Generalized Network Psychometrics: Combining Network and
Latent Variable Models

arXiv:1605.09288v4 [math.ST] 11 Sep 2017

Sacha Epskamp,1 Mijke T. Rhemtulla,2 and Denny Borsboom1
1. University of Amsterdam: Department of Psychological Methods
2. University of California, Davis: Department of Psychology

Abstract
We introduce the network model as a formal psychometric model, conceptualizing the covariance between psychometric indicators as resulting from
pairwise interactions between observable variables in a network structure.
This contrasts with standard psychometric models, in which the covariance
between test items arises from the influence of one or more common latent
variables. Here, we present two generalizations of the network model that
encompass latent variable structures, establishing network modeling as parts
of the more general framework of Structural Equation Modeling (SEM). In
the first generalization, we model the covariance structure of latent variables
as a network. We term this framework Latent Network Modeling (LNM) and
show that, with LNM, a unique structure of conditional independence relationships between latent variables can be obtained in an explorative manner. In the second generalization, the residual variance-covariance structure
of indicators is modeled as a network. We term this generalization Residual
Network Modeling (RNM) and show that, within this framework, identifiable
models can be obtained in which local independence is structurally violated.
These generalizations allow for a general modeling framework that can be
used to fit, and compare, SEM models, network models, and the RNM and
LNM generalizations. This methodology has been implemented in the freeto-use software package lvnet, which contains confirmatory model testing as
well as two exploratory search algorithms: stepwise search algorithms for
low-dimensional datasets and penalized maximum likelihood estimation for
larger datasets. We show in simulation studies that these search algorithms
performs adequately in identifying the structure of the relevant residual or
latent networks. We further demonstrate the utility of these generalizations
in an empirical example on a personality inventory dataset.

Manuscript accepted for publication in Psychometrika.

GENERALIZED NETWORK PSYCHOMETRICS

2

Introduction
Recent years have seen an emergence of network modeling in psychometrics (Borsboom, 2008; Schmittmann et al., 2013), with applications in clinical psychology (e.g., van
Borkulo et al. 2015; McNally et al. 2015; Fried et al. 2015), psychiatry (e.g., Isvoranu et al.
2016b,a), health sciences (e.g., Kossakowski et al. 2015), social psychology (e.g., Dalege
et al. 2016; Cramer et al. 2012), and other fields (see for a review of recent literature Fried
and van Borkulo 2016). This line of literature stems from the network perspective of psychology (Cramer et al., 2010), in which noticeable macroscopic behavior—the co-occurrence
of aspects of psychology such as cognitive abilities, psychopathological symptoms, or a set
of behaviors—is hypothesized to not be due to the influence of unobserved common causes,
such as general intelligence, psychopathological disorders, or personality traits, but rather
to emergent behavior in a network of interacting psychological, sociological, and biological
components. Network models are used to gain insight into this potentially high-dimensional
interplay (Epskamp et al., 2016). In practice, network models can be used as a sparse representation of the joint distribution of observed indicators, and as such these models show
great promise in psychometrics by providing a perspective that complements latent variable modeling. Network modeling highlights variance that is unique to pairs of variables,
whereas latent variable modeling focuses on variance that is shared across all variables
(Costantini et al., 2015). As a result, network modeling and latent variable modeling can
complement—rather than exclude—one-another.
In this paper, we introduce the reader to this field of network psychometrics (Epskamp
et al., ress) and formalize the network model for multivariate normal data, the Gaussian
Graphical Model (GGM; Lauritzen 1996), as a formal psychometric model. We contrast
the GGM to the Structural Equation Model (SEM; Wright 1921; Kaplan 2000) and show
that the GGM can be seen as another way to approach modeling covariance structures as is
typically done in psychometrics. In particular, rather than modeling the covariance matrix,
the GGM models the inverse of a covariance matrix. The GGM and SEM are thus very
closely related: every GGM model and every SEM model imply a constrained covariance
structure. We make use of this relationship to show that, through a reparameterization of
the SEM model, the GGM model can be obtained in two different ways: first, as a network
structure that relates a number of latent variables to each other, and second, as a network
between residuals that remain given a fitted latent variable model. As such, the GGM
can be modeled and estimated in SEM, which allows for network modeling of psychometric
data to be carried out in a framework familiar to psychometricians and methodologists. In
addition, this allows for one to assess the fit of a GGM, compare GGMs to one-another and
compare a GGM to a SEM model.
However, the combination of GGM and SEM allows for more than fitting network
models. As we will show, the strength of one framework can help overcome shortcomings
of the other framework. In particular, SEM falls short in that exploratory estimation is
complicated and there is a strong reliance on local independence, whereas the GGM falls
short in that it assumes no latent variables. In this paper, we introduce network models
for latent covariances and for residual covariances as two distinct generalized frameworks of
both the SEM and GGM. The first framework, Latent Network Modeling (LNM), formulates a network among latent variables. This framework allows researchers to exploratively

GENERALIZED NETWORK PSYCHOMETRICS

3

estimate conditional independence relationships between latent variables through model
search algorithms; this estimation is difficult in the SEM framework due to the presence
of equivalent models (MacCallum et al., 1993). The second framework, which we denote
Residual Network Modeling (RNM), formulates a network structure on the residuals of a
SEM model. With this framework, researchers can circumvent critical assumptions of both
SEM and the GGM: SEM typically relies on the assumption of local independence, whereas
network modeling typically relies on the assumption that the covariance structure among a
set of the items is not due to latent variables at all. The RNM framework allows researchers
to estimate SEM models without the assumption of local independence (all residuals can
be correlated, albeit due to a constrained structure on the inverse residual covariance matrix) as well as to estimate a network structure, while taking into account the fact that the
covariance between items may be partly due to latent factors.
While the powerful combination of SEM and GGM allows for confirmative testing of
network structures both with and without latent variables, we recognize that few researchers
have yet formulated strict confirmatory hypotheses in the relatively new field of network
psychometrics. Often, researchers are more interested in exploratively searching a plausible
network structure. To this end, we present two exploratory search algorithms. The first is a
step-wise model search algorithm that adds and removes edges of a network as long as fit is
improved, and the second uses penalized maximum likelihood estimation (Tibshirani, 1996)
to estimate a sparse model. We evaluate the performance of these search methods in four
simulation studies. Finally, the proposed methods have been implemented in a free-to-use
R package, lvnet, which we illustrate in an empirical example on personality inventory
items (Revelle, 2010).
Modeling Multivariate Gaussian Data
Let y be the response vector of a random subject on P items1 . We assume y is
centered and follows a multivariate Gaussian density:
y ∼ NP (00, Σ ) ,

Σ.
In which Σ is a P × P variance–covariance matrix, estimated by some model-implied Σ̂
Σ is often done through some form of maximum likelihood estimation. If we
Estimating Σ̂
measure N independent samples of y we can formulate the N × P matrix Y containing
realization y >
i as its ith row. Let S represent the sample variance–covariance matrix of Y :
S=

1
Y >Y .
N −1

In maximum likelihood estimation, we use S to compute and minimize −2 times the logΣ (Lawley, 1940; Jöreskog, 1967; Jacobucci et al., 2016):
likelihood function to find Σ̂
h

 



Σ + Trace S Σ̂
Σ
min log det Σ̂
Σ
Σ̂

1

−1 

 

i

S −P .
− log det Ŝ

(1)

Throughout this paper, vectors will be represented with lowercase boldfaced letters and matrices will be
denoted by capital boldfaced letters. Roman letters will be used to denote observed variables and parameters
(such as the number of nodes) and Greek letters will be used to denote latent variables and parameters that
need to be estimated. The subscript i will be used to denote the realized response vector of subject i and
omission of this subscript will be used to denote the response of a random subject.

GENERALIZED NETWORK PSYCHOMETRICS

4

Σ should be estimated as closely as possible to S and perfect
To optimize this expression, Σ̂
Σ = S . A properly identified model with the same number of parameters
fit is obtained if Σ̂
Σ as there are unique elements in S (P (P + 1)/2 parameters) will lead to
(K) used to form Σ̂
Σ = S and therefore a saturated model. The goal of modeling multivariate Gaussian data
Σ̂
Σ with positive degrees of freedom, K < P (P + 1)/2, in which
is to obtain some model for Σ̂
Σ resembles S closely.
Σ̂
Structural Equation Modeling
In Confirmatory Factor Analysis (CFA), Y is typically assumed to be a causal linear
effect of a set of M centered latent variables, η , and independent residuals or error, ε :
y = Λη + ε .
Here, Λ represents a P × M matrix of factor loadings. This model implies the following
Σ:
model for Σ̂
Σ = Λ ΨΛ > + Θ ,
Σ̂
(2)
in which Ψ = Var (ηη ) and Θ = Var (εε). In Structural Equation Modeling (SEM), Var (ηη )
can further be modeled by adding structural linear relations between the latent variables2 :
η = Bη + ζ ,
in which ζ is a vector of residuals and B is an M × M matrix of regression coefficients.
Σ can be more extensively modeled as:
Now, Σ̂
Σ = Λ (II − B )−1 Ψ (II − B )−1> Λ > + Θ ,
Σ̂

(3)

in which now Ψ = Var (ζζ ). This framework can be used to model direct causal effects
between observed variables by setting Λ = I and Θ = O , which is often called path analysis
(Wright, 1934).
Σ and S , a P × P matrix; if Θ is fully estimated—contains no
The Θ matrix is, like Σ̂
restricted elements—then Θ alone constitutes a saturated model. Therefore, to make either
(2) or (3) identifiable, Θ must be strongly restricted. Typically, Θ is set to be diagonal,
a restriction often termed local independence (Lord et al., 1968; Holland and Rosenbaum,
1986) because indicators are independent of each other after conditioning on the set of latent
variables. To improve fit, select off-diagonal elements of Θ can be estimated, but systematic
violations of local independence—many nonzero elements in Θ —are not possible as that will
quickly make (2) and (3) saturated or even over-identified. More precisely, Θ can not be
fully-populated—some elements of Θ must be set to equal zero—when latent variables are
used. An element of Θ being fixed to zero indicates that two variables are locally independent
after conditioning on the set of latent variables. As such, local independence is a critical
assumption in both CFA and SEM; if local independence is systematically violated, CFA
and SEM will never result in correct models.
The assumption of local independence has led to critiques of the factor model and
its usage in psychology; local independence appears to be frequently violated due to direct
2

We make use here of the convenient all-y notation and do not distinguish between exogenous and
endogenous latent variables (Hayduk, 1987).

GENERALIZED NETWORK PSYCHOMETRICS

5

causal effects, semantic overlap, or reciprocal interactions between putative indicators of
a latent variable (Borsboom, 2008; Cramer et al., 2010; Borsboom et al., 2011; Cramer
et al., 2012; Schmittmann et al., 2013). In psychopathology research, local independence of
symptoms given a person’s level of a latent mental disorder has been questioned (Borsboom
and Cramer, 2013). For example, three problems associated with depression are “fatigue”,
“concentration problems” and “rumination”. It is plausible that a person who suffers from
fatigue will also concentrate more poorly, as a direct result of being fatigued and regardless
of his or her level of depression. Similarly, rumination might lead to poor concentration.
In another example, Kossakowski et al. (2015) describe the often-used SF-36 questionnaire
(Ware Jr and Sherbourne, 1992) designed to measure health related quality of life. The
SF-36 contains items such as “can you walk for more than one kilometer” and “can you
walk a few hundred meters”. Clearly, these items can never be locally independent after
conditioning on any latent trait, as one item (the ability to walk a few hundred meters)
is a prerequisite for the other (walking more than a kilometer). In typical applications,
the excessive covariance between items of this type is typically left unmodeled, and treated
instead by combining items into a subscale or total score that is subsequently subjected to
factor analysis; of course, however, this is tantamount to ignoring the relevant psychometric
problem rather than solving it.
Given the many theoretically expected violations of local independence in psychometric applications, many elements of Θ in both (2) and (3) should ordinarily be freely
estimated. Especially when violations of local independence are expected to be due to
causal effects of partial overlap, residual correlations should not be constrained to zero;
in addition, a chain of causal relationships between indicators can lead to all residuals to
become correlated. Thus, even when latent factors cause much of the covariation between
measured items, fitting a latent variable model that involves local independence may not
fully account for correlation structure between measured items. Of course, in practice,
many psychometricians are aware of this problem, which is typically addressed by freeing
up correlations between residuals to improve model fit. However, this is usually done in
an ad-hoc fashion, on the basis of inspection of modification indices and freeing up error
covariances one by one, which is post hoc, suboptimal, and involves an uncontrolled journey
through the model space. As a result, it is often difficult to impossible to tell how exactly
authors arrived at their final reported models. As we will show later in this paper, this
process can be optimized and systematized using network models to connect residuals on
top of a latent variable structure.
Network Modeling
Recent authors have suggested that the potential presence of causal relationships between measured variables may allow the explanation of the covariance structure without the
need to invoke any latent variables (Borsboom, 2008; Cramer et al., 2010; Borsboom et al.,
2011; Schmittmann et al., 2013). The interactions between indicators can instead be modeled as a network, in which indicators are represented as nodes that are connected by edges
representing pairwise interactions. Such interactions indicate the presence of covariances
that cannot be explained by any other variable in the model and can represent—possibly
reciprocal—causal relationships. Estimating a network structure on psychometric data is
termed network psychometrics (Epskamp et al., ress). Such a network of interacting compo-

GENERALIZED NETWORK PSYCHOMETRICS

6

Y2

Y1

Y3

Figure 1 . Example of a pairwise Markov Random Field model. Edges in this model indicate
pairwise interactions, and are drawn using undirected edges to distinguish from (bidirectional) covariances. Rather than a model for marginal associations (such as a network
indicating covariances), this is a model for conditional associations. The network above
encodes that Y1 and Y3 are independent after conditioning on Y2 . Such a model allows all
three variables to correlate while retaining one degree of freedom (the model only has two
parameters)
nents can generate data that fit factor models well, as is commonly the case in psychology.
van der Maas et al. (2006) showed that the positive manifold of intelligence—which is commonly explained with the general factor for intelligence, g—can emerge from a network
of mutually benefiting cognitive abilities. Borsboom et al. (2011) showed that a network
of psychopathological symptoms, in which disorders are modeled as clusters of symptoms,
could explain comorbidity between disorders. Furthermore, Epskamp et al. (ress) showed
that the Ising model for ferromagnetism (Ising, 1925), which models magnetism as a network
of particles, is equivalent to multidimensional item response theory (Reckase, 2009).
In network psychometrics, psychometric data are modeled through directed or undirected networks. Directed networks are equivalent to path analysis models. For modeling
undirected networks, pairwise Markov Random Fields (Lauritzen, 1996; Murphy, 2012) are
used. In these models, each variable is represented by a node, and nodes are connected
by a set of edges. If two nodes, yj and yk , are not connected by an edge, then this means
they are independent after conditioning on the set of all other nodes, y −(j,k) . Whenever two
nodes cannot be rendered independent conditional on the other nodes in the system, they
are said to feature in a pairwise interaction, which is represented by an undirected edge—an
edge with no arrows—to contrast such an effect from covariances typically represented in
the SEM literature with bidirectional edges. Figure 1 represents such a network model,
in which nodes y1 and y3 are independent after conditioning on node y2 . Such a model
can readily arise from direct interactions between the nodes. For example, this conditional

GENERALIZED NETWORK PSYCHOMETRICS

7

independence structure would emerge if y2 is a common cause of y1 and y3 , or if y2 is the
mediator in a causal path between y1 and y3 . In general, it is important to note that pairwise interactions are not mere correlations; two variables may be strongly correlated but
unconnected (e.g., when both are caused by another variable in the system) and they may
be uncorrelated but strongly connected in the network (e.g., when they have a common
effect in the system). For instance, in the present example the model does not indicate that
y1 and y3 are uncorrelated, but merely indicates that any correlation between y1 and y3 is
due to their mutual interaction with y2 ; a network model in which either directly or indirectly connected paths exist between all pairs of nodes typically implies a fully populated
(no zero elements) variance–covariance matrix.
In the case of multivariate Gaussian data, this model is termed the Gaussian Graphical
Model (GGM; Lauritzen 1996). Here, the partial correlation coefficient is sufficient to test
the degree of conditional independence of two variables after conditioning on all other
variables; if the partial correlation coefficient is zero, there is conditional independence and
hence no edge in the network. As such, partial correlation coefficients can directly be used
in the network as edge weights; the strength of connection between two nodes3 . Such a
network is typically encoded in a symmetrical and real valued p × p weight matrix, Ω , in
which element ωjk represents the edge weight between node j and node k:




Cor yj , yk | y −(j,k) = ωjk = ωkj .
The partial correlation coefficients can be directly obtained from the inverse of variance–
Σ, also termed the precision matrix K̂
K (Lauritzen, 1996):
covariance matrix Σ̂




Cor yj , yk | y −(j,k) = − √

κjk
√ .
κkk κjj

Thus, element κjk of the precision matrix is proportional to to the partial correlation coefficient of variables yj and yk after conditioning on all other variables. Since this process
simply involves standardizing the precision matrix, we propose the following model4 :
Σ = K̂
K
Σ̂

−1

= ∆ (II − Ω )−1 ∆ ,

(4)

−1

in which ∆ is a diagonal matrix with δjj = κjj 2 and Ω has zeroes on the diagonal. This model
allows for confirmative testing of the GGM structures on psychometric data. Furthermore,
the model can be compared to a saturated model (fully populated off-diagonal values of Ω )
Ω = O ), allowing one to obtain χ2 fit statistics as well as fit
and the independence model (Ω
indices such as the RMSEA (Browne and Cudeck, 1992) and CFI (Bentler, 1990). Such
methods of assessing model fit have not yet been used in network psychometrics.
Similar to CFA and SEM, the GGM relies on a critical assumption; namely, that
covariances between observed variables are not caused by any latent or unobserved variable.
If we estimate a GGM in a case where, in fact, a latent factor model was the true data
3

A saturated GGM is also called a partial correlation network because it contains the sample partial
correlation coefficients as edge weights.
4
To our knowledge, the GGM has not yet been framed in this form. We chose this form because it allows
for clear modeling and interpretation of the network parameters.

GENERALIZED NETWORK PSYCHOMETRICS

8

generating structure, then generally we would expect the GGM to be saturated—i.e., there
would be no missing edges in the GGM (Chandrasekaran et al., 2012). A missing edge in
the GGM indicates the presence of conditional independence between two indicators given
all other indicators; we do not expect indicators to become independent given subsets of
other indicators (see also Ellis and Junker 1997; Holland and Rosenbaum 1986). Again, this
critical assumption might not be plausible. While variables such as “Am indifferent to the
feelings of others” and “Inquire about others’ well-being” quite probably interact with each
other, it might be far-fetched to assume that no unobserved variable, such as a personality
trait, in part also causes some of the variance in responses on these items.
Generalizing Factor Analysis and Network Modeling
We propose two generalizations of both SEM and the GGM that both allow the
modeling of network structures in SEM. In the first generalization, we adopt the CFA5
decomposition in (2) and model the variance–covariance matrix of latent variables as a
GGM:
Ψ = ∆Ψ (II − ΩΨ )−1 ∆Ψ .
This framework can be seen as modeling conditional independencies between latent variables
not by directed effects (as in SEM) but as an undirected network. As such, we term this
framework latent network modeling (LNM).
In the second generalization, we adopt the SEM decomposition of the variance–
covariance matrix in (3) and allow the residual variance–covariance matrix Θ to be modeled
as a GGM:
Θ = ∆Θ (II − ΩΘ )−1 ∆Θ .
Because this framework conceptualizes associations between residuals as pairwise interactions, rather than correlations, we term this framework Residual Network Modeling (RNM).
Using this framework allows—as will be described below—for a powerful way of fitting a
confirmatory factor structure even though local independence is systematically violated and
all residuals are correlated.
Figure 2 shows four different examples of possible models that are attainable under
the SEM, LNM and RNM frameworks. Panel (a) shows a typical SEM model in which
one latent variable functions as a common cause of two others. Panel (b) shows a network
model which can be estimated using both the RNM and the LNM frameworks. Panel (c)
shows a completely equivalent LNM model to the SEM model of Panel (a) in which the
direction of effect between latent variables is not modeled. Finally, panel (d) shows a model
in which three exogenous latent variables underlie a set of indicators of which the residuals
form a network. The remainder of this section will describe RNM and LNM in more detail
and will outline the class of situations in which using these models is advantageous over
CFA or SEM.
5
We use the CFA framework instead of the SEM framework here as the main application of this framework
is in exploratively estimating relationships between latent variables.

GENERALIZED NETWORK PSYCHOMETRICS

(a) Structural Equation Modeling

(b) Network Modeling

(c) Latent Network Modeling

(d) Residual Network Modeling

9

Figure 2 . Examples of possible models under four different modeling frameworks. Circular nodes indicate latent variables, square nodes indicate manifest variables and gray
nodes indicate residuals. Directed edges indicate factor loadings or regression parameters
and undirected edges indicate pairwise interactions. Note that such undirected edges do
not indicate covariances, which are typically denoted with bidirectional edges. Replacing
covariances with interactions is where the network models differ from typical SEM.

GENERALIZED NETWORK PSYCHOMETRICS

10

Latent Network Modeling
The LNM framework models the latent variance–covariance matrix of a CFA model
as a GGM:
Σ = Λ∆Ψ (II − ΩΨ )−1 ∆Ψ Λ > + Θ .
Σ̂
(5)

This allows researchers to model conditional independence relationships between latent
variables without making the implicit assumptions of directionality or acyclicness. In SEM,
B is typically modeled as a directed acyclic graph (DAG), meaning that elements of B
can be represented by directed edges, and following along the path of these edges it is not
possible to return to any node (latent variable). The edges in such a DAG can be interpreted
as causal, and in general they imply a specific set of conditional independence relationships
between the nodes (Pearl, 2000).
While modeling conditional independence relationships between latent variables as a
DAG is a powerful tool for testing strictly confirmatory hypotheses, it is less useful for more
exploratory estimation. Though there have been recent advances in exploratory estimation
of DAGs within an SEM framework (e.g., Gates and Molenaar 2012; Rosa et al. 2012),
many equivalent DAGs can imply the same conditional independence relationships, and
thus fit the data equally well even though their causal interpretation can be strikingly
different (MacCallum et al., 1993). Furthermore, the assumption that the generating model
is acyclic—which, in practice, often is made on purely pragmatic grounds to identify a
model—is problematic in that much psychological behavior can be assumed to have at
least some cyclic and complex behavior and feedback (Schmittmann et al., 2013). Thus,
the true conditional independence relationships in a dataset can lead to many equivalent
compositions of B , and possibly none of them are the true model.
In psychometrics and SEM the GGM representation has not been very prominent,
even though it has some manifest benefits over the attempt to identify DAGs directly. For
example, by modeling conditional independence relationships between latent variables as
a GGM, many relationships can be modeled in a simpler way as compared to a DAG. In
addition, in the GGM each set of conditional independence relations only corresponds to one
model: there are no equivalent GGMs with the same nodes. Figure 3 shows a comparison of
several conditional independence relations that can be modeled equivalently or not by using
a GGM or by using a DAG. Panel (a) and Panel (b) show two DAGs that represent the
same conditional independence relations, y1 ⊥
⊥ y3 | y2 , which can both be represented by
the same GGM shown in Panel (e) and Panel (f). There are some conditional independence
relations that a GGM cannot represent in the same number of parameters as a DAG;
Panel (c) shows a collider structure that cannot be exactly represented by a GGM (the best
fitting GGM would feature three edges instead of two). On the other hand, there are also
conditional independence relationships that a GGM can represent and a DAG cannot; the
cycle of Panel (h) cannot be represented by a DAG. Further equivalences and differences
between GGMs and DAGs are beyond the scope of this paper, but haven been well described
in the literature (e.g., chapter 3 of Lauritzen 1996; Koller and Friedman 2009; Kolaczyk
2009). In sum, the GGM offers a natural middle ground between zero-order correlations and
DAGs: every set of zero-order correlations implies exactly one GGM, and every DAG implies
exactly one GGM. In a sense, the road from correlations to DAGs (including hierarchical
factor models) thus always must pass through the realm of GGMs, which acts as a bridge

GENERALIZED NETWORK PSYCHOMETRICS

Y2

Y2

11

Y2
No equivalent GGM with
the same # of parameters

Y3

Y1

(a)

Y3

Y1

(b)

Y3

Y1

(c)

Y2

(d)

Y2

Y2
No equivalent GGM with
the same # of parameters

Y1

(e)

Y3

Y1

(f)

Y3

Y1

Y3

Y4

(g)

(h)

Figure 3 . Equivalent models between directed acyclic graphs (DAG; top) and Gaussian
Graphical Models (GGM; bottom). Each column of graphs show two models that are
equivalent. Panels (a), (b), (e) and (f) all represent the same conditional independence
structure: Y1 and Y3 are independent after conditioning on Y2 . Panel (c) represents that
Y1 and Y3 are marginally independent even though they both are correlated with Y2 , a
structure that cannot be represented in a GGM with only two edges. Panel (h) shows
that Y1 and Y3 are independent after conditioning on the set Y2 and Y4 , which cannot be
represented with a DAG of four edges.
between the correlational and causal worlds.
Because there are no equivalent undirected models possible, LNM offers a powerful
tool for exploratory estimation of relationships between latent variables. For example,
suppose one encounters data generated by the SEM model in Figure 2, Panel (a). Without
prior theory on the relations between latent variables, exploratory estimation on this dataset
would lead to three completely equivalent models: the one shown in Figure 2, Panel (c)
and two models in which the common cause instead is the middle node in a causal chain.
As the number of latent variables increases, the potential number of equivalent models that
encode the same conditional independence relationships grows without bound. The LNM
model in Panel (c) of Figure 2 portrays the same conditional independence relationship as
the SEM model in Panel (a) of Figure 2, while having no equivalent model. Exploratory
estimation could easily find this model, and portrays the retrieved relationship in a clear
and unambiguous way.
A final benefit of using LNM models is that they allow network analysts to construct a network while taking measurement error into account. So far, networks have been
constructed based on single indicators only and no attempt has been made to remediate
measurement error. By forming a network on graspable small concepts measured by a few
indicators, the LNM framework can be used to control for measurement error.
Residual Network Modeling
In the RNM framework the residual structure of SEM is modeled as a GGM:
Σ = Λ (II − B )−1 Ψ (II − B )−1> Λ > + ∆Θ (II − ΩΘ )−1 ∆Θ .
Σ̂

(6)

GENERALIZED NETWORK PSYCHOMETRICS

12

This modeling framework conceptualizes latent variable and network modeling as two sides
of the same coin, and offers immediate benefits to both. In latent variable modeling, RNM
allows for the estimation of a factor structure (possibly including structural relationships
between the latent variables), while having no uncorrelated errors and thus no local independence. The error-correlations, however, are still highly structured due to the residual
network structure. This can be seen as a compromise between the ideas of network analysis
and factor modeling; while we agree that local independence is plausibly violated in many
psychometric tests, we think the assumption of no underlying latent traits and therefore a
sparse GGM may often be too strict. For network modeling, RNM allows a researcher to
estimate a sparse network structure while taking into account that some of the covariation
between items was caused by a set of latent variables. Not taking this into account would
lead to a saturated model (Chandrasekaran et al., 2012), whereas the residual network
structure can be sparse.
To avoid confusion between residual correlations, we will denote edges in ΩΘ residual
interactions. Residual interactions can be understood as pairwise linear effects, possibly due
to some causal influence or partial overlap between indicators that is left after controlling
for the latent structure. Consider again the indicators for agreeableness “Am indifferent to
the feelings of others” and “Inquire about others’ well-being”. It seems clear that we would
not expect these indicators to be locally independent after conditioning on agreeableness;
being indifferent to the feelings of others will cause one to not inquire about other’s wellbeing. Thus, we could expect these indicators to feature a residual interaction; some degree
of correlation between these indicators is expected to remain, even after conditioning on
the latent variable and all other indicators in the model.
The RNM framework in particular offers a new way of improving the fit of confirmatory factor models. In contrast to increasingly popular methods such as exploratory
SEM (ESEM; Marsh et al. 2014) or LASSO regularized SEM models (Jacobucci et al.,
2016), the RNM framework improves the fit by adding residual interactions rather than
allowing for more cross-loadings. The factor structure is kept exactly intact as specified
in the confirmatory model. Importantly, therefore, the interpretation of the latent factor
does not change. This can be highly valuable in the presence of a strong theory on the
latent variables structure underlying a dataset even in the presence of violations of local
independence.
Exploratory Network Estimation
Both the LNM and RNM modeling frameworks allow for confirmative testing of network structures. Confirmatory estimation is straightforward and similar to estimating SEM
models, with the exception that instead of modeling Ψ or Θ now the latent network ΩΨ or
ΩΘ is modeled. Furthermore, both modeling frameworks allow for the confirmatory fit of a
network model. In LNM, a confirmatory network structure can be tested by setting Λ = I
and Θ = O ; in RNM, a confirmatory network model can be tested by omitting any latent
variables. We have developed the R package lvnet6 , which utilizes OpenMx (Neale et al.,
2016) for confirmative testing of RNM and LNM models (as well as a combination of the
6

Developmental version: http://www.github.com/sachaepskamp/lvnet; stable version:
https://cran.r-project.org/package=lvnet

GENERALIZED NETWORK PSYCHOMETRICS

13

two). The lvnet function can be used for this purpose by specifying the fixed and the free
elements of model matrices. The package returns model fit indices (e.g., the RMSEA, CFI
and χ2 value), parameter estimates, and allows for model comparison tests.
Often the network structure, either at the residual or the latent level, is unknown
and needs to be estimated. To this end, the package includes two exploratory search
algorithms described below: step-wise model search and penalized maximum likelihood
estimation. For both model frameworks and both search algorithms, we present simulation
studies to investigate the performance of these procedures. As is typical in simulation
studies investigating the performance of network estimation techniques, we investigated
the sensitivity and specificity (van Borkulo et al., 2014). These measures investigate the
estimated edges versus the edges in the true model, with a ‘positive’ indicating an estimated
edge and a ‘negative’ indicating an edge that is estimated to be zero. Sensitivity, also termed
the true positive rate, gives the ratio of the number of true edges that were detected in the
estimation versus the total number of edges in the true model:
sensitivity =

# true positives
# true positives + # of false negatives

Specificity, also termed the true negative rate, gives the ratio of true missing edges detected
in the estimation versus the total number of absent edges in the true model:
specificity =

# true negatives
# true negatives + # false positives

The specificity can be seen as a function of the number of false positives: a high specificity
indicates that there were not many edges detected to be nonzero that are zero in the true
model. To favor degrees of freedom, model sparsity and interpretability, specificity should be
high all-around—estimation techniques should not result in many false positives—whereas
sensitivity should increase as a function of the sample size.
Simulating Gaussian Graphical models
In all simulation studies reported here, networks were constructed in the same way as
done by Yin and Li (2011) in order to obtain a positive definite inverse-covariance matrix
K . First, a network structure was generated without weights. Next, weights were drawn
randomly from a uniform distribution between 0.5 and 1, and made negative with 50%
probability. The diagonal elements of K were then set to 1.5 times the sum of all absolute
values in the corresponding row, or 1 if this sum was zero. Next, all values in each row were
divided by the diagonal value, ensuring that the diagonal values become 1. Finally, the
matrix was made symmetric by averaging the lower and upper triangular elements. In the
chain graphs used in the following simulations, this algorithm created networks in which
the non-zero partial correlations had a mean of 0.33 and a standard deviation of 0.04.
Stepwise Model Search
In exploratory search, we are interested in recovering the network structure of either
ΩΨ in LNM or ΩΘ in RNM. This can be done through a step-wise model search, either based

GENERALIZED NETWORK PSYCHOMETRICS

14

on χ2 difference tests (Algorithm 1) or on minimization of some information criterion (Algorithm 2) such as the Akaike information criterion (AIC), Bayesian information criterion
(BIC) or the extended Bayesian information criterion (EBIC; Chen and Chen 2008) which
is now often used in network estimation (van Borkulo et al., 2014; Foygel and Drton, 2010).
In LNM, removing edges from ΩΨ cannot improve the fit beyond that of an already fitting
CFA model. Hence, model search for ΩΨ should start at a fully populated initial setup for
ΩΨ . In RNM, on the other hand, a densely populated ΩΘ would lead to an over-identified
model, and hence the step-wise model search should start at an empty network ΩΘ = O .
The function lvnetSearch in the lvnet package can be used for both search algorithms.
Algorithm 1 Stepwise network estimation by χ2 difference testing.
Start with initial setup for Ω
repeat
for all Unique elements of Ω do
Remove edge if present or add edge if absent
Fit model with changed edge
end for
if Adding an edge significantly improves fit (α = 0.05) then
Add edge that improves fit the most
else if Removing an edge does not significantly worsen fit (α = 0.05) then
Remove edge that worsens fit the least
end if
until No added edge significantly improves fit and removing any edge significantly worsens
fit
Algorithm 2 Stepwise network estimation by AIC, BIC or EBIC optimization.
Start with initial setup for Ω
repeat
for all Unique elements of Ω do
Remove edge if present or add edge if absent
Fit model with changed edge
end for
if Any changed edge improved AIC, BIC or EBIC then
Change edge that improved AIC, BIC or EBIC the most
end if
until No changed edge improves AIC, BIC or EBIC
Simulation Study 1: Latent Network Modeling. We performed a simulation
study to assess the performance of the above mentioned step-wise search algorithms in LNM
models. Figure 4 shows the LNM model under which we simulated data. In this model,
four latent factors with three indicators each were connected in a latent network. The latent
network was a chain network, leading all latent variables to be correlated according to a
structure that cannot be represented in SEM. Factor loadings and residual variances were
set to 1, and the network weights were simulated as described in the section “Simulating

GENERALIZED NETWORK PSYCHOMETRICS

15

Figure 4 . Model used in simulation study 1: step-wise model search in latent network
modeling. Four latent variables were each measured by three items. Latent variables
covary due to the structure of a latent Gaussian graphical model in which edges indicate
partial correlation coefficients. This model has the form of a chain graph, which cannot be
represented in a structural equation model. Factor loadings, residual variances and latent
variances were set to 1 and the latent partial correlations had an average of 0.33 with a
standard deviation of 0.04.
Gaussian Graphical models”. The simulation study followed a 5 × 4 design: the sample
size was varied between 50, 100, 250, 500 and 1 000 to represent typical sample sizes in
psychological research, and the stepwise evaluation criterion was either χ2 difference testing,
AIC, BIC or EBIC (using a tuning parameter of 0.5). Each condition was simulated 1 000
times, resulting in 20 000 total simulated datasets.
Figure 5 shows the results of the simulation study. Data is represented in standard
boxplots (McGill et al., 1978): the box shows the 25th, 50th (median) 75th quantiles, the
whiskers range from the largest values in 1.5 times the inter-quantile range (75th - 25th
quantile) and points indicate outliers outside that range. In each condition, we investigated
the sensitivity and specificity. The top panel shows that sensitivity improves with sample
size, with AIC performing best and EBIC worst. From sample sizes of 500 and higher all
estimation criterion performed well in retrieving the edges. The bottom panel shows that
specificity is generally very high, with EBIC performing best and AIC worst. These results
indicate that the step-wise procedure is conservative and prefers simpler models to more
complex models; missing edges are adequately detected but present edges in the true model
might go unnoticed except in larger samples. With sample sizes over 500, all four estimation

GENERALIZED NETWORK PSYCHOMETRICS

16

Sensitivity (true positive rate)
AIC

Chi−square

BIC

EBIC

Sensitivity

1.00

●

0.75

●

●

0.50

●

●

0.25

0.00

●

●

●

●

●

●

●

●

●

●

50

●

●

100 250 500 1000

50

100 250 500 1000

50

●

100 250 500 1000

50

100 250 500 1000

Sample Size

Specificity (true negative rate)
AIC

Chi−square

BIC

EBIC

1.00

Specificity

0.75

0.50

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

0.25

0.00
50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

Sample Size

Figure 5 . Simulation results of simulation study 1: step-wise model search in latent network
modeling. Each condition was replicated 1 000 times, leading to 20 000 total simulated
datasets. High sensitivity indicates that the method is able to detect edges in the true
model, and high specificity indicates that the method does not detect edges that are zero
in the true model.

GENERALIZED NETWORK PSYCHOMETRICS

17

Figure 6 . Model used in simulation study 2: step-wise model search in residual network
modeling. Two latent variables were each measured by five items; a Gaussian graphical
model, in which edges indicate partial correlation coefficients, leads to all residuals to be
correlated due to a chain graph between residuals, which cannot be represented in a structural equation model. Factor loadings, residual variances and latent variances were set to
1, the factor covariance was set to 0.25 and the latent partial correlations had an average
of 0.33 with a standard deviation of 0.04.
methods show both a high sensitivity and specificity.
Simulation Study 2: Residual Network Modeling. We conducted a second
simulation study to assess the performance of step-wise model selection in RNM models.
Figure 7 shows the model under which data were simulated: two latent variables with 5
indicators each. The residual network was constructed to be a chain graph linking a residual
of an indicator of one latent variable to two indicators of the other latent variable. This
structure cannot be represented by a DAG and causes all residuals to be connected, so
that Θ is fully populated. Factor loadings and residual variances were set to 1, the factor
covariance was set to 0.25, and the network weights were simulated as described in the
section “Simulating Gaussian Graphical models”.
The simulation study followed a 5×4 design; sample size was again varied between 50,
100, 250, 500 and 1 000, and models were estimated using either χ2 significance testing, AIC,
BIC or EBIC. Factor loadings and factor variances were set to 1 and the factor correlation
was set to 0.25. The weights in ΩΘ were chosen as described in the section “Simulating
Gaussian Graphical models”. Each condition was simulated 1, 000 times, leading to 20 000
total datasets.
Figure 7 shows the results of the simulation study. The top panel shows that sensitivity increases with sample size and performs best when using AIC as the criterion. BIC
performed comparably in sensitivity to χ2 testing and EBIC performed the worst. The
bottom panel shows that specificity was very high for all sample sizes and all criteria, with
EBIC performing best and AIC worst. These results indicate that the number of false pos-

GENERALIZED NETWORK PSYCHOMETRICS

18

Sensitivity (true positive rate)
AIC

Chi−square

BIC

1.00

EBIC

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

Sensitivity

0.75

0.50

0.25
●

●

●

●

●

●

0.00

●

50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

50

●

100 250 500 1000

Sample Size

Specificity (true negative rate)
AIC
1.00

Specificity

0.75

●

●
●
●
●

Chi−square

●

●
●
●

BIC

EBIC

●

●
●
●
●

●
●
●
●
●
●

●
●
●
●

●

●

●
●
●

●
●
●
●
●

●
●
●
●
●
●
●
●
●

●
●
●

●
●
●

●
●

●
●
●
●
●

●
●
●
●
●
●

●
●
●
●
●
●
●
●
●

●
●
●

●

●
●
●
●
●
●

●
●

●
●
●

●
●
●
●
●
●
●

●
●
●
●
●
●
●
●
●

●

●
●

0.50
●

●
●

0.25

0.00
50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

Sample Size

Figure 7 . Simulation results of simulation study 2: step-wise model search in residual
network modeling. Each condition was replicated 1 000 times, leading to 20 000 total simulated datasets. High sensitivity indicates that the method is able to detect edges in the
true model, and high specificity indicates that the method does not detect edges that are
zero in the true model.

GENERALIZED NETWORK PSYCHOMETRICS

19

itives is very low and that the method is on average well capable of discovering true edges
for sample size larger than 250. In sum, all four criteria perform well with EBIC erring on
the side of caution and AIC erring on the side of discovery.
LASSO Regularization
While the step-wise model selection algorithms perform well in retrieving the correct
network structure, they are very slow when the number of nodes in the network increases
(e.g., more than 10 nodes). This is particularly important in the context of RNM, in which
the number of indicators can be larger than 10 even in small models. A popular method
for fast estimation of high-dimensional network structures is by applying the least absolute
shrinkage and selection operator (LASSO; Tibshirani 1996). LASSO regularization has
also recently been introduced in the SEM literature (Jacobucci et al., 2016) as a method for
obtaining sparser structures of Λ and B . In the LASSO, instead of optimizing the likelihood
function as described in (1), a penalized likelihood is optimized (Jacobucci et al., 2016):
h

 



Σ + Trace S Σ̂
Σ
min log det Σ̂
Σ
Σ̂

−1 

 

i

S − P + νPenalty ,
− log det Ŝ

(7)

in which ν denotes a tuning parameter controlling the level of penalization. The penalty
here is taken to be the sum of absolute parameters:
Penalty =

X

<i,j>

|ωij |,

in which ωij denotes an element from either ΩΨ or ΩΘ . Other penalty functions may be
used as well—such as summing the squares of parameter estimates (ridge regression; Hoerl
and Kennard 1970) or combining both absolute and squared values (elastic net; Zou and
Hastie 2005)—but these are not currently implemented in lvnet. The benefit of the LASSO
is that it returns models that perform better in cross-validation. In addition, the LASSO
yields sparse models in which many relationships are estimated to be zero.
Algorithm 3 LASSO estimation for exploratory network search.
for all Sequence of tuning parameters ν1 , ν2 , . . . do
Estimate LASSO regularized model using given tuning parameter
Count the number of parameters for which the absolute estimate is larger than 
Determine information criterion (AIC, BIC or EBIC) given fit and number of parameters
end for
Select model with best AIC, BIC or EBIC
Refit this model without LASSO in which absolute parameters smaller than  are fixed
to zero
ΩΘ ,
The lvnet function allows for LASSO regularization for a given model matrix (Ω
ΩΨ , Θ , Ψ , Λ or B ) and a given value for the tuning parameter ν. The optimizer used
in lvnet does not return exact zeroes. To circumvent this issue, any absolute parameter
below some small value  (by default  = 0.0001) is treated as zero in counting the number of

GENERALIZED NETWORK PSYCHOMETRICS

20

parameters and degrees of freedom (Zou et al., 2007). The lvnetLasso function implements
the search algorithm described in Algorithm 3 to automatically choose an appropriate tuning
parameter, use that for model selection and rerun the model to obtain a comparable fit to
non-regularized models. In this algorithm, a sequence of tuning parameters is tested, which
is set by default to a logarithmically spaced sequence of 20 values between 0.01 and 1.
Simulation Study 3: Latent Network Modeling. We studied the performance
of LASSO penalization in estimating the latent network structure in a similar simulation
study to the study of the step-wise procedure described above. Data were simulated under
a similar model to the one shown in Figure 4, except that now 8 latent variables were
used leading to a total of 24 observed variables. All parameter values were the same as in
simulation study 1. The simulation followed a 5 × 3 design. Sample size was varied between
100, 250, 500, 1 000 and 2 500, and for each sample size 1 000 datasets were simulated leading
to a total of 5 000 generated datasets. On these datasets the best model was selected using
either AIC, BIC or EBIC, leading to 15 000 total replications. In each replication, sensitivity
and specificity were computed. Figure 8 shows that AIC had a relatively poor specificity
all-around, but a high sensitivity. EBIC performed well with sample sizes of 500 and higher.
Simulation Study 4: Residual Network Modeling. To assess the performance
of LASSO in estimating the residual network structure we simulated data as in Figure 6,
except that in this case four latent variables were used, each with 5 indicators, the residuals
of which were linked via a chain graph. All parameter values were the same as in simulation
study 2. The design was the same as in simulation study 3, leading to 5 000 generated
datasets on which AIC, BIC or EBIC were used to select the best model. While Figure 9
shows good performance of the LASSO in retrieving the residual network structure and
similar results as before: AIC performs the worst in specificity and EBIC the best.
Empirical Example: Personality Inventory
In this section, we demonstrate LNM and RNM models by confirmatively testing a
model and exploratively searching a residual and a latent network structure. We use the
lvnet package, which can be installed in R as follows:
> install.packages("lvnet")
> library("lvnet")
To exemplify the method, we will use a dataset from the psych package (Revelle, 2010)
on the Big 5 personality traits (Benet-Martinez and John, 1998; Digman, 1989; Goldberg,
1990, 1993; McCrae and Costa, 1997). This dataset consists of 2800 observations of 25 items
designed to measure the 5 central personality traits with 5 items per trait. We estimated
the CFA model on the BFI dataset. Next, we used LASSO estimation to the RNM model
using 100 different tuning parameters and using EBIC as criterion to maximize specificity
and search for a sparse residual network. The fully correlated covariance matrix of latent
variables is equivalent to a fully connected latent network structure. Thus, after fitting a
RNM model, we can again apply LASSO to estimate a latent network in the resulting model,
which we abbreviate here to an RNM+LNM model. The R code used for this analysis can
be found in the supplementary materials.

GENERALIZED NETWORK PSYCHOMETRICS

21

Sensitivity (true positive rate)
AIC

BIC

EBIC

1.00

●

●

●

Sensitivity

0.75

●

●

●

●

0.50

●

●

●

0.25

●

●

0.00

●

100

250

500

1000

2500

100

250

500

1000

2500

100

250

500

1000

2500

Sample Size

Specificity (true negative rate)
AIC

BIC

1.00

EBIC

●
●
●
●
●

0.75

●

Specificity

●

●
●

0.50

●

●

●

●

●

●

●

●

●

●

●

●

●

●
●

●

0.25
●

●

●

●

●

●

●

●

●

●

●

250

500

1000

●

0.00

●

100

2500

100

250

500

1000

2500

100

250

500

1000

2500

Sample Size

Figure 8 . Simulation results of simulation study 3: model selection via penalized maximum
likelihood estimation in latent network modeling. The same model as in Figure 4 was used
except now with 4 latent variables leading to 24 observed variables. For each sample size
1 000 datasets were generated, leading to 5 000 total simulated datasets on which AIC, BIC
or EBIC was used to select the best model. High sensitivity indicates that the method is
able to detect edges in the true model, and high specificity indicates that the method does
not detect edges that are zero in the true model.

GENERALIZED NETWORK PSYCHOMETRICS

22

Sensitivity (true positive rate)
AIC

BIC

EBIC

1.00

●
●

●

●

●

●

●
●
●

Sensitivity

0.75

●
●

●

●

●

●

●

●

0.50

●

●

●

●

●

●

●

●

●

0.25

●

●

●

●

●

●

●

●

●

0.00
100

250

500

1000

2500

100

250

500

1000

2500

100

250

500

1000

2500

Sample Size

Specificity (true negative rate)
AIC

BIC

EBIC

1.00
●
●

Specificity

0.75

●

●
●
●
●
●
●
●
●
●
●
●

●
●
●
●

●
●
●

●
●
●
●

●
●
●
●
●
●

●

●
●
●
●
●
●

●
●
●
●

250

500

1000

●
●
●
●

●
●
●
●
●
●
●
●
●
●
●
●
●
●

0.50

●
●
●
●
●
●

●
●
●
●
●
●
●
●

●

●
●
●
●
●

●

0.25

0.00
100

250

500

1000

2500

100

250

500

1000

2500

100

2500

Sample Size

Figure 9 . Simulation results of simulation study 4: model selection via penalized maximum
likelihood estimation in residual network modeling. The same model as in Figure 6 was
used except now with 4 latent variables leading to 20 observed variables. For each sample
size 1 000 datasets were generated, leading to 5 000 total simulated datasets on which AIC,
BIC or EBIC was used to select the best model. High sensitivity indicates that the method
is able to detect edges in the true model, and high specificity indicates that the method
does not detect edges that are zero in the true model.

GENERALIZED NETWORK PSYCHOMETRICS
χ2
4713.94
806.63
843.18

df
265
172
176

CFA
RNM
RNM+LNM

AIC
183233.7
179511.0
179539.5

BIC
183589.9
180419.4
180424.2

23
EBIC
184542.4
182848.2
182789.5

RMSEA
0.08
0.04
0.04

TLI
0.75
0.94
0.94

CFI
0.78
0.97
0.97

Table 1
Fit measures for three models estimated on the BFI dataset in the psych R package. CFA is
the correlated five-factor model. RNM is the same model as the CFA model with a residual
network. RNM+LNM denotes the same model as the RNM model in which edges of the
latent network have been removed.
Factor structure & latent network

Residual network

●
●

A1

O5

●

C2

A2

●

A3

O4

●

O3

E1

A4

A4
●

E2

O2

A5

O2

●
●

O4
O

A

●

E4

O1

C1

C1
N5

N5

N4

●

C2

N3
A5

N

N4

C

C3

●
●
●
●

C5

C4
O3

N3

●

O5

N2

O1

N1

●

C4

●

E

●

E5
N2

C5

●

A1

●

Agreeableness
A1: Am indifferent to the feelings of others.
A2: Inquire about others' well−being.
A3: Know how to comfort others.
A4: Love children.
A5: Make people feel at ease.
Conscientiousness
C1: Am exacting in my work.
C2: Continue until everything is perfect.
C3: Do things according to a plan.
C4: Do things in a half−way manner.
C5: Waste my time.
Extraversion
E1: Don't talk a lot.
E2: Find it difficult to approach others.
E3: Know how to captivate people.
E4: Make friends easily.
E5: Take charge.
Neuroticism
N1: Get angry easily.
N2: Get irritated easily.
N3: Have frequent mood swings.
N4: Often feel blue.
N5: Panic easily.

E3
N1

E1
E5

A3

E2
E4

E3

A2

C3

●
●
●
●
●

Openness
O1: Am full of ideas.
O2: Avoid difficult reading material.
O3: Carry the conversation to a higher level.
O4: Spend time reflecting on things.
O5: Will not probe deeply into a subject.

Figure 10 . Visualization of the factor structure and latent network (left) and the residual
network (right) of the BFI personality dataset from the psych package in R. LASSO estimation with 100 different tuning parameters in combination with EBIC model selection
was used to first estimate the residual network structure and following the latent network
structure.
Table 1 shows the fit of the three models. The CFA model fits poorly. The RNM
model has substantively improved fit and resulted in good fit indices overall. The estimated RNM+LNM model showed that 5 edges could be removed from the latent network
after taking residual interactions into account. Figure 10 shows the factor structure and
residual network of the final RNM+LNM model. It can be seen that Agreeableness is now
only connected to extraversion: after taking into account someone’s level of extraversion
agreeableness is independent of the other three personality traits. Extraversion is the most
central node in this network and the only trait that is directly linked to all other traits. The
residual network shows many meaningful connections. While seemingly densely connected,
this network only has 30% of all possible edges in a network of that size, leading the model
to have 176 degrees of freedom. The corresponding residual covariance structure is fully
populated with no zero elements.
It should be noted that the procedures used in this example are highly explorative.
While Table 1 shows that the RNM fits better than the CFA model, the CFA model is
solely based on theory whereas the RNM model was found through high-dimensional model
search. As a result, to substantially interpret the structure found in Figure 10 it should
first be replicated in independent samples.

GENERALIZED NETWORK PSYCHOMETRICS

24

Conclusion
In this paper we introduced a formal psychometric model for network modeling of multivariate normal data. We contrasted this model with latent variable models as commonly
used in CFA and SEM. Furthermore, using the CFA and SEM frameworks, we proposed
two generalizations of the network model to encompass latent variable structures within
the network paradigm. In the first generalization, LNM, we construct a network among the
latent variables, whereas in the second generalization, RNM, a network is formed among
residuals of indicators. Both frameworks offer powerful benefits over both latent variable
and network modeling. From the perspective of latent variable modeling, the LNM framework allows one to exploratively search for conditional independence relationships between
latent variables without the need for prior theory, and the RNM framework allows one to
model latent common causes without assuming local independence. From the perspective
of network modeling, the LNM framework allows one to model network structures while
taking measurement error into account, and the RNM framework allows one to estimate a
network structure, even when all nodes are in part caused by unobserved or latent variables.
In addition, both frameworks allow for network models to be fitted and compared to SEM
models. The discussed methodology has been implemented in the freely available R package
lvnet.
Simulation studies showed that step-wise search and penalized maximum likelihood estimation of the residual or latent network structures resulted in high specificity all around—
the methods did not result often in false positives—and rapidly increasing sensitivity as a
function of the sample size; the higher the sample size, the more true edges were detected in
the algorithm. These numbers are comparable to state-of-the-art network estimation techniques in sample and model sizes that are plausible in psychological settings (van Borkulo
et al., 2014). In all four simulation studies, using AIC as the model selection criterion led
to the best sensitivity and using EBIC led to the best specificity. However, it is important
to note that the choice of a particular information criterion cannot be argued by these
numbers alone, and depends on the relative importance one assigns to the side of discovery
(optimizing sensitivity) or the side of caution (optimizing specificity; Dziak et al. 2012).
Furthermore, it should be noted that these simulation results are specific for the particular
model setup and sample sizes used; results might be different for other kinds of models or
sample size ranges.
In addition to the LNM and RNM frameworks, other combinations of CFA, SEM
and network modeling are possible as well. For example, a framework can be constructed
which contains both a latent and a residual network (as shown in our empirical example),
or directed regression paths as in the SEM model can be added to the LNM model. While
these models are all estimable in the lvnet software, in the current paper we chose to focus
on the distinct benefits that modeling a residual or latent network presents. Thus, in this
paper, we only described the modeling of multivariate normal data. More advanced models
are possible, but not yet implemented in the lvnet software. In the case of binary variables,
the appropriate model to use is the Ising model, which has been shown to be equivalent
to multivariate item response models (Epskamp et al., ress; Marsman et al., 2015). Future
research could aim at constructing Ising models among binary latent variables in latent
class analysis, or constructing residual networks in models with binary indicators. Finally,

GENERALIZED NETWORK PSYCHOMETRICS

25

the expressions optimized in Equations (1) and (7) are based on summary statistics and
therefore only truly applicable to complete data. With incomplete data, the appropriate
estimation method is to use full-information maximum likelihood (FIML; Arbuckle et al.
1996); however, FIML has not yet been implemented in the lvnet software.
In our view, the presented modeling framework is a versatile and promising addition
to the spectrum of psychometric models. The GGM, which has a central place in this
modeling framework, acts as a natural interface between correlation and causality, and
we think this representation should receive more attention in psychometrics. From the
point of view afforded by the current paper, the typical attempt to determine directed
SEMs from correlation structures in fact appears somewhat haphazard in psychology, a
historical accident in a field that has been prematurely directed to hypothesis testing at
the expense of systematic exploration. Perhaps, psychometrics as a field should consider
taking a step back to focus on the consistent identification of GGMs, instead of wanting
to jump to the causal conclusion immediately. In this regard, the fact that GGMs do
not have equivalent models would appear to be a major benefit, as they allow us to focus
on charting connections between variables systematically, without being forced to adhere
to one particular causal interpretation or another. In addition, because the GGM does
not specify the nature or direction of interactions between variables, it appears a natural
model for research situations where no temporal information or experimental interventions
are present, so that associations may arise for a multitude of reasons: the GGM can be
consistently interpreted regardless of whether associations arise from direct causal relations,
reciprocal causation, latent common causes, semantic overlap between items, or homeostatic
couplings of parameters. While this can be seen as a downside of the GGM—the lack of
directionality leads to fewer falsifiable hypotheses—it can also be a major asset in a field
like psychology, where strong causal theory is sparse and the identification of DAGs often
appears a bridge too far.
References
Arbuckle, J. L., Marcoulides, G. A., and Schumacker, R. E. (1996). Full information estimation in the presence of incomplete data. Advanced structural equation modeling: Issues
and techniques, pages 243–277.
Benet-Martinez, V. and John, O. (1998). Los Cinco Grandes across Cultures and Ethnic
Groups: Multitrait Multimethod Analyses of the Big Five in Spanish and English. Journal
of Personality and Social Psychology, 75:729–750.
Bentler, P. M. (1990). Comparative fit indexes in structural models. Psychological bulletin,
107(2):238–346.
Borsboom, D. (2008). Psychometric Perspectives on Diagnostic Systems. Journal of clinical
psychology, 64(9):1089–1108.
Borsboom, D. and Cramer, A. O. J. (2013). Network analysis: An integrative approach to
the structure of psychopathology. Annual Review of Clinical Psychology, 9:91–121.
Borsboom, D., Cramer, A. O. J., Schmittmann, V. D., Epskamp, S., and Waldorp, L. J.
(2011). The small world of psychopathology. PloS one, 6(11):e27407.
Browne, M. W. and Cudeck, R. (1992). Alternative ways of assessing model fit. Sociological
Methods & Research, 21(2):230–258.

GENERALIZED NETWORK PSYCHOMETRICS

26

Chandrasekaran, V., Parrilo, P. A., and Willsky, A. S. (2012). Latent variable graphical
model selection via convex optimization (with discussion). The Annals of Statistics,
40(4):1935–1967.
Chen, J. and Chen, Z. (2008). Extended bayesian information criteria for model selection
with large model spaces. Biometrika, 95(3):759–771.
Costantini, G., Epskamp, S., Borsboom, D., Perugini, M., Mõttus, R., Waldorp, L. J., and
Cramer, A. O. J. (2015). State of the aRt personality research: A tutorial on network
analysis of personality data in R. Journal of Research in Personality, 54:13–29.
Cramer, A. O. J., Sluis, S., Noordhof, A., Wichers, M., Geschwind, N., Aggen, S. H.,
Kendler, K. S., and Borsboom, D. (2012). Dimensions of normal personality as networks
in search of equilibrium: You can’t like parties if you don’t like people. European Journal
of Personality, 26(4):414–431.
Cramer, A. O. J., Waldorp, L., van der Maas, H., and Borsboom, D. (2010). Comorbidity:
A Network Perspective. Behavioral and Brain Sciences, 33(2-3):137–150.
Dalege, J., Borsboom, D., van Harreveld, F., van den Berg, H., Conner, M., and van der
Maas, H. L. J. (2016). Toward a formalized account of attitudes: The Causal Attitude
Network (CAN) model. Psychological review, 123(1):2–22.
Digman, J. (1989). Five Robust Trait Dimensions: Development, Stability, and Utility.
Journal of Personality, 57(2):195–214.
Dziak, J. J., Coffman, D. L., Lanza, S. T., and Li, R. (2012). Sensitivity and specificity of
information criteria. The Methodology Center and Department of Statistics, Penn State,
The Pennsylvania State University.
Ellis, J. L. and Junker, B. W. (1997). Tail-measurability in monotone latent variable models.
Psychometrika, 62(4):495–523.
Epskamp, S., Maris, G., Waldorp, L., and Borsboom, D. (in press). Network psychometrics.
In Irwing, P., Hughes, D., and Booth, T., editors, Handbook of Psychometrics. Wiley, New
York, NY, USA.
Epskamp, S., Waldorp, L. J., Mõttus, R., and Borsboom, D. (2016). Discovering psychological dynamics in time-series data. arXiv preprint, page arXiv:1609.04156.
Foygel, R. and Drton, M. (2010). Extended Bayesian information criteria for Gaussian
graphical models. Advances in Neural Information Processing Systems, 23:2020–2028.
Fried, E. I., Bockting, C., Arjadi, R., Borsboom, D., Amshoff, M., Cramer, O. J., Epskamp, S., Tuerlinckx, F., Carr, D., and Stroebe, M. (2015). From loss to loneliness:
The relationship between bereavement and depressive symptoms. Journal of abnormal
psychology, 124(2):256–265.
Fried, E. I. and van Borkulo, C. (2016). Mental disorders as networks of problems: a review
of recent insights.
Gates, K. M. and Molenaar, P. C. (2012). Group search algorithm recovers effective connectivity maps for individuals in homogeneous and heterogeneous samples. NeuroImage,
63(1):310–319.
Goldberg, L. (1993). The Structure of Phenotypic Personality Traits. American Psychologist, 48(1):26–34.
Goldberg, L. R. (1990). An alternative “description of personality”: the big-five factor
structure. Journal of personality and social psychology, 59(6):1216–1229.

GENERALIZED NETWORK PSYCHOMETRICS

27

Hayduk, L. A. (1987). Structural equation modeling with LISREL: Essentials and advances.
Johns Hopkins University Press, Baltimore, MD, USA.
Hoerl, A. E. and Kennard, R. W. (1970). Ridge regression: Biased estimation for nonorthogonal problems. Technometrics, 12(1):55–67.
Holland, P. W. and Rosenbaum, P. R. (1986). Conditional association and unidimensionality
in monotone latent variable models. The Annals of Statistics, 14:1523–1543.
Ising, E. (1925). Beitrag zur theorie des ferromagnetismus. Zeitschrift für Physik A Hadrons
and Nuclei, 31(1):253–258.
Isvoranu, A. M., Borsboom, D., van Os, J., and Guloksuz, S. (2016a). A Network Approach to Environmental Impact in Psychotic Disorders: Brief Theoretical Framework.
Schizophrenia Bulletin, 42(4):870–873.
Isvoranu, A. M., van Borkulo, C. D., Boyette, L., Wigman, J. T. W., Vinkers, C. H.,
Borsboom, D., and GROUP Investigators (2016b). A Network Approach to Psychosis:
Pathways between Childhood Trauma and Psychotic Symptoms. Schizophrenia Bulletin.
Advance Access published May 10, 2016.
Jacobucci, R., Grimm, K. J., and McArdle, J. J. (2016). Regularized structural equation
modeling. Structural Equation Modeling: A Multidisciplinary Journal, 23(4):555–566.
Jöreskog, K. G. (1967). A general approach to confirmatory maximum likelihood factor
analysis. ETS Research Bulletin Series, 1967(2):183–202.
Kaplan, D. (2000). Structural equation modeling: Foundations and extensions. Sage, Thousand Oaks, CA, USA.
Kolaczyk, E. D. (2009). Statistical analysis of network data. Springer, New York, NY, USA.
Koller, D. and Friedman, N. (2009). Probabilistic graphical models: Principles and techniques. MIT press, Cambridge, MA, USA.
Kossakowski, J. J., Epskamp, S., Kieffer, J. M., van Borkulo, C. D., Rhemtulla, M., and
Borsboom, D. (2015). The application of a network approach to health-related quality of
life (HRQoL): Introducing a new method for assessing hrqol in healthy adults and cancer
patient. Quality of Life Research, 25:781–92.
Lauritzen, S. L. (1996). Graphical models. Clarendon Press, Oxford, UK.
Lawley, D. N. (1940). VI.—the estimation of factor loadings by the method of maximum
likelihood. Proceedings of the Royal Society of Edinburgh, 60(01):64–82.
Lord, F. M., Novick, M. R., and Birnbaum, A. (1968). Statistical theories of mental test
scores. Addison-Wesley, Oxford, UK.
MacCallum, R. C., Wegener, D. T., Uchino, B. N., and Fabrigar, L. R. (1993). The problem of equivalent models in applications of covariance structure analysis. Psychological
Bulletin, 114(1):185–199.
Marsh, H. W., Morin, A. J., Parker, P. D., and Kaur, G. (2014). Exploratory structural
equation modeling: An integration of the best features of exploratory and confirmatory
factor analysis. Annual Review of Clinical Psychology, 10:85–110.
Marsman, M., Maris, G., Bechger, T., and Glas, C. (2015). Bayesian inference for low-rank
ising networks. Scientific reports, 5(9050):1–7.
McCrae, R. R. and Costa, P. T. (1997). Personality trait structure as a human universal.
American Psychologist, 52(5):509–516.
McGill, R., Tukey, J. W., and Larsen, W. A. (1978). Variations of box plots. The American
Statistician, 32(1):12–16.

GENERALIZED NETWORK PSYCHOMETRICS

28

McNally, R. J., Robinaugh, D. J., Wu, G. W., Wang, L., Deserno, M. K., and Borsboom, D.
(2015). Mental disorders as causal systems a network approach to posttraumatic stress
disorder. Clinical Psychological Science, 3(6):836–849.
Murphy, K. P. (2012). Machine learning: A probabilistic perspective. MIT press, Cambridge,
MA, USA.
Neale, M. C., Hunter, M. D., Pritikin, J. N., Zahery, M., Brick, T. R., Kirkpatrick, R. M.,
Estabrook, R., Bates, T. C., Maes, H. H., and Boker, S. M. (2016). Openmx 2.0: Extended
structural equation and statistical modeling. Psychometrika, 81(2):535–549.
Pearl, J. (2000). Causality: Models, Reasoning, and Inference. Cambridge University Press,
New York, NY.
Reckase, M. D. (2009). Multidimensional item response theory. Springer, New York, NY,
USA.
Revelle, W. (2010). psych: Procedures for Psychological, Psychometric, and Personality
Research (R package version 1.0-93). Northwestern University, Evanston, Illinois.
Rosa, M., Friston, K., and Penny, W. (2012). Post-hoc selection of dynamic causal models.
Journal of neuroscience methods, 208(1):66–78.
Schmittmann, V. D., Cramer, A. O. J., Waldorp, L. J., Epskamp, S., Kievit, R. A., and
Borsboom, D. (2013). Deconstructing the construct: A network perspective on psychological phenomena. New Ideas in Psychology, 31(1):43–53.
Tibshirani, R. (1996). Regression shrinkage and selection via the lasso. Journal of the Royal
Statistical Society. Series B (Methodological), 58:267–288.
van Borkulo, C. D., Borsboom, D., Epskamp, S., Blanken, T. F., Boschloo, L., Schoevers,
R. A., and Waldorp, L. J. (2014). A new method for constructing networks from binary
data. Scientific Reports, 4(5918):1–10.
van Borkulo, C. D., Boschloo, L., Borsboom, D., Penninx, B. W. J. H., Waldorp, L. J., and
Schoevers, R. A. (2015). Association of symptom network structure with the course of
depression. JAMA Psychiatry, 72(12):1219–1226.
van der Maas, H. L., Dolan, C. V., Grasman, R. P., Wicherts, J. M., Huizenga, H. M.,
and Raijmakers, M. E. (2006). A dynamical model of general intelligence: The positive
manifold of intelligence by mutualism. Psychological Review, 113(4):842–861.
Ware Jr, J. E. and Sherbourne, C. D. (1992). The MOS 36-item short-form health survey
(SF-36): I. Conceptual framework and item selection. Medical Care, 30:473–483.
Wright, S. (1921). Correlation and causation. Journal of agricultural research, 20(7):557–
585.
Wright, S. (1934). The method of path coefficients. The Annals of Mathematical Statistics,
5(3):161–215.
Yin, J. and Li, H. (2011). A sparse conditional gaussian graphical model for analysis of
genetical genomics data. The Annals of Applied Statistics, 5(4):2630–2650.
Zou, H. and Hastie, T. (2005). Regularization and variable selection via the elastic net.
Journal of the Royal Statistical Society: Series B (Statistical Methodology), 67(2):301–
320.
Zou, H., Hastie, T., Tibshirani, R., et al. (2007). On the “degrees of freedom” of the lasso.
The Annals of Statistics, 35(5):2173–2192.

GENERALIZED NETWORK PSYCHOMETRICS: COMBINING
NETWORK AND LATENT VARIABLE MODELS

arXiv:1605.09288v4 [math.ST] 11 Sep 2017

Sacha Epskamp, Mijke Rhemtulla & Denny Borsboom
university of amsterdam

September 13, 2017

Manuscript submitted to “Psychometrika”
Correspondence should be sent to
E-Mail: sacha.epskamp@gmail.com
Address: POBox 15906, 1001 NK Amsterdam
Website: www.sachaepskamp.com

Psychometrika Submission

September 13, 2017

2

GENERALIZED NETWORK PSYCHOMETRICS: COMBINING NETWORK AND
LATENT VARIABLE MODELS
Abstract
We introduce the network model as a formal psychometric model,
conceptualizing the covariance between psychometric indicators as resulting
from pairwise interactions between observable variables in a network
structure. This contrasts with standard psychometric models, in which the
covariance between test items arises from the influence of one or more
common latent variables. Here, we present two generalizations of the
network model that encompass latent variable structures, establishing
network modeling as parts of the more general framework of Structural
Equation Modeling (SEM). In the first generalization, we model the
covariance structure of latent variables as a network. We term this
framework Latent Network Modeling (LNM) and show that, with LNM, a
unique structure of conditional independence relationships between latent
variables can be obtained in an explorative manner. In the second
generalization, the residual variance-covariance structure of indicators is
modeled as a network. We term this generalization Residual Network
Modeling (RNM) and show that, within this framework, identifiable models
can be obtained in which local independence is structurally violated. Both
generalizations have been implemented in the free-to-use software package
lvnet, which contains confirmatory model testing as well as two exploratory
search algorithms: stepwise search algorithms for low-dimensional datasets
and penalized maximum likelihood estimation for larger datasets. We show

Psychometrika Submission

September 13, 2017

in simulation studies that these search algorithms performs adequately in
identifying the structure of the relevant residual or latent networks. We
further demonstrate the utility of these generalizations in an empirical
example on a personality inventory dataset.
Key words: networks, structural equation modeling, simulation study

3

Psychometrika Submission

September 13, 2017

4

Introduction
Recent years have seen an emergence of network modeling in psychometrics
(Borsboom, 2008; Borsboom et al., 2011; Cramer et al., 2012; Schmittmann et al.,
2013). In this paper, we introduce the reader to this field of network psychometrics
(Epskamp et al., in press) and formalize the network model for multivariate normal
data, the Gaussian Graphical Model (GGM; Lauritzen 1996), as a formal psychometric
model. We contrast the GGM to the Structural Equation Model (SEM; Wright 1921;
Kaplan 2000) and show that the GGM can be seen as another way to approach
modeling covariance structures as is typically done in psychometrics. In particular,
rather than modeling the covariance matrix, the GGM models the inverse of a
covariance matrix. The GGM and SEM are thus very closely related: every GGM
model and every SEM model imply a constrained covariance structure. We make use of
this relationship to show that, through a reparameterization of the SEM model, the
GGM model can be obtained in two different ways: first, as a network structure that
relates a number of latent variables to each other, and second, as a network between
residuals that remain given a fitted latent variable model. As such, the GGM can be
modeled and estimated in SEM, which allows for network modeling of psychometric
data to be carried out in a framework familiar to psychometricians and methodologists.
However, the combination of GGM and SEM allows for more than fitting network
models. As we will show, the strength of one framework can help the other framework
and vice versa. In this paper, we introduce network models for latent covariances and
for residual covariances as two distinct generalized frameworks of both the SEM and
GGM. The first framework, Latent Network Modeling (LNM), formulates a network
among latent variables. This framework allows researchers to exploratively estimate
conditional independence relationships between latent variables through model search
algorithms; this estimation is difficult in the SEM framework due to the presence of
equivalent models (MacCallum et al., 1993). The second framework, which we denote
Residual Network Modeling (RNM), formulates a network structure on the residuals of
a SEM model. With this framework, researchers can circumvent critical assumptions of

Psychometrika Submission

September 13, 2017

5

both SEM and the GGM: SEM typically relies on the assumption of local
independence, whereas network modeling typically relies on the assumption that the
covariance structure among a set of the items is not due to latent variables at all. The
RNM framework allows researchers to estimate SEM models without the assumption of
local independence (all residuals can be correlated, albeit due to a constrained
structure on the inverse residual covariance matrix) as well as to estimate a network
structure, while taking into account the fact that the covariance between items may be
partly due to latent factors.
While the powerful combination of SEM and GGM allows for confirmative testing
of network structures both with and without latent variables, we recognize that few
researchers have yet formulated strict confirmatory hypotheses in the relatively new
field of network psychometrics. Often, researchers are more interested in exploratively
searching a plausible network structure. To this end, we present two exploratory search
algorithms. The first is a step-wise model search algorithm that adds and removes
edges of a network as long as fit is improved, and the second uses penalized maximum
likelihood estimation (Tibshirani, 1996) to estimate a sparse model. We evaluate the
performance of these search methods in four simulation studies. Finally, the proposed
methods have been implemented in a free-to-use R package, lvnet, which we illustrate
in an empirical example on personality inventory items (Revelle, 2015).
Modeling Multivariate Gaussian Data
Let y be the response vector of a random subject on P items1 . We assume y is
centered and follows a multivariate Gaussian density:
y ∼ NP (00, Σ ) ,
1

Throughout this paper, vectors will be represented with lowercase boldfaced letters and matrices

will be denoted by capital boldfaced letters. Roman letters will be used to denote observed variables
and parameters (such as the number of nodes) and Greek letters will be used to denote latent variables
and parameters that need to be estimated. The subscript i will be used to denote the realized response
vector of subject i and omission of this subscript will be used to denote the response of a random subject.

Psychometrika Submission

September 13, 2017

6

Σ.
In which Σ is a P × P variance-covariance matrix, estimated by some model-implied Σ̂
Σ is often done through some form of maximum likelihood estimation. If we
Estimating Σ̂
measure N independent samples of y we can formulate the N × P matrix Y containing
realization y >
i as its ith row. Let S represent the sample variance-covariance matrix of

Y:
S=

1
Y >Y .
N −1

In maximum likelihood estimation, we use S to compute and minimize −2 times the the
Σ (Lawley, 1940; Jöreskog, 1967; Jacobucci et al., 2016):
log-likelihood function to find Σ̂
h
 
 −1 
 
i
Σ + Trace S Σ̂
Σ
S −P .
min log det Σ̂
− log det Ŝ

(1)

Σ
Σ̂

Σ should be estimated as closely as possible to S and
To optimize this expression, Σ̂
Σ = S . A properly identified model with the same number of
perfect fit is obtained if Σ̂
Σ as there are unique elements in S (P (P + 1)/2
parameters (K) used to form Σ̂
Σ = S and therefore a saturated model. The goal of modeling
parameters) will lead to Σ̂
Σ with positive degrees of
multivariate Gaussian data is to obtain some model for Σ̂
Σ resembles S closely.
freedom, K < P (P + 1)/2, in which Σ̂
Structural Equation Modeling
In Confirmatory Factor Analysis (CFA), Y is typically assumed to be a causal
linear effect of a set of M centered latent variables, η , and independent residuals or
error, ε :
y = Λη + ε .
Here, Λ represents a P × M matrix of factor loadings. This model implies the following
Σ:
model for Σ̂
Σ = Λ ΨΛ > + Θ ,
Σ̂

(2)

in which Ψ = Var (ηη ) and Θ = Var (εε). In Structural Equation Modeling (SEM),
Var (ηη ) can further be modeled by adding structural linear relations between the latent

Psychometrika Submission

September 13, 2017

7

variables2 :
η = Bη + ζ ,
in which ζ is a vector of residuals and B is an M × M matrix of regression coefficients.
Σ can be more extensively modeled as:
Now, Σ̂
Σ = Λ (II − B )−1 Ψ (II − B )−1> Λ> + Θ,
Σ̂

(3)

in which now Ψ = Var (ζζ ). This framework can be used to model direct causal effects
between observed variables by setting Λ = I and Θ = O , which is often called path
analysis (Wright, 1934).
Σ and S , a P × P matrix; if Θ is fully estimated—contains
The Θ matrix is, like Σ̂
no restricted elements—then Θ alone constitutes a saturated model. Therefore, to
make either (2) or (3) identifiable, Θ must be strongly restricted. Typically, Θ is set to
be diagonal, a restriction often termed local independence (Lord et al., 1968; Holland
and Rosenbaum, 1986) because indicators are independent of each other after
conditioning on the set of latent variables. To improve fit, select off-diagonal elements
of Θ can be estimated, but systematic violations of local independence—many nonzero
elements in Θ —are not possible as that will quickly make (2) and (3) saturated or even
over-identified. More precisely, Θ can not be fully-populated—some elements of Θ
must be set to equal zero—when latent variables are used. An element of Θ being fixed
to zero indicates that two variables are locally independent after conditioning on the
set of latent variables. As such, local independence is a critical assumption in both
CFA and SEM; if local independence is systematically violated, CFA and SEM will
never result in correct models.
The assumption of local independence has led to critiques of the factor model and
its usage in psychology; local independence appears to be frequently violated due to
direct causal effects, semantic overlap, or reciprocal interactions between putative
indicators of a latent variable (Borsboom, 2008; Cramer et al., 2010; Borsboom et al.,
2

We make use here of the convenient all-y notation and do not distinguish between exogenous and

endogenous latent variables (Hayduk, 1987).

Psychometrika Submission

September 13, 2017

8

2011; Cramer et al., 2012; Schmittmann et al., 2013). For example, in psychopathology
research, local independence of symptoms given a person’s level of a latent mental
disorder has been questioned (Borsboom and Cramer, 2013). For example, three
problems associated with depression are “fatigue”, “concentration problems” and
“rumination”. It is plausible that a person who suffers from fatigue will also
concentrate more poorly, as a direct result of being fatigued and regardless of his or her
level of depression. Similarly, rumination might lead to poor concentration. In another
example, Kossakowski et al. (2015) describe the often-used SF-36 questionnaire (Ware
Jr and Sherbourne, 1992) designed to measure health related quality of life. The SF-36
contains items such as “can you walk for more than one kilometer” and “can you walk
a few hundred meters”. Clearly, these items can never be locally independent after
conditioning on any latent trait, as one item (the ability to walk a few hundred meters)
is a prerequisite for the other (walking more than a kilometer). In typical applications,
the excessive covariance between items of this type is typically left unmodeled, and
treated instead by combining items into a subscale or total score that is subsequently
subjected to factor analysis; of course, however, this is tantamount to ignoring the
relevant psychometric problem rather than solving it.
Given the many theoretically expected violations of local independence in
psychometric applications, many elements of Θ in both (2) and (3) should ordinarily
be freely estimated. Especially when violations of local independence are expected to
be due to causal effects of partial overlap, residual correlations should not be
constrained to zero; in addition, a chain of causal relationships between indicators can
lead to all residuals to become correlated. Thus, even when latent factors cause much
of the covariation between measured items, fitting a latent variable model that involves
local independence may not fully account for correlation structure between measured
items. Of course, in practice, many psychometricians are aware of this problem, which
is typically addressed by freeing up correlations between residuals to improve model fit.
However, this is usually done in an ad-hoc fashion, on the basis of inspection of
modification indices and freeing up error covariances one by one, which is post hoc,

Psychometrika Submission

September 13, 2017

9

suboptimal, and involves an uncontrolled journey through the model space. As a
result, it is often difficult to impossible to tell how exactly authors arrived at their final
reported models. As we will show later in this paper, this process can be optimized
and systematized using network models to connect residuals on top of a latent variable
structure.
Network Modeling
Recent authors have suggested that the potential presence of causal relationships
between measured variables may allow the explanation of the covariance structure
without the need to invoke any latent variables (Borsboom, 2008; Cramer et al., 2010;
Borsboom et al., 2011; Schmittmann et al., 2013). The interactions between indicators
can instead be modeled as a network, in which indicators are represented as nodes that
are connected by edges representing pairwise interactions. Such interactions indicate
the presence of covariances that cannot be explained by any other variable in the model
and can represent—possibly reciprocal—causal relationships. Estimating a network
structure on psychometric data is termed network psychometrics (Epskamp et al., in
press). Such a network of interacting components can generate data that fit factor
models well, as is commonly the case in psychology. Van der Maas et al. (2006) showed
that the positive manifold of intelligence—which is commonly explained with the
general factor for intelligence, g—can emerge from a network of mutually benefiting
cognitive abilities. Borsboom et al. (2011) showed that a network of psychopathological
symptoms, in which disorders are modeled as clusters of symptoms, could explain
comorbidity between disorders. Furthermore, Epskamp et al. (in press) showed that the
Ising model for ferromagnetism (Ising, 1925), which models magnetism as a network of
particles, is equivalent to multidimensional item response theory (Reckase, 2009).
In network psychometrics, psychometric data are modeled through directed or
undirected networks. Directed networks are equivalent to path analysis models. For
modeling undirected networks, pairwise Markov Random Fields (Lauritzen, 1996;
Murphy, 2012) are used. In these models, each variable is represented by a node, and

Psychometrika Submission

September 13, 2017

10

Y2

Y1

Y3

Figure 1: Example of a pairwise Markov Random Field model. Edges in this model
indicate pairwise interactions, and are drawn using undirected edges to distinguish from
(bidirectional) covariances. Rather than a model for marginal associations (such as a
network indicating covariances), this is a model for conditional associations. The network
above encodes that Y1 and Y3 are independent after conditioning on Y2 . Such a model
allows all three variables to correlate while retaining one degree of freedom (the model
only has two parameters)

Psychometrika Submission

September 13, 2017

11

nodes are connected by a set of edges. If two nodes, yj and yk , are not connected by an
edge, then this means they are independent after conditioning on the set of all other
nodes, y −(j,k) . Whenever two nodes cannot be rendered independent conditional on the
other nodes in the system, they are said to feature in a pairwise interaction, which is
represented by an undirected edge—an edge with no arrows—to contrast such an effect
from covariances typically represented in the SEM literature with bidirectional edges.
Figure 1 represents such a network model, in which nodes y1 and y3 are independent
after conditioning on node y2 . Such a model can readily arise from direct interactions
between the nodes. For example, this conditional independence structure would emerge
if y2 is a common cause of y1 and y3 , or if y2 is the mediator in a causal path between
y1 and y3 . In general, it is important to note that pairwise interactions are not mere
correlations; two variables may be strongly correlated but unconnected (e.g., when
both are caused by another variable in the system) and they may be uncorrelated but
strongly connected in the network (e.g., when they have a common effect in the
system). For instance, in the present example the model does not indicate that y1 and
y3 are uncorrelated, but merely indicates that any correlation between y1 and y3 is due
to their mutual interaction with y2 ; a network model in which either directly or
indirectly connected paths exist between all pairs of nodes typically implies a fully
populated (no zero elements) variance-covariance matrix.
In the case of multivariate Gaussian data this model is termed the Gaussian
Graphical Model (GGM; Lauritzen 1996). In the case of multivariate normality, the
partial correlation coefficient is sufficient to test the degree of conditional independence
of two variables after conditioning on all other variables; if the partial correlation
coefficient is zero, there is conditional independence and hence no edge in the network.
As such, partial correlation coefficients can directly be used in the network as edge
weights; the strength of connection between two nodes3 . Such a network is typically
encoded in a symmetrical and real valued p × p weight matrix, Ω , in which element ωjk
3

A saturated GGM is also called a partial correlation network because it contains the sample partial

correlation coefficients as edge weights.

Psychometrika Submission

September 13, 2017

12

represents the edge weight between node j and node k:


Cor yj , yk | y −(j,k) = ωjk = ωkj .
The partial correlation coefficients can be directly obtained from the inverse of
Σ, also termed the precision matrix K̂
K (Lauritzen, 1996):
variance-covariance matrix Σ̂


κjk
Cor yj , yk | y −(j,k) = − √ √ .
κkk κjj
Thus, element κjk of the precision matrix is proportional to the partial correlation
coefficient of variables yj and yk after conditioning on all other variables. Since this
process simply involves standardizing the precision matrix, we propose the following
model4 :
Σ = K̂
K
Σ̂

−1

= ∆ (II − Ω )−1 ∆ ,

(4)

−1

in which ∆ is a diagonal matrix with δjj = κjj 2 and Ω has zeroes on the diagonal. This
model allows for confirmative testing of the GGM structures on psychometric data.
Furthermore, the model can be compared to a saturated model (fully populated
Ω = O ), allowing one to obtain
off-diagonal values of Ω ) and the independence model (Ω
χ2 fit statistics as well as fit indices such as the RMSEA (Browne and Cudeck, 1992)
and CFI (Bentler, 1990). Such methods of assessing model fit have not yet been used
in network psychometrics.
Similar to CFA and SEM, the GGM relies on a critical assumption; namely, that
covariances between observed variables are not caused by any latent or unobserved
variable. If we estimate a GGM in a case where, in fact, a latent factor model was the
true data generating structure, then generally we would expect the GGM to be
saturated—i.e., there would be no missing edges in the GGM (Chandrasekaran et al.,
2010). A missing edge in the GGM indicates the presence of conditional independence
between two indicators given all other indicators; we do not expect indicators to
become independent given subsets of other indicators (see also Ellis and Junker 1997;
4

To our knowledge, the GGM has not yet been framed in this form. We chose this form because it

allows for clear modeling and interpretation of the network parameters.

Psychometrika Submission

September 13, 2017

13

Holland and Rosenbaum 1986). Again, this critical assumption might not be plausible.
While variables such as “Am indifferent to the feelings of others” and “Inquire about
others’ well-being” quite probably interact with each other, it might be far-fetched to
assume that no unobserved variable, such as a personality trait, in part also causes
some of the variance in responses on these items.
Generalizing Factor Analysis and Network Modeling
We propose two generalizations of both SEM and the GGM that both allow the
modeling of network structures in SEM. In the first generalization, we adopt the CFA5
decomposition in (2) and model the variance-covariance matrix of latent variables as a
GGM:
Ψ = ∆Ψ (II − ΩΨ )−1 ∆Ψ .
This framework can be seen as modeling conditional independencies between latent
variables not by directed effects (as in SEM) but as an undirected network. As such,
we term this framework latent network modeling (LNM).
In the second generalization, we adopt the SEM decomposition of the
variance-covariance matrix in (3) and allow the residual variance-covariance matrix Θ
to be modeled as a GGM:
Θ = ∆Θ (II − ΩΘ )−1 ∆Θ .
Because this framework conceptualizes associations between residuals as pairwise
interactions, rather than correlations, we term this framework Residual Network
Modeling (RNM). Using this framework allows—as will be described below—for a
powerful way of fitting a confirmatory factor structure even though local independence
is systematically violated and all residuals are correlated.
5

We use the CFA framework instead of the SEM framework here as the main application of this

framework is in exploratively estimating relationships between latent variables.

Psychometrika Submission

September 13, 2017

(a) Structural Equation Modeling

(b) Network Modeling

(c) Latent Network Modeling

(d) Residual Network Modeling

14

Figure 2: Examples of possible models under four different modeling frameworks. Circular nodes indicate latent variables, square nodes indicate manifest variables and gray
nodes indicate residuals. Directed edges indicate factor loadings or regression parameters and undirected edges indicate pairwise interactions. Note that such undirected
edges do not indicate covariances, which are typically denoted with bidirectional edges.
Replacing covariances with interactions is where the network models differ from typical
SEM.

Psychometrika Submission

September 13, 2017

15

Figure 2 shows four different examples of possible models that are attainable
under the SEM, LNM and RNM frameworks. Panel 2a shows a typical SEM model in
which one latent variable functions as a common cause of two others. Panel 2b shows a
network model that can be estimated using both the RNM and the LNM frameworks.
Panel 2c shows a completely equivalent LNM model to the SEM model of Panel 2a in
which the direction of effect between latent variables is not modeled. Finally, panel 2d
shows a model in which three exogenous latent variables underlie a set of indicators of
which the residuals form a network. The remainder of this section will describe RNM
and LNM in more detail and will outline the class of situations in which using these
models is advantageous over CFA or SEM.
Latent Network Modeling
The LNM framework models the latent variance-covariance matrix of a CFA
model as a GGM:
Σ = Λ∆Ψ (II − ΩΨ )−1 ∆Ψ Λ > + Θ .
Σ̂

(5)

This allows researchers to model conditional independence relationships between latent
variables without making the implicit assumptions of directionality or acyclicness. In
SEM, B is typically modeled as a directed acyclic graph (DAG), meaning that
elements of B can be represented by directed edges, and following along the path of
these edges it is not possible to return to any node (latent variable). The edges in such
a DAG can be interpreted as causal, and in general they imply a specific set of
conditional independence relationships between the nodes (Pearl, 2000).
While modeling conditional independence relationships between latent variables as
a DAG is a powerful tool for testing strictly confirmatory hypotheses, it is less useful
for more exploratory estimation. Though there have been recent advances in
exploratory estimation of DAGs within an SEM framework (e.g., Gates and Molenaar
2012; Rosa et al. 2012), many equivalent DAGs can imply the same conditional
independence relationships, and thus fit the data equally well even though their causal
interpretation can be strikingly different (MacCallum et al., 1993). Furthermore, the

Psychometrika Submission

Y2

September 13, 2017

Y2

16

Y2
No equivalent GGM with
the same # of parameters

Y3

Y1

Y3

Y1

(a)

(b)

Y2

Y2

Y3

Y1

(c)

Y2
No equivalent GGM with
the same # of parameters

Y3

Y1

(e)

Y3

Y1

(d)

(f)

Y3

Y1
Y4

(g)

(h)

Figure 3: Equivalent models between directed acyclic graphs (DAG; top) and Gaussian
Graphical Models (GGM; bottom). Each column of graphs show two models that are
equivalent. Panels (a), (b), (e) and (f) all represent the same conditional independence
structure: Y1 and Y3 are independent after conditioning on Y2 . Panel (c) represents that
Y1 and Y3 are marginally independent even though they both are correlated with Y2 , a
structure that cannot be represented in a GGM with only two edges. Panel (h) shows
that Y1 and Y3 are independent after conditioning on the set Y2 and Y4 , which cannot
be represented with a DAG of four edges.

Psychometrika Submission

September 13, 2017

17

assumption that the generating model is acyclic—which, in practice, often is made on
purely pragmatic grounds to identify a model—is problematic in that much
psychological behavior can be assumed to have at least some cyclic and complex
behavior and feedback (Schmittmann et al., 2013). Thus, the true conditional
independence relationships in a dataset can lead to many equivalent compositions of B ,
and possibly none of them are the true model.
In psychometrics and SEM the GGM representation has not been very prominent,
even though it has some manifest benefits over the attempt to identify DAGs directly.
For example, by modeling conditional independence relationships between latent
variables as a GGM, many relationships can be modeled in a simpler way as compared
to a DAG. In addition, in the GGM each set of conditional independence relations only
corresponds to one model: there are no equivalent GGMs with the same nodes. Figure
3 shows a comparison of several conditional independence relations that can be
modeled equivalently or not by using a GGM or by using a DAG. Figures 3a and 3b
show two DAGs that represent the same conditional independence relations,
y1 ⊥⊥ y3 | y2 , which can both be represented by the same GGM shown in Figure 3e.
There are some conditional independence relations that a GGM cannot represent in
the same number of parameters as a DAG; Figure 3c shows a collider structure that
cannot be exactly represented by a GGM (the best fitting GGM would feature three
edges instead of two). On the other hand, there are also conditional independence
relationships that a GGM can represent and a DAG cannot; the cycle of Figure 3d
cannot be represented by a DAG. Further equivalences and differences between GGMs
and DAGs are beyond the scope of this paper, but haven been well described in the
literature (e.g., chapter 3 of Lauritzen 1996; Koller and Friedman 2009; Kolaczyk and
Csárdi 2014). In sum, the GGM offers a natural middle ground between zero-order
correlations and DAGs: every set of zero-order correlations implies exactly one GGM,
and every DAG implies exactly one GGM. In a sense, the road from correlations to
DAGs (including hierarchical factor models) thus always must pass through the realm
of GGMs, which acts as a bridge between the correlational and causal worlds.

Psychometrika Submission

September 13, 2017

18

Because there are no equivalent undirected models possible, LNM offers a powerful
tool for exploratory estimation of relationships between latent variables. For example,
suppose one encounters data generated by the SEM model in Figure 2a. Without prior
theory on the relations between latent variables, exploratory estimation on this dataset
would lead to three completely equivalent models: the one shown in Figure 2c and two
models in which the common cause instead is the middle node in a causal chain. As
the number of latent variables increases, the potential number of equivalent models
that encode the same conditional independence relationships grows without bound.
The LNM model in Figure 2c portrays the same conditional independence relationship
as the SEM model in 2a while having no equivalent model. Exploratory estimation
could easily find this model, and portrays the retrieved relationship in a clear and
unambiguous way.
A final benefit of using LNM models is that they allow network analysts to
construct a network while taking measurement error into account. So far, networks
have been constructed based on single indicators only and no attempt has been made
to remediate measurement error. By forming a network on graspable small concepts
measured by a few indicators, the LNM framework can be used to control for
measurement error.
Residual Network Modeling
In the RNM framework the residual structure of SEM is modeled as a GGM:
Σ = Λ (II − B )−1 Ψ (II − B )−1> Λ > + ∆Θ (II − ΩΘ )−1 ∆Θ .
Σ̂

(6)

This modeling framework conceptualizes latent variable and network modeling as two
sides of the same coin, and offers immediate benefits to both. In latent variable
modeling, RNM allows for the estimation of a factor structure (possibly including
structural relationships between the latent variables), while having no uncorrelated
errors and thus no local independence. The error-correlations, however, are still highly
structured due to the residual network structure. This can be seen as a compromise

Psychometrika Submission

September 13, 2017

19

between the ideas of network analysis and factor modeling; while we agree that local
independence is plausibly violated in many psychometric tests, we think the
assumption of no underlying latent traits and therefore a sparse GGM may often be
too strict. For network modeling, RNM allows a researcher to estimate a sparse
network structure while taking into account that some of the covariation between items
was caused by a set of latent variables. Not taking this into account would lead to a
saturated model (Chandrasekaran et al., 2010), whereas the residual network structure
can be sparse.
To avoid confusion between residual correlations, we will denote edges in ΩΘ
residual interactions. Residual interactions can be understood as pairwise linear effects,
possibly due to some causal influence or partial overlap between indicators that is left
after controlling for the latent structure. Consider again the indicators for
agreeableness “Am indifferent to the feelings of others” and “Inquire about others’
well-being”. It seems clear that we would not expect these indicators to be locally
independent after conditioning on agreeableness; being indifferent to the feelings of
others will cause one to not inquire about other’s well-being. Thus, we could expect
these indicators to feature a residual interaction; some degree of correlation between
these indicators is expected to remain, even after conditioning on the latent variable
and all other indicators in the model.
The RNM framework in particular offers a new way of improving the fit of
confirmatory factor models. In contrast to increasingly popular methods such as
exploratory SEM (ESEM; Marsh et al. 2014) or LASSO regularized SEM models
(Jacobucci et al., 2016), the RNM framework improves the fit by adding residual
interactions rather than allowing for more cross-loadings. The factor structure is kept
exactly intact as specified in the confirmatory model. Importantly, therefore, the
interpretation of the latent factor does not change. This can be highly valuable in the
presence of a strong theory on the latent variables structure underlying a dataset even
in the presence of violations of local independence.

Psychometrika Submission

September 13, 2017

20

Exploratory Network Estimation
Both the LNM and RNM modeling frameworks allow for confirmative testing of
network structures. Confirmatory estimation is straightforward and similar to
estimating SEM models, with the exception that instead of modeling Ψ or Θ now the
latent network ΩΨ or ΩΘ is modeled. Furthermore, both modeling frameworks allow
for the confirmatory fit of a network model. In LNM, a confirmatory network structure
can be tested by setting Λ = I and Θ = O ; in RNM, a confirmatory network model can
be tested by omitting any latent variables. We have developed the R package lvnet6 ,
which utilizes OpenMx (Neale et al., 2015) for confirmative testing of RNM and LNM
models (as well as a combination of the two). The lvnet function can be used for this
purpose by specifying the fixed and the free elements of model matrices. The package
returns model fit indices (e.g., the RMSEA, CFI and χ2 value), parameter estimates,
and allows for model comparison tests.
Often the network structure, either at the residual or the latent level, is unknown
and needs to be estimated. To this end, the package includes two exploratory search
algorithms described below: step-wise model search and penalized maximum likelihood
estimation. For both model frameworks and both search algorithms, we present
simulation studies to investigate the performance of these procedures. As is typical in
simulation studies investigating the performance of network estimation techniques, we
investigated the sensitivity and specificity (van Borkulo et al., 2014). These measures
investigate the estimated edges versus the edges in the true model, with a ‘positive’
indicating an estimated edge and a ‘negative’ indicating an edge that is estimated to
be zero. Sensitivity, also termed the true positive rate, gives the ratio of the number of
true edges that were detected in the estimation versus the total number of edges in the
true model:
sensitivity =

# true positives
# true positives + # of false negatives

Specificity, also termed the true negative rate, gives the ratio of true missing edges
6

github.com/sachaepskamp/lvnet

Psychometrika Submission

September 13, 2017

21

detected in the estimation versus the total number of absent edges in the true model:
specificity =

# true negatives
# true negatives + # false positives

The specificity can be seen as a function of the number of false positives: a high
specificity indicates that there were not many edges detected to be nonzero that are
zero in the true model. To favor degrees of freedom, model sparsity and interpretability,
specificity should be high all-around—estimation techniques should not result in many
false positives—whereas sensitivity should increase as a function of the sample size.
Simulating Gaussian Graphical models
In all simulation studies reported here, networks were constructed in the same way
as done by Yin and Li (2011) in order to obtain a positive definite inverse-covariance
matrix K . First, a network structure was generated without weights. Next, weights
were drawn randomly from a uniform distribution between 0.5 and 1, and made
negative with 50% probability. The diagonal elements of K were then set to 1.5 times
the sum of all absolute values in the corresponding row, or 1 if this sum was zero.
Next, all values in each row were divided by the diagonal value, ensuring that the
diagonal values become 1. Finally, the matrix was made symmetric by averaging the
lower and upper triangular elements. In the chain graphs used in the following
simulations, this algorithm created networks in which the non-zero partial correlations
had a mean of 0.33 and a standard deviation of 0.04.
Stepwise Model Search
In exploratory search, we are interested in recovering the network structure of
either ΩΨ in LNM or ΩΘ in RNM. This can be done through a step-wise model search,
either based on χ2 difference tests (Algorithm 1) or on minimization of some
information criterion (Algorithm 2) such as the Akaike information criterion (AIC),
Bayesian information criterion (BIC) or the extended Bayesian information criterion
(EBIC; Chen and Chen 2008) which is now often used in network estimation (van

Psychometrika Submission

September 13, 2017

22

Borkulo et al., 2014; Foygel and Drton, 2010). In LNM, removing edges from ΩΨ
cannot improve the fit beyond that of an already fitting CFA model. Hence, model
search for ΩΨ should start at a fully populated initial setup for ΩΨ . In RNM, on the
other hand, a densely populated ΩΘ would lead to an over-identified model, and hence
the step-wise model search should start at an empty network ΩΘ = O . The function
lvnetSearch in the lvnet package can be used for both search algorithms.
Algorithm 1 Stepwise network estimation by χ2 difference testing.
Start with initial setup for Ω
repeat
for all Unique elements of Ω do
Remove edge if present or add edge if absent
Fit model with changed edge
end for
if Adding an edge significantly improves fit (α = 0.05) then
Add edge that improves fit the most
else if Removing an edge does not significantly worsen fit (α = 0.05) then
Remove edge that worsens fit the least
end if
until No added edge significantly improves fit and removing any edge significantly
worsens fit

Psychometrika Submission

September 13, 2017

23

Algorithm 2 Stepwise network estimation by AIC, BIC or EBIC optimization.
Start with initial setup for Ω
repeat
for all Unique elements of Ω do
Remove edge if present or add edge if absent
Fit model with changed edge
end for
if Any changed edge improved AIC, BIC or EBIC then
Change edge that improved AIC, BIC or EBIC the most
end if
until No changed edge improves AIC, BIC or EBIC

Simulation Study 1: Latent Network Modeling
We performed a simulation study to assess the performance of the above
mentioned step-wise search algorithms in LNM models. Figure 4 shows the LNM
model under which we simulated data. In this model, four latent factors with three
indicators each were connected in a latent network. The latent network was a chain
network, leading all latent variables to be correlated according to a structure that
cannot be represented in SEM. Factor loadings and residual variances were set to 1,
and the network weights were simulated as described in the section “Simulating
Gaussian Graphical models”. The simulation study followed a 5 × 4 design: the sample
size was varied between 50, 100, 250, 500 and 1 000 to represent typical sample sizes in
psychological research, and the stepwise evaluation criterion was either χ2 difference
testing, AIC, BIC or EBIC (using a tuning parameter of 0.5). Each condition was
simulated 1 000 times, resulting in 20 000 total simulated datasets.
Figure 5 shows the results of the simulation study. Data is represented in standard
boxplots (McGill et al., 1978): the box shows the 25th, 50th (median) 75th quantiles,
the whiskers range from the largest values in 1.5 times the inter-quantile range (75th 25th quantile) and points indicate outliers outside that range. In each condition, we

Psychometrika Submission

September 13, 2017

24

Figure 4: Model used in simulation study 1: step-wise model search in latent network
modeling. Four latent variables were each measured by three items. Latent variables
covary due to the structure of a latent Gaussian graphical model in which edges indicate
partial correlation coefficients. This model has the form of a chain graph, which cannot
be represented in a structural equation model. Factor loadings, residual variances and
latent variances were set to 1 and the latent partial correlations had an average of 0.33
with a standard deviation of 0.04.

Psychometrika Submission

September 13, 2017

25

Sensitivity (true positive rate)
AIC

Chi−square

BIC

EBIC

Sensitivity

1.00

●

0.75

●

●

0.50

●

●

0.25

0.00

●

●

●

●

●

●

●

●

●

●

50

●

●

100 250 500 1000

50

100 250 500 1000

50

●

100 250 500 1000

50

100 250 500 1000

Sample Size

Specificity (true negative rate)
AIC

Chi−square

BIC

EBIC

1.00

Specificity

0.75

0.50

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

0.25

0.00
50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

Sample Size

Figure 5: Simulation results of simulation study 1: step-wise model search in latent
network modeling. Each condition was replicated 1 000 times, leading to 20 000 total
simulated datasets. High sensitivity indicates that the method is able to detect edges
in the true model, and high specificity indicates that the method does not detect edges
that are zero in the true model.

Psychometrika Submission

September 13, 2017

26

investigated the sensitivity and specificity. The top panel shows that sensitivity
improves with sample size, with AIC performing best and EBIC worst. From sample
sizes of 500 and higher all estimation criterion performed well in retrieving the edges.
The bottom panel shows that specificity is generally very high, with EBIC performing
best and AIC worst. These results indicate that the step-wise procedure is conservative
and prefers simpler models to more complex models; missing edges are adequately
detected but present edges in the true model might go unnoticed except in larger
samples. With sample sizes over 500, all four estimation methods show both a high
sensitivity and specificity.
Simulation Study 2: Residual Network Modeling
We conducted a second simulation study to assess the performance of step-wise
model selection in RNM models. Figure 7 shows the model under which data were
simulated: two latent variables with 5 indicators each. The residual network was
constructed to be a chain graph linking a residual of an indicator of one latent variable
to two indicators of the other latent variable. This structure cannot be represented by
a DAG and causes all residuals to be connected, so that Θ is fully populated. Factor
loadings and residual variances were set to 1, the factor covariance was set to 0.25, and
the network weights were simulated as described in the section “Simulating Gaussian
Graphical models”.
The simulation study followed a 5 × 4 design; sample size was again varied

between 50, 100, 250, 500 and 1 000, and models were estimated using either χ2
significance testing, AIC, BIC or EBIC. Factor loadings and factor variances were set
to 1 and the factor correlation was set to 0.25. The weights in ΩΘ were chosen as
described in the section “Simulating Gaussian Graphical models”. Each condition was
simulated 1, 000 times, leading to 20 000 total datasets.
Figure 7 shows the results of the simulation study. The top panel shows that
sensitivity increases with sample size and performs best when using AIC as the
criterion. BIC performed comparably in sensitivity to χ2 testing and EBIC performed

Psychometrika Submission

September 13, 2017

27

Figure 6: Model used in simulation study 2: step-wise model search in residual network
modeling. Two latent variables were each measured by five items; a Gaussian graphical
model, in which edges indicate partial correlation coefficients, leads to all residuals to
be correlated due to a chain graph between residuals, which cannot be represented in a
structural equation model. Factor loadings, residual variances and latent variances were
set to 1, the factor covariance was set to 0.25 and the latent partial correlations had an
average of 0.33 with a standard deviation of 0.04.

Psychometrika Submission

September 13, 2017

28

Sensitivity (true positive rate)
AIC

Chi−square

BIC

1.00

EBIC

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

●

Sensitivity

0.75

0.50

0.25
●

●

●

●

●

●

0.00

●

50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

50

●

100 250 500 1000

Sample Size

Specificity (true negative rate)
AIC
1.00

Specificity

0.75

●

●
●
●
●

Chi−square

●

●
●
●

BIC

EBIC

●

●
●
●
●

●
●
●
●
●
●

●
●
●
●

●

●

●
●
●

●
●
●
●
●

●
●
●
●
●
●
●
●
●

●
●
●

●
●
●

●
●

●
●
●
●
●

●
●
●
●
●
●

●
●
●
●
●
●
●
●
●

●
●
●

●

●
●
●
●
●
●

●
●

●
●
●

●
●
●
●
●
●
●

●
●
●
●
●
●
●
●
●

●

●
●

0.50
●

●
●

0.25

0.00
50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

50

100 250 500 1000

Sample Size

Figure 7: Simulation results of simulation study 2: step-wise model search in residual
network modeling. Each condition was replicated 1 000 times, leading to 20 000 total
simulated datasets. High sensitivity indicates that the method is able to detect edges
in the true model, and high specificity indicates that the method does not detect edges
that are zero in the true model.

Psychometrika Submission

September 13, 2017

29

the worst. The bottom panel shows that specificity was very high for all sample sizes
and all criteria, with EBIC performing best and AIC worst. These results indicate that
the number of false positives is very low and that the method is on average well
capable of discovering true edges for sample size larger than 250. In sum, all four
criteria perform well with EBIC erring on the side of caution and AIC erring on the
side of discovery.
LASSO regularization
While the step-wise model selection algorithms perform well in retrieving the
correct network structure, they are very slow when the number of nodes in the network
increases (e.g., more than 10 nodes). This is particularly important in the context of
RNM, in which the number of indicators can be larger than 10 even in small models. A
popular method for fast estimation of high-dimensional network structures is by
applying the least absolute shrinkage and selection operator (LASSO; Tibshirani 1996).
LASSO regularization has also recently been introduced in the SEM literature
(Jacobucci et al., 2016) as a method for obtaining sparser structures of Λ and B . In
the LASSO, instead of optimizing the likelihood function as described in (1), a
penalized likelihood is optimized (Jacobucci et al., 2016):
h
 
 −1 
 
i
Σ + Trace S Σ̂
Σ
S − P + νPenalty ,
min log det Σ̂
− log det Ŝ

(7)

Σ
Σ̂

in which ν denotes a tuning parameter controlling the level of penalization. The
penalty here is taken to be the sum of absolute parameters:
Penalty =

X

|ωij |,

<i,j>

in which ωij denotes an element from either ΩΨ or ΩΘ . Other penalty functions may
be used as well—such as summing the squares of parameter estimates (ridge regression;
Hoerl and Kennard 1970) or combining both absolute and squared values (elastic net;
Zou and Hastie 2005)—but these are not currently implemented in lvnet. The benefit
of the LASSO is that it returns models that perform better in cross-validation. In

Psychometrika Submission

September 13, 2017

30

addition, the LASSO yields sparse models in which many relationships are estimated
to be zero.
Algorithm 3 LASSO estimation for exploratory network search.
for all Sequence of tuning parameters ν1 , ν2 , . . . do
Estimate LASSO regularized model using given tuning parameter
Count the number of parameters for which the absolute estimate is larger than 
Determine information criterion AIC or BIC given fit and number of parameters
end for
Select model with best AIC, BIC or EBIC
Refit this model without LASSO in which absolute parameters smaller than  are fixed
to zero
The lvnet function allows for LASSO regularization for a given model matrix
ΩΘ , ΩΨ , Θ , Ψ , Λ or B ) and a given value for the tuning parameter ν. The optimizer
(Ω
used in lvnet does not return exact zeroes. To circumvent this issue, any absolute
parameter below some small value  (by default  = 0.0001) is treated as zero in
counting the number of parameters and degrees of freedom (Zou et al., 2007). The
lvnetLasso function implements the search algorithm described in Algorithm 3 to
automatically choose an appropriate tuning parameter, use that for model selection
and rerun the model to obtain a comparable fit to non-regularized models. In this
algorithm, a sequence of tuning parameters is tested, which is set by default to a
logorithmically spaced sequence of 20 values between 0.01 and 1.
Simulation Study 3: Latent Network Modeling
We studied the performance of LASSO penalization in estimating the latent
network structure in a similar simulation study to the study of the step-wise procedure
described above. Data were simulated under a similar model to the one shown in
Figure 4, except that now 8 latent variables were used leading to a total of 24 observed
variables. All parameter values were the same as in simulation study 1. The simulation

Psychometrika Submission

September 13, 2017

31

followed a 5 × 3 design. Sample size was varied between 100, 250, 500, 1 000 and 2 500,
and for each sample size 1 000 datasets were simulated leading to a total of 5 000
generated datasets. On these datasets the best model was selected using either AIC,
BIC or EBIC, leading to 15 000 total replications. In each replication, sensitivity and
specificity were computed. Figure 8 shows that AIC had a relatively poor specificity
all-around, but a high sensitivity. EBIC performed well with sample sizes of 500 and
higher.
Simulation Study 4: Residual Network Modeling
To assess the performance of LASSO in estimating the residual network structure
we simulated data as in Figure 6, except that in this case four latent variables were
used, each with 5 indicators, the residuals of which were linked via a chain graph. All
parameter values were the same as in simulation study 2. The design was the same as
in simulation study 3, leading to 5 000 generated datasets on which AIC, BIC or EBIC
were used to select the best model. While Figure 9 shows good performance of the
LASSO in retrieving the residual network structure and similar results as before: AIC
performs the worst in specificity and EBIC the best.
Empirical Example: Personality Inventory
In this section, we demonstrate LNM and RNM models by confirmative testing a
model and exploratively searching a residual and a latent network structure. We use
the lvnet package, which can be installed directly from Github using the devtools
package in R:
> library("devtools")
> install_github("sachaepskamp/lvnet")
> library("lvnet")
To exemplify the method, we will use a dataset from the psych package (Revelle,
2015) on the Big 5 personality traits (Benet-Martinez and John, 1998; Digman, 1989;

Psychometrika Submission

September 13, 2017

32

Sensitivity (true positive rate)
AIC

BIC

EBIC

1.00

●

●

●

Sensitivity

0.75

●

●

●

●

0.50

●

●

●

0.25

●

●

0.00

●

100

250

500

1000

2500

100

250

500

1000

2500

100

250

500

1000

2500

Sample Size

Specificity (true negative rate)
AIC

BIC

1.00

EBIC

●
●
●
●
●

0.75

●

Specificity

●

●
●

0.50

●

●

●

●

●

●

●

●

●

●

●

●

●

●
●

●

0.25
●

●

●

●

●

●

●

●

●

●

●

250

500

1000

●

0.00

●

100

2500

100

250

500

1000

2500

100

250

500

1000

2500

Sample Size

Figure 8: Simulation results of simulation study 3: model selection via penalized maximum likelihood estimation in latent network modeling. The same model as in Figure 4
was used except now with 4 latent variables leading to 24 observed variables. For each
sample size 1 000 datasets were generated, leading to 5 000 total simulated datasets on
which AIC, BIC or EBIC was used to select the best model. High sensitivity indicates
that the method is able to detect edges in the true model, and high specificity indicates
that the method does not detect edges that are zero in the true model.

Psychometrika Submission

September 13, 2017

33

Sensitivity (true positive rate)
AIC

BIC

EBIC

1.00

●
●

●

●

●

●

●
●
●

0.75

●

Sensitivity

●

●

●

●

●

●

●

0.50

●

●

●

●

●

●

●

●

●

0.25

●

●

●

●

●

●

●

●

●

0.00
100

250

500

1000

2500

100

250

500

1000

2500

100

250

500

1000

2500

Sample Size

Specificity (true negative rate)
AIC

BIC

EBIC

1.00
●
●

Specificity

0.75

●

●
●
●
●
●
●
●
●
●
●
●

●
●
●
●

●
●
●

●
●
●
●

●
●
●
●
●
●

●

●
●
●
●
●
●

●
●
●
●

250

500

1000

●
●
●
●

●
●
●
●
●
●
●
●
●
●
●
●
●
●

0.50

●
●
●
●
●
●

●
●
●
●
●
●
●

●
●
●
●
●

●

●

●

0.25

0.00
100

250

500

1000

2500

100

250

500

1000

2500

100

2500

Sample Size

Figure 9: Simulation results of simulation study 4: model selection via penalized maximum likelihood estimation in residual network modeling. The same model as in Figure
6 was used except now with 4 latent variables leading to 20 observed variables. For each
sample size 1 000 datasets were generated, leading to 5 000 total simulated datasets on
which AIC, BIC or EBIC was used to select the best model. High sensitivity indicates
that the method is able to detect edges in the true model, and high specificity indicates
that the method does not detect edges that are zero in the true model.

Psychometrika Submission

September 13, 2017

34

Goldberg, 2001, 1998; McCrae and Costa, 1997). This dataset consists of 2800
observations of 25 items designed to measure the 5 central personality traits with 5
items per trait. We estimated the CFA model on the BFI dataset. Next, we used
LASSO estimation to the RNM model using 100 different tuning parameters and using
EBIC as criterion to maximize specificity and search for a sparse residual network. The
fully correlated covariance matrix of latent variables is equivalent to a fully connected
latent network structure. Thus, after fitting a RNM model, we can again apply LASSO
to estimate a latent network in the resulting model, which we term an RNM+LNM
model. The R code used for this analysis can be found in the supplementary materials.
df

χ2

AIC

BIC

EBIC

RMSEA

TLI

CFI

CFA

265

4713.94

183233.67

183589.91

184542.39

0.08

0.75

0.78

RNM

172

806.63

179510.97

180419.39

182848.22

0.04

0.94

0.97

RNM+LNM

176

843.18

179539.52

180424.19

182789.53

0.04

0.94

0.97

Table 1: Fit measures for three models estimated on the BFI dataset in the psych R
package. CFA is the correlated five-factor model, RNM the same model but with a
residual network, and RNM+LNM the RNM model in which edges of the latent network
have been removed.
Table 1 shows the fit of the three models. The CFA model fits poorly. The RNM
model has substantively improved fit and resulted in good fit indices overall. The
estimated RNM+LNM showed that 5 edges could be removed from the latent network
after taking residual interactions into account. Figure 10 shows the factor structure
and residual network of the final RNM+LNM model. It can be seen that Agreeableness
is now only connected to extraversion: after taking into account someone’s level of
extraversion agreeableness is independent of the other three personality traits.
Extraversion is the most central node in this network and the only trait that is directly
linked to all other traits. The residual network shows many meaningful connections.
While seemingly densely connected, this network only has 30% of all possible edges in

Psychometrika Submission
Factor structure & latent network

September 13, 2017
Residual network

35

●
●

A1

O5

●

C2

A2

●

A3

O4

●

O3

E1

A4

A4
●

E2

O2

A5

O2

●
●

O4
O

A

●

E4

O1

C1

C1
N5

N5

N4

●

C2

N3
A5

N

N4

C

C3

●
●
●
●

C5

C4
O3

N3

●

O5

N2

O1

N1

●

C4

●

E

●

E5
N2

C5

●

A1

●

Agreeableness
A1: Am indifferent to the feelings of others.
A2: Inquire about others' well−being.
A3: Know how to comfort others.
A4: Love children.
A5: Make people feel at ease.
Conscientiousness
C1: Am exacting in my work.
C2: Continue until everything is perfect.
C3: Do things according to a plan.
C4: Do things in a half−way manner.
C5: Waste my time.
Extraversion
E1: Don't talk a lot.
E2: Find it difficult to approach others.
E3: Know how to captivate people.
E4: Make friends easily.
E5: Take charge.
Neuroticism
N1: Get angry easily.
N2: Get irritated easily.
N3: Have frequent mood swings.
N4: Often feel blue.
N5: Panic easily.

E3
N1

E1
E5

A3

E2
E4

E3

A2

C3

●
●
●
●
●

Openness
O1: Am full of ideas.
O2: Avoid difficult reading material.
O3: Carry the conversation to a higher level.
O4: Spend time reflecting on things.
O5: Will not probe deeply into a subject.

Figure 10: Visualization of the factor structure and latent network (left) and the residual network (right) of the BFI personality dataset from the psych package in R. LASSO
estimation with 100 different tuning parameters in combination with EBIC model selection was used to first estimate the residual network structure and following the latent
network structure.
a network of that size, leading the model to have 176 degrees of freedom. The
corresponding residual covariance structure is fully populated with no zero elements.
Conclusion
In this paper we introduced a formal psychometric model for network modeling of
multivariate normal data. We contrasted this model with latent variable models as
commonly used in CFA and SEM. Furthermore, using the CFA and SEM frameworks,
we proposed two generalizations of the network model to encompass latent variable
structures within the network paradigm. In the first generalization, LNM, we construct
a network among the latent variables, whereas in the second generalization, RNM, a
network is formed among residuals of indicators. Both frameworks offer powerful
benefits over both latent variable and network modeling. From the perspective of
latent variable modeling, the LNM framework allows one to exploratively search for
conditional independence relationships between latent variables without the need for

Psychometrika Submission

September 13, 2017

36

prior theory, and the RNM framework allows one to model latent common causes
without assuming local independence. From the perspective of network modeling, the
LNM framework allows one to model network structures while taking measurement
error into account, and the RNM framework allows one to estimate a network
structure, even when all nodes are in part caused by unobserved or latent variables.
Both frameworks have been implemented in the freely available R package lvnet.
Simulation studies showed that step-wise search and penalized maximum
likelihood estimation of the residual or latent network structures resulted in high
specificity all around—the methods did not result often in false positives—and rapidly
increasing sensitivity as a function of the sample size; the higher the sample size, the
more true edges were detected in the algorithm. These numbers are comparable to
state-of-the-art network estimation techniques in sample and model sizes that are
plausible in psychological settings (van Borkulo et al., 2014). In all four simulation
studies, using AIC as the model selection criterion led to the best sensitivity and using
EBIC led to the best specificity. However, it is important to note that the choice of a
particular information criterion cannot be argued by these numbers alone, and depends
on the relative importance one assigns to the side of discovery (optimizing sensitivity)
or the side of caution (optimizing specificity; Dziak et al. 2012). Furthermore, it should
be noted that these simulation results are specific for the particular model setup and
sample sizes used; results might be different for other kinds of models or sample size
ranges.
In addition to the LNM and RNM frameworks, other combinations of CFA, SEM
and network modeling are possible as well. For example, a framework can be
constructed which contains both a latent and a residual network (as shown in our
empirical example), or directed regression paths as in the SEM model can be added to
the LNM model. While these models are all estimable in the lvnet software, in the
current paper we chose to focus on the distinct benefits that modeling a residual or
latent network presents. Thus, in this paper, we only described the modeling of
multivariate normal data. More advanced models are possible, but not yet implemented

Psychometrika Submission

September 13, 2017

37

in the lvnet software. In the case of binary variables, the appropriate model to use is
the Ising model, which has been shown to be equivalent to multivariate item response
models (see Epskamp et al. in press). Future research could aim at constructing Ising
models among binary latent variables in latent class analysis, or constructing residual
networks in models with binary indicators. Finally, the expressions optimized in
Equations (1) and (7) are based on summary statistics and therefore only truly
applicable to complete data. With incomplete data, the appropriate estimation method
is to use full-information maximum likelihood (FIML; Arbuckle et al. 1996); however,
FIML has not yet been implemented in the lvnet software.
In our view, the presented modelling framework is a versatile and promising
addition to the spectrum of psychometric models. The GGM, which has a central place
in this modelling framework, acts as a natural interface between correlation and
causality, and we think this representation should receive more attention in
psychometrics. From the point of view afforded by the current paper, the typical
attempt to determine directed SEMs from correlation structures in fact appears
somewhat haphazard in psychology, a historical accident in a field that has been
prematurely directed to hypothesis testing at the expense of systematic exploration.
Perhaps, psychometrics as a field should consider taking a step back to focus on the
consistent identification of GGMs, instead of wanting to jump to the causal conclusion
immediately. In this regard, the fact that GGMs do not have equivalent models would
appear to be a major benefit, as they allow us to focus on charting connections
between variables systematically, without being forced to adhere to one particular
causal interpretation or another. In addition, because the GGM does not specify the
nature or direction of interactions between variables, it appears a natural model for
research situations where no temporal information or experimental interventions are
present, so that associations may arise for a multitude of reasons: the GGM can be
consistently interpreted regardless of whether associations arise from direct causal
relations, reciprocal causation, latent common causes, semantic overlap between items,
or homeostatic couplings of parameters. This appears to be a major asset in a field like

Psychometrika Submission

September 13, 2017

38

psychology, where strong causal theory is sparse and the identification of DAGs often
appears a bridge too far.

Psychometrika Submission

September 13, 2017

39

References
Arbuckle, J. L., Marcoulides, G. A., and Schumacker, R. E. (1996). Full information
estimation in the presence of incomplete data. Advanced structural equation
modeling: Issues and techniques, pages 243–277.
Benet-Martinez, V. and John, O. (1998). Los Cinco Grandes across cultures and ethnic
groups: Multitrait multimethod analyses of the Big Five in Spanish and English.
Journal of Personality and Social Psychology, 75:729–750.
Bentler, P. M. (1990). Comparative fit indexes in structural models. Psychological
bulletin, 107(2):238–346.
Borsboom, D. (2008). Psychometric perspectives on diagnostic systems. Journal of
clinical psychology, 64(9):1089–1108.
Borsboom, D. and Cramer, A. O. J. (2013). Network analysis: an integrative approach
to the structure of psychopathology. Annual review of clinical psychology, 9:91–121.
Borsboom, D., Cramer, A. O. J., Schmittmann, V. D., Epskamp, S., and Waldorp,
L. J. (2011). The small world of psychopathology. PloS one, 6(11):e27407.
Browne, M. W. and Cudeck, R. (1992). Alternative ways of assessing model fit.
Sociological Methods & Research, 21(2):230–258.
Chandrasekaran, V., Parrilo, P., Willsky, A. S., et al. (2010). Latent variable graphical
model selection via convex optimization. In Communication, Control, and
Computing (Allerton), 2010 48th Annual Allerton Conference on, pages 1610–1613.
IEEE.
Chen, J. and Chen, Z. (2008). Extended bayesian information criteria for model
selection with large model spaces. Biometrika, 95(3):759–771.
Cramer, A. O. J., Sluis, S., Noordhof, A., Wichers, M., Geschwind, N., Aggen, S. H.,
Kendler, K. S., and Borsboom, D. (2012). Dimensions of normal personality as

Psychometrika Submission

September 13, 2017

40

networks in search of equilibrium: You can’t like parties if you don’t like people.
European Journal of Personality, 26(4):414–431.
Cramer, A. O. J., Waldorp, L. J., van der Maas, H. L. J., and Borsboom, D. (2010).
Comorbidity: a network perspective. Behavioral and Brain Sciences,
33(2-3):137–150.
Digman, J. (1989). Five robust trait dimensions: Development, stability, and utility.
Journal of Personality, 57(2):195–214.
Dziak, J. J., Coffman, D. L., Lanza, S. T., and Li, R. (2012). Sensitivity and specificity
of information criteria. The Methodology Center and Department of Statistics, Penn
State, The Pennsylvania State University.
Ellis, J. L. and Junker, B. W. (1997). Tail-measurability in monotone latent variable
models. Psychometrika, 62(4):495–523.
Epskamp, S., Maris, G., Waldorp, L., and Borsboom, D. (in press). Network
psychometrics. In Irwing, P., Hughes, D., and Booth, T., editors, Handbook of
Psychometrics. Wiley, New York, NY, USA.
Foygel, R. and Drton, M. (2010). Extended Bayesian information criteria for Gaussian
graphical models. In Advances in Neural Information Processing Systems 23, pages
604–612.
Gates, K. M. and Molenaar, P. C. (2012). Group search algorithm recovers effective
connectivity maps for individuals in homogeneous and heterogeneous samples.
NeuroImage, 63(1):310–319.
Goldberg, L. (1998). The structure of phenotypic personality traits. Personality:
critical concepts in psychology, 48:26–34.
Goldberg, L. (2001). An alternative” description of personality”: The Big-Five factor
structure. The Science of Mental Health: Personality and personality disorder,
page 34.

Psychometrika Submission

September 13, 2017

41

Hayduk, L. A. (1987). Structural equation modeling with LISREL: Essentials and
advances. Johns Hopkins University Press, Baltimore, MD, USA.
Hoerl, A. E. and Kennard, R. W. (1970). Ridge regression: Biased estimation for
nonorthogonal problems. Technometrics, 12(1):55–67.
Holland, P. W. and Rosenbaum, P. R. (1986). Conditional association and
unidimensionality in monotone latent variable models. The Annals of Statistics,
14:1523–1543.
Ising, E. (1925). Beitrag zur theorie des ferromagnetismus. Zeitschrift für Physik A
Hadrons and Nuclei, 31(1):253–258.
Jacobucci, R., Grimm, K. J., and McArdle, J. J. (2016). Regularized structural
equation modeling. Structural Equation Modeling: A Multidisciplinary Journal,
pages 1–12.
Jöreskog, K. G. (1967). A general approach to confirmatory maximum likelihood factor
analysis. ETS Research Bulletin Series, 1967(2):183–202.
Kaplan, D. (2000). Structural equation modeling: Foundations and extensions. Sage,
Thousand Oaks, CA, USA.
Kolaczyk, E. D. and Csárdi, G. (2014). Statistical analysis of network data with R.
Springer, New York, NY, USA.
Koller, D. and Friedman, N. (2009). Probabilistic graphical models: principles and
techniques. MIT press, Cambridge, MA, USA.
Kossakowski, J. J., Epskamp, S., Kieffer, J. M., van Borkulo, C. D., Rhemtulla, M.,
and Borsboom, D. (2015). The application of a network approach to health-related
quality of life (HRQoL): Introducing a new method for assessing hrqol in healthy
adults and cancer patient. Quality of Life Research, 0962-9343:1–12.
Lauritzen, S. L. (1996). Graphical models. Oxford University Press, Oxford, UK.

Psychometrika Submission

September 13, 2017

42

Lawley, D. N. (1940). Vi.the estimation of factor loadings by the method of maximum
likelihood. Proceedings of the Royal Society of Edinburgh, 60(01):64–82.
Lord, F. M., Novick, M. R., and Birnbaum, A. (1968). Statistical theories of mental
test scores. Addison-Wesley, Oxford, UK.
MacCallum, R. C., Wegener, D. T., Uchino, B. N., and Fabrigar, L. R. (1993). The
problem of equivalent models in applications of covariance structure analysis.
Psychological bulletin, 114(1):185 – 199.
Marsh, H. W., Morin, A. J., Parker, P. D., and Kaur, G. (2014). Exploratory
structural equation modeling: An integration of the best features of exploratory and
confirmatory factor analysis. Annual Review of Clinical Psychology, 10:85–110.
McCrae, R. and Costa, P. (1997). Personality trait structure as a human universal.
American psychologist, 52(5):509–516.
McGill, R., Tukey, J. W., and Larsen, W. A. (1978). Variations of box plots. The
American Statistician, 32(1):12–16.
Murphy, K. P. (2012). Machine learning: a probabilistic perspective. MIT press,
Cambridge, MA, USA.
Neale, M. C., Hunter, M. D., Pritikin, J. N., Zahery, M., Brick, T. R., Kirkpatrick,
R. M., Estabrook, R., Bates, T. C., Maes, H. H., and Boker, S. M. (2015). Openmx
2.0: Extended structural equation and statistical modeling. Psychometrika, pages
1–15.
Pearl, J. (2000). Causality: models, reasoning and inference, volume 29. Cambridge
Univ Press, Cambridge, UK.
Reckase, M. D. (2009). Multidimensional item response theory. Springer, New York,
NY, USA.

Psychometrika Submission

September 13, 2017

43

Revelle, W. (2015). psych: Procedures for Psychological, Psychometric, and Personality
Research. Northwestern University, Evanston, Illinois. R package version 1.5.8.
Rosa, M., Friston, K., and Penny, W. (2012). Post-hoc selection of dynamic causal
models. Journal of neuroscience methods, 208(1):66–78.
Schmittmann, V. D., Cramer, A. O. J., Waldorp, L. J., Epskamp, S., Kievit, R. A.,
and Borsboom, D. (2013). Deconstructing the construct: A network perspective on
psychological phenomena. New ideas in psychology, 31(1):43–53.
Tibshirani, R. (1996). Regression shrinkage and selection via the lasso. Journal of the
Royal Statistical Society. Series B (Methodological), pages 267–288.
van Borkulo, C. D., Borsboom, D., Epskamp, S., Blanken, T. F., Boschloo, L.,
Schoevers, R. A., and Waldorp, L. J. (2014). A new method for constructing
networks from binary data. Scientific Reports, 4:Article number: 5918.
Van der Maas, H. L. J., Dolan, C. V., Grasman, R. P. P. P., Wicherts, J. M., Huizenga,
H. M., and Raijmakers, M. E. J. (2006). A dynamical model of general intelligence:
the positive manifold of intelligence by mutualism. Psychological review, 113:842–861.
Ware Jr, J. E. and Sherbourne, C. D. (1992). The MOS 36-item short-form health
survey (SF-36): I. Conceptual framework and item selection. Medical Care,
30:473–483.
Wright, S. (1921). Correlation and causation. Journal of agricultural research,
20(7):557–585.
Wright, S. (1934). The method of path coefficients. The Annals of Mathematical
Statistics, 5(3):161–215.
Yin, J. and Li, H. (2011). A sparse conditional gaussian graphical model for analysis of
genetical genomics data. The annals of applied statistics, 5(4):26302650.

Psychometrika Submission

September 13, 2017

44

Zou, H. and Hastie, T. (2005). Regularization and variable selection via the elastic net.
Journal of the Royal Statistical Society: Series B (Statistical Methodology),
67(2):301–320.
Zou, H., Hastie, T., Tibshirani, R., et al. (2007). On the degrees of freedom of the
lasso. The Annals of Statistics, 35(5):2173–2192.

