Rpackage_lavaan
2023-06-28
摘要
结构方程建模（SEM）已被社会和行为科学中的许多应用研究者广泛使用。这些年来，许多用于结构方程建模的软件包已经被开发出来，其中既有免费的也有商用的。R包lavaan为应用研究人员、教师和统计学家提供了一个免费的、完全开源的，但拥有商业质量的潜变量建模包。本文解释了开发lavaan包背后的目的，概述了它最重要的特性，并提供了一些例子来说明lavaan包在实践中是如何工作的。
介绍
本文介绍了R包lavaan。这是一个在R系统中实现的结构方程建模软件包。该软件包可从http://CRAN.R-project的综合R档案网络（CRAN）获得，也可以在http://lavaan.org/网站查询相关信息。lavaan是潜变量分析的首字母缩写，它的名字揭示了我们的目标：提供一套工具用来探索、估计和理解一个潜变量模型，其中的功能包括因素分析、构建结构方程。然而，lavaan包的构建才刚刚开始，要实现这一雄心勃勃的目标，还有许多工作要做。迄今为止，lavaan包的重点主要集中在具有连续观测变量的结构方程建模（SEM）上（Bollen 1989）。结构方程模型包含了广泛的多元统计技术。 在R环境中，有两种方法来估计结构方程模型。第一种方法是将R与外部的商用SEM程序连接起来。这在模拟研究中通常很有用（例如，参见van de Schoot、Hoijtink和Dekovi‘c （2010）。第二种方法是使用一个专用的R包来进行结构方程的建模。在撰写本文时，除了lavaan之外，还有两种备选方案可供选择。由约翰·福克斯开发的sem软件包从2001年就开始存在了，在很长一段时间里，它是R语言中唯一的SEM包。第二个包是OpenMx(Boker,Neale, Maes, Wilde, Spiegel, Brick, Spies, Estabrook, Kenny, Bates, Mehta, and Fox 2011)，可从http://openmx.psyc.virginia.edu/获得。正如软件包的名称所暗示的那样，OpenMx是对Mx程序的重写，它包括三个部分：一个用R语言编写的前端，一个用C语言编写的后端，以及一个第三方商业优化器（NPSOL）。OpenMx除了NPSOL优化器外的所有部分都是开源的。 本文结构如下。首先，我们简单介绍了作者开发lavaan包的原因与过程。接下来，我们介绍了lavaan的模型语法，并用两个来自SEM文献的著名的例子（一个CFA例子和一个SEM例子）来说明lavaan包在实践中的使用。最后，我们将讨论多组模型的应用。
我们为什么需要lavaan？
如上所述，许多SEM软件包都是可用的，包括免费的和商用的，其中还包括可以在R环境中运行的两个软件包。那么，为什么还需要另一个SEM软件包呢？这个问题的答案有三方面： 1. lavaan的目标受众是大量的应用研究人员，他们需要SEM软件来回答他们的问题。许多应用研究人员以前没有使用过R语言，并且习惯了商用SEM软件。应用研究人员通常喜欢直观且有丰富的建模功能的软件，而lavaan试图实现这两个目标。 2. lavaan旨在吸引那些教授SEM课程或指导SEM研讨会的人；理想情况下，教师应该有一个易于使用且完整的SEM程序，并且要容易在计算机教室里安装。 3. lavaan的目标是吸引在SEM领域工作的统计学家。为了实现更多功能，制造一个能够直接访问代码的开源SEM程序是有利的。 第一点可以说是最难实现的目标。如果作者希望商用SEM程序的用户使用lavaan，就必须有足够吸引人的点。lavaan是免费的，这通常是无关紧要的。对许多应用研究人员来说，最重要的是(1)该软件使用简单直观，(2)该软件拥有他们想要的所有功能，(3)lavaan的结果与他们目前的商用软件非常接近。为了确保软件易于使用和直观，作者开发了“lavaan模型语法”，它提供了一种拟合结构方程模型的简洁方法。许多应用研究人员经常要求的两个功能是支持非正常（但连续的）数据，和处理缺失的数据。这两个功能在lavaan包中都得到了实现。最后，为了确保lavaan报告的结果与商用软件的输出相似，lavaan中的所有拟合函数都包含一个模拟选项。如果mimic= “Mplus”，lavaan将努力产生类似Mplus的输出。如果mimic=“EQS”，lavaan产生的输出将接近于EQS的输出。在未来发布的lavaan中，我们计划添加mimic=“LISREL”和mimic=“AMOS”（但这些程序的用户目前可以使用mimic=“EQS”进行替代）。 第二点是针对我们这些在课堂上或研讨会上教授SEM技术的人。对教师来说，lavaan是免费的这个事实很重要。如果该软件是免费的，那么就不再需要安装有限制的“学生版本”来教授SEM课程。当然，教师也会喜欢一个简单和直观的软件，这样他们就可以花更多的时间来讨论和解释SEM分析的实质性结果，而不是花时间来解释一个软件的笨拙的模型语法。最后，模拟选项可以使lavaana平稳过渡到其他的商用软件。第三点是针对在结构方程建模领域工作的专业统计学家。长期以来，该领域一直由闭源代码的商用软件所主导。在实践中，这意味着该领域的许多技术贡献都是由那些能够访问源代码的研究小组（及其合作者）实现的。不幸的是，由于缺乏可以使用的开源软件，其他研究人员很难实现新的计算和统计方面的进步。这与统计遗传学或神经成像等其他领域形成了鲜明的对比，这些领域几乎取得了指数级的进展，部分原因是这两个领域都有开源软件包来促进发展。因此，作者使lavaan包完全开源，而不依赖于商用组件。
安装与调用
install.packages("lavaan")library(lavaan)例1：验证性因子分析(CFA)
lavaan包提供了一个内置数据集叫做HolzingerSwineford，输入以下命令查看数据集描述：
??HolzingerSwineford数据形式可以通过一下命令获取：
data <- HolzingerSwineford1939
data此数据集包含了来自两个学校的七、八年级孩子的智力能力测验分数。在我们的版本里，只包含原有26个测试中的9个，这9个测试分数作为9个测量变量分别对应3个潜变量：
视觉因子(visual)对应x1，x2，x3
文本因子(textual)对应x4，x5，x6
速度因子(speed)对应x7，x8，x9
模型如下图所示：
在接下来的内容中，我们将把这个3因素模型称为“H&S模型”，可用图1表示。注意，图中的路径图是简化过的，它不表示观测变量的残差，也不表示外生潜在变量的方差。尽管如此，它还是抓住了模型的精髓。在讨论这个模型的lavaan语法之前，有必要先确定这个模型中的自由参数。该模型中有三个潜变量(因子)，每个变量有三个指标，因此需要估计九个因子的负荷。潜变量之间还存在三个协方差——也就是另外三个参数。这12个参数在路径图中分别用单箭头和双箭头表示。此外，我们还需要估计9个观测变量的残差方差和潜在变量的方差，从而得到12个额外的自由参数。如上所说，总共有24个参数。但是模型还没有被定义好，因为我们还需要度量潜变量。通常有两种方法可以做到这一点:
对于每个潜变量，将其中一个指标(通常是第一个)的因子负荷固定为常数(通常为1.0)；
标准化潜变量的方差。
不管怎样，我们固定了其中的三个参数，剩下的21个参数仍然是自由的。表2是由parTable()方法生成的，包含了该模型的所有相关参数的概述，包括三个固定的因子负荷。表中的每一行对应一个参数。“rhs”、“op”和”lhs”列定义了模型的参数。带有”=~“操作符的所有参数都是因子负荷，而带有”~~“操作符的所有参数都是方差或协方差。“free”列中的非零元素是模型的自由参数。“free”列中的零元素对应于固定参数，其值在”ustart”列中找到。
简单来说，例1是一个定义模型、拟合模型、提取结果的过程。 建立模型代码如下：
HS.model <- 'visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9'
# 然后拟合cfa函数，第一个参数是模型，第二个参数是数据集
fit <- cfa(HS.model, data = HolzingerSwineford1939)
# 再通过summary函数给出结果
summary(fit, fit.measure = TRUE)## lavaan 0.6.15 ended normally after 35 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 21
##
## Number of observations 301
##
## Model Test User Model:
##
## Test statistic 85.306
## Degrees of freedom 24
## P-value (Chi-square) 0.000
##
## Model Test Baseline Model:
##
## Test statistic 918.852
## Degrees of freedom 36
## P-value 0.000
##
## User Model versus Baseline Model:
##
## Comparative Fit Index (CFI) 0.931
## Tucker-Lewis Index (TLI) 0.896
##
## Loglikelihood and Information Criteria:
##
## Loglikelihood user model (H0) -3737.745
## Loglikelihood unrestricted model (H1) -3695.092
##
## Akaike (AIC) 7517.490
## Bayesian (BIC) 7595.339
## Sample-size adjusted Bayesian (SABIC) 7528.739
##
## Root Mean Square Error of Approximation:
##
## RMSEA 0.092
## 90 Percent confidence interval - lower 0.071
## 90 Percent confidence interval - upper 0.114
## P-value H_0: RMSEA <= 0.050 0.001
## P-value H_0: RMSEA >= 0.080 0.840
##
## Standardized Root Mean Square Residual:
##
## SRMR 0.065
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 0.554 0.100 5.554 0.000
## x3 0.729 0.109 6.685 0.000
## textual =~
## x4 1.000
## x5 1.113 0.065 17.014 0.000
## x6 0.926 0.055 16.703 0.000
## speed =~
## x7 1.000
## x8 1.180 0.165 7.152 0.000
## x9 1.082 0.151 7.155 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.408 0.074 5.552 0.000
## speed 0.262 0.056 4.660 0.000
## textual ~~
## speed 0.173 0.049 3.518 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.549 0.114 4.833 0.000
## .x2 1.134 0.102 11.146 0.000
## .x3 0.844 0.091 9.317 0.000
## .x4 0.371 0.048 7.779 0.000
## .x5 0.446 0.058 7.642 0.000
## .x6 0.356 0.043 8.277 0.000
## .x7 0.799 0.081 9.823 0.000
## .x8 0.488 0.074 6.573 0.000
## .x9 0.566 0.071 8.003 0.000
## visual 0.809 0.145 5.564 0.000
## textual 0.979 0.112 8.737 0.000
## speed 0.384 0.086 4.451 0.000输出由三部分组成。前九行称为标题，包含以下内容信息:
lavaan版本号
优化是否正常结束，需要多少次迭代
使用的估计器(这里:ML，表示最大似然)
用于为该估计器找到最佳拟合参数值的优化器(这里:NLMINB)
模型参数的数量(这里为21)
分析中有效使用的观察值的数量(这里是301)
一个叫做“模型测试用户模型”的部分:它提供了一个测试统计量、自由度和一个用户指定的模型的p值。
下一节包含了其他的度量值，之所以显示它，是因为我们使用了可选参数fit.measures = TRUE。它从Model Test Baseline Model开始，并以SRMR的值做结尾。最后一部分包含参数估计。它包含了标准误差的信息。然后，它将模型中包含的所有自由(和固定)参数制成表格。通常，首先是显示潜变量，然后是协方差和(残差)方差。第一列(Estimate)包含每个模型参数的(估计或固定)参数值;第二列(Std.err)包含每个估计参数的标准误差;第三列(z值)包含Wald统计量(即简单地通过将参数值除以其标准误差获得)，最后一列(P(>|z|))包含检验总体中参数值为零的零假设的p值。 注意，在Variances部分中，观察到的变量名称前面有一个点。这是因为他们是因变量(或内生变量)(由潜变量预测)，因此，输出的方差是残差方差的估计值。相比之下，在潜变量名之前没有点，因为它们是这个模型中的外生变量(没有单箭头指向它们)。这里的方差是潜变量的估计总方差。
为了总结第一个例子，我们总结了这个三因素模型所需的完整代码:
# load the lavaan package (only needed once per session)
library(lavaan)
# specify the model
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9'
# fit the model
fit <- cfa(HS.model, data = HolzingerSwineford1939)
# display summary output
summary(fit, fit.measures = TRUE)简单地复制这段代码并粘贴到R中应该可以工作。典型的语法有:
使用lavaan模型语法指定您的模型。在本例中，只有潜变量定义被使用。在下面的示例中，将使用其他公式类型。
拟合模型。这需要一个包含观察到的变量(或者样本)的数据协集。在本例中，我们使用了cfa()函数。lavaan包中的其他函数是sem()和growth()，它们用于拟合全结构方程模型和生长曲线模型。
从拟合模型中提取信息。
例2：结构方程(SEM)
在第二个示例中，我们将使用内置的PoliticalDemocracy数据集。这是Bollen在他1989年关于结构方程建模的书中(以及其他地方)使用的数据集。要了解有关该数据集的更多信息，请参阅其帮助页面和其中的参考资料。我们想要拟合的模型如图3所示，各变量的含义如下图4所示。
PoliticalDemocracy中的数据形式如下：
data <- PoliticalDemocracy
data
该模型对应的lavaan语法如下:
model <- '# measurement model 测量模型
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem65 =~ y5 + y6 + y7 + y8
# regressions 回归
dem60 ~ ind60
dem65 ~ ind60 + dem60
# residual correlations 残余相关
y1 ~~ y5
y2 ~~ y4 + y6
y3 ~~ y7
y4 ~~ y8
y6 ~~ y8'在本例中，我们使用三种不同的公式类型:潜在变量定义(使用=~操作符)、回归公式(使用~操作符)和(协)方差公式(使用~~操作符)。回归公式与R语言中的普通公式类似，(协)方差公式通常有如下形式: 变量 ~~ 变量
变量可以是观察变量，也可以是潜在变量。如果两个变量名称相同，则该表达式指的是该变量的方差(或残差方差)。如果两个变量的名称不同，则表示这两个变量之间的(残差)协方差。lavaan包会自动区分方差和残差方差。
在我们的例子中，表达式y1 ~~ y5允许两个观测变量的残差相关。如果我们认为y1和y5这两个变量有一些相同之处，而这些共同点没有被潜在变量捕捉到，我们就用这条式子来表示。在这个例子中，这两个变量指的是在两个不同的年份(分别是1960年和1965年)测量但有相同含义的分数。请注意，两个表达式y2 ~~ y4和y2 ~~ y6，可以组合成表达式y2 ~~ y4 + y6，这是因为~~运算符左边的变量（y2）是相同的。这只是一个速记符号。
为了拟合模型并查看结果，我们可以输入:
# 拟合SEM
fit <- sem(model, data = PoliticalDemocracy)
# 提取结果
summary(fit, standardized = TRUE) ## lavaan 0.6.15 ended normally after 68 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 31
##
## Number of observations 75
##
## Model Test User Model:
##
## Test statistic 38.125
## Degrees of freedom 35
## P-value (Chi-square) 0.329
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## ind60 =~
## x1 1.000 0.670 0.920
## x2 2.180 0.139 15.742 0.000 1.460 0.973
## x3 1.819 0.152 11.967 0.000 1.218 0.872
## dem60 =~
## y1 1.000 2.223 0.850
## y2 1.257 0.182 6.889 0.000 2.794 0.717
## y3 1.058 0.151 6.987 0.000 2.351 0.722
## y4 1.265 0.145 8.722 0.000 2.812 0.846
## dem65 =~
## y5 1.000 2.103 0.808
## y6 1.186 0.169 7.024 0.000 2.493 0.746
## y7 1.280 0.160 8.002 0.000 2.691 0.824
## y8 1.266 0.158 8.007 0.000 2.662 0.828
##
## Regressions:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## dem60 ~
## ind60 1.483 0.399 3.715 0.000 0.447 0.447
## dem65 ~
## ind60 0.572 0.221 2.586 0.010 0.182 0.182
## dem60 0.837 0.098 8.514 0.000 0.885 0.885
##
## Covariances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .y1 ~~
## .y5 0.624 0.358 1.741 0.082 0.624 0.296
## .y2 ~~
## .y4 1.313 0.702 1.871 0.061 1.313 0.273
## .y6 2.153 0.734 2.934 0.003 2.153 0.356
## .y3 ~~
## .y7 0.795 0.608 1.308 0.191 0.795 0.191
## .y4 ~~
## .y8 0.348 0.442 0.787 0.431 0.348 0.109
## .y6 ~~
## .y8 1.356 0.568 2.386 0.017 1.356 0.338
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .x1 0.082 0.019 4.184 0.000 0.082 0.154
## .x2 0.120 0.070 1.718 0.086 0.120 0.053
## .x3 0.467 0.090 5.177 0.000 0.467 0.239
## .y1 1.891 0.444 4.256 0.000 1.891 0.277
## .y2 7.373 1.374 5.366 0.000 7.373 0.486
## .y3 5.067 0.952 5.324 0.000 5.067 0.478
## .y4 3.148 0.739 4.261 0.000 3.148 0.285
## .y5 2.351 0.480 4.895 0.000 2.351 0.347
## .y6 4.954 0.914 5.419 0.000 4.954 0.443
## .y7 3.431 0.713 4.814 0.000 3.431 0.322
## .y8 3.254 0.695 4.685 0.000 3.254 0.315
## ind60 0.448 0.087 5.173 0.000 1.000 1.000
## .dem60 3.956 0.921 4.295 0.000 0.800 0.800
## .dem65 0.172 0.215 0.803 0.422 0.039 0.039与上例不同，这里我们忽略了参数fit.measure = TRUE，用standardized = TRUE来标准化参数值）,结果如上。 可以看出sem()和cfa()这两个函数很相似，实际上，这两个函数现在几乎是一样的，但在将来可能发生变化。在summary()函数中，我们省略了参数fit.measure = TRUE，因此，您只能得到基本的卡方检验统计量。我们改用参数standardized = TRUE来输出标准化参数值。在添加了standardized参数后，出现了Std.lv， Std.all两列，前者只有潜变量被标准化了，后者为所有变量都标准化了，也被称为“完全标准化解”。
这个案例的完整代码如下所示：
model <- '# measurement model 测量模型
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem65 =~ y5 + y6 + y7 + y8
# regressions 回归
dem60 ~ ind60
dem65 ~ ind60 + dem60
# residual correlations 残余相关
y1 ~~ y5
y2 ~~ y4 + y6
y3 ~~ y7
y4 ~~ y8
y6 ~~ y8'
# 拟合SEM
fit <- sem(model, data = PoliticalDemocracy)
# 提取结果
summary(fit, standardized = TRUE)
#与上例不同，这里我们忽略了参数fit.measure = TRUE，用standardized = TRUE来标准化参数值）例3：结构方程的多组分析
lavaan包完全支持多组分析。为了进行多组分析，你需要在数据组中加上组变量的名称。默认情况下，所有组会拟合相同的模型。在以下示例中，我们针对两个学校(Pasteur和Grant-White，可以看到数据集中有关学校的变量名为school)拟合例1的模型：
HS.model <- 'visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9'
fit <- cfa(HS.model,
data = HolzingerSwineford1939,
group = "school")
summary(fit)## lavaan 0.6.15 ended normally after 57 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 60
##
## Number of observations per group:
## Pasteur 156
## Grant-White 145
##
## Model Test User Model:
##
## Test statistic 115.851
## Degrees of freedom 48
## P-value (Chi-square) 0.000
## Test statistic for each group:
## Pasteur 64.309
## Grant-White 51.542
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
##
## Group 1 [Pasteur]:
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 0.394 0.122 3.220 0.001
## x3 0.570 0.140 4.076 0.000
## textual =~
## x4 1.000
## x5 1.183 0.102 11.613 0.000
## x6 0.875 0.077 11.421 0.000
## speed =~
## x7 1.000
## x8 1.125 0.277 4.057 0.000
## x9 0.922 0.225 4.104 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.479 0.106 4.531 0.000
## speed 0.185 0.077 2.397 0.017
## textual ~~
## speed 0.182 0.069 2.628 0.009
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .x1 4.941 0.095 52.249 0.000
## .x2 5.984 0.098 60.949 0.000
## .x3 2.487 0.093 26.778 0.000
## .x4 2.823 0.092 30.689 0.000
## .x5 3.995 0.105 38.183 0.000
## .x6 1.922 0.079 24.321 0.000
## .x7 4.432 0.087 51.181 0.000
## .x8 5.563 0.078 71.214 0.000
## .x9 5.418 0.079 68.440 0.000
## visual 0.000
## textual 0.000
## speed 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.298 0.232 1.286 0.198
## .x2 1.334 0.158 8.464 0.000
## .x3 0.989 0.136 7.271 0.000
## .x4 0.425 0.069 6.138 0.000
## .x5 0.456 0.086 5.292 0.000
## .x6 0.290 0.050 5.780 0.000
## .x7 0.820 0.125 6.580 0.000
## .x8 0.510 0.116 4.406 0.000
## .x9 0.680 0.104 6.516 0.000
## visual 1.097 0.276 3.967 0.000
## textual 0.894 0.150 5.963 0.000
## speed 0.350 0.126 2.778 0.005
##
##
## Group 2 [Grant-White]:
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 0.736 0.155 4.760 0.000
## x3 0.925 0.166 5.583 0.000
## textual =~
## x4 1.000
## x5 0.990 0.087 11.418 0.000
## x6 0.963 0.085 11.377 0.000
## speed =~
## x7 1.000
## x8 1.226 0.187 6.569 0.000
## x9 1.058 0.165 6.429 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.408 0.098 4.153 0.000
## speed 0.276 0.076 3.639 0.000
## textual ~~
## speed 0.222 0.073 3.022 0.003
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .x1 4.930 0.095 51.696 0.000
## .x2 6.200 0.092 67.416 0.000
## .x3 1.996 0.086 23.195 0.000
## .x4 3.317 0.093 35.625 0.000
## .x5 4.712 0.096 48.986 0.000
## .x6 2.469 0.094 26.277 0.000
## .x7 3.921 0.086 45.819 0.000
## .x8 5.488 0.087 63.174 0.000
## .x9 5.327 0.085 62.571 0.000
## visual 0.000
## textual 0.000
## speed 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.715 0.126 5.676 0.000
## .x2 0.899 0.123 7.339 0.000
## .x3 0.557 0.103 5.409 0.000
## .x4 0.315 0.065 4.870 0.000
## .x5 0.419 0.072 5.812 0.000
## .x6 0.406 0.069 5.880 0.000
## .x7 0.600 0.091 6.584 0.000
## .x8 0.401 0.094 4.249 0.000
## .x9 0.535 0.089 6.010 0.000
## visual 0.604 0.160 3.762 0.000
## textual 0.942 0.152 6.177 0.000
## speed 0.461 0.118 3.910 0.000如果你想要固定参数，或者自己设定初值，可以使用相同的预乘法技术，但是此时单个参数将被参数向量取代，该参数向量每组都应该有一个。如果您使用单个参数而不是向量(不推荐使用)，则该参数将应用于所有组。如果您指定单个标签，则会出现警告，因为这将表示跨组的等式约束。比如：
HS.model <- 'visual =~ x1 + 0.5*x2 + c(0.6, 0.8)*x3
textual =~ x4 + start(c(1.2, 0.6))*x5 + a*x6
speed =~ x7 + x8 + x9'在潜在因子“visual”的定义中，我们将指标x3在第一组中的因子负荷设定为“0.6”，在第二组中的因子负荷设定为“0.8”，指标x2在两组中的的因子负荷均固定为“0.5”。在“textual”因子的定义中，我们为x5指标提供了两个不同的初值（每组一个）。另外，我们将指标x6在第一组中的因子负荷标记为a1，在第二组中标记为a2。直接写a*x6可能很有诱惑力。但是，在多组设置中使用单个标签具有双重效果:它给了两个组的x6的因子负载相同的标签a，因此，这两个参数现在被约束为相等。因为这可能不是使用者的本意，lavaan将对此发出警告。如果确实想这样做，最好使用向量标签:c(a, a)*x6。
为了检验修饰符的效果，我们对模型进行了一些修改:
fit <- cfa(HS.model,
data = HolzingerSwineford1939,
group = "school")
summary(fit)## lavaan 0.6.15 ended normally after 57 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 60
##
## Number of observations per group:
## Pasteur 156
## Grant-White 145
##
## Model Test User Model:
##
## Test statistic 115.851
## Degrees of freedom 48
## P-value (Chi-square) 0.000
## Test statistic for each group:
## Pasteur 64.309
## Grant-White 51.542
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
##
## Group 1 [Pasteur]:
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 0.394 0.122 3.220 0.001
## x3 0.570 0.140 4.076 0.000
## textual =~
## x4 1.000
## x5 1.183 0.102 11.613 0.000
## x6 0.875 0.077 11.421 0.000
## speed =~
## x7 1.000
## x8 1.125 0.277 4.057 0.000
## x9 0.922 0.225 4.104 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.479 0.106 4.531 0.000
## speed 0.185 0.077 2.397 0.017
## textual ~~
## speed 0.182 0.069 2.628 0.009
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .x1 4.941 0.095 52.249 0.000
## .x2 5.984 0.098 60.949 0.000
## .x3 2.487 0.093 26.778 0.000
## .x4 2.823 0.092 30.689 0.000
## .x5 3.995 0.105 38.183 0.000
## .x6 1.922 0.079 24.321 0.000
## .x7 4.432 0.087 51.181 0.000
## .x8 5.563 0.078 71.214 0.000
## .x9 5.418 0.079 68.440 0.000
## visual 0.000
## textual 0.000
## speed 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.298 0.232 1.286 0.198
## .x2 1.334 0.158 8.464 0.000
## .x3 0.989 0.136 7.271 0.000
## .x4 0.425 0.069 6.138 0.000
## .x5 0.456 0.086 5.292 0.000
## .x6 0.290 0.050 5.780 0.000
## .x7 0.820 0.125 6.580 0.000
## .x8 0.510 0.116 4.406 0.000
## .x9 0.680 0.104 6.516 0.000
## visual 1.097 0.276 3.967 0.000
## textual 0.894 0.150 5.963 0.000
## speed 0.350 0.126 2.778 0.005
##
##
## Group 2 [Grant-White]:
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 0.736 0.155 4.760 0.000
## x3 0.925 0.166 5.583 0.000
## textual =~
## x4 1.000
## x5 0.990 0.087 11.418 0.000
## x6 0.963 0.085 11.377 0.000
## speed =~
## x7 1.000
## x8 1.226 0.187 6.569 0.000
## x9 1.058 0.165 6.429 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.408 0.098 4.153 0.000
## speed 0.276 0.076 3.639 0.000
## textual ~~
## speed 0.222 0.073 3.022 0.003
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .x1 4.930 0.095 51.696 0.000
## .x2 6.200 0.092 67.416 0.000
## .x3 1.996 0.086 23.195 0.000
## .x4 3.317 0.093 35.625 0.000
## .x5 4.712 0.096 48.986 0.000
## .x6 2.469 0.094 26.277 0.000
## .x7 3.921 0.086 45.819 0.000
## .x8 5.488 0.087 63.174 0.000
## .x9 5.327 0.085 62.571 0.000
## visual 0.000
## textual 0.000
## speed 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.715 0.126 5.676 0.000
## .x2 0.899 0.123 7.339 0.000
## .x3 0.557 0.103 5.409 0.000
## .x4 0.315 0.065 4.870 0.000
## .x5 0.419 0.072 5.812 0.000
## .x6 0.406 0.069 5.880 0.000
## .x7 0.600 0.091 6.584 0.000
## .x8 0.401 0.094 4.249 0.000
## .x9 0.535 0.089 6.010 0.000
## visual 0.604 0.160 3.762 0.000
## textual 0.942 0.152 6.177 0.000
## speed 0.461 0.118 3.910 0.0003.1 在部分组中固定参数
有时，我们希望在除了一个特定的组之外的所有组中固定一个参数的值。在特定组中，我们希望自由地估计该参数的值。这个参数的修饰符仍然是一个向量，它包含每个组对这个参数的设定值，但是我们可以使用NA来强制让一个参数在一个(或多个)组中是自由的。假设我们有四组。我们用三个指标来定义一个潜变量(比如f)。我们希望在除第二组之外的所有组中，将指标item2的因子负荷固定为1.0。我们可以这样写：
f =~ item1 + c(1,NA,1,1)*item2 + item33.2 约束一个参数使其在各组中相等
如果你想要将一个或多个参数约束为跨组相等，则需要为它们提供相同的标签。例如，若要约束指标x3的因子负载在两组中相等，可以这样写:
HS.model <- 'visual =~ x1 + x2 + c(v3,v3)*x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9'和之前一样，无论是在组内还是跨组，相同的标签都意味着相同的参数。
3.3 约束一组参数使其在各组中相等
虽然提供相同的标签是为几个参数指定相等约束的一种非常灵活的方法，但有一种更方便的方法可以对整个参数集(例如:所有因子负荷或所有截距)施加相等约束。我们称这种类型的约束为“组间相等约束”，它们可以通过通过group.equal()来实现。例如，要约束(所有)因子负载在组间相等，您可以按照以下步骤进行:
HS.model <- 'visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9'
fit <- cfa(HS.model,
data = HolzingerSwineford1939,
group = "school",
group.equal = c("loadings"))
summary(fit)## lavaan 0.6.15 ended normally after 42 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 60
## Number of equality constraints 6
##
## Number of observations per group:
## Pasteur 156
## Grant-White 145
##
## Model Test User Model:
##
## Test statistic 124.044
## Degrees of freedom 54
## P-value (Chi-square) 0.000
## Test statistic for each group:
## Pasteur 68.825
## Grant-White 55.219
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
##
## Group 1 [Pasteur]:
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 (.p2.) 0.599 0.100 5.979 0.000
## x3 (.p3.) 0.784 0.108 7.267 0.000
## textual =~
## x4 1.000
## x5 (.p5.) 1.083 0.067 16.049 0.000
## x6 (.p6.) 0.912 0.058 15.785 0.000
## speed =~
## x7 1.000
## x8 (.p8.) 1.201 0.155 7.738 0.000
## x9 (.p9.) 1.038 0.136 7.629 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.416 0.097 4.271 0.000
## speed 0.169 0.064 2.643 0.008
## textual ~~
## speed 0.176 0.061 2.882 0.004
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .x1 4.941 0.093 52.991 0.000
## .x2 5.984 0.100 60.096 0.000
## .x3 2.487 0.094 26.465 0.000
## .x4 2.823 0.093 30.371 0.000
## .x5 3.995 0.101 39.714 0.000
## .x6 1.922 0.081 23.711 0.000
## .x7 4.432 0.086 51.540 0.000
## .x8 5.563 0.078 71.087 0.000
## .x9 5.418 0.079 68.153 0.000
## visual 0.000
## textual 0.000
## speed 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.551 0.137 4.010 0.000
## .x2 1.258 0.155 8.117 0.000
## .x3 0.882 0.128 6.884 0.000
## .x4 0.434 0.070 6.238 0.000
## .x5 0.508 0.082 6.229 0.000
## .x6 0.266 0.050 5.294 0.000
## .x7 0.849 0.114 7.468 0.000
## .x8 0.515 0.095 5.409 0.000
## .x9 0.658 0.096 6.865 0.000
## visual 0.805 0.171 4.714 0.000
## textual 0.913 0.137 6.651 0.000
## speed 0.305 0.078 3.920 0.000
##
##
## Group 2 [Grant-White]:
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 (.p2.) 0.599 0.100 5.979 0.000
## x3 (.p3.) 0.784 0.108 7.267 0.000
## textual =~
## x4 1.000
## x5 (.p5.) 1.083 0.067 16.049 0.000
## x6 (.p6.) 0.912 0.058 15.785 0.000
## speed =~
## x7 1.000
## x8 (.p8.) 1.201 0.155 7.738 0.000
## x9 (.p9.) 1.038 0.136 7.629 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.437 0.099 4.423 0.000
## speed 0.314 0.079 3.958 0.000
## textual ~~
## speed 0.226 0.072 3.144 0.002
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .x1 4.930 0.097 50.763 0.000
## .x2 6.200 0.091 68.379 0.000
## .x3 1.996 0.085 23.455 0.000
## .x4 3.317 0.092 35.950 0.000
## .x5 4.712 0.100 47.173 0.000
## .x6 2.469 0.091 27.248 0.000
## .x7 3.921 0.086 45.555 0.000
## .x8 5.488 0.087 63.257 0.000
## .x9 5.327 0.085 62.786 0.000
## visual 0.000
## textual 0.000
## speed 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.645 0.127 5.084 0.000
## .x2 0.933 0.121 7.732 0.000
## .x3 0.605 0.096 6.282 0.000
## .x4 0.329 0.062 5.279 0.000
## .x5 0.384 0.073 5.270 0.000
## .x6 0.437 0.067 6.576 0.000
## .x7 0.599 0.090 6.651 0.000
## .x8 0.406 0.089 4.541 0.000
## .x9 0.532 0.086 6.202 0.000
## visual 0.722 0.161 4.490 0.000
## textual 0.906 0.136 6.646 0.000
## speed 0.475 0.109 4.347 0.000用.p2., .p3.， .p5.，…标记的输出结果在已自动施加相等约束。可以添加更多的“组间相等约束”。除了因子负荷外，group.equal还可以用来施加其他的相等约束。
intercepts: 观测变量的截距
means: 在变量的截距/均值
residuals: 观测变量的残差方差
residual.covariances: 观测变量的残差协方差
lv.variances: 潜在变量的(残差)方差
lv.covariances: 潜在变量的(残差)协方差
regressions:模型中的所有回归系数
如果省略group.equal，每组的所有参数都是自由估计的(但模型结构是相同的)。
但是，如果您希望除了一两个参数是自由的以外，每一组的所有参数都是相同的(比如所有因子加载和截距)，该怎么办呢?在这个情况下，您可以使用参数group.partial，其中要包含那些需要保持自由的参数的名称。例如:
fit <- cfa(HS.model,
data = HolzingerSwineford1939,
group = "school",
group.equal = c("loadings", "intercepts"),
group.partial = c("visual=~x2", "x7~1"))3.4 测量不变性
在我们比较各个组的潜在变量均值之前，我们首先需要测量不变性。因子分析的测量不变性(Measurement Invariance)是指测验的观察变量与潜变量之间的关系在不同组别或时间点上保持不变。测量不变性的成立是进行有意义组间比较的重要前提。当数据连续时，测量不变性测试涉及固定序列的模型比较测试。一个典型的序列包含三个模型:
模型1:配置不变性。所有组因素结构都是相同的。
模型2:弱不变性。各组之间的因子负荷被限制为相等。
模型3:强不变性。因子载荷和截距在组间被限制为相等。
在lavaan中，我们可以进行如下操作:
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9 '
# configural invariance
fit1 <- cfa(HS.model, data = HolzingerSwineford1939, group = "school")
# weak invariance
fit2 <- cfa(HS.model, data = HolzingerSwineford1939, group = "school",
group.equal = "loadings")
# strong invariance
fit3 <- cfa(HS.model, data = HolzingerSwineford1939, group = "school",
group.equal = c("intercepts", "loadings"))
# model comparison tests
lavTestLRT(fit1, fit2, fit3)##
## Chi-Squared Difference Test
##
## Df AIC BIC Chisq Chisq diff RMSEA Df diff Pr(>Chisq)
## fit1 48 7484.4 7706.8 115.85
## fit2 54 7480.6 7680.8 124.04 8.192 0.049272 6 0.2244
## fit3 60 7508.6 7686.6 164.10 40.059 0.194211 6 4.435e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1lavTestLRT()函数可用于模型比较测试。因为我们提供了三个模型拟合，它将产生两个测试:第一个测试将第一个模型与第二个模型进行比较，而第二个测试将第二个模型与第三个模型进行比较。由于第一个p值不显著，我们可以得出结论，该数据集支持弱不变性(等因子负载)。然而，由于第二个p值是显著的，所以不具有强不变性。因此，直接比较两组的潜变量均值是不明智的。
我们也可以用semTools包的measurementInvariance()函数来按特定次序执行多组分析，每个模型都与基线模型比较并且对先前的模型进行卡方差异性检验，并且展示拟合度量差别。不过这个方法仍旧比较原始。
library(semTools)
measurementInvariance(HS.model,
data = HolzingerSwineford1939,
group = "school")例3：增长曲线模型
(一) 相关函数
growth()
(1) 简要描述：
用于拟合一个增长曲线模型，但只有当模型中全部潜在变量都是增长因子时才有用，如果要解决更复杂的模型，使用lavaan()会更好
(2) 用法：
growth(model = NULL, data = NULL, ordered = NULL, sampling.weights = NULL,sample.cov = NULL, sample.mean = NULL, sample.th = NULL,sample.nobs = NULL, group = NULL, cluster = NULL, constraints = ““, WLS.V = NULL, NACOV = NULL, ov.order =”model”,…)
(3) 常用参数介绍：
model: 用于介绍用户设定好的结构方程模型，并且这个模型需要使用lavaan模型语法来表示。
data：表示输入的数据，包含模型中使用的观察变量。如果某些变量被声明为有序因子，lavaan将把它们视为有序变量。
ordered：表示指定变量是否为有序的，默认为FALSE
(4) 查看结果：
运行growth函数后，可以查看拟合优度、参数估计值、标准误和置信区间等模型统计信息。拟合优度可以反映模型对实际数据的拟合程度，一般使用常见的拟合度量指标，如CFI、RMSEA、SRMR等。
(二) 应用例子：含时间变化协变量的线性增长模型
第一步 调用lavaan包：
library("lavaan")第二步定义线性增长模型
model.syntax <- '
# 定义两个系数固定的含时间的隐变量。其中，s具有一个随时间递增的线性趋势。
i =~ 1*t1 + 1*t2 + 1*t3 + 1*t4
s =~ 0*t1 + 1*t2 + 2*t3 + 3*t4
# 定义隐变量i和s来探讨对隐变量增长趋势的影响。
i ~ x1 + x2
s ~ x1 + x2
# 定义时间变化协变量t1、t2、t3、t4，并分别对指定的外部变量进行回归预测
t1 ~ c1
t2 ~ c2
t3 ~ c3
t4 ~ c4
'代码详细解释：
i ~ x1 + x2
s ~ x1 + x2
i ~ x1 + x2表示隐变量i的取值随着外部变量x1和x2的取值的变化而发生变化，这里的~符号表示”受影响于”的意思。在回归模型中，x1和x2分别是自变量，可以用来预测隐变量i的取值。系数可以被理解为不同自变量对于隐变量影响的权重。同理，s ~ x1 + x2表示隐变量s的取值随着外部变量x1和x2的取值而发生变化，而x1和x2的系数则表示这些变量对隐变量s增长趋势的影响权重。
对外部变量 x1 和 x2 进行回归的主要目的是探究它们是否对隐变量 i 和 s 的增长趋势产生影响。在对 i 和 s 进行建模时，我们希望：（1）对隐变量增长趋势进行描述和预测，（2）确定哪些因素可能会影响这些增长趋势，（3）量化这些影响。
简单地说，对外部变量进行回归分析可以帮助我们建立一种关系模型，用来描述隐变量和外部变量之间的关系。回归模型可以提供一些关键指标，如变量之间的关联程度，影响因素的权重系数以及误差等，这些指标可以帮助我们更好地了解隐变量如何随着时间和其他影响因素的变化而变化。
特别地，在这个模型中，隐变量i和s分别表示随着时间的推移呈线性增长的变量，而x1和x2则代表各自的影响因素，那么通过对i和s分别对x1和x2进行回归，可以帮助我们定量描述这些影响因素对隐变量增长趋势的实际影响，从而更好地理解隐变量的行为模式和演变规律。
t1 ~ c1
t2 ~ c2
t3 ~ c3
t4 ~ c4让时间变量 t1、t2、t3、t4 分别对外部变量 c1、c2、c3、c4 进行回归的主要目的是探究随着时间推移，外部变量的取值变化是否对隐变量的增长趋势产生影响。回归模型可以通过对时间变量和外部变量之间的关系进行建模，帮助我们研究变量之间的动态演化特征。
在上述的模型中，时间变量 t1、t2、t3、t4 分别与外部变量 c1、c2、c3、c4 进行回归，就是要探究时间和外部变量之间是否存在影响，以及这种影响可能具有的动态特性。例如，如果外部变量以某种规律与时间相关，那么外部变量的取值可能会对隐变量的增长趋势产生影响；而且随着时间推移，这种影响还可能会发生变化。通过对时间变量和外部变量的回归分析，我们可以得到对隐变量增长趋势最为关键的外部影响因素及其随时间变化的权重系数。
第三步 进行潜在变量增长模型建模
fit <- growth(model.syntax, data = Demo.growth)
summary(fit)代码详细解释：
fit <- growth(model.syntax, data = Demo.growth)这行代码的执行过程是，把潜在变量增长模型的结构定义为model.syntax，用到的数据为 Demo.growth。函数 growth 会根据这两个参数进行拟合，并把结果保存在 fit 中。
summary(fit)这行代码的执行过程是，对模型拟合结果进行汇总统计，并打印出来。summary() 这个函数可以用来统计和汇总不同类型的拟合结果（比如回归分析、ANOVA分析等等），并且打印这些结果的详细信息。对于隐变量增长模型来说，打印出来的内容通常包括模型的参数估计值、标准差、拟合优度指标等等。这些信息对于评估模型的质量和调整模型的参数都非常重要。
第四步 结果的读取与解释：
## lavaan 0.6.15 ended normally after 31 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 17
##
## Number of observations 400
##
## Model Test User Model:
##
## Test statistic 26.059
## Degrees of freedom 21
## P-value (Chi-square) 0.204
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## i =~
## t1 1.000
## t2 1.000
## t3 1.000
## t4 1.000
## s =~
## t1 0.000
## t2 1.000
## t3 2.000
## t4 3.000
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## i ~
## x1 0.608 0.060 10.134 0.000
## x2 0.604 0.064 9.412 0.000
## s ~
## x1 0.262 0.029 9.198 0.000
## x2 0.522 0.031 17.083 0.000
## t1 ~
## c1 0.143 0.050 2.883 0.004
## t2 ~
## c2 0.289 0.046 6.295 0.000
## t3 ~
## c3 0.328 0.044 7.361 0.000
## t4 ~
## c4 0.330 0.058 5.655 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .i ~~
## .s 0.075 0.040 1.855 0.064
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|)
## .t1 0.000
## .t2 0.000
## .t3 0.000
## .t4 0.000
## .i 0.580 0.062 9.368 0.000
## .s 0.958 0.029 32.552 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .t1 0.580 0.080 7.230 0.000
## .t2 0.596 0.054 10.969 0.000
## .t3 0.481 0.055 8.745 0.000
## .t4 0.535 0.098 5.466 0.000
## .i 1.079 0.112 9.609 0.000
## .s 0.224 0.027 8.429 0.000结果中，Regressions部分报告了外部变量对隐变量以及时间变量的回归系数；Covariances部分报告了两个隐变量之间的协方差；Intercepts部分报告了每个变量的截距参数，从结果中可以看到，s的截距值非常高，意味着该变量再t1时刻的基础上就有很高的取值；Variances部分报告的是每个变量的方差值，方差值越高则说明该变量的取值在样本中的变异性就越大，其中该结果显示这些变量的P值都为0，说明这些方差值都是显著的。
例4：将协方差矩阵作为输入
输入协方差矩阵，并通过包中的getcov函数计算得到完整的协方差矩阵
首先定义一个协方差矩阵数据的数值:
lower <- '
11.834
6.947 9.364
6.819 5.091 12.532
4.783 5.028 7.495 9.986
-3.839 -3.889 -3.841 -3.625 9.610
-21.899 -18.831 -21.748 -18.775 35.522 450.288'利用gercov函数建立完整的协方差矩阵，输出完整的协方差矩阵
wheaton.cov <- getCov(lower, names = c("anomia67", "powerless67",
"anomia71", "powerless71",
"education", "sei"))其中getcov函数里面的names参数的作用是给矩阵的行列进行命名 查看输出的矩阵数据
## anomia67 powerless67 anomia71 powerless71 education sei
## anomia67 11.834 6.947 6.819 4.783 -3.839 -21.899
## powerless67 6.947 9.364 5.091 5.028 -3.889 -18.831
## anomia71 6.819 5.091 12.532 7.495 -3.841 -21.748
## powerless71 4.783 5.028 7.495 9.986 -3.625 -18.775
## education -3.839 -3.889 -3.841 -3.625 9.610 35.522
## sei -21.899 -18.831 -21.748 -18.775 35.522 450.288下面将建立模型
wheaton.model <- '
# latent variables
ses =~ education + sei
alien67 =~ anomia67 + powerless67
alien71 =~ anomia71 + powerless71
# regressions
alien71 ~ alien67 + ses
alien67 ~ ses
# correlated residuals
anomia67 ~~ anomia71
powerless67 ~~ powerless71'通过以上代码建立一个模型，模型中第一部分是定义潜在变量，通过符号=来实现，第二部分是回归模型，通过符号来展示变量间的关系，第三部分是对残差部分进行定义。
模型建立完整后，需要通过sem函数来估计方程中的一些参数，其中参数sample.cov是调用上文中计算得到的完整方差协方差矩阵。sample.nobs指的是观测变量的个数。拟合模型后，通过summary函数显示出整体模型的估计效果。
summary(fit, standardized = TRUE)## lavaan 0.6.15 ended normally after 31 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 17
##
## Number of observations 400
##
## Model Test User Model:
##
## Test statistic 26.059
## Degrees of freedom 21
## P-value (Chi-square) 0.204
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## i =~
## t1 1.000 1.386 0.875
## t2 1.000 1.386 0.660
## t3 1.000 1.386 0.507
## t4 1.000 1.386 0.412
## s =~
## t1 0.000 0.000 0.000
## t2 1.000 0.768 0.366
## t3 2.000 1.536 0.562
## t4 3.000 2.304 0.684
##
## Regressions:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## i ~
## x1 0.608 0.060 10.134 0.000 0.439 0.451
## x2 0.604 0.064 9.412 0.000 0.436 0.419
## s ~
## x1 0.262 0.029 9.198 0.000 0.341 0.351
## x2 0.522 0.031 17.083 0.000 0.679 0.653
## t1 ~
## c1 0.143 0.050 2.883 0.004 0.143 0.089
## t2 ~
## c2 0.289 0.046 6.295 0.000 0.289 0.131
## t3 ~
## c3 0.328 0.044 7.361 0.000 0.328 0.112
## t4 ~
## c4 0.330 0.058 5.655 0.000 0.330 0.091
##
## Covariances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .i ~~
## .s 0.075 0.040 1.855 0.064 0.152 0.152
##
## Intercepts:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .t1 0.000 0.000 0.000
## .t2 0.000 0.000 0.000
## .t3 0.000 0.000 0.000
## .t4 0.000 0.000 0.000
## .i 0.580 0.062 9.368 0.000 0.419 0.419
## .s 0.958 0.029 32.552 0.000 1.247 1.247
##
## Variances:
## Estimate Std.Err z-value P(>|z|) Std.lv Std.all
## .t1 0.580 0.080 7.230 0.000 0.580 0.231
## .t2 0.596 0.054 10.969 0.000 0.596 0.135
## .t3 0.481 0.055 8.745 0.000 0.481 0.064
## .t4 0.535 0.098 5.466 0.000 0.535 0.047
## .i 1.079 0.112 9.609 0.000 0.562 0.562
## .s 0.224 0.027 8.429 0.000 0.379 0.379上述显示出结果包含有参数估计部分，协方差以及房差的参数估计部分。第一列是变量名，第二列是估计得到的参数估计值，第三列是标准误的估计值。展示出不同变量间的结构关系及影响作用。
例5：估计方法：极大似然估计法
模型估计的第一步是构建需要估计的模型
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9'接下来通过cfa函数对模型进行拟合
fit <- cfa(HS.model,data = HolzingerSwineford1939)函数内第一个参数是上述构建的模型，第二个参数data是需要计算的数据集
同理，可以通过summary函数查看估计的结果
summary(fit)## lavaan 0.6.15 ended normally after 35 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 21
##
## Number of observations 301
##
## Model Test User Model:
##
## Test statistic 85.306
## Degrees of freedom 24
## P-value (Chi-square) 0.000
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Latent Variables:
## Estimate Std.Err z-value P(>|z|)
## visual =~
## x1 1.000
## x2 0.554 0.100 5.554 0.000
## x3 0.729 0.109 6.685 0.000
## textual =~
## x4 1.000
## x5 1.113 0.065 17.014 0.000
## x6 0.926 0.055 16.703 0.000
## speed =~
## x7 1.000
## x8 1.180 0.165 7.152 0.000
## x9 1.082 0.151 7.155 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## visual ~~
## textual 0.408 0.074 5.552 0.000
## speed 0.262 0.056 4.660 0.000
## textual ~~
## speed 0.173 0.049 3.518 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .x1 0.549 0.114 4.833 0.000
## .x2 1.134 0.102 11.146 0.000
## .x3 0.844 0.091 9.317 0.000
## .x4 0.371 0.048 7.779 0.000
## .x5 0.446 0.058 7.642 0.000
## .x6 0.356 0.043 8.277 0.000
## .x7 0.799 0.081 9.823 0.000
## .x8 0.488 0.074 6.573 0.000
## .x9 0.566 0.071 8.003 0.000
## visual 0.809 0.145 5.564 0.000
## textual 0.979 0.112 8.737 0.000
## speed 0.384 0.086 4.451 0.000上述估计结果的描述性统计可以清楚看到估计的方法为极大似然估计，模型的参数量、观测量潜在变量的估计结果以及协方差、方差等结果的参数。
例6: 间接效应和中介效应
假定存在这样一个中介分析X->M->Y
下面通过生成随机数的方式来展示该包用于处理间接效应和中介效应
首先通过set.seed函数设立随机种子
set.seed(1234)随机生成变量X，假定变量X服从与标准正态分布，生成100个数值
X <- rnorm(100)根据生成的X计算得到M以及Y的随机数值
M <- 0.6*X + rnorm(100)
Y <- 0.68*M + rnorm(100)将三个变量纳入到数据框中，命名数据框为Data
Data <- data.frame(X = X, Y = Y, M = M)接下通过建立间接效应和中介效应的模型来进行计算和估计
model <- '# direct effect
Y ~ c*X
# mediator
M ~ a*X
Y ~ b*M
# indirect effect (a*b)
ab := a*b
# total effect
total := c + (a*b)'通过sem函数对上述构建的模型进行估计，函数内第一个参数是模型，第二个参数是上述纳入多变量的数据框Data
fit <- sem(model, data = Data)估计结果后，通过summary函数来查看整体估计的结果
summary(fit) ## lavaan 0.6.15 ended normally after 1 iteration
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 5
##
## Number of observations 100
##
## Model Test User Model:
##
## Test statistic 0.000
## Degrees of freedom 0
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## Y ~
## X (c) 0.028 0.109 0.253 0.800
## M ~
## X (a) 0.574 0.103 5.586 0.000
## Y ~
## M (b) 0.768 0.092 8.323 0.000
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .Y 0.898 0.127 7.071 0.000
## .M 1.054 0.149 7.071 0.000
##
## Defined Parameters:
## Estimate Std.Err z-value P(>|z|)
## ab 0.441 0.095 4.638 0.000
## total 0.468 0.123 3.797 0.000根据上述估计结果的描述性统计可以看到，估计方法为极大似然法，参数为5，观测量为100，regression为参数的回归估计结果，Variances为方差，回归结果展示了不同变量间的关系。
例7：修正指标
在修正指标的过程中可以通过在summary函数中加入加参数modindices = TRUE，来得出修正指标。
fit <- cfa(HS.model,
data = HolzingerSwineford1939)
mi <- modindices(fit)
mi[mi$op == "=~",]## lhs op rhs mi epc sepc.lv sepc.all sepc.nox
## 25 visual =~ x4 1.211 0.077 0.069 0.059 0.059
## 26 visual =~ x5 7.441 -0.210 -0.189 -0.147 -0.147
## 27 visual =~ x6 2.843 0.111 0.100 0.092 0.092
## 28 visual =~ x7 18.631 -0.422 -0.380 -0.349 -0.349
## 29 visual =~ x8 4.295 -0.210 -0.189 -0.187 -0.187
## 30 visual =~ x9 36.411 0.577 0.519 0.515 0.515
## 31 textual =~ x1 8.903 0.350 0.347 0.297 0.297
## 32 textual =~ x2 0.017 -0.011 -0.011 -0.010 -0.010
## 33 textual =~ x3 9.151 -0.272 -0.269 -0.238 -0.238
## 34 textual =~ x7 0.098 -0.021 -0.021 -0.019 -0.019
## 35 textual =~ x8 3.359 -0.121 -0.120 -0.118 -0.118
## 36 textual =~ x9 4.796 0.138 0.137 0.136 0.136
## 37 speed =~ x1 0.014 0.024 0.015 0.013 0.013
## 38 speed =~ x2 1.580 -0.198 -0.123 -0.105 -0.105
## 39 speed =~ x3 0.716 0.136 0.084 0.075 0.075
## 40 speed =~ x4 0.003 -0.005 -0.003 -0.003 -0.003
## 41 speed =~ x5 0.201 -0.044 -0.027 -0.021 -0.021
## 42 speed =~ x6 0.273 0.044 0.027 0.025 0.025总结与改进
lavaan 包是一个功能强大的R包，用于结构方程模型 (SEM)分析，但它也存在一些缺点或可改进的方面。以下是一些可能的缺点或改进点: 1. 缺乏一些高级功能：虽然lavaan提供了许多常见的结构方程模型分析功能，但在某些特定领域或复杂模型的分析中可能缺乏一些高级功能，例如，对于长期数据（如时间序列数据）或高级统计方法（如贝叶斯方法）的支持相对有限。 2. 计算效率：lavaan在处理大规模数据集或复杂模型时可能会面临一些计算效率的挑战。一些改进的算法或优化策略可以提高计算速度和效率，以适应更大规模的数据和复杂的模型。