A first look at structured equation models using the Lavaan package
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
# A FIRST LOOK AT LAVAAN | |
## By Grace Charles, presented at Davis R Users' Group on May 15, 2015 | |
## adapted from Jim Grace's SEM workshop and Lavaan tutorials | |
# Set your working directory | |
setwd("~/Desktop/DAVIS/sem workshop") | |
###Load Libraries | |
library(lavaan) | |
library(qgraph) | |
#Built in dataset | |
data(PoliticalDemocracy) | |
PD<-PoliticalDemocracy | |
head(PD) | |
## here is an example model from one of Lavaan's build-in datasets | |
## the measurement model equations are "latent" and represented by =~ | |
## regressions are indicated by ~ | |
## residual correlations (in this case because they represent different years of the same measurement) are represented by ~~ | |
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 | |
y2~~ y6 | |
y3 ~~ y7 | |
y4 ~~ y8 | |
y6 ~~ y8 | |
' | |
#fit your SEM | |
fit <- sem(model, data = PD) | |
#summarize results | |
summary(fit, standardized = TRUE, rsq = T) | |
##plot results using semPaths function in qgraph | |
semPaths(fit, "std", edge.label.cex = 0.5, curvePivot = TRUE, layout = "tree") | |
##check to see if you missed anything. High mi values suggest that there is a path that you missed. | |
modindices(fit) | |
## looks good | |
##can also look at variance tables | |
vartable(fit) | |
## sometimes you get warnings about the scale of your variables | |
#Warning message: | |
# In getDataFull(data = data, group = group, group.label = group.label, : | |
# lavaan WARNING: some observed variances are (at least) a factor 100 times larger than others; please rescale | |
# in this case, all you have to do to make this error go away is rescale variables | |
#model comparison | |
#you can compare alternative pathway models using AIC, BIC, etc: | |
#create second alternative model | |
names(PD) | |
model2 <- ' | |
# measurement model | |
ind60 =~ x1 + x2 + x3 | |
dem60 =~ y1 + y2 + y3 + y4 | |
dem65 =~ y5 + y6 + y7 + y8 | |
# regressions | |
dem60 ~ ind60 | |
dem65 ~ dem60 | |
#took out ind60 from regression | |
# residual correlations | |
y1 ~~ y5 | |
y2 ~~ y4 + y6 | |
y3 ~~ y7 | |
y4 ~~ y8 | |
y6 ~~ y8 | |
' | |
fit2 <- sem(model2, data = PD) | |
summary(fit2) | |
AIC(fit, fit2) | |
## what about nonlinear data? | |
# Set your working directory | |
("~/Desktop/sem workshop") | |
#commands in bold | |
# Load data and name file ?k.dat? | |
k.dat<-read.csv("./Keeley_rawdata_select4.csv") | |
# Examine contents of keeley data file | |
names(k.dat) | |
head(k.dat) | |
# Write lavaan code for this single equation model | |
mod <- ' | |
rich ~ cover | |
cover ~ firesev | |
' | |
k.dat$cov2<-k.dat$cover^2 | |
mod2<- ' | |
rich ~ cover + cov2 | |
cover ~ firesev | |
cover ~~ cov2 | |
cov2 ~~ firesev | |
' | |
# Fit the model (i.e. est. parameters) | |
mod1.fit <- sem(mod, data=k.dat) | |
mod2.fit<- sem(mod2, data=k.dat,fixed.x=FALSE) | |
#need to rescale data. | |
vartable(mod1.fit) | |
k.dat$rich<-k.dat$rich/100 | |
# Output a summary of the computed results - summary of mod2 suggests that both cover and cover squared can impact | |
summary(mod1.fit, rsq=T) # rsq=T means output the r-sqr | |
summary(mod2.fit, rsq=T) | |
semPaths(mod1.fit, "std", edge.label.cex = 0.5, curvePivot = TRUE, layout = "tree") | |
semPaths(mod2.fit, "std", edge.label.cex = 0.5, curvePivot = TRUE, layout = "tree") | |
Thanks, but where is the "Keeley_rawdata_select4.csv" file. Could we use this file too please?
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
It happens to me too when I try to load the semPlot package. Is the problem solved?