mgsem: Estimation of Multiple-Group Structural Equation Models

In sirt: Supplementary Item Response Theory Models

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. \Sexpr[results=rd]{tools:::Rd_expr_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)

sirt documentation built on May 29, 2024, 8:43 a.m.