今天收到了北京大学老师打来的电话，问我如果没有被数据科学方向的导师录取，愿不愿意去读生物统计的博士。

我婉拒了，些许遗憾，但不后悔，原因全是个人选择，读博挺好的，但是我决定换一种环境，去工作了。

从去年11月开始申请，到一系列的纠结，到现在做下决定，确实是释然了很多。读书很自由，但读书并不是适合每一个人的选择，或者说不是适合一个人特定时期的选择，也许工作不顺意又想去读也说不定。

真的是越长大越体会到人生重要决定时的艰难迷茫。怎么选都会遗憾，只有勇敢的走下去吧。

感慨一下哈，今天继续给大家写潜增长模型。

随机效应

要理解随机效应，首先还是要理解嵌套，比如同一个学校的学生，或者同一个人的纵向测量数据，这些数据我们不能想当然地认为它们之间是独立的，但是线性回归的前提假设便是观测的独立性，所以对于嵌套数据一般不能用普通的回归分析法。

对于嵌套数据我们要用混合效应模型，混合效应模型之前的文章写了很多，大家可以自己翻翻哦，为啥叫混合效应模型呢？

They are mixed, because there is a mixture of fixed effects and random effects.

因为这个模型可以同时估计固定效应和随机效应，所以就叫它混合效应模型，而其中的固定效应就和我们普通回归分析的回归系数一样，是我们主要关心的东西。

加上随机效应就是为了让我们能够把嵌套数据的变异分解的更加的清晰：我们可以让每个学校有每个学校自己的特征（随机截距或随即斜率），也可以让纵向测量中的每个人有自己的特征（随机截距或随即斜率），之所以要加随机效应，是我们想让嵌套数据的变异分解地更好从而使得我们的固定效应估计更准确。

混合效应模型的构成

给大家写写混合效应模型构成的实例，比如现在我要研究y和x的关系，但是数据是从学生中收集的，学生来自不同的地区，那么这就是一个嵌套，学生嵌套在不同的地区。我们规定学生用i表示，地区用c表示，拟合的模型如下：

上面这个模型中有一个随机截距b0c，就是说不同的地区学生的起始x允许不一样：

怎么不一样呢？

就是这个随机截距b0c为整体截距b0和不同地区截距uc的和（地区水平的变异可大可小，但我们要考虑它，这就是随机效应的作用），如下式：

我们会假定所有的uc都是服从正态分布的。

此时我们把随机截距模型改写一下就成了下面这个样子：

或者：

看第二个式子你就会明白随机截距只是在常规回归模型上多了一个截距项uc，所以叫做随机截距。

把随机截距放在结构方程模型中理解

刚刚写了随机截距是服从正态分布的，回想我们在结构方程中潜变量也是服从正态分布的，有没有什么联想呢？

是不是可以把随机截距当作潜变量处理呢？

我们来试试，首先模拟数据，我们模拟一个纵向数据，共500个观测，4个时间点：

set.seed(1234)
n = 500
timepoints = 4
time = rep(0:3, times=n)
subject = rep(1:n, each=4)

我们会有固定效应，还有随机效应，对于本例来讲，测量是嵌套在每个观测上的（每个观测测了4次），所以我们的随机效应便是个体水平的截距或斜率，我们让固定效应的截距和斜率分别为 0.5 和 0.25，随机截距和随机斜率的相关为0.2

intercept = .5
slope = .25
randomEffectsCorr = matrix(c(1,.2,.2,1), ncol=2)
randomEffects = MASS::mvrnorm(n, mu=c(0,0), Sigma = randomEffectsCorr, empirical=T) %>% data.frame()
colnames(randomEffects) = c('Int', 'Slope')

于是此时我们的数据长这样：

我们接着模拟自变量和应变量，我们的思路就是把随机效应加在固定效应上，就是把随机截距和固定截距相加，把随机斜率和固定斜率相加，自变量我们就用时间time就行，应变量的模拟代码如下：

sigma = .5
y1 = (intercept + randomEffects$Int[subject]) + # 随机截距(slope + randomEffects$Slope[subject])*time + #随机斜率rnorm(n*timepoints, mean=0, sd=sigma)d = data.frame(subject, time, y1)

