Instrumental Variable (MIIV-2SLS) Estimation
estimator_iv.RdEstimate a structural equation model using model-implied instrumental
variables and two-stage least squares (MIIV-2SLS), a non-iterative,
equation-by-equation alternative to maximum likelihood (Bollen, 1996).
This page documents how to invoke the estimator, the options that can be
passed via the estimator = list() argument, how to specify
instruments manually with the |~ operator, and the lavInspect()
fields that are relevant for IV estimation.
Details
Instead of minimizing a single discrepancy function for the whole model, MIIV-2SLS estimates each measurement and structural equation separately by two-stage least squares (2SLS). For each equation, the instruments are other observed variables in the model that are (under the hypothesized model) uncorrelated with the equation's composite error term: the model-implied instrumental variables (MIIVs). Because the equations are estimated separately, a misspecification in one equation does not necessarily propagate to the others, the method is non-iterative, and no starting values are required.
The estimator is selected by setting estimator = "IV" (the alias
"MIIV" is also accepted; names are case-insensitive), for example
fit <- sem(model, data = Data, estimator = "IV")Estimation proceeds in two stages. In the first stage, the directed effects (factor loadings and regression coefficients) are estimated equation-by-equation by 2SLS. In the second stage, the undirected effects (residual variances and covariances) are estimated by a weighted least-squares step applied to the residual moments. Standard errors are available for both stages (see the options below).
Both continuous and ordered-categorical data are supported (for ordered data,
the polychoric-based PIV estimator of Bollen & Maydeu-Olivares, 2007 is used).
Simple equality constraints (e.g. equal factor loadings) and general linear
equality constraints (e.g. a == 2*b) are honored. Multiple groups are
supported, with standard errors (for continuous, two-stage-missing and
categorical data), including the mean structure and cross-group equality
constraints for measurement invariance testing (configural, metric and scalar,
e.g. via group.equal = c("loadings", "intercepts")). For scalar
invariance the latent-mean (and intercept) estimates are obtained by a joint
GLS solve over the mean structure (see iv_mean_structure), which matches
the maximum likelihood mean estimates given the model-implied covariance
matrix. The IV estimator is restricted to single-level models.
Estimator options
Options specific to the IV estimator are passed as named elements of a list
to the estimator argument, alongside the estimator name:
fit <- sem(model, data = Data,
estimator = list(estimator = "IV",
iv_varcov_method = "2RLS",
iv_sargan_adjust = "BH"))The option names use snake_case. For backward compatibility, dot.case names
(e.g. iv.varcov.method) are also accepted and silently converted.
iv_method:Character. The instrumental-variable method. Currently only
"2SLS"(two-stage least squares) is available. Default is"2SLS".iv_samplestats:Logical. If
TRUE(the default), the equations are estimated from the sample moments (covariances and means) rather than from the raw data. This is required for categorical data and when only sample moments are provided as input.iv_varcov_method:Character. The second-stage method used to estimate the residual variances and covariances. One of
"ULS","GLS","2RLS","RLS"or"NONE". Default is"RLS"(reweighted least squares). For categorical data,"ULS"is used. If"NONE", the variance/covariance parameters are not estimated.iv_sargan:Logical. If
TRUE(the default), overidentification tests are computed for each overidentified equation. Two tests are reported side by side: the classical Sargan test (Sargan, 1958), which assumes normality (continuous data) and is not valid for polychoric correlations; and a robust, residual-based test that is valid more generally. The robust test is Browne's (1984) residual statistic applied to the single equation, computed with the asymptotic covariance matrix (ACOV) of the sample statistics: the distribution-free (ADF) ACOV of the sample covariances for continuous data (so it is robust to non-normality, and is then equivalent to Hansen's overidentification test), or the polychoric ACOV for categorical data. Both tests share the same degrees of freedom (number of instruments minus number of regressors). For continuous data with missing values the robust test is not available. The results can be retrieved withlavInspect(fit, "sargan")(or, equivalently,lavInspect(fit, "hansen")).iv_sargan_adjust:Character. A multiple-comparison adjustment applied to the per-equation overidentification p-values, using any method accepted by
p.adjust("holm","hochberg","hommel","bonferroni","BH","BY","fdr", or"none"). Default is"none". When a method other than"none"is requested, extrasargan.pval.adjandbrowne.pval.adjcolumns are added to the table returned bylavInspect(fit, "sargan").iv_weak:Character. How to respond to weak instruments, diagnosed per equation by the first-stage F-statistic (Staiger & Stock, 1997). One of:
"warn":(the default) report the affected equations and their first-stage F, and suggest supplying stronger instruments manually; the estimates are not changed.
"prune":greedily drop the weakest excess instruments from weak, overidentified equations (never below just-identified) and report what was dropped.
"none":skip the weak-instrument check.
iv_weak_threshold:Numeric. The first-stage F-statistic below which instruments are deemed weak. Default is
10.iv_vcov_stage1:Character. The standard-error method for the first-stage (directed) coefficients. One of
"lm.vcov","lm.vcov.dfres","gamma"or"none". Default is"lm.vcov.dfres"for continuous data and"gamma"for categorical data. The"gamma"method requiresiv_samplestats = TRUE. When the model contains simple equality constraints among the directed coefficients (e.g. equal factor loadings), the"lm.vcov"/"lm.vcov.dfres"methods use a variance-weighted restricted-2SLS covariance for the constrained coefficients (as in the MIIVsem package); setiv_vcov_stage1 = "gamma"to obtain the delta-method (moment-Jacobian) standard errors instead.iv_vcov_stage2:Character. The standard-error method for the second-stage (undirected) parameters. One of
"h2","delta"or"none". Default is"h2". The"h2"method requiresiv_samplestats = TRUE.iv_vcov_gamma_modelbased:Logical. If
TRUE(the default), the normal-theory weight matrix (gamma) used for the standard errors is based on the model-implied covariance matrix; ifFALSE, it is based on the unrestricted (H1) sample covariance matrix.iv_mean_structure:Character. How the free mean-structure parameters (observed intercepts and latent means) are estimated. With
"wls"(the default) they are obtained by a joint GLS solve that fits the model-implied means to the sample means across all groups, pooling parameters that are constrained equal across groups (e.g. intercepts for scalar invariance); this matches the maximum likelihood mean estimates (given the model-implied covariance matrix). With"moments"the intercepts are recomputed per equation and shared intercepts are pooled by a (simpler) nobs-weighted average. The two agree unless intercepts are constrained across groups (scalar invariance).iv_vcov_jack_numerical,iv_vcov_jaca_numerical,iv_vcov_jacb_numerical:Logical. If
TRUE, the corresponding Jacobians used in the standard-error computation are obtained numerically rather than from analytic expressions. The default isFALSE; these are mainly intended for checking the analytic derivatives.iv_mimic_ml:Logical. If
FALSE(the default), the IV estimator uses unbiased divisors: the sample covariance matrix is divided by \(N-1\) (sample.cov.rescale = FALSE) and the per-equation residual variances (the RSS terms entering the standard errors and the Sargan test) use the \(N-P\) divisor. IfTRUE, the maximum likelihood divisor \(N\) is used throughout: the sample covariance is rescaled to the ML divisor \(N\) and the residual variances use the \(N\) divisor as well. As a result, for just-identified models the IV point estimates and standard errors match the ML solution exactly.
Specifying instruments manually (the |~ operator)
By default lavaan determines the model-implied instruments for each equation
automatically. The |~ operator overrides this choice for a specific
equation: the left-hand side is the dependent variable of the equation (for a
latent variable, its scaling indicator), and the right-hand side lists the
instruments to use.
model <- '
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem60 ~ ind60
# use only x2 and x3 as instruments for the
# dem60 ~ ind60 equation
y1 |~ x2 + x3
'
fit <- sem(model, data = PoliticalDemocracy, estimator = "IV")This is useful when the automatically selected instruments are weak (see
iv_weak above) or when subject-matter knowledge suggests a different
set of instruments.
External instruments. The instruments listed on the right-hand side of
|~ need not be part of the model. If an instrument is an observed
variable that does not otherwise appear in the model (it is not an indicator,
a predictor, an outcome, or otherwise mentioned), it is treated as an
external instrument: it is read from the data so that its covariances
with the model variables are available to the 2SLS estimator, but it never
enters the model-implied summary statistics. The model degrees of freedom, the
fit measures, and the model-implied moments are therefore exactly the same as
if the instrument had not been mentioned; only the equations for which it is
used are affected. External instruments require raw data (a data=
argument); they are not available when only sample moments are supplied. For
categorical data the external instruments must themselves be ordered
(categorical); the augmented polychoric correlations and the asymptotic
covariance of the sample statistics are then used for the estimates and the
standard errors. Single-group and multiple-group categorical models are
supported, but (for now) not equality constraints among the directed
coefficients (e.g. measurement invariance) together with external instruments.
## z is not part of the model; it is used only as an
## instrument for the (endogenous) regression y ~ x
model <- '
y ~ x
y |~ z
'
fit <- sem(model, data = Data, estimator = "IV")For external instruments the point estimates and the directed-coefficient (loading and regression) standard errors are the usual 2SLS quantities. The second-stage (residual variance and covariance) standard errors also account for the sampling variability of the instrument moments: the directed coefficients depend on the instrument-model covariances, and this uncertainty is propagated into the variance/covariance parameters through the full (augmented) moment Jacobian.
Missing data
For continuous data with missing values, the IV estimator defaults to a
two-stage (FIML) approach: the saturated mean vector and covariance matrix are
estimated by the EM algorithm, the model is then estimated from these moments,
and the standard errors are corrected for the additional uncertainty stemming
from the EM estimates (Savalei & Bentler, 2009). This default is used when the
data contain missing values and the user does not set the missing
argument. Set missing = "robust.two.stage" for sandwich-corrected
standard errors (Savalei & Falk, 2014), or missing = "listwise" to
delete incomplete cases instead. Two-stage missing-data estimation is also
available for multiple-group models.
Multiple groups and measurement invariance
Multiple-group models are supported, with point estimates and standard errors
for continuous, two-stage-missing and categorical data. Cross-group equality
constraints are imposed in the usual way, either through shared parameter
labels in the model syntax or through the group.equal argument, so the
standard measurement-invariance sequence is available: configural
(group.equal = NULL), metric (group.equal = "loadings") and
scalar (group.equal = c("loadings", "intercepts")).
A parameter that is constrained equal across groups is estimated by pooling the
information from all groups; for a configural model the per-group estimates (and
their standard errors) reproduce the separate single-group fits. The mean
structure (observed intercepts and latent means) is estimated by a joint GLS
solve over all groups (see iv_mean_structure), which reproduces the
maximum likelihood mean estimates given the model-implied covariance matrix.
Inspecting the instruments and the Sargan test
Several lavInspect fields are specific to IV estimation:
lavInspect(fit, "iv"):(aliases
"miiv","instruments") a per-equation table with the dependent variable (lhs), the regressors (rhs), the observed variables actually used (lhs_new,rhs_new), the equation type (type:"miiv","ols"or"user"), and the instruments (iv) for that equation.lavInspect(fit, "sargan"):a per-equation table with the overidentification tests: the degrees of freedom (
df), the classical Sargan statistic (sargan.stat) and its p-value (sargan.pval), and the robust Browne residual-based statistic (browne.stat) and its p-value (browne.pval); seeiv_sarganabove. For categorical data the classical Sargan p-value isNA(not valid), so use thebrowne.*columns. Adjusted p-values (sargan.pval.adjandbrowne.pval.adj) are added wheniv_sargan_adjustis set.lavInspect(fit, "hansen")is an alias.lavInspect(fit, "eqs"):the underlying per-equation information.
References
Bollen, K. A. (1996). An alternative two stage least squares (2SLS) estimator for latent variable equations. Psychometrika, 61(1), 109-121.
Bollen, K. A., & Maydeu-Olivares, A. (2007). A polychoric instrumental variable (PIV) estimator for structural equation models with categorical variables. Psychometrika, 72(3), 309-326.
Browne, M. W. (1984). Asymptotically distribution-free methods for the analysis of covariance structures. British Journal of Mathematical and Statistical Psychology, 37(1), 62-83.
Sargan, J. D. (1958). The estimation of economic relationships using instrumental variables. Econometrica, 26(3), 393-415.
Savalei, V., & Bentler, P. M. (2009). A two-stage approach to missing data: Theory and application to auxiliary variables. Structural Equation Modeling, 16(3), 477-497.
Savalei, V., & Falk, C. F. (2014). Robust two-stage approach outperforms robust full information maximum likelihood with incomplete nonnormal data. Structural Equation Modeling, 21(2), 280-302.
Staiger, D., & Stock, J. H. (1997). Instrumental variables regression with weak instruments. Econometrica, 65(3), 557-586.
Examples
## The classic Bollen (1989) Political Democracy example
model <- '
# measurement model
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem65 =~ y5 + y6 + y7 + y8
# regressions
dem60 ~ ind60
dem65 ~ ind60 + dem60
'
## fit using MIIV-2SLS instead of maximum likelihood
fit <- sem(model, data = PoliticalDemocracy, estimator = "IV")
summary(fit)
#> lavaan 0.7-1.3097 ended normally after 1 iteration
#>
#> Estimator IV
#> Optimization method NONITER
#> Number of model parameters 25
#>
#> Number of observations 75
#>
#> Model Test User Model:
#>
#> Browne's residual-based (NT) test
#> Test statistic 55.032
#> Degrees of freedom 41
#> P-value (Chi-square) 0.070
#>
#> Parameter Estimates:
#>
#> Standard errors Two-stage
#> Standard errors stage 1 lm.vcov.dfres
#> Standard errors stage 2 h2
#> Gamma matrix NT
#>
#> Latent Variables:
#> Estimate Std.Err z-value P(>|z|)
#> ind60 =~
#> x1 1.000
#> x2 2.078 0.130 15.954 0.000
#> x3 1.751 0.151 11.623 0.000
#> dem60 =~
#> y1 1.000
#> y2 1.227 0.171 7.187 0.000
#> y3 0.986 0.131 7.515 0.000
#> y4 1.178 0.133 8.850 0.000
#> dem65 =~
#> y5 1.000
#> y6 1.086 0.159 6.821 0.000
#> y7 1.132 0.139 8.161 0.000
#> y8 1.149 0.146 7.854 0.000
#>
#> Regressions:
#> Estimate Std.Err z-value P(>|z|)
#> dem60 ~
#> ind60 1.261 0.431 2.923 0.003
#> dem65 ~
#> ind60 1.123 0.319 3.525 0.000
#> dem60 0.724 0.104 6.996 0.000
#>
#> Variances:
#> Estimate Std.Err z-value P(>|z|)
#> .x1 0.077 0.020 3.944 0.000
#> .x2 0.152 0.068 2.242 0.025
#> .x3 0.470 0.090 5.243 0.000
#> .y1 1.767 0.424 4.168 0.000
#> .y2 6.597 1.224 5.390 0.000
#> .y3 5.204 0.945 5.504 0.000
#> .y4 3.101 0.670 4.628 0.000
#> .y5 2.228 0.443 5.034 0.000
#> .y6 4.715 0.855 5.515 0.000
#> .y7 3.702 0.702 5.274 0.000
#> .y8 3.145 0.620 5.073 0.000
#> ind60 0.480 0.090 5.321 0.000
#> .dem60 4.891 1.125 4.346 0.000
#> .dem65 0.363 0.257 1.414 0.157
#>
## the model-implied instruments used for each equation
lavInspect(fit, "iv")
#> lhs rhs lhs.new rhs.new type instruments
#> 1 x2 ind60 x2 x1 miiv x3, y1, y2, y3, y4, y5, y6, y7, y8
#> 2 x3 ind60 x3 x1 miiv x2, y1, y2, y3, y4, y5, y6, y7, y8
#> 3 y2 dem60 y2 y1 miiv x1, x2, x3, y3, y4, y5, y6, y7, y8
#> 4 y3 dem60 y3 y1 miiv x1, x2, x3, y2, y4, y5, y6, y7, y8
#> 5 y4 dem60 y4 y1 miiv x1, x2, x3, y2, y3, y5, y6, y7, y8
#> 6 y6 dem65 y6 y5 miiv x1, x2, x3, y1, y2, y3, y4, y7, y8
#> 7 y7 dem65 y7 y5 miiv x1, x2, x3, y1, y2, y3, y4, y6, y8
#> 8 y8 dem65 y8 y5 miiv x1, x2, x3, y1, y2, y3, y4, y6, y7
#> 9 dem60 ind60 y1 x1 miiv x2, x3
#> 10 dem65 ind60 + dem60 y5 x1 + y1 miiv x2, x3, y2, y3, y4
## the per-equation Sargan overidentification tests
lavInspect(fit, "sargan")
#> lhs rhs df sargan.stat sargan.pval
#> 1 x2 ind60 8 8.080 0.426
#> 2 x3 ind60 8 8.505 0.386
#> 3 y2 dem60 8 18.360 0.019
#> 4 y3 dem60 8 9.884 0.273
#> 5 y4 dem60 8 14.487 0.070
#> 6 y6 dem65 8 20.021 0.010
#> 7 y7 dem65 8 13.400 0.099
#> 8 y8 dem65 8 14.893 0.061
#> 9 dem60 ind60 1 0.489 0.484
#> 10 dem65 ind60 + dem60 3 0.769 0.857
## pass estimator options via estimator = list(...)
fit2 <- sem(model, data = PoliticalDemocracy,
estimator = list(estimator = "IV",
iv_varcov_method = "2RLS",
iv_sargan_adjust = "BH",
iv_weak = "warn"))
lavInspect(fit2, "sargan")
#> lhs rhs df sargan.stat sargan.pval sargan.pval.adj
#> 1 x2 ind60 8 8.080 0.426 0.532
#> 2 x3 ind60 8 8.505 0.386 0.532
#> 3 y2 dem60 8 18.360 0.019 0.093
#> 4 y3 dem60 8 9.884 0.273 0.455
#> 5 y4 dem60 8 14.487 0.070 0.175
#> 6 y6 dem65 8 20.021 0.010 0.093
#> 7 y7 dem65 8 13.400 0.099 0.198
#> 8 y8 dem65 8 14.893 0.061 0.175
#> 9 dem60 ind60 1 0.489 0.484 0.538
#> 10 dem65 ind60 + dem60 3 0.769 0.857 0.857
## specify instruments manually with the |~ operator:
## use only x2, x3 as instruments for the dem60 ~ ind60 equation
model.iv <- '
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + y2 + y3 + y4
dem60 ~ ind60
y1 |~ x2 + x3
'
fit3 <- sem(model.iv, data = PoliticalDemocracy, estimator = "IV")
lavInspect(fit3, "iv")
#> lhs rhs lhs.new rhs.new type instruments
#> 1 x2 ind60 x2 x1 miiv x3, y1, y2, y3, y4
#> 2 x3 ind60 x3 x1 miiv x2, y1, y2, y3, y4
#> 3 y2 dem60 y2 y1 miiv x1, x2, x3, y3, y4
#> 4 y3 dem60 y3 y1 miiv x1, x2, x3, y2, y4
#> 5 y4 dem60 y4 y1 miiv x1, x2, x3, y2, y3
#> 6 dem60 ind60 y1 x1 user x2, x3
if (FALSE) { # \dontrun{
## equality constraints are honored (here: equal loadings across factors)
model.eq <- '
ind60 =~ x1 + x2 + x3
dem60 =~ y1 + a*y2 + b*y3 + c*y4
dem65 =~ y5 + a*y6 + b*y7 + c*y8
dem60 ~ ind60
dem65 ~ ind60 + dem60
'
fit4 <- sem(model.eq, data = PoliticalDemocracy, estimator = "IV")
coef(fit4)
## multiple groups: metric measurement invariance (equal loadings)
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9 '
fit.metric <- sem(HS.model, data = HolzingerSwineford1939, estimator = "IV",
group = "school", group.equal = "loadings")
summary(fit.metric)
} # }