предсказать () в регрессии lmer, но мне нужно только 2 категории

Я пытаюсь оценить многоуровневую модель. Мой код:

fullModel2 <- lmer(pharmexp_2001 ~ gdp_1000_gm + health_exp_per_cap_1000_gm + life_exp +
                   labour_cost_1000_gm + (year_gm|lowerID), data=adat, REML=F)

что приводит к следующей модели:

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: pharmexp_2001 ~ gdp_1000_gm + health_exp_per_cap_1000_gm + life_exp +      
         labour_cost_1000_gm + (year_gm | lowerID)
   Data: adat

     AIC      BIC   logLik deviance df.resid 
  1830.2   1859.9   -906.1   1812.2      191 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.5360 -0.6853 -0.0842  0.4923  4.0051 

Random effects:
 Groups   Name        Variance Std.Dev. Corr 
 lowerID  (Intercept) 134.6851 11.6054       
          year_gm       0.4214  0.6492  -1.00
 Residual             487.5324 22.0801       
Number of obs: 200, groups: lowerID, 2

Fixed effects:
                            Estimate Std. Error t value
(Intercept)                -563.7924    75.4125  -7.476
gdp_1000_gm                  -0.9050     0.2051  -4.413
health_exp_per_cap_1000_gm   37.5394     6.3943   5.871
life_exp                      8.8571     0.9498   9.326
labour_cost_1000_gm          -1.3573     0.4684  -2.898

Correlation of Fixed Effects:
            (Intr) g_1000 h____1 lif_xp
gdp_1000_gm -0.068                     
hl____1000_  0.374 -0.254              
life_exp    -0.996  0.072 -0.393       
lbr_c_1000_ -0.133 -0.139 -0.802  0.142

Я знаю, что проблема в том, что корреляция равна -1 из-за случайных эффектов, но у меня проблема посерьезнее. Мне нужно построить свои результаты, но мне нужны только 2 строки: когда lowerID=0 и когда lowerID=1. Итак, я хочу построить pharmaexp_2001 по оси Y против year по оси X, но мне нужно только 2 линии (по lowerID). Я знаю, что мне нужно использовать predict.merMod, но как я могу построить эти результаты, отобразив только эти две линии? В настоящее время на моем графике 21 линия (потому что я анализирую расходы на лекарства в 21 стране).


r data-visualization lmer
person Eszter Takács    schedule 25.03.2014    source источник
comment
Не могли бы вы объяснить значение year_gm в части случайных эффектов?   -  person Michael M    schedule 25.03.2014

Ответы (1)


arrow_upward
3
arrow_downward

Добро пожаловать на сайт, @Eszter Takács!

Вам нужно только указать два идентификатора в newdata. Вот пример, основанный на данных sleepstudy в R. Я предполагаю, что вы хотите отобразить прогнозируемые значения по оси Y. Просто замените код своими данными и переменными, и вы получите предсказанные значения для lowerID==0 и lowerID==1. Затем вы можете использовать свой код для построения двух линий для двух идентификаторов.

> (fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy, REML=F))
Linear mixed model fit by maximum likelihood ['lmerMod']
Formula: Reaction ~ Days + (Days | Subject) 
   Data: sleepstudy 
      AIC       BIC    logLik  deviance 
1763.9393 1783.0971 -875.9697 1751.9393 
Random effects:
 Groups   Name        Std.Dev. Corr
 Subject  (Intercept) 23.781       
          Days         5.717   0.08
 Residual             25.592       
Number of obs: 180, groups: Subject, 18
Fixed Effects:
(Intercept)         Days  
     251.41        10.47  

> newdata = sleepstudy[sleepstudy$Subject==308 | sleepstudy$Subject==333,]
> str(p <- predict(fm1,newdata)) # new data, all RE
 Named num [1:20] 254 274 293 313 332 ...
 - attr(*, "names")= chr [1:20] "1" "2" "3" "4" ...
person Randel    schedule 25.03.2014