R: Estimation of Multiple-Group Structural Equation Models

mgsem {sirt}

R Documentation

Estimation of Multiple-Group Structural Equation Models

Description

Estimates a multiple-group structural equation model. The function allows arbitrary prior distributions on model parameters and regularized estimation with the SCAD and the LASSO penalty. Moreover, it can also conduct robust moment estimation using the L_p loss function \rho(x)=|x|^p for p \ge 0. See Robitzsch (2023) for more details.

Usage

mgsem(suffstat, model, data=NULL, group=NULL, weights=NULL, estimator="ML",
     p_me=2, p_pen=1, pen_type="scad", diffpar_pen=NULL, pen_sample_size=TRUE,
     a_scad=3.7, eps_approx=0.001, comp_se=TRUE, se_delta_formula=FALSE,
     prior_list=NULL, hessian=TRUE, fixed_parms=FALSE, cd=FALSE,
     cd_control=list(maxiter=20, tol=5*1e-04, interval_length=0.05, method="exact"),
     partable_start=NULL, num_approx=FALSE, technical=NULL, control=list())

Arguments

`suffstat`	List containing sufficient statistics
`model`	Model specification, see examples. Can have entries `est`, `index`, `lower`, `upper`, `prior`, `pen_l2`, `pen_lp`, `pen_difflp`. Each entry can be defined for model matrices `ALPHA`, `NU`, `LAM`, `PHI`, and `PSI`.
`data`	Optional data frame
`group`	Optional vector of group identifiers
`weights`	Optional vector of sampling weights
`estimator`	Character. Can be either `"ML"` for maximum likelihood fitting function or `"ME"` for robust moment estimation.
`p_me`	Power in $L_p$ loss function for robust moment estimation
`p_pen`	Power for penalty in regularized estimation. For regular LASSO and SCAD penalty functions, it is $p=1$.
`pen_type`	Penalty type. Can be either `"scad"` or `"lasso"`.
`diffpar_pen`	List containing values of regularization parameters in fused lasso estimation
`pen_sample_size`	List containing values for sample sizes for regularized estimation
`a_scad`	Parameter $a$ used in SCAD penalty
`eps_approx`	Approximation value for nondifferentiable robust moment fitting function or penalty function
`comp_se`	Logical indicating whether standard errors should be computed
`se_delta_formula`	Logical indicating whether standard errors should be computed according to the delta formula
`prior_list`	List containing specifications of the prior distributions
`hessian`	Logical indicating whether the Hessian matrix should be computed
`fixed_parms`	Logical indicating whether all model parameters should be fixed
`cd`	Logical indicating whether coordinate descent should be used for estimation
`cd_control`	Control parameters for coordinate descent estimation
`partable_start`	Starting values for parameter estimation
`num_approx`	Logical indicating whether derivatives should be computed based on numerical differentiation
`technical`	Parameters used for optimization in `sirt_optimizer`
`control`	Control paramaters for optimization

Details

[MORE INFORMATION TO BE ADDED]

Value

A list with following values

`coef`	Coeffients
`vcov`	Variance matrix
`se`	Vector of standard errors
`partable`	Parameter table
`model`	Specified model
`opt_res`	Result from optimization
`opt_value`	Value of fitting function
`residuals`	Residuals of sufficient statistics
`ic`	Information criteria
`technical`	Specifications of optimizer
`suffstat_vcov`	Variance matrix of sufficient statistics
`me_delta_method`	Input and output matrices for delta method if `estimator="ME"`
`data_proc`	Processed data
`case_ll`	Casewise log-likelihood function
`...`	Further values

References

Robitzsch, A. (2023). Model-robust estimation of multiple-group structural equation models. Algorithms, 16(4), 210. doi:10.3390/a16040210

Examples

## Not run: 
#############################################################################
# EXAMPLE 1: Noninvariant item intercepts in a multiple-group SEM
#############################################################################

#---- simulate data
set.seed(65)
G <- 3  # number of groups
I <- 5  # number of items
# define lambda and nu parameters
lambda <- matrix(1, nrow=G, ncol=I)
nu <- matrix(0, nrow=G, ncol=I)
err_var <- matrix(1, nrow=G, ncol=I)

