R语言实现潜变量模型构建
结构方程模型是基于变量的协方差矩阵来分析变量之间关系的一种统计方法。作为多元数据分析的重要工具。其可以分成两部分:测量模型和结构模型。所谓测量模型主要是研究潜变量(因子)和显变量(测量指标)的关系;结构模型指的是研究潜变量(因子)之间的关系。今天我们给大家介绍一个集合各种潜变量分析的模型的R包lavaan。
首先,我们看下需要安装的包:
install.packages("lavaan")
install.packages("semPlot")
install.packages("semTools")
接下来我们看下在此包中回归方程的几种表示形式:
y ~ f1 + f2 + x1 + x2 #回归方程~左边为因变量
f1 =~ y1 + y2 + y3 #潜变量表示形式
y1 ~~ y1 # 方差相关性
y1 ~~ y2 # 残差相关性
y1 ~ 1 #只有截距项的回归方程
其中y,y1,y2,y3为因变量,f1为潜变量,x1,x2为显变量。那么将上面的模型展现形式综合在一起就是下面的格式:
model <- '
# 定义潜变量
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + a*y2 + b*y3 + c*y4
dem65 =~ y5 + a*y6 + b*y7 + c*y8
# 回归方程
dem60 ~ ind60
dem65 ~ ind60 + dem60
# 残差相关性
y1 ~~ y5
y2 ~~ y4 + y6
y3 ~~ y7
y4 ~~ y8
y6 ~~ y8
# 截距项
1
f1 ~ 1
'
当然,如果模型的方程式太多也可以保存在txt的文件中,并且后缀名改为.lav。接下来就借助包内的函数readLines()读取模型的内容。另外,我们也可以直接以表格形式查看模型的构成,需要用到函数lavaanify(),我们直接看实例:
model <- '
#two latent variables with fixed loadings
f1 =~ 1*y1a + 1*y1b + 1*y1c
f2 =~ 1*y2a + 1*y2b + 1*y2c
#intercepts constrained to be equal
#using the default names
y1a ~ 1
y1b ~ equal("y1a~1") * 1
y1c ~ equal("y1a~1") * 1
#intercepts constrained to be equal
#using a custom label
y2a ~ int2*1
y2b ~ int2*1
y2c ~ int2*1
'
接下来我们看下包中自带的三个主要模型:
1. 验证性因子模型(Confirmatory FactorAnalysis (CFA) model)。我们首先看下其主要函数:
我们一般会把默认的参数作为最好的设置,所以我们不一一展开参数解释了,我们直接看实例:
首先我们看下数据源:
## The famous Holzinger and Swineford(1939) example
<- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9 '
fit <- cfa(HS.model, data =HolzingerSwineford1939)
fit.measures = TRUE)
检查模型的拟合度,使用χ2/df(小于3.0)、AGFI (>0.90)、CFI(>0.90)、 RMSEA (<0.07)与SRMR(<0.08)等指标。
接下来进行信效度的分析需要用到函数reliability()。所谓信度(Reliability)即可靠性,是指采用同一方法对同一对象进行调查时,问卷调查结果的稳定性和一致性,即测量工具(问卷或量表)能否稳定地测量所测的事物或变量;效度(Validity)即有效性,它是指测量工具或手段能够准确测出所需测量的事物的程度。信度与效度二者的关系:信度低,效度不可能高;效度高,信度也必然高。
其中alpha>0.8,内部一致性极好;0.6-0.8较好;<0.6一致性较差。
在上述结果中omega值为CR(信度系数),其中alpha值以及AVE值都已列举出来,可以看出AVE并未都大于0.5,alpha值并未都在0.8以上,CR值也并未都在0.8以上,说明模型信效度不是很良好。
最后我们看下我们构建的模型的结构:
semPaths(fit)
2. 结构方程模型(StructuralEquation Model (SEM))。其参数与功能和cfa基本一致,我们直接看下实例:
model <- '
#latent variable definitions
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + a*y2 + b*y3 + c*y4
dem65 =~ y5 + a*y6 + b*y7 + c*y8
#regressions
dem60 ~ ind60
dem65 ~ ind60 + dem60
#residual correlations
y1 ~~ y5
y2 ~~ y4 + y6
y3 ~~ y7
y4 ~~ y8
y6 ~~ y8
'
fit <- sem(model, data =PoliticalDemocracy)
fit.measures = TRUE)
semPaths(fit)
3. 增长曲线模型(FitGrowth Curve Models)。一般来说,事物总是经过发生、发展、成熟三个阶段,而每一个阶段的发展速度各不相同。通常在发生阶段,变化速度较为缓慢;在发展阶段,变化速度加快;在成熟阶段,变化速度又趋缓慢,按上述三个阶段发展规律得到的变化曲线称为生长曲线。我们直接看实例:
## linear growth model with a time-varyingcovariate
<- '
#intercept and slope with fixed coefficients
i=~ 1*t1 + 1*t2 + 1*t3 + 1*t4
s=~ 0*t1 + 1*t2 + 2*t3 + 3*t4
#regressions
x1 + x2
x1 + x2
#time-varying covariates
t1 ~ c1
t2 ~ c2
t3 ~ c3
t4 ~ c4
'
fit <- growth(model.syntax, data =Demo.growth)
fit.measures = TRUE)
semPaths(fit)
以上个模型都是基于lavaan函数的相关参数进行的默认运行模型,如果想更加细致的修改各参数,需要直接调用lavaan函数。同时,此包还设置了一些对summary结构化的输出函数:
1. parameterEstimates 直接提取各方程式的评估结果。我们直接看下实例:
parameterEstimates(fit)
2. parTable() 模型参数列表,输出为data.frame格式。
parTable(fit)
当我们在模型中加入等式约束时,我们可以评估等式约束的可靠性需要用到函数lavTestScore。我们直接看下包的实例:
<- '
visual =~ x1 + b1*x2 + x3
textual =~ x4 + b2*x5 + x6
speed =~ x7 + b3*x8 + x9
b1 == b2
b2 == b3
'
fit <- cfa(HS.model,data=HolzingerSwineford1939)
# test 1: release both two equality constraints
cumulative = TRUE)
与此同时此包还添加了模型之间的对比函数lavTestLRT,可以对比两个模型是否有差异。实例:
<- '
visual =~ x1 + b1*x2 + x3
textual =~ x4 + b2*x5 + x6
speed =~ x7 + b3*x8 + x9
'
fit1 <- cfa(HS.model, data =HolzingerSwineford1939)
fit0 <- cfa(HS.model, data =HolzingerSwineford1939,
orthogonal = TRUE)
fit0)
至此就是目前此包的重点内容了,欢迎大家互相学习交流!