Posterior_survfit () не использует nd - PullRequest
0 голосов
/ 03 августа 2020

У меня проблемы с созданием апостериорных прогнозов с помощью posterior_survfit (). Я пытаюсь использовать новый фрейм данных, но он не использует новый фрейм данных, а вместо этого использует значения из набора данных, которые я использовал для соответствия модели. Подбираемые переменные в модели: New.Treatment (6 обработок = категориальные), Openness (непрерывный световой индекс min = 2,22, средний = 6,903221 и max = 10,54), subplot_by_site (категориальные - 720 сайтов), New.Species.name ( категориальный - 165 видов). В моем новом фрейме данных 94 строки, а posterior_survfit () дает мне 3017800 строк. Помогите, пожалуйста!

head(nd)
      New.Treatment Openness
1          BE                  5
2          BE                  6
3          BE                  7
4          BE                  8
5          BE                  9
6          BE                 10


fit= stan_surv(formula = Surv(days, Status_surv) ~ New.Treatment*Openness + (1 |subplot_by_site)+(1|New.Species.name),
  data = dataset,
  basehaz = "weibull",
  chains=4,
  iter = 2000,
  cores =4 )

Post=posterior_survfit(fit, type="surv",
                        newdata=nd5)

head(Post)
  id cond_time    time median  ci_lb  ci_ub
1  1        NA 62.0000 0.9626 0.9623 1.0000
2  1        NA 69.1313 0.9603 0.9600 0.9997
3  1        NA 76.2626 0.9581 0.9579 0.9696
4  1        NA 83.3939 0.9561 0.9557 0.9665
5  1        NA 90.5253 0.9541 0.9537 0.9545
6  1        NA 97.6566 0.9522 0.9517 0.9526

##Here some reproducible code to explain my problem:

library(rstanarm)

data_NHN<- expand.grid(New.Treatment = c("A","B","C"), Openness = c(seq(2, 11, by=0.15)))
data_NHN$subplot_by_site=c(rep("P1",63),rep("P2",60),rep("P3",60))
data_NHN$Status_surv=sample(0:1,183, replace=TRUE) 
data_NHN$New.Species.name=c(rep("sp1",10),rep("sp2",40),rep("sp1",80),rep("sp2",20),rep("sp1",33))
data_NHN$days=sample(10, size = nrow(data_NHN), replace = TRUE)

nd_t<- expand.grid(New.Treatment = c("A","B","C"), Openness = c(seq(2, 11, by=1)))


mod= stan_surv(formula = Surv(days, Status_surv) ~ New.Treatment+Openness + (1 |subplot_by_site)+(1|New.Species.name),
                  data =data_NHN,
                  basehaz = "weibull",
                  chains=4,
                  iter = 30,
                  cores =4)

summary(mod)
pos=posterior_survfit(mod, type="surv",
                        newdataEvent=nd_t,
                      times = 0)
head(pos)

 #I am interested in predicting values for specific Openess values  
 #(nd_t=20 rows)but I am getting instead values for each point in time 
 #(pos=18300rows)

Операционная система: Ma c OS Catalina 10.15.6 Версия R: 4.0 Версия rstan: 2.21.2 Версия rstanarm: rstanarm_2.21.2 Любые предложения, почему она не работает. неясно, как создать какой-то график влияния одной переменной во взаимодействии по мере изменения другой и связанной с этим неопределенности (например, график предельных эффектов). В моем примере я заинтересован в получении значений в определенных c "значениях открытости", а не в каждом заданном c времени, как показано в апостериорных результатах. TIA.

1 Ответ

0 голосов
/ 04 августа 2020

Я могу дать вам лишь частичный ответ. Вы должны знать, что это довольно актуально; хотя в документе ArXiv (от февраля 2020 г.) говорится:

Надеюсь, к тому времени, когда вы это прочтете, функциональность будет доступна в стабильном выпуске в Comprehensive R Archive Network (CRAN)

но пока это даже не так; его даже нет в основной ветке GitHub, поэтому я использовал remotes::install_github("stan-dev/rstanarm@feature/survival"), чтобы установить его из источника.

Ближайшая проблема заключается в том, что вы должны указать новый фрейм данных как newdata, а не newdataEvent. Кажется, есть много несоответствий между мастером и этой ветвью, а также между документами и кодом ... newdataEvent используется в старом методе для моделей stanjm, но не для моделей stansurv. Вы можете посмотреть код здесь или использовать formals(rstanarm:::posterior_survfit.stansurv). К сожалению, поскольку этот метод имеет аргумент (неиспользуемый, не отмеченный) ..., это означает, что любые неверно названные аргументы будут игнорироваться.

Следующая проблема заключается в том, что если вы укажете новые данные в своем примере как newdata вы получите

Ошибка: в данных отсутствуют следующие переменные: subplot_by_site, New.Species.name

То есть не похоже чтобы быть очевидным способом создания апостериорного прогноза на уровне популяции. (Установка переменных группировки случайных эффектов на NA не допускается.) Если вы хотите сделать это, вы можете:

  • развернуть newdata, чтобы включить все комбинации группирующих переменных в вашем наборе данных и самостоятельно усредняйте результаты по уровням;
  • опубликуйте проблему на GitHub или свяжитесь с сопровождающими ...
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...