# define extent of noninvariance
dif_int <- .5

#- 1st group: N(0,1)
nu[1,4] <- dif_int
#- 2nd group: N(0.3,1.5)
gg <- 2 ;
nu[gg,1] <- -dif_int
#- 3nd group: N(.8,1.2)
gg <- 3
nu[gg,2] <- -dif_int
#- define distributions of groups
mu <- c(0,.3,.8)
sigma <- sqrt(c(1,1.5,1.2))
N <- rep(1000,3) # sample sizes per group

exact <- FALSE
suffstat <- sirt::invariance_alignment_simulate(nu, lambda, err_var, mu, sigma, N,
                output="suffstat", groupwise=TRUE, exact=exact)

#---- model specification

# model specifications joint group
est <- list(
        ALPHA=matrix( c(0), ncol=1),
        NU=matrix( 0, nrow=I, ncol=1),
        LAM=matrix(1, nrow=I, ncol=1),
        PHI=matrix(0,nrow=1,ncol=1),
        PSI=diag(rep(1,I))
        )

# parameter index
index <- list(
        ALPHA=0*est$ALPHA,
        NU=1+0*est$NU,
        LAM=1+0*est$LAM,
        PHI=0*est$PHI,
        PSI=diag(1,I)
        )

# lower bounds
lower <- list(
        PSI=diag(rep(0.01,I)), PHI=matrix(0.01,1,1)
        )

#*** joint parameters
group0 <- list(est=est, index=index, lower=lower)

#*** group1
est <- list(
        ALPHA=matrix( c(0), ncol=1),
        NU=matrix( 0, nrow=I, ncol=1),
        LAM=matrix(0, nrow=I, ncol=1),
        PHI=matrix(1,nrow=1,ncol=1)
            )

# parameter index
index <- list(
        ALPHA=0*est$ALPHA,
        NU=0*est$NU,
        LAM=1*est$LAM,
        PHI=0*est$PHI
        )

group1 <- list(est=est, index=index, lower=lower)

#*** group 2 and group 3

# modify parameter index
index$ALPHA <- 1+0*est$ALPHA
index$PHI <- 1+0*est$PHI
group3 <- group2 <- list(est=est, index=index, lower=lower)

#*** define model
model <- list(group0=group0, group1=group1, group2=group2, group3=group3)

#-- estimate model with ML
res1 <- sirt::mgsem( suffstat=suffstat, model=model2, eps_approx=1e-4, estimator="ML",
                    technical=list(maxiter=500, optimizer="optim"),
                    hessian=FALSE, comp_se=FALSE, control=list(trace=1) )
str(res1)

#-- robust moment estimation with p=0.5

optimizer <- "optim"
technical <- list(maxiter=500, optimizer=optimizer)
eps_approx <- 1e-3

res2 <- sirt::mgsem( suffstat=suffstat, model=res1$model, p_me=0.5,
                    eps_approx=eps_approx, estimator="ME", technical=technical,
                    hessian=FALSE, comp_se=FALSE, control=list(trace=1) )

#---- regularized estimation

nu_lam <- 0.1    # regularization parameter

# redefine model
define_model <- function(model, nu_lam)
{
    pen_lp <- list( NU=nu_lam+0*model$group1$est$NU)
    ee <- "group1"
    for (ee in c("group1","group2","group3"))
    {
        model[[ee]]$index$NU <- 1+0*index$NU
        model[[ee]]$pen_lp <- pen_lp
    }
    return(model)
}

model3 <- define_model(model=model, nu_lam=nu_lam)
pen_type <- "scad"

res3 <- sirt::mgsem( suffstat=suffstat, model=model3, p_pen=1, pen_type=pen_type,
                    eps_approx=eps_approx, estimator="ML",
                    technical=list(maxiter=500, optimizer="optim"),
                    hessian=FALSE, comp_se=FALSE, control=list(trace=1) )
str(res3)

## End(Not run)

[Package sirt version 4.1-15 Index]