好了，现在这个d便是我们要用的数据集，它长这样：

因为这个数据集是我们亲自模拟出来的，所以我们知道它既有随机斜率又有随机截距，而且我们还知道固定效应和随机效应大小分别是多少，我们现在就来逆操作，看看我们拟合一个混合模型能不能得到我们预想的系数：

library(lme4)
mixedModel = lmer(y1 ~ time + (1 + time|subject), data=d)  # 带有随机截距和随机斜率的混合模型
## summary(mixedModel)

运行上面的代码，就可以见证奇迹了：

我们可以看到，我们设定的固定斜率为0.25，图中结果输出为0.261，我们设定的固定截距为0.5，图中结果输出为0.487，因为我们本身还加了随机误差，除过这个随机误差之后可以说是完全吻合。

跑潜增长模型

上面依然是给大家回忆和混合效应模型，希望大家对混合模型已经没问题了，我们接着看潜增长模型：

在潜增长模型中我们会把随机截距和随机斜率当作潜变量处理，相应地我们会固定住因子载荷，看下面这个示意图，我们有随机截距和随机斜率，我们还把载荷都给固定了，注意，我们把不同时间点的截距的载荷都固定为1，斜率载荷的固定一定要反映出我们数据的实际间隔，通常是从0开始。

我们做潜增长模型的时候每个时间点的测量都是以变量纳入的，所以我们的数据得是宽型数据，转换方法如下：

dWide = d %>%  spread(time, y1) %>% rename_at(vars(-subject), function(x) paste0('y', x))
head(dWide)

我们的数据就成了每个观测分别在不同的时间点时的y值，此时时间点变成了变量而非变量的水平（看不懂这句话的话去看看tidydata）。

数据处理好之后我们就可以来跑增长模型了，用到的是lavaan包中的growth方法：

model = "# 将随机效应作为潜变量i =~ 1*y0 + 1*y1 + 1*y2 + 1*y3s =~ 0*y0 + 1*y1 + 2*y2 + 3*y3
"
growthCurveModel = growth(model, data=dWide)
summary(growthCurveModel)

运行上面的代码便可以得到增长模型的结果，如下图：

大家注意，结果很多地方都是空的，我们需要看的结果就是随机效应的截距项：

结果中有我们的潜变量i和s的截距，这个东西啥意思呢？

回忆一下，我们现在是在跑随机效应，我们是将随机截距和随机斜率当作了潜变量，那么这个i和s的截距就是整体的均值，大家可能又糊涂了，就是说，我们的随机效应是服从正态分布的，它的均值就是固定效应，随机效应会在固定效应周围，而这个潜变量i和s的截距就是均值，所以你就可以把i和s的截距认为是固定效应，希望我写明白了，还不懂的话私信我吧。

还有点怀疑？

哈哈，没事我们来验证下嘛。记住上图的截距分别是0.487和0.261。

我们跑潜增长模型的数据和混合效应的数据完全是一样的，我模拟数据的时候是规定截距和斜率分别是0.5和0.25的，我们i和s的截距是0.487和0.261，考虑我加了随机误差，所以这些值是完全一致的，我们还可以用以下的代码调出来固定效应：

fixef(mixedModel)

看到没，依然是0.487和0.261，和我们的增长模型中的i和s的截距是一模一样的。所以说，大家记住，增长模型是跑的随机效应，所以我们要关心截距输出，截距输出就是整体均值，也就是固定效应。

增长模型的结果输出中还有一部分结果是方差的估计：

这个该怎么看呢？这个不重要，但是大家可以了解一下

潜增长模型给每个时间点都会估计一个方差，但是在混合模型中认为不同时点的方差是一样的，所以你会看到在混合模型中std.Dev为0.48，而在增长模型中每个时点的variance都不一样，平均起来也是0.48.这个就是混合模型和增长模型的区别之一。

为了大家更好地理解上面的叙述，再给大家写个例子：

我们用同样的数据我们做增长模型的时候把方差固定住：

#将增长模型的变量方差固定就和混合模型结果一样了
model = "# intercept and slope with fixed coefficientsi =~ 1*y0 + 1*y1 + 1*y2 + 1*y3s =~ 0*y0 + 1*y1 + 2*y2 + 3*y3y0 ~~ resvar*y0    y1 ~~ resvar*y1y2 ~~ resvar*y2y3 ~~ resvar*y3
"growthCurveModel = growth(model, data=dWide)
summary(growthCurveModel)

结果如下：

还是同样的数据，我们再对比其与混合模型的结果：

可以看出来，增长模型固定方差后就和混合模型的结果是一样的。也就是说大家记住，增长模型的和混合模型的区别就是增长模型认为不同时间变量的变异是不一样的。

小结

今天从混合效应模型出发，进一步给大家解释了混合效应模型和潜增长模型的区别，感谢大家耐心看完，自己的文章都写的很细，代码都在原文中，希望大家都可以自己做一做，请关注后私信回复“数据链接”获取所有数据和本人收集的学习资料。如果对您有用请先收藏，再点赞转发。

也欢迎大家的意见和建议，大家想了解什么统计方法都可以在文章下留言，说不定我看见了就会给你写教程哦。

R数据分析：潜增长模型LGM的做法和解释，及其与混合模型对比相关推荐

R数据分析：工具变量回归的做法和解释，实例解析
什么是工具变量,以及什么是孟德尔随机化,以及孟德尔随机化怎么实现都给大家写了(大家去翻翻之前的文章呀),因为孟德尔随机化的工具变量是基因变量,所以我们会用专门的R包去做,普通的工具变量研究,我们要用的 ...
mplus数据分析：增长模型潜增长模型与增长混合模型再解释
Python微信订餐小程序课程视频 https://edu.csdn.net/course/detail/36074 Python实战量化交易理财系统 https://edu.csdn.net/cou ...
R数据分析：潜在转化分析LTA的做法和解释（一）
之前给大家写了很多潜在类别分析的教程Mplus教程:如何做潜在类别分析LCA R数据分析:用R语言做潜类别分析LCA Mplus数据分析:潜在类别分析(LCA)流程(详细版) R数据分析:再写潜在类别 ...
CDA数据分析——AARRR增长模型的介绍、使用
AARRR增长模型:将产品的营收路径拆分为:激活→注册→留存→下单→传播. 下面是对 AARRR增长模型中各渠道的定义及运营方式做详细讲解: Acquisition(用户获取)激活: 这是流量来源的 ...
R数据分析：生存分析与有竞争事件的生存分析的做法和解释
今天被粉丝发的文章给难住了,又偷偷去学习了一下竞争风险模型,想起之前写的关于竞争风险模型的做法,真的都是皮毛哟,大家见笑了.想着就顺便把所有的生存分析的知识和R语言的做法和论文报告方法都给大家梳理一遍 ...
R数据分析：如何做数据的非线性关系，多项式回归的做法和解释
线性关系其实是最常见也是最有效,同时还是最好解释的,不过变量间复杂的关系我们用多项式回归做出来可能会更加的准确.刚好有位粉丝的数据需要用到多项式回归,今天就给大家写写. 要理解非线性关系,首先我们看看 ...
Mplus数据分析:随机截距交叉之后的做法和如何加协变量，写给粉丝
记得之前有写过如何用R做随机截距交叉滞后,有些粉丝完全是R小白,还是希望我用mplus做,今天就给大家写写如何用mplus做随机截距交叉滞后. 做之前我们需要知道一些Mplus的默认的设定: obse ...
R数据分析：交叉滞后模型基础与实例解析
最近问纵向数据分析的同学贼多,像潜增长,GEE,多水平,之前都有写,今天偷空出个简易的交叉滞后教程哈,大家只要遇到像causal models,cross- lagged panel models,l ...
R数据分析：生存分析的做法和结果解释
今天给大家写写生存分析: Survival analysis corresponds to a set of statistical approaches used to investigate th ...

R数据分析：潜增长模型LGM的做法和解释，及其与混合模型对比

随机效应

混合效应模型的构成

把随机截距放在结构方程模型中理解

跑潜增长模型

小结

R数据分析：潜增长模型LGM的做法和解释，及其与混合模型对比相关推荐

最新文章

热门文章