Skip to contents

Introduction

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.

Package example

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/graphviz_and_mermaid.html#node-attributes https://rich-iannone.github.io/DiagrammeR/graphviz_and_mermaid.html#edge-attributes

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.

## This is lavaan 0.6-17
## lavaan is FREE software! Please report any bugs.
library(lavaanPlot)

model <- 'mpg ~ cyl + disp + hp
          qsec ~ disp + hp + wt'

fit <- sem(model, data = mtcars)
summary(fit)
## lavaan 0.6.17 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.000

Then using that model fit object, simply call the lavaanPlot function, specifying your desired graph parameters.

lavaanPlot(model = fit, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = FALSE)
plot cyl cyl mpg mpg cyl->mpg disp disp disp->mpg qsec qsec disp->qsec hp hp hp->mpg hp->qsec wt wt wt->qsec

Plot customization options

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)
plot cyl Cylinders mpg Miles Per Gallon cyl->mpg disp Displacement disp->mpg qsec Speed disp->qsec hp Horsepower hp->mpg hp->qsec wt Weight wt->qsec

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"))
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual visual visual->x1 visual->x2 visual->x3 textual textual textual->x4 textual->x5 textual->x6 speed speed speed->x7 speed->x8 speed->x9

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.17 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.000
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE)
plot cyl Cylinders mpg Miles Per Gallon cyl->mpg -0.99 disp Displacement disp->mpg -0.02 qsec Speed disp->qsec -0.01 hp Horsepower hp->mpg -0.02 hp->qsec -0.02 wt Weight wt->qsec 1.69

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)
plot cyl Cylinders mpg Miles Per Gallon cyl->mpg disp Displacement disp->mpg -0.02 qsec Speed disp->qsec -0.01 hp Horsepower hp->mpg hp->qsec -0.02 wt Weight wt->qsec 1.69

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)
plot cyl Cylinders mpg Miles Per Gallon cyl->mpg -0.29 disp Displacement disp->mpg -0.44 qsec Speed disp->qsec -0.56 hp Horsepower hp->mpg -0.19 hp->qsec -0.85 wt Weight wt->qsec 0.91

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)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.55 visual->x3 0.73 textual Textual Ability textual->x4 1 textual->x5 1.11 textual->x6 0.93 speed Speed Ability speed->x7 1 speed->x8 1.18 speed->x9 1.08
# Significant paths
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, sig = .05)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.55 visual->x3 0.73 textual Textual Ability textual->x4 1 textual->x5 1.11 textual->x6 0.93 speed Speed Ability speed->x7 1 speed->x8 1.18 speed->x9 1.08
# All paths standardized
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, stand = TRUE)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 0.77 visual->x2 0.42 visual->x3 0.58 textual Textual Ability textual->x4 0.85 textual->x5 0.86 textual->x6 0.84 speed Speed Ability speed->x7 0.57 speed->x8 0.72 speed->x9 0.67

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)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.55 visual->x3 0.73 textual Textual Ability textual->x4 1 textual->x5 1.11 textual->x6 0.93 textual->visual 0.41 speed Speed Ability speed->x7 1 speed->x8 1.18 speed->x9 1.08 speed->visual 0.26 speed->textual 0.17

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")
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.55 visual->x3 0.73 textual Textual Ability textual->x4 1 textual->x5 1.11 textual->x6 0.93 textual->visual 0.41*** speed Speed Ability speed->x7 1 speed->x8 1.18 speed->x9 1.08 speed->visual 0.26*** speed->textual 0.17***
lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = "latent")
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1*** visual->x2 0.55*** visual->x3 0.73*** textual Textual Ability textual->x4 1*** textual->x5 1.11*** textual->x6 0.93*** textual->visual 0.41 speed Speed Ability speed->x7 1*** speed->x8 1.18*** speed->x9 1.08*** speed->visual 0.26 speed->textual 0.17

You can change the number of decimal places in the coefficient value labels as well with the digits option.

lavaanPlot(model = fit, labels = labels, node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2

Graph Options

Orientations with rankdir

You 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)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2
lavaanPlot(model = fit, labels = labels, graph_options = list(rankdir = "RL"), node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2
lavaanPlot(model = fit, labels = labels, graph_options = list(rankdir = "BT"), node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2

Layouts with layout

The 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)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2
lavaanPlot(model = fit, labels = labels, graph_options = list(layout = "circo"), node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2
lavaanPlot(model = fit, labels = labels, graph_options = list(layout = "twopi"), node_options = list(shape = "box", fontname = "Helvetica"), edge_options = list(color = "grey"), coefs = TRUE, covs = TRUE, stars = TRUE, digits = 1)
plot x1 x1 x2 x2 x3 x3 x4 x4 x5 x5 x6 x6 x7 x7 x8 x8 x9 x9 visual Visual Ability visual->x1 1 visual->x2 0.6 visual->x3 0.7 textual Textual Ability textual->x4 1 textual->x5 1.1 textual->x6 0.9 textual->visual 0.4 speed Speed Ability speed->x7 1 speed->x8 1.2 speed->x9 1.1 speed->visual 0.3 speed->textual 0.2