R：在混合模型效果中绘制RANEF

Question

因此，我对R相当陌生。我正在研究一个项目，其中有一个数据集用于评估鸟类年龄。我们有来自95个人的> 400,000个观测值。我的任务是这样的：

“如果图上有一些数据点，那么这个图会更有说服力，你可以从lme4的getcall中获取ranef的数据点，因此，运行模型时不需要年龄，然后，抽出随机效应项（每个个体的值），然后绘制它们（y =随机效应，x =年龄）。

这是有问题的图（它是在Excel中使用我们的模型制作）：Graph

所以我跑这

> lmbdna<-lmer(Gs ~ (1 | Individual) + Bin + year + Mass) 
> ranef(lmbdna) 
>$Individual 

> plot(Age, ranef(lmbdna)) 

而且我得到了在XY这个

错误.coords（x，y，xlabel，ylabel，log）： 'x'和'y'长度不同

我真的迷失了，我不是如何绘制时代的ranefs。有没有办法将一个年龄与个人联系起来以摆脱这个错误？

下面是我的一些数据：

Indiv Age Mass MaxDepth Depth  Gs  PDBA Bin year 
69903 12 1015 3.806 3.025 0.1854302 92.7151 A N 
52712 20 957.5 3.806 3.025 0.204678 102.339 A T 
55969 19 1002.5 3.806 3.025 0.222338 111.169 A T 
64442 15 1040 3.806 3.025 0.1872954 93.6477 A T 
76252 11 940  3.806 3.025 0.223136 111.568 A T 
53391 21 1022.5 3.806 3.025 0.234452 117.226 A E 
53391 21 1022.5 3.806 3.025 0.299438 149.719 A E 
60117 18 937.5 3.806 3.025 0.1469442 73.4721 A E 
60151 18 970  3.806 3.025 0.1941052 97.0526 A E 
52712 20 957.5 3.855 3.025 0.1812926 90.6463 A T 
52712 20 957.5 3.855 3.025 0.25101  125.505 A T 
64442 15 1040 3.855 3.025 0.1850976 92.5488 A T 
64442 15 1040 3.855 3.025 0.1026478 51.3239 A T 
76252 11 940  3.855 3.025 0.235822 117.911 A T 
78712 10 880  3.855 3.025 0.1638106 81.9053 A T 
87819 7 1000 3.855 3.025 0.166391 83.1955 A T 
90281 6 957.5 3.855 3.025 0.1493948 74.6974 A T 
60151 18 970  3.855 3.025 0.1904232 95.2116 A E 
69944 12 915  3.904 3.025 0.256504 128.252 A N 
3260 24 960  3.904 3.025 0.168019 84.0095 A T 
52712 20 957.5 3.904 3.025 0.270704 135.352 A T 
64442 15 1040 3.904 3.025 0.1507102 75.3551 A T 
71432 12 970  3.904 3.025 0.1238154 61.9077 A T 
80538 15 917.5 3.904 3.025 0.236976 118.488 A E 
76583 14 870  3.952 3.025 0.295982 147.991 A N 
84420 7 1005 3.952 3.025 0.1861876 93.0938 A N 
87819 7 1000 3.952 3.025 0.178179 89.0895 A T 
53391 21 1022.5 3.952 3.025 0.1917954 95.8977 A E 
53391 21 1022.5 3.952 3.025 0.1482036 74.1018 A E 
53391 21 1022.5 3.952 3.025 0.1999868 99.9934 A E 
53391 21 1022.5 3.952 3.025 0.276334 138.167 A E 
60151 18 970  3.952 3.025 0.1776108 88.8054 A E 
80538 15 917.5 3.952 3.025 0.188733 94.3665 A E 
69944 12 915  4.001 3.025 0.2596  129.8 A N 
3260 24 960  4.001 3.025 0.1824546 91.2273 A T 

任何帮助表示赞赏。谢谢

Answer 1

主要问题是你提取的随机效应是每个人的个人，而你试图绘制的年龄数据是针对每个观察。您需要将其汇总到个人级别（例如，在所有观察结果中为每个人获取最大值，以获得您要查找的结果。下面的示例使用您在上面提供的文本）

不完全确定是否这样是你所追求的，但主要的问题是，你的向量的长度不同。希望下面的代码会给你一些想法。本来也极大地得到更大的数据样本。

require(lme4) 

# the data snapshot supplied above only contains one level of the factor bin, let's add another 
dta[dta$Indiv < 7000, "Bin"] <- "B" 

# estimate model 
lmbdna <- lmer(Gs ~ (1 | Indiv) + Bin + year + Mass, dta) 

# random effects are for indivuals, so need to aggregate individuals' ages 
plt_dta <- aggregate(dta$Age, by = list(dta$Indiv), max) 

# create one dataset by merging with random effects 
plt_dta <- merge(plt_dta, ranef(lmbdna)$Indiv, by.x = 1, by.y = 0) 

# plot data 
plot(plt_dta[,2], plt_dta[,3], xlab = "age", ylab = "ranef") 

Answer 2

我不要以为我们可以根据你分享的数据得到一个有意义的模型，所以让我们使用类似?ranef中的示例：

library("lme4") 
fm1 = lmer(Reaction ~ Days + (1 | Subject), data = sleepstudy) 

my_ranef = ranef(fm1) 
str(my_ranef) 
# List of 1 
# $ Subject:'data.frame': 18 obs. of 1 variable: 
# ..$ (Intercept): num [1:18] 40.78 -77.85 -63.11 4.41 10.22 ... 
# - attr(*, "class")= chr "ranef.mer" 

str()用于查看对象的结构。在这里我们可以看到my_ranef是1项目的列表，并且该项目是具有18行和1列（名为(Intercept)）的数据帧（名称为Subject）。让我们在数据帧仔细一看：

head(my_ranef$Subject) 
#  (Intercept) 
# 308 40.783710 
# 309 -77.849554 
# 310 -63.108567 
# 330 4.406442 
# 331 10.216189 
# 332 8.221238 

所以数据帧具有对应于每个个体对象，然后随机效应截距是在（截距）栏列名。因此，我们可以画出这样的随机效应：

plot(row.names(my_ranef$Subject), my_ranef$Subject[, "(Intercept)"]) 

你有问题，你是给无论是整个列表或一个数据帧作为y参数的情节，而你需要提取刚刚矢量。

我们也可以提取随机效应

est_ranef = my_ranef$Subject[, "(Intercept)"] 

（注意括号使"(Intercept)"其增添了些许困难，提取其非标准的列名，一个简单的$是行不通一拉my_ranef$Subject$(Intercept)，但我们可以用[和引用列名作为以上，或反引号可以$像这样使用：

my_ranef$Subject$`(Intercept)` 

如果您在ranef()查看更复杂的模型，您会看到为什么使用数据帧结构列表。