我正在尝试为预测值生成标准误差。我使用下面的代码来生成预测值,但它也没有给出标准错误。
ord6 <- veg$ord1-2
laimod.group = lmer(log(lai+0.000019) ~ ord6*plant_growth_form +
(1|plot.code) +
(1|species.code),
data=veg,
REML=FALSE)
summary(laimod.group)
new.ord6 <- c(-1,0,1,2,3,4,5,6,7)
new.plant_growth_form <- c("fern", "grass", "herb","herbaceous climber",
"herbaceous shrub", "moss", "tree sapling",
"undet", "woody climber", "woody shrub")
newdat <- expand.grid(
ord6=new.ord6,plant_growth_form=new.plant_growth_form)
newdat$pred <- predict(laimod.group,newdat, se.fit=TRUE, re.form=NA)
newdat
评论 1:laimod.group = 使用 lmer 比较五个模型后选择的最终模型(包 lme4)
评论 2:predictSE.mer 需要包 AICcmodavg
我确实尝试了以下代码作为替代,但继续收到以下错误消息: fam.link.mer(mod) 中的错误:找不到对象'out.link'
newdat$pred <- predictSE.mer(laimod.group, newdat, se.fit = TRUE, type = "response",
level = 0, print.matrix = FALSE)
请查看我的数据的可重现子集:
structure(list(plot.code = structure(c(1L, 2L, 3L, 4L, 5L, 5L,
5L, 5L, 5L, 5L, 6L, 7L, 7L, 8L, 8L, 8L, 8L, 8L, 8L, 8L, 9L, 9L,
9L, 9L, 10L, 11L, 11L, 11L, 11L, 11L, 12L, 13L, 14L, 15L, 15L,
15L, 15L, 15L, 15L, 15L, 16L, 17L, 18L, 19L, 19L, 19L, 20L, 21L,
22L, 23L, 24L, 25L, 26L, 27L, 27L, 28L, 28L, 28L, 28L, 28L, 29L,
29L, 30L, 30L, 30L, 30L, 30L, 30L, 30L, 31L, 32L, 34L, 34L, 34L,
34L, 34L, 34L, 34L, 35L, 36L, 36L, 36L, 37L, 38L, 39L, 39L, 39L,
40L, 40L, 33L, 33L), .Label = c("a100f1r", "a100m562r", "a10m562r",
"a1f56r", "a1m5r", "b100f177r", "b100m17r", "b100m5r", "c100f17r",
"c100f1r", "c100f5r", "d100m56r", "d100m5r", "d10f1r", "d10f5r",
"e100m17r(old)", "e100m1r", "e100m5r", "e10f177r", "e10f17r(old)",
"e10f5r(old)", "e1f17r", "e1f5r", "f100m177r", "f10f177r", "f10f17r",
"f1m177r", "f1m56r", "lf1f1r", "lf1f5r", "lf1m1r", "og100f5r",
"og10f1r", "og10m1r", "og10m5r", "op100f562r", "op100m177r",
"op10f1r", "op10f5r", "op10m562r"), class = "factor"), species.code = structure(c(69L,
59L, 67L, 69L, 20L, 44L, 28L, 32L, 31L, 7L, 13L, 63L, 69L, 52L,
69L, 14L, 54L, 57L, 42L, 9L, 62L, 10L, 22L, 69L, 35L, 49L, 38L,
11L, 41L, 39L, 16L, 40L, 69L, 32L, 33L, 41L, 22L, 69L, 43L, 4L,
68L, 48L, 6L, 34L, 53L, 3L, 15L, 30L, 13L, 31L, 66L, 64L, 38L,
46L, 61L, 29L, 61L, 27L, 8L, 41L, 55L, 58L, 23L, 25L, 18L, 45L,
26L, 13L, 65L, 12L, 51L, 50L, 60L, 47L, 17L, 5L, 19L, 61L, 1L,
37L, 13L, 36L, 13L, 2L, 11L, 24L, 44L, 13L, 49L, 56L, 21L), .Label = c("agetri",
"alb214", "annunk", "arimin", "baudip", "beg032", "blurip", "buc009",
"cal079", "calplu", "chrodo", "cishas", "clihir", "cos049", "cycari",
"cypunk", "cyr075", "cyrped1", "dae205", "dalpin", "diapla1",
"dio063", "diosum", "emison", "ery046", "eryborb", "fic119",
"ficmeg", "friacu", "graunk", "indunk", "jactom", "lauunk", "leeind",
"luvsar", "lyccer", "mac068", "melmal", "mergra", "miccra1",
"mikcor", "mitken", "nep127", "nepbis", "paldas", "palunk", "panunk",
"penlax", "poaunk", "pol019", "pop246", "ptecog", "ptesub1",
"rubcle", "ryphul", "scamac", "scl051", "sclsum", "selcup", "selfro",
"spa098", "sphste1", "stitrut", "tet055", "tetdie", "tetdie1",
"tetkor", "xanfla", "zinunk"), class = "factor"), plant_growth_form = structure(c(3L,
6L, 9L, 3L, 7L, 1L, 7L, 4L, 8L, 5L, 5L, 1L, 3L, 7L, 3L, 3L, 9L,
2L, 9L, 7L, 9L, 7L, 7L, 3L, 9L, 2L, 10L, 4L, 4L, 9L, 2L, 7L,
3L, 4L, 7L, 4L, 7L, 3L, 1L, 4L, 7L, 7L, 3L, 10L, 7L, 7L, 1L,
2L, 5L, 8L, 9L, 9L, 10L, 7L, 9L, 9L, 9L, 7L, 7L, 4L, 7L, 2L,
7L, 10L, 3L, 7L, 10L, 5L, 9L, 9L, 7L, 7L, 6L, 7L, 3L, 9L, 9L,
9L, 9L, 7L, 5L, 6L, 5L, 9L, 4L, 3L, 1L, 5L, 2L, 7L, 7L), .Label = c("fern",
"grass", "herb", "herbaceous climber", "herbaceous shrub", "moss",
"tree sapling", "undet", "woody climber", "woody shrub"), class = "factor"),
ord1 = c(9L, 5L, 7L, 9L, 4L, 4L, 5L, 5L, 5L, 2L, 9L, 5L,
4L, 6L, 8L, 6L, 3L, 3L, 5L, 3L, 4L, 5L, 3L, 5L, 3L, 9L, 6L,
4L, 4L, 6L, 2L, 5L, 5L, 9L, 3L, 4L, 3L, 5L, 3L, 4L, 1L, 8L,
1L, 5L, 7L, 6L, 9L, 1L, 9L, 1L, 4L, 4L, 2L, 5L, 2L, 3L, 5L,
1L, 3L, 3L, 3L, 2L, 6L, 5L, 2L, 6L, 5L, 2L, 5L, 3L, 6L, 5L,
6L, 3L, 3L, 4L, 7L, 4L, 6L, 1L, 2L, 2L, 4L, 3L, 3L, 3L, 3L,
4L, 4L, 3L, 3L), lai = c(4.525068022, 0.325399379, 0.229222148,
4.076350538, 0.006889889, 0.003279268, 0.037268428, 0.056032134,
0.013573973, 0.001304667, 0.696949844, 1.256477431, 0.122569437,
0.191398415, 1.606070777, 0.425381508, 0.03013251, 0.00181661,
0.017317993, 0.014455456, 0.102704752, 0.031065374, 0.000923601,
0.453384679, 0.017859983, 7.765697214, 0.127071322, 0.102178413,
0.049099766, 0.427983019, 4.22e-05, 0.229034333, 0.694745347,
0.068069112, 0.218354525, 0.05883256, 0.032252145, 0.304812298,
0.009320025, 0.036424481, 0, 0.326, 0.000201724, 0.286106787,
0.556249444, 0.274764132, 4.21, 0, 0.695663959, 0.000213763,
0.00476907, 0.000205017, 3.77e-05, 0.134661951, 0.005631489,
0.0971, 0.172154618, 5.91e-05, 0.000371101, 0.000145266,
0.013382779, 0.00025348, 0.11016712, 0.0616302, 0.018011524,
0.107619537, 0.189926726, 0.000857257, 0.041252452, 0, 0.00475341,
0.077329281, 0.633865958, 0.038182437, 0.015560589, 0.010375148,
1.515423445, 0.008559863, 0.003636564, 0.000424537, 0.002786085,
0.091458876, 0.014216177, 0.165042816, 0.009187705, 0.00115711,
0.000920496, 0.009072635, 0.001443384, 0.001595447, 0.023263507
)), .Names = c("plot.code", "species.code", "plant_growth_form",
"ord1", "lai"), class = "data.frame", row.names = c(NA, -91L))