I am attempting to estimate a multilevel model. My code is:
我正在尝试估计一个多层次模型。我的代码是:
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)
which results in the following model:
其结果为:
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
I know that it's a problem that the correlation is -1 by random effects, but I have a bigger problem. I have to plot my results, but only I need 2 lines: when lowerID=0 and when lowerID=1. So I want to plot pharmaexp_2001 on the y-axis against year on the x-axis, but I need only 2 lines (by lowerID). I know that I have to use predict.merMod, but how can I plot these results, plotting only these two lines? Currently my plot has 21 lines (because I analyse pharmaceutical expenditure in 21 countries).
我知道随机效应的相关性是-1是一个问题，但我有一个更大的问题。我需要画出结果，但只需要两行:当小写id =0和当小写id =1。我想在y轴上绘制pharmaexp_2001，在x轴上绘制年份，但我只需要2行(用小写id)。我知道我必须使用预测。merMod，但是我怎么画这些结果，只画这两条线呢?目前我的地块有21条线(因为我分析了21个国家的医药支出)。
1 个解决方案
#1
3
Welcome to the site, @Eszter Takács!
欢迎来到网站@Eszter Takacs !
You only need to specify the two IDs in newdata. Here is an example based on sleepstudy data in R. I assume you want to plot the predicted values on the y-axis. Just replace the code with your data and variables, you will obtain the predicted values for lowerID==0 and lowerID==1. Then you can use your code to plot the two lines for the two IDs.
您只需在newdata中指定这两个id。这里有一个基于r的睡眠研究数据的例子，我假设你想在y轴上绘制预测值。只需用您的数据和变量替换代码，您将获得小写id =0和小写id =1的预测值。然后您可以使用您的代码为这两个id绘制两行代码。
> (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" ...
#1
3
Welcome to the site, @Eszter Takács!
欢迎来到网站@Eszter Takacs !
You only need to specify the two IDs in newdata. Here is an example based on sleepstudy data in R. I assume you want to plot the predicted values on the y-axis. Just replace the code with your data and variables, you will obtain the predicted values for lowerID==0 and lowerID==1. Then you can use your code to plot the two lines for the two IDs.
您只需在newdata中指定这两个id。这里有一个基于r的睡眠研究数据的例子，我假设你想在y轴上绘制预测值。只需用您的数据和变量替换代码，您将获得小写id =0和小写id =1的预测值。然后您可以使用您的代码为这两个id绘制两行代码。
> (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" ...