Predict the values of latent variables (and their indicators)
lavPredict.RdThe main purpose of the lavPredict() function is to compute (or
`predict') individual scores for the latent variables in the model
(`factor scores'). NOTE: the goal of this
function is NOT to predict future values of dependent variables as in the
regression framework! (For models with only continuous observed variables, the function lavPredictY() supports this.)
Usage
lavPredict(object, newdata = NULL, type = "lv", method = "EBM",
transform = FALSE, se = "none", acov = "none",
label = TRUE, fsm = FALSE, mdist = FALSE, rel = FALSE,
append_data = FALSE, assemble = FALSE,
level = 1L, optim_method = "bfgs", eta = NULL,
parallel = c("auto", "no", "multicore", "snow"),
ncpus = NULL, cl = NULL,
drop_list_single_group = TRUE,
mdist_draws = 2000L,
...)Arguments
- object
An object of class
lavaan.- newdata
An optional data.frame, containing the same variables as the data.frame used when fitting the model in object. For multilevel models, the cluster variable(s) must be included; the cluster structure and (for
missing = "ml") the missing-data patterns are rebuilt from newdata, so the predictions are conditional on the newdata only (the cluster labels need not match those of the original data).- type
A character string. If
"lv", estimated values for the latent variables in the model are computed. If"ov", model predicted values for the indicators of the latent variables in the model are computed. If"yhat", the estimated values for the observed indicators are computed, given the user-specified values for the latent variables provided by theetaargument. If"fy", densities (or probabilities) are computed for each observed indicator, given the user-specified values for the latent variables provided by theetaargument.- method
A character string. In the linear case (when the indicators are continuous), the possible options are
"regression"or"Bartlett". In the categorical case, the two options are"EBM"for the Empirical Bayes Modal approach, and"ML"for the maximum likelihood approach. For higher-order factors (which have no observed indicators of their own), the"Bartlett"scores are computed from the collapsed measurement model (regressing the observed indicators of the lower-order factors on the higher-order factors), as insam.- transform
Logical. If
TRUE, transform the factor scores (per group) so that their mean and variance-covariance matrix matches the model-implied mean and variance-covariance matrix. This may be useful if the individual factor scores will be used in a follow-up (regression) analysis. Note: the standard errors (if requested) are not transformed (yet). The resulting factor scores are often called correlation-preserving factor scores.- se
Character. If
"none", no standard errors are computed. If"standard", standard errors are computed (assuming the parameters of the measurement model are known). The standard errors are returned as an attribute. For continuous data, the (naive) standard errors are the same for every observation. For categorical data, where the factor scores are obtained by numerical optimization (method = "EBM"ormethod = "ML"), the standard errors are based on the curvature of the objective at the optimum (the inverse of its Hessian, a Laplace/observed -information approximation) and therefore differ from one response pattern to the next; the result is an (nobs x nfactor) matrix per group.- acov
Similar to the
"se"argument, but optionally returns the full sampling covariance matrix of the factor scores as an attribute. For continuous data this is a single matrix per group; for categorical data it is a list with one (nfactor x nfactor) matrix per observation.- label
Logical. If TRUE, the columns in the output are labeled.
- fsm
Logical. If TRUE, return the factor score matrix as an attribute. Only for numeric data.
- mdist
Logical. If TRUE, the (squared) Mahalanobis distances of the factor scores (if
type = "lv") or the casewise residuals (iftype = "resid") are returned as an attribute. If (some of) the data is ordered categorical, the distances are generalized as in Mansolf and Reise (2017): the continuous distance is replaced by its expected value over the region of the (multivariate normally distributed) latent response vector that is consistent with the observed response pattern; this expectation is approximated by Monte Carlo integration (usingmdist_drawsrandom draws per distinct response pattern; the result therefore depends on the state of the random number generator, andset.seed()can be used for reproducibility). Observed continuous variables are conditioned on, and missing values are integrated out. In this categorical setting,type = "resid"is allowed, and returns the expected casewise residuals of the latent responses (a by-product of the same Monte Carlo integration).- rel
Logical. Only used if
type = "lv". If TRUE, the factor reliabilities are returned as an attribute. (The squared values are often called the factor determinacies.)- append_data
Logical. Only used when
type = "lv". If TRUE, the original data (or the data provided in the newdata argument) is appended to the factor scores.- assemble
Logical. If TRUE, the separate groups are reassembled into a single data.frame with a group column, having the same dimensions as the original (or newdata) dataset.
- level
Integer. Only used in a multilevel SEM. If
level = 1, only factor scores for latent variable defined at the first (within) level are computed; iflevel = 2, only factor scores for latent variables defined at the second (between) level are computed. For two-level models withmissing = "ml", the posterior means of the cluster-level components and of the missing values (given all the observed data; computed as in the E-step of the EM algorithm) are used, so that the factor scores are the exact posterior means given all the observed data (new in 0.7-1).- optim_method
Character string. Only used in the categorical case. If
"nlminb"(the default in 0.5), the"nlminb()"function is used for the optimization. If"bfgs"or"BFGS"(the default in 0.6), the"optim()"function is used with the BFGS method.- eta
An optional matrix or list, containing latent variable values for each observation. Used for computations when
type = "ov".- parallel
Character. Only used in the categorical case, where factor scores are computed by a per-observation numerical optimization. The options are
"no"(serial),"multicore"(fork-based, not available on Windows) and"snow"(a PSOCK cluster). The default,"auto", uses"multicore"automatically (on non-Windows platforms) when a large number of distinct optimizations is required, and runs serially otherwise. Because the computation is deterministic, the results do not depend on the value of this argument.- ncpus
Integer. The number of processes to use in parallel computation. The default is to use all available cores, minus two.
- cl
An optional parallel or snow cluster for use when
parallel = "snow". If not supplied, a cluster on the local machine is created for the duration of the call.- drop_list_single_group
Logical. If
FALSE, the results are returned as a list, where each element corresponds to a group (even if there is only a single group). IfTRUE, the list will be unlisted if there is only a single group.- mdist_draws
Integer. Only used if
mdist = TRUEand (some of) the data is ordered categorical: the number of Monte Carlo draws per distinct response pattern used to approximate the expected Mahalanobis distances.- ...
To support old argument names.
Details
The predict() function calls the lavPredict() function
with its default options.
If there are no latent variables in the model, type = "ov" will
simply return the values of the observed variables. Note that this function
can not be used to `predict' values of dependent variables, given the
values of independent variables (in the regression sense). In other words,
the structural component is completely ignored (for now).
References
For an overview (and evaluation) of the various factor score methods, see:
Grice, J. W. (2001). Computing and evaluating factor scores. Psychological Methods, 6(4), 430-450. doi:10.1037/1082-989X.6.4.430
For the (continuous) regression and Bartlett methods, see:
Bartlett, M. S. (1937). The statistical conception of mental factors. British Journal of Psychology, 28, 97-104. doi:10.1111/j.2044-8295.1937.tb00863.x
Bentler, P. M., & Yuan, K.-H. (1997). Optimal conditionally unbiased equivariant factor score estimators. In M. Berkane (Ed.), Latent variable modeling and applications to causality (pp. 259-281). New York: Springer-Verlag. doi:10.1007/978-1-4612-1842-5_14
For the (categorical) Empirical Bayes Modal (EBM) and Maximum Likelihood (ML) methods, see:
Skrondal, A., & Rabe-Hesketh, S. (2004). Generalized latent variable modeling: Multilevel, longitudinal, and structural equation models. Boca Raton, FL: Chapman & Hall/CRC.
For the Mahalanobis distances with ordered categorical data
(mdist = TRUE), see:
Mansolf, M., & Reise, S. P. (2017). Case diagnostics for factor analysis of ordered categorical data with applications to person-fit measurement. Structural Equation Modeling: A Multidisciplinary Journal, 25(1), 86-100. doi:10.1080/10705511.2017.1367926
For the correlation-preserving factor scores (transform = TRUE), see:
ten Berge, J. M. F., Krijnen, W. P., Wansbeek, T., & Shapiro, A. (1999). Some new results on correlation-preserving factor scores prediction methods. Linear Algebra and its Applications, 289(1-3), 311-318. doi:10.1016/S0024-3795(97)10007-6
See also
lavPredictY to predict y-variables given x-variables.
Examples
data(HolzingerSwineford1939)
## fit model
HS.model <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6
speed =~ x7 + x8 + x9 '
fit <- cfa(HS.model, data = HolzingerSwineford1939)
head(lavPredict(fit))
#> visual textual speed
#> [1,] -0.81767524 -0.13754501 0.06150726
#> [2,] 0.04951940 -1.01272402 0.62549360
#> [3,] -0.76139670 -1.87228634 -0.84057276
#> [4,] 0.41934153 0.01848569 -0.27133710
#> [5,] -0.41590481 -0.12225009 0.19432951
#> [6,] 0.02325632 -1.32981727 0.70885348
head(lavPredict(fit, type = "ov"))
#> x1 x2 x3 x4 x5 x6 x7 x8
#> [1,] 4.118094 5.635456 1.654027 2.923363 4.187433 2.0581851 4.247409 5.599652
#> [2,] 4.985289 6.115449 2.286533 2.048184 3.213292 1.2476414 4.811396 6.265128
#> [3,] 4.174373 5.666607 1.695075 1.188622 2.256533 0.4515610 3.345329 4.535242
#> [4,] 5.355111 6.320146 2.556271 3.079394 4.361108 2.2026924 3.914565 5.206912
#> [5,] 4.519865 5.857836 1.947067 2.938658 4.204458 2.0723504 4.380232 5.756376
#> [6,] 4.959026 6.100912 2.267378 1.731091 2.860343 0.9539666 4.894756 6.363489
#> x9
#> [1,] 5.440645
#> [2,] 6.050613
#> [3,] 4.465019
#> [4,] 5.080664
#> [5,] 5.584297
#> [6,] 6.140770
## ---------------------------------------------
## standard errors for the factor scores (se =)
## ---------------------------------------------
## the standard errors are returned as the "se" attribute (a list, one
## (nobs x nfactor) matrix per group)
## for continuous indicators, the (naive) standard errors are the same
## for every observation
fscores <- lavPredict(fit, se = "standard")
attr(fscores, "se")[[1]]
#> visual textual speed
#> [1,] 0.4747366 0.3316382 0.3276657
## for categorical indicators, the factor scores are obtained by numerical
## optimization, and their standard errors differ from one response pattern
## (observation) to the next
## (a two-factor model is used here, as the numerical optimization and the
## Monte Carlo integration below both get expensive as the number of
## ordered indicators grows)
HS.model.ord <- ' visual =~ x1 + x2 + x3
textual =~ x4 + x5 + x6 '
HS.ord <- HolzingerSwineford1939
HS.ord[ , paste0("x", 1:6)] <-
lapply(HS.ord[ , paste0("x", 1:6)], function(x) ordered(cut(x, 3)))
fit.ord <- cfa(HS.model.ord, data = HS.ord, ordered = paste0("x", 1:6))
fscores <- lavPredict(fit.ord, se = "standard")
head(attr(fscores, "se")[[1]])
#> visual textual
#> [1,] 0.5256836 0.2973992
#> [2,] 0.5429494 0.3913308
#> [3,] 0.5429494 0.3913308
#> [4,] 0.5025922 0.3350009
#> [5,] 0.5529813 0.3422373
#> [6,] 0.5429494 0.3913308
## -------------------------------------------------
## casewise Mahalanobis distances (case diagnostics)
## -------------------------------------------------
## continuous data: the squared distances are exact
fscores <- lavPredict(fit, mdist = TRUE)
head(attr(fscores, "mdist")[[1]])
#> [1] 2.270885 3.971682 5.494760 1.489545 1.067721 5.892432
## ordered categorical data: expected (squared) distances of the latent
## responses, obtained by Monte Carlo integration (Mansolf & Reise, 2017);
## set the seed for reproducibility
set.seed(123)
resid.ord <- lavPredict(fit.ord, type = "resid", mdist = TRUE)
## expected casewise residuals of the latent responses:
head(resid.ord)
#> x1 x2 x3 x4 x5 x6
#> [1,] -0.11373505 1.2144354 -0.9378889 -0.14164761 0.3768976 -0.587162724
#> [2,] 0.03387823 -0.1080701 0.3427065 -0.23396444 -0.1214889 -0.005852471
#> [3,] 0.03387823 -0.1080701 0.3427065 -0.23396444 -0.1214889 -0.005852471
#> [4,] -0.18803703 1.1277673 -0.0332436 -0.04067443 -0.1407833 0.438579839
#> [5,] 0.02544435 -0.1134666 -0.7694199 -0.02062462 -0.1141321 0.477697843
#> [6,] 0.03387823 -0.1080701 0.3427065 -0.23396444 -0.1214889 -0.005852471
## expected squared residual-based distances (person fit):
head(attr(resid.ord, "mdist")[[1]])
#> [1] 6.644701 2.420138 2.420138 4.043837 3.559749 2.420138
## ------------------------------------------
## merge factor scores to original data.frame
## ------------------------------------------
idx <- lavInspect(fit, "case.idx")
fscores <- lavPredict(fit)
## loop over factors
for (fs in colnames(fscores)) {
HolzingerSwineford1939[idx, fs] <- fscores[ , fs]
}
head(HolzingerSwineford1939)
#> id sex ageyr agemo school grade x1 x2 x3 x4 x5 x6
#> 1 1 1 13 1 Pasteur 7 3.333333 7.75 0.375 2.333333 5.75 1.2857143
#> 2 2 2 13 7 Pasteur 7 5.333333 5.25 2.125 1.666667 3.00 1.2857143
#> 3 3 2 13 1 Pasteur 7 4.500000 5.25 1.875 1.000000 1.75 0.4285714
#> 4 4 1 13 2 Pasteur 7 5.333333 7.75 3.000 2.666667 4.50 2.4285714
#> 5 5 2 12 2 Pasteur 7 4.833333 4.75 0.875 2.666667 4.00 2.5714286
#> 6 6 2 14 1 Pasteur 7 5.333333 5.00 2.250 1.000000 3.00 0.8571429
#> x7 x8 x9 visual textual speed
#> 1 3.391304 5.75 6.361111 -0.81767524 -0.13754501 0.06150726
#> 2 3.782609 6.25 7.916667 0.04951940 -1.01272402 0.62549360
#> 3 3.260870 3.90 4.416667 -0.76139670 -1.87228634 -0.84057276
#> 4 3.000000 5.30 4.861111 0.41934153 0.01848569 -0.27133710
#> 5 3.695652 6.30 5.916667 -0.41590481 -0.12225009 0.19432951
#> 6 4.347826 6.65 7.500000 0.02325632 -1.32981727 0.70885348
## multigroup models return a list of factor scores (one per group)
data(HolzingerSwineford1939)
mgfit <- update(fit, group = "school", group.equal = c("loadings","intercepts"))
#> Warning: lavaan->lav_step11_estoptim():
#> Model estimation FAILED! Returning starting values.
#> Error in lav_mvn_loglik_samp(sample_mean = lavsamplestats@mean[[g]], sample_cov = lavsamplestats@cov[[g]], sample_nobs = lavsamplestats@nobs[[g]], mu = mu, sigma_1 = lavimplied$cov[[g]], x_idx = [email protected][[g]], x_mean = [email protected][[g]], x_cov = [email protected][[g]], sinv_method = "eigen", sigma_inv = NULL): non-conformable arguments
idx <- lavInspect(mgfit, "case.idx") # list: 1 vector per group
#> Error: object 'mgfit' not found
fscores <- lavPredict(mgfit) # list: 1 matrix per group
#> Error: object 'mgfit' not found
## loop over groups and factors
for (g in seq_along(fscores)) {
for (fs in colnames(fscores[[g]])) {
HolzingerSwineford1939[ idx[[g]], fs] <- fscores[[g]][ , fs]
}
}
head(HolzingerSwineford1939)
#> id sex ageyr agemo school grade x1 x2 x3 x4 x5 x6
#> 1 1 1 13 1 Pasteur 7 3.333333 7.75 0.375 2.333333 5.75 1.2857143
#> 2 2 2 13 7 Pasteur 7 5.333333 5.25 2.125 1.666667 3.00 1.2857143
#> 3 3 2 13 1 Pasteur 7 4.500000 5.25 1.875 1.000000 1.75 0.4285714
#> 4 4 1 13 2 Pasteur 7 5.333333 7.75 3.000 2.666667 4.50 2.4285714
#> 5 5 2 12 2 Pasteur 7 4.833333 4.75 0.875 2.666667 4.00 2.5714286
#> 6 6 2 14 1 Pasteur 7 5.333333 5.00 2.250 1.000000 3.00 0.8571429
#> x7 x8 x9 visual textual speed
#> 1 3.391304 5.75 6.361111 -0.81767524 -0.13754501 0.06150726
#> 2 3.782609 6.25 7.916667 0.04951940 -1.01272402 0.62549360
#> 3 3.260870 3.90 4.416667 -0.76139670 -1.87228634 -0.84057276
#> 4 3.000000 5.30 4.861111 0.41934153 0.01848569 -0.27133710
#> 5 3.695652 6.30 5.916667 -0.41590481 -0.12225009 0.19432951
#> 6 4.347826 6.65 7.500000 0.02325632 -1.32981727 0.70885348
## -------------------------------------
## Use factor scores in subsequent models
## -------------------------------------
## see Examples in semTools package: ?plausibleValues