The lavaan package is an excellent package for structural equation models, and the DiagrammeR package is an excellent package for producing nice looking graph diagrams. As of right now, the lavaan package has no built in plotting functions for models, and the available options from external packages don’t look as nice and aren’t as easy to use as DiagrammeR, in my opinion. Of course, you can use DiagrammeR to build path diagrams for your models, but it requires you to build the diagram specification manually. This package exists to streamline that process, allowing you to plot your lavaan models directly, without having to translate them into the DOT language specification that DiagrammeR uses.
The package is very straightforward to use, simply call the
lavaanPlot function with your lavaan model, adding whatever
graph, node and edge attributes you want as a named list (graph
attributes are specified as a standard default value that shows you what
the other attribute lists should look like). For your reference, the
available attributes can be found here:
https://rich-iannone.github.io/DiagrammeR/reference/node_aes.html https://rich-iannone.github.io/DiagrammeR/reference/edge_aes.html
Here’s a quick example using the mtcars data set.
First fit your lavaan model. The package supports plotting lavaan regression relationships and latent variable - indicator relationships.
library(lavaan)
library(lavaanPlot)
model <- 'mpg ~ cyl + disp + hp
qsec ~ disp + hp + wt'
fit <- sem(model, data = mtcars)
summary(fit)## lavaan 0.6.15 ended normally after 32 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 9
##
## Number of observations 32
##
## Model Test User Model:
##
## Test statistic 18.266
## Degrees of freedom 2
## P-value (Chi-square) 0.000
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## mpg ~
## cyl -0.987 0.738 -1.337 0.181
## disp -0.021 0.010 -2.178 0.029
## hp -0.017 0.014 -1.218 0.223
## qsec ~
## disp -0.008 0.004 -2.122 0.034
## hp -0.023 0.004 -5.229 0.000
## wt 1.695 0.398 4.256 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .mpg ~~
## .qsec 0.447 0.511 0.874 0.382
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .mpg 8.194 2.049 4.000 0.000
## .qsec 0.996 0.249 4.000 0.000Then using that model fit object, simply call the
lavaanPlot function, specifying your desired graph
parameters.
You can also specify different variable labels using the a named list
and the labels argument.
labels <- list(mpg = "Miles Per Gallon", cyl = "Cylinders", disp = "Displacement", hp = "Horsepower", qsec = "Speed", wt = "Weight")
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = FALSE)An example showing latent variable relationships:
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
'
fit <- cfa(HS.model, data=HolzingerSwineford1939)
lavaanPlot(model = fit, edge_options = list(color = "grey"))You can label the plot edges with the coefficient values using
coefs = TRUE.
model <- 'mpg ~ cyl + disp + hp
qsec ~ disp + hp + wt'
fit <- sem(model, data = mtcars)
summary(fit)## lavaan 0.6.15 ended normally after 32 iterations
##
## Estimator ML
## Optimization method NLMINB
## Number of model parameters 9
##
## Number of observations 32
##
## Model Test User Model:
##
## Test statistic 18.266
## Degrees of freedom 2
## P-value (Chi-square) 0.000
##
## Parameter Estimates:
##
## Standard errors Standard
## Information Expected
## Information saturated (h1) model Structured
##
## Regressions:
## Estimate Std.Err z-value P(>|z|)
## mpg ~
## cyl -0.987 0.738 -1.337 0.181
## disp -0.021 0.010 -2.178 0.029
## hp -0.017 0.014 -1.218 0.223
## qsec ~
## disp -0.008 0.004 -2.122 0.034
## hp -0.023 0.004 -5.229 0.000
## wt 1.695 0.398 4.256 0.000
##
## Covariances:
## Estimate Std.Err z-value P(>|z|)
## .mpg ~~
## .qsec 0.447 0.511 0.874 0.382
##
## Variances:
## Estimate Std.Err z-value P(>|z|)
## .mpg 8.194 2.049 4.000 0.000
## .qsec 0.996 0.249 4.000 0.000lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE)By default it will show all paths, but you can also specify whatever
significance level you want using the sig argument to only
show significant paths.
# significant standardized paths only
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, sig = .05)You can also show standardized paths only with the stand
argument.
# All paths unstandardized
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, stand = TRUE)Same works for latent variable loadings
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
'
fit <- cfa(HS.model, data=HolzingerSwineford1939)
labels = list(visual = "Visual Ability", textual = "Textual Ability", speed = "Speed Ability")
# Show coefs
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE)# Significant paths
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, sig = .05)# All paths standardized
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, stand = TRUE)You can also include double-sided edges to represent model covariances if you want:
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9
'
fit <- cfa(HS.model, data=HolzingerSwineford1939)
labels = list(visual = "Visual Ability", textual = "Textual Ability", speed = "Speed Ability")
# significant standardized paths only
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE)You can include significance stars as well using the
stars option. You can choose stars for regression paths,
latent paths, or covariances. Specify which of the 3 you want
(“regress”, “latent”, “covs”).
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = "covs")lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = "latent")You can change the number of decimal places in the coefficient value
labels as well with the digits option.
rankdirYou can also do additional options to modify the layout of your plots in certain ways.
The graph option rankdir allows you to change the
overall orientation of the plot between TB (top-bottom; default), BT
(bottom-top), RL (right-left) and LR (left-right):
lavaanPlot(model = fit, labels = labels, graph_options = list(rankdir = "LR"), node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)layoutThe graph option layout allows you to change the layout
engine used to generate the plot between dot (the default), neato,
circo, twopi, plus some others, see: https://graphviz.org/docs/layouts/. Your mileage may
vary as to whether any of these work for a given graph.
lavaanPlot(model = fit, labels = labels, graph_options = list(layout = "neato"), node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)