Title:	Likelihood-Based Confidence Interval in Structural Equation Models
Version:	0.11.3
Description:	Forms likelihood-based confidence intervals (LBCIs) for parameters in structural equation modeling, introduced in Cheung and Pesigan (2023) <doi:10.1080/10705511.2023.2183860>. Currently implements the algorithm illustrated by Pek and Wu (2018) <doi:10.1037/met0000163>, and supports the robust LBCI proposed by Falk (2018) <doi:10.1080/10705511.2017.1367254>.
URL:	https://sfcheung.github.io/semlbci/
BugReports:	https://github.com/sfcheung/semlbci/issues
License:	GPL-3
Encoding:	UTF-8
LazyData:	true
RoxygenNote:	7.3.2
Suggests:	testthat (≥ 3.0.0), knitr, rmarkdown
Depends:	R (≥ 4.0.0)
Imports:	lavaan (≥ 0.6.13), nloptr, stats, utils, MASS, ggplot2, ggrepel, rlang, pbapply, callr, methods
VignetteBuilder:	knitr
Config/testthat/parallel:	true
Config/testthat/edition:	3
Config/testthat/start-first:	semlbci_wn_mg_*, ur_ci_bound_ur*
NeedsCompilation:	no
Packaged:	2025-01-25 02:30:12 UTC; shufa
Author:	Shu Fai Cheung [aut, cre], Ivan Jacob Agaloos Pesigan [ctb]
Maintainer:	Shu Fai Cheung <shufai.cheung@gmail.com>
Repository:	CRAN
Date/Publication:	2025-01-25 07:30:02 UTC

semlbci: Likelihood-Based Confidence Interval in Structural Equation Models

Description

Forms likelihood-based confidence intervals (LBCIs) for parameters in structural equation modeling, introduced in Cheung and Pesigan (2023) doi:10.1080/10705511.2023.2183860. Currently implements the algorithm illustrated by Pek and Wu (2018) doi:10.1037/met0000163, and supports the robust LBCI proposed by Falk (2018) doi:10.1080/10705511.2017.1367254.

Author(s)

Maintainer: Shu Fai Cheung shufai.cheung@gmail.com (ORCID)

Other contributors:

• Ivan Jacob Agaloos Pesigan r.jeksterslab@gmail.com (ORCID) [contributor]

See Also

Useful links:

• https://sfcheung.github.io/semlbci/

• Report bugs at https://github.com/sfcheung/semlbci/issues

Dataset (CFA, Two Factors, One Standardized Error Variance Close to Zero)

Description

Generated from a two-factor model, with one standardized error variance close to zero.

Usage

cfa_evar_near_zero

Format

A data frame with 120 rows and six variables, x1 to x6

Details

This model is used for examples like this one:

# If fitted by the following model, the standardized
# error variance of `x3` is close to zero.
# Consequently, the R-square of `x3` is close to one:

library(lavaan)
mod <- "f1 =~ x1 + x2 + x3
        f2 =~ x4 + x5 + x6"
fit <- cfa(mod, cfa_evar_near_zero)
summary(fit, standardized = TRUE, rsquare = TRUE)

Examples

print(head(cfa_evar_near_zero), digits = 3)
nrow(cfa_evar_near_zero)

Dataset (CFA, Two Factors, Six Variables)

Description

Generated from a two-factor model with six variables, n = 500

Usage

cfa_two_factors

Format

A data frame with 500 rows and six variables, x1 to x6.

Details

This model is used for examples like this one:

library(lavaan)
mod <- "f1 =~ x1 + x2 + x3
        f2 =~ x4 + x5 + x6"
fit <- cfa(mod, cfa_two_factors)
summary(fit)

Examples

print(head(cfa_two_factors), digits = 3)
nrow(cfa_two_factors)

Dataset (CFA, Two Factors, Six Variables, Two Groups)

Description

Generated from a two-factor model with six variables, n = 500, two groups, n = 250 each.

Usage

cfa_two_factors_mg

Format

A data frame with 500 rows, one grouping variable, gp, six variables, x1 to x6.

Details

This model is used for examples like this one:

library(lavaan)
mod <- "f1 =~ x1 + x2 + x3
        f2 =~ x4 + x5 + x6"
fit <- cfa(mod, cfa_two_factors_mg, group = "gp")
summary(fit)

Examples

print(head(cfa_two_factors_mg), digits = 3)
nrow(cfa_two_factors_mg)
table(cfa_two_factors_mg$gp)

Pre-analysis Check For 'semlbci'

Description

Check the output passed to semlbci()

Usage

check_sem_out(
  sem_out,
  robust = c("none", "satorra.2000"),
  multigroup_ok = TRUE
)

Arguments

Details

It checks whether the model and the estimation method in the sem_out object passed to semlbci() are supported by the current version of semlbci(). This function is to be used by semlbci() but is exported such that the compatibility of an SEM output can be checked directly.

Estimation methods (estimator in lavaan::lavaan()) currently supported:

• Maximum likelihood (ML) and its variants (e.g., MLM, MLR). For methods with robust test statistics (e.g., MLR), only robust LBCIs (robust = "satorra.2000" in calling semlbci()) can be requested.

Estimation methods not yet supported:

• Generalized least squares (GLS).

• Weighted least squares (a.k.a. asymptotically distribution free) (WLS) and its variants (e.g., WLSMV).

• Unweighted least squares (ULS).

• Diagonally weighted least squares (DWLS).

• Other methods not listed.

Models supported:

• Single-group models with continuous variables.

• Multiple-group models with continuous variables.

Models not tested:

• Models with categorical variables.

Models not yet supported:

• Models with formative factors.

• Multilevel models.

Value

A numeric vector of one element. If 0, the model and estimation method are officially supported. If larger than zero, then the model and method are not officially supported but users can still try to use semlbci() on it at their own risks. If less than zero, then the model and/or the method are officially not supported.

The attributes info contains the reason for a value other than zero.

See Also

semlbci(), ci_i_one()

Examples

library(lavaan)
data(cfa_two_factors)
mod <-
"
f1 =~ x1 + x2 + x3
f2 =~ x4 + x5 + x6
"

fit <- sem(mod, cfa_two_factors)

# Should be 0
check_sem_out(fit)

fit2 <- sem(mod, cfa_two_factors, estimator = "DWLS")

# Should be negative because DWLS is officially not supported
check_sem_out(fit2)

fit3 <- sem(mod, cfa_two_factors, estimator = "MLR")

# Should be negative because MLR is supported only if
# robust is set to "satorra.2000"
check_sem_out(fit3)

# Should be zero because robust is set to "satorra.2000"
check_sem_out(fit3, robust = "satorra.2000")

Find a Likelihood-Based Confidence Bound By Root Finding

Description

Find the lower or upper bound of the likelihood-based confidence interval (LBCI) for one parameter in a structural equation model fitted in lavaan::lavaan() using uniroot().

Usage

ci_bound_ur(
  sem_out,
  func,
  ...,
  level = 0.95,
  which = c("lbound", "ubound"),
  interval = NULL,
  progress = FALSE,
  method = "uniroot",
  lrt_method = "default",
  tol = 5e-04,
  root_target = c("chisq", "pvalue"),
  d = 5,
  uniroot_extendInt = switch(which, lbound = "downX", ubound = "upX"),
  uniroot_trace = 0,
  uniroot_maxiter = 1000,
  use_callr = TRUE,
  rs = NULL
)

gen_est_i(i, sem_out, standardized = FALSE)

Arguments

Details

This function is called xby ci_bound_ur_i(). This function is exported because it is a stand-alone function that can be used directly for any function that receives a lavaan object and returns a scalar.

The function ci_bound_ur_i() is a wrapper of this function, with an interface similar to that of ci_bound_wn_i() and returns a cibound-class object. The user-parameter function is generated internally by ci_bound_wn_i().

This function, on the other hand, requires users to supply the function directly through the func argument. This provides the flexibility to find the bound for any function of the model parameter, even one that cannot be easily coded in lavaan model syntax.

Value

The function ci_bound_ur() returns a list with the following elements:

• bound: The bound found.

• optimize_out: THe output of the root finding function, uniroot() for now. (Called optimize_out because an earlier version of this function also uses optimize()).

• sem_out_bound: The lavaan model with the user-defined parameter fixed to the bound.

• lrt: The output of lavaan::lavTestLRT() comparing sem_out and sem_out_bound.

• bound_start: The Wald or delta method confidence bound returned when determining the interval internally.

• user_est: The estimate of the user-defined parameter when determining the interval internally.

The function gen_est_i() returns a special function can inspects the Model slot (and implied slot if necessary) of a modified lavaan object and return the parameter estimate. This function is to be used by ci_bound_ur() or gen_sem_out_userp().

Examples

library(lavaan)
data(simple_med)
dat <- simple_med
mod <-
"
m ~ x
y ~ m
"
fit_med <- lavaan::sem(mod, simple_med, fixed.x = FALSE)
parameterTable(fit_med)
# Create a function to get the second parameter
est_i <- gen_est_i(i = 2, sem_out = fit_med)
# Find the lower bound of the likelihood-based confidence interval
# of the second parameter.
# user_callr should be TRUE or omitted in read research.
# Remove interval in read research. It is added to speed up the example.
out1l <- ci_bound_ur(sem_out = fit_med,
                     func = est_i,
                     which = "lbound",
                     use_callr = FALSE,
                     interval = c(.39070, .39075))
out1l

Likelihood-Based Confidence Bound By Root Finding

Description

Using root finding to find the lower or upper bound of the likelihood-based confidence interval (LBCI) for one parameter in a structural equation model fitted in lavaan::lavaan().

Usage

ci_bound_ur_i(
  i = NULL,
  npar = NULL,
  sem_out = NULL,
  f_constr = NULL,
  which = NULL,
  history = FALSE,
  perturbation_factor = 0.9,
  lb_var = -Inf,
  standardized = FALSE,
  wald_ci_start = !standardized,
  opts = list(),
  ciperc = 0.95,
  ci_limit_ratio_tol = 1.5,
  verbose = FALSE,
  sf = 1,
  sf2 = 0,
  p_tol = 5e-04,
  std_method = "internal",
  bounds = "none",
  xtol_rel_factor = 1,
  ftol_rel_factor = 1,
  lb_prop = 0.05,
  lb_se_k = 3,
  d = 5,
  ...
)

Arguments

Details

Important Notice

This function is not supposed to be used directly by users in typical scenarios. Its interface is user-unfriendly because it should be used through semlbci(). It is exported such that interested users can examine how a confidence bound is found, or use it for experiments or simulations.

Usage

This function is the lowest level function used by semlbci(). semlbci() calls this function once for each bound of each parameter.

For consistency in the interface, most of the arguments in ci_bound_wn_i() are also included in this function, even those not used internally.

Algorithm

This function, unlike ci_bound_wn_i(), use a simple root finding algorithm. Basically, it tries fixing the target parameter to different values until the likelihood ratio test p-value, or the corresponding chi-square difference, is equal to the value corresponding to the desired level of confidence. (Internally, the difference between the p-value and the target p-value, that for the chi-square difference, is the function value.)

For finding the bound, this algorithm can be inefficient compared to the one proposed by Wu and Neale (2012). The difference can be less than one second versus 10 seconds. It is included as a backup algorithm for parameters which are difficult for the method by Wu and Neale.

Internally, it uses uniroot() to find the root.

Limitation(s)

This function does not handle an estimate close to an attainable bound using the method proposed by Wu and Neale (2012). Use it for such parameters with cautions.

Value

A cibound-class object which is a list with three elements:

• bound: A single number. The value of the bound located. NA is the search failed for various reasons.

• diag: A list of diagnostic information.

• call: The original call.

A detailed and organized output can be printed by the default print method (print.cibound()).

References

Wu, H., & Neale, M. C. (2012). Adjusted confidence intervals for a bounded parameter. Behavior Genetics, 42(6), 886-898. doi:10.1007/s10519-012-9560-z

See Also

print.cibound(), semlbci(), ci_i_one(); see ci_bound_wn_i() on the version for the method by Wu and Neale (2012).

Examples

library(lavaan)
data(simple_med)
dat <- simple_med
mod <-
"
m ~ x
y ~ m
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
# Remove `opts` in real cases.
# The options are added just to speed up the example
out1l <- ci_bound_ur_i(i = 1,
                       sem_out = fit_med,
                       which = "lbound",
                       opts = list(use_callr = FALSE,
                                   interval = c(0.8277, 0.8278)))
out1l

Likelihood-based Confidence Bound By Wu-Neale-2012

Description

User the method proposed by Wu and Neale (2012) to find the lower or upper bound of the likelihood-based confidence interval (LBCI) for one parameter in a structural equation model fitted in lavaan::lavaan().

Usage

ci_bound_wn_i(
  i = NULL,
  npar = NULL,
  sem_out = NULL,
  f_constr = NULL,
  which = NULL,
  history = FALSE,
  perturbation_factor = 0.9,
  lb_var = -Inf,
  standardized = FALSE,
  wald_ci_start = !standardized,
  opts = list(),
  ciperc = 0.95,
  ci_limit_ratio_tol = 1.5,
  verbose = FALSE,
  sf = 1,
  sf2 = 0,
  p_tol = 5e-04,
  std_method = "internal",
  bounds = "none",
  xtol_rel_factor = 1,
  ftol_rel_factor = 1,
  lb_prop = 0.05,
  lb_se_k = 3,
  try_harder = 0,
  fit_lb = -Inf,
  fit_ub = +Inf,
  timeout = 300,
  ...
)

Arguments

Details

Important Notice

This function is not supposed to be used directly by users in typical scenarios. Its interface is user-unfriendly because it should be used through semlbci(). It is exported such that interested users can examine how a confidence bound is found, or use it for experiments or simulations.

Usage

This function is the lowest level function used by semlbci(). semlbci() calls this function once for each bound of each parameter. To use it, set_constraint() needs to be called first to create the equality constraint required by the algorithm proposed by Wu and Neale (2012).

Algorithm

This function implements the algorithm presented in Wu and Neale (2012; see also Pek & Wu, 2015, Equation 12) that estimates all free parameters in the optimization.

Limitation(s)

This function does not yet implement the method by Wu and Neale (2012) for an estimate close to an attainable bound.

Value

A cibound-class object which is a list with three elements:

• bound: A single number. The value of the bound located. NA is the search failed for various reasons.

• diag: A list of diagnostic information.

• call: The original call.

A detailed and organized output can be printed by the default print method (print.cibound()).

References

Pek, J., & Wu, H. (2015). Profile likelihood-based confidence intervals and regions for structural equation models. Psychometrika, 80(4), 1123-1145. doi:10.1007/s11336-015-9461-1

Wu, H., & Neale, M. C. (2012). Adjusted confidence intervals for a bounded parameter. Behavior Genetics, 42(6), 886-898. doi:10.1007/s10519-012-9560-z

See Also

print.cibound(), semlbci(), ci_i_one()

Examples

data(simple_med)
dat <- simple_med

mod <-
"
m ~ x
y ~ m
"

fit_med <- lavaan::sem(mod, simple_med, fixed.x = FALSE)

fn_constr0 <- set_constraint(fit_med)

out1l <- ci_bound_wn_i(i = 1,
                       npar = 5,
                       sem_out = fit_med,
                       f_constr = fn_constr0,
                       which = "lbound")
out1l

Likelihood-Based Confidence Bound for One Parameter

Description

Find the likelihood-based confidence bound for one parameter.

Usage

ci_i_one(
  i,
  which = NULL,
  sem_out,
  method = c("wn", "ur"),
  standardized = FALSE,
  robust = "none",
  sf_full = NA,
  sf_args = list(),
  sem_out_name = NULL,
  try_k_more_times = 0,
  ...
)

Arguments

Details

Important Notice

This function is not supposed to be used directly by users in typical scenarios. Its interface is user-unfriendly because it should be used through semlbci(). It is exported such that interested users can examine how a confidence bound is found, or use it for experiments or simulations.

Usage

ci_i_one() is the link between semlbci() and the lowest level function (currently ci_bound_wn_i()). When called by semlbci() to find the bound of a parameter, ci_i_one() calls a function (ci_bound_wn_i() by default) one or more times to find the bound (limit) for a likelihood-based confidence interval.

Value

A list of the following elements.

• bound: The bound located. NA if the search failed.

• diags: Diagnostic information.

• method: Method used. Currently only "wn" is the only possible value.

• times: Total time used in the search.

• sf_full: The scaling and shift factors used.

• ci_bound_i_out: The original output from ci_bound_wn_i().

• attempt_lb_var: How many attempts used to reduce the lower bounds of free variances.

• attempt_more_times: How many additional attempts used to search for the bounds. Controlled by try_k_more_times.

See Also

semlbci(), ci_bound_wn_i()

Examples

data(simple_med)

library(lavaan)
mod <-
"
m ~ x
y ~ m
"
fit_med <- lavaan::sem(mod, simple_med, fixed.x = FALSE)

parameterTable(fit_med)

# Find the LBCI for the first parameter
# The method "wn" needs the constraint function.
# Use set_constraint() to generate this function:
fn_constr0 <- set_constraint(fit_med)

# Call ci_i to find the bound, the lower bound in this example.
# The constraint function, assigned to f_constr, is passed
# to ci_bound_wn_i().
# npar is an argument for ci_bound_wn_i().
out <- ci_i_one(i = 1,
                which = "lbound",
                sem_out = fit_med,
                npar = 5,
                f_constr = fn_constr0)
out$bounds

Check The Order of Bounds in a List of semlbci Objects

Description

Check whether the LBCIs in a list of semlbci-class of objects are consistent with their levels of confidence.

Usage

ci_order(semlbci_list)

## S3 method for class 'ci_order'
print(x, digits = 3, ...)

Arguments

Value

A ci_order-class object with a print method print.ci_order(). The number of rows is equal to the number of parameters in semlbci_list, and the columns stores the confidence limits from the list, ordered according to the level of confidence.

x is returned invisibly. Called for its side effect.

Methods (by generic)

• print(ci_order): The print method of the output of ci_order().

Author(s)

Shu Fai Cheung https://orcid.org/0000-0002-9871-9448

See Also

nearby_levels(), semlbci()

Examples

library(lavaan)
mod <-
"
m ~ x
y ~ m
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
lbci_fit <- semlbci(fit_med)
lbci_fit_nb <- nearby_levels(lbci_fit,
                             ciperc_levels = c(-.050, .050))

# Check the order of the confidence bounds.
# A confidence interval with a higher level of confidence
# should enclose a confidence interval with
# a lower level of confidence.
ci_order(lbci_fit_nb)

Confidence Intervals for a 'smelbci' Object

Description

Return the confidence intervals of the parameters in the output of semlbci().

Usage

## S3 method for class 'semlbci'
confint(object, parm, level = 0.95, ...)

Arguments

Details

It returns the likelihood-based confidence intervals in the output of semlbci().

Value

A two-column matrix of the confidence intervals.

Author(s)

Shu Fai Cheung https://orcid.org/0000-0002-9871-9448

See Also

semlbci()

Examples

library(lavaan)
mod <-
"
m ~ a*x
y ~ b*m
ab := a * b
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
p_table <- parameterTable(fit_med)
p_table
lbci_med <- semlbci(fit_med,
                    pars = "ab :=")
lbci_med

confint(lbci_med)

Create a Wrapper To Be Used in 'lavaan' Models

Description

Make a function on lavaan object usable in a lavaan model syntax.

Usage

gen_userp(func, sem_out)

gen_sem_out_userp(
  userp,
  sem_out,
  userp_name = "semlbciuserp1234",
  fix = TRUE,
  control_args = list(),
  iter.max = 10000,
  max_attempts = 5
)

Arguments

Details

gen_userp

There are cases in which we want to create a user parameter which is a function of other free parameters, computed by a function. However such a function may work only on a lavaan object.

If the target function works by extracting parameter estimates stored in the Model slot and/or the implied slot, then gen_userp() can be used to convert it to a function that retrieves the parameter estimates when being called by lavaan::lavaan() or its wrappers, modifies the stored lavaan object using lavaan::lav_model_set_parameters() and lavaan::lav_model_implied() to change the estimates, and call the target function.

Note that this is an unconventional way to define a user parameter and the generated function should always be checked to see whether it works as expected.

As shown in the examples, the parameter computed this may not have standard error nor p-value.

The main purpose is for the point estimate, for searching the likelihood-based confidence bound using ci_bound_ur() and ci_bound_ur_i().

Note that the target function specified in func should work directly on the parameter estimates stored in the Model slot and then get the estimates using lavaan::lav_model_get_parameters(). Functions that work on the unmodified output generated by lavaan::lavaan() usually do not work.

Users are not recommended to use gen_userp() and gen_sem_out_userp() directly because they require unconventional way to extract parameter estimates from a lavaan model. However, developers may use them to include functions they wrote in a lavaan model. This is the technique used by ci_bound_ur_i() to constrain any parameter in a model to an arbitrary value.

gen_sem_out_userp

The function gen_sem_out_userp() is to be used internally for generating a function for searching a likelihood-based confidence bound. It is exported because it needs to be run in an fresh external R process, usually created by callr in other internal functions.

Value

gen_userp

It returns a function that accepts a numeric vector of length equals to the number of free parameters in sem_out, and returns a scalar which is the output of func. If this vector is not supplied, it will try to find it in the parent.frame(). This is how it works inside a lavaan model.

gen_sem_out_userp

If fix is TRUE, it returns a function with these arguments:

• target: The value to which the user-defined parameter will be fixed to.

• verbose: If TRUE, additional information will be printed when fitting the model.

• control: The values to be passed as a list to the argument of the same name in lavaan::lavaan().

• seed: Numeric. If supplied, it will be used in set.seed() to initialize the random number generator. Necessary to reproduce some results because random numbers are used in some steps in lavaan. If NULL, the default, set.seed() will not be called.

If fix is 'FALSE, then it returns a function with optional arguments that will be ignored, Calling it will simply fit the modified model to the data. Useful for getting the value of the user-defined parameter.

Examples

library(lavaan)

data(simple_med)
dat <- simple_med
mod <-
"
m ~ a*x
y ~ b*m
ab := a*b
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
parameterEstimates(fit_med)

# A trivial example for verifying the results
my_ab <- function(object) {
    # Need to use lav_model_get_parameters()
    # because the object is only a modified
    # lavaan-object, not one directly
    # generated by lavaan function
    est <- lavaan::lav_model_get_parameters(object@Model, type = "user")
    unname(est[1] * est[2])
  }

# Check the function
my_ab(fit_med)
coef(fit_med, type = "user")["ab"]

# Create the function
my_userp <- gen_userp(func = my_ab,
                      sem_out = fit_med)
# Try it on the vector of free parameters
my_userp(coef(fit_med))

# Generate a modified lavaan model
fit_userp <- gen_sem_out_userp(userp = my_userp,
                               userp_name = "my_userp",
                               sem_out = fit_med)

# This function can then be used in the model syntax.

# Note that the following example only work when called inside the
# workspace or inside other functions such as ci_bound_ur()`
# and `ci_bound_ur_i()` because `lavaan::sem()` will
# search `my_userp()` in the global environment.

# Therefore, the following lines are commented out.
# They should be run only in a "TRUE" interactive
# session.

# mod2 <-
# "
# m ~ x
# y ~ m
# ab := my_userp()
# "
# fit_med2 <- sem(mod2, simple_med, fixed.x = FALSE)
# parameterEstimates(fit_med2)
#
# # Fit the model with the output of the function, a*b
# # fixed to .50
#
# fit_new <- fit_userp(.50)
#
# # Check if the parameter ab is fixed to .50
# parameterEstimates(fit_new)

A 'cibound' Output From a 'semlbci' Object

Description

Get the cibound output of a bound from a semlbci object, the output of semlbci().

Usage

get_cibound(x, row_id, which = c("lbound", "ubound"))

get_cibound_status_not_0(x)

Arguments

Details

The function get_cibound() returns the original output of ci_bound_wn_i() for a bound. Usually for diagnosis.

The function get_cibound_status_not_0() checks the status code of each bound, and returns the cibound outputs of bounds with status code not equal to zero (i.e., something wrong in the search). Printing it can print the diagnostic information for all bounds that failed in the search.

Value

get_cibound() returns a cibound-class object. See ci_bound_wn_i() for details. get_cibound_status_not_0() returns a list of cibound-class objects with status not equal to zero. If all bounds have status equal to zero, it returns an empty list.

Author(s)

Shu Fai Cheung https://orcid.org/0000-0002-9871-9448

See Also

semlbci()

Examples

library(lavaan)
mod <-
"
m ~ a*x
y ~ b*m
ab := a * b
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
p_table <- parameterTable(fit_med)
p_table
lbci_med <- semlbci(fit_med,
                    pars = c("ab :="))
lbci_med

# Get the output of ci_bound_wn_i() of the lower
# bound of the LBCI for the indirect effect:
get_cibound(lbci_med, row_id = 6, which = "lbound")

# Get the output of ci_bound_wn_i() of the upper
# bound of the LBCI for the indirect effect:
get_cibound(lbci_med, row_id = 6, which = "ubound")

Log Profile likelihood of a Parameter

Description

These functions compute the log profile likelihood of a parameter when it is fixed to a value or a range of values

Usage

loglike_compare(
  sem_out,
  semlbci_out = NULL,
  par_i,
  confidence = 0.95,
  n_points = 21,
  start = "default",
  try_k_more = 5,
  parallel = FALSE,
  ncpus = parallel::detectCores(logical = FALSE) - 1,
  use_pbapply = TRUE
)

loglike_range(
  sem_out,
  par_i,
  confidence = 0.95,
  n_points = 21,
  interval = NULL,
  verbose = FALSE,
  start = "default",
  try_k_more = 5,
  parallel = FALSE,
  ncpus = parallel::detectCores(logical = FALSE) - 1,
  use_pbapply = TRUE
)

loglike_point(
  theta0,
  sem_out,
  par_i,
  verbose = FALSE,
  start = "default",
  try_k_more = 5
)

loglike_quad_range(
  sem_out,
  par_i,
  confidence = 0.95,
  n_points = 21,
  interval = NULL,
  parallel = FALSE,
  ncpus = parallel::detectCores(logical = FALSE) - 1,
  use_pbapply = TRUE,
  try_k_more = 5,
  start = "default"
)

loglike_quad_point(theta0, sem_out, par_i)

loglike_compare_ur(
  sem_out,
  semlbci_out = NULL,
  par_i,
  confidence = 0.95,
  n_points = 21,
  standardized = FALSE,
  parallel = FALSE,
  ncpus = parallel::detectCores(logical = FALSE) - 1,
  use_pbapply = TRUE,
  loadbalancing = TRUE
)

loglike_range_ur(
  sem_out,
  par_i,
  standardized = FALSE,
  confidence = 0.95,
  n_points = 21,
  interval = NULL,
  verbose = FALSE,
  parallel = FALSE,
  ncpus = parallel::detectCores(logical = FALSE) - 1,
  use_pbapply = TRUE,
  loadbalancing = TRUE
)

loglike_point_ur(theta0, sem_out, par_i, standardized = FALSE, verbose = FALSE)

loglike_quad_range_ur(
  sem_out,
  par_i,
  confidence = 0.95,
  standardized = FALSE,
  n_points = 21,
  interval = NULL,
  parallel = FALSE,
  ncpus = parallel::detectCores(logical = FALSE) - 1,
  use_pbapply = TRUE,
  loadbalancing = TRUE
)

loglike_quad_point_ur(theta0, sem_out, par_i, standardized = FALSE)

Arguments

Details

It uses the methods presented in Pawitan (2013) to compute and visualize the log profile likelihood of a parameter in a structural equation model when this parameter is fixed to a value or a range of values. loglike_range() and loglike_point() compute the so-called "true" log profile likelihood, while loglike_quad_range() and loglike_quad_point() approximate the log profile likelihood by a quadratic function.

These functions are for creating illustrative examples and learning only, not for research use. Therefore, they are not as versatile as semlbci() in the types of models and parameters supported. They can be used for free parameters and user-defined parameters not involved in any constraints. Only a model fitted by maximum likelihood is supported.

They will not check whether the computation is appropriate for a model. It is the responsibility of the users to ensure that the computation is appropriate for the model and parameter.

In version 0.11.2.1, added variants of the function, with suffix ⁠_ur⁠, which use root finding ("ur" in semlbci()). These variants are slower to run but can be used for parameter in the standardized solution. Therefore, they can used to generate plots for parameters such as standardized regression paths and correlations.

Value

loglike_compare() calls loglike_range() and loglike_quad_range() and returns their results in a loglike_compare-class object, a list with these elements:

• quadratic: The output of loglike_quad_range().

• loglikelihood: The output of loglike_range().

• pvalue_quadratic: The likelihood ratio test p-values at the quadratic approximation confidence bounds.

• pvalue_loglikelihood: The likelihood ratio test p-values at the likelihood-based confidence bounds.

• est: The point estimate of the parameter in sem_out.

loglike_compare-class object has a plot method (plot.loglike_compare()) that can be used to plot the log profile likelihood.

loglike_point() returns a list with these elements:

• loglike: The log profile likelihood of the parameter when it is fixed to theta0.

• pvalue: The p-values based on the likelihood ratio difference test between the original model and the model with the parameter fixed to theta0.

• fit: A lavaan::lavaan object. The original model with the parameter fixed to theta0.

• lrt: The output of lavaan::lavTestLRT(), comparing the original model to the model with the parameter fixed to theta0.

loglike_quad_range() returns a data frame with these columns:

• theta: The values to which the parameter is fixed to.

• loglike: The log profile likelihood values of the parameter using quadratic approximation.

• pvalue: The p-values based on the likelihood ratio difference test between the original model and the model with the parameter fixed to theta.

loglike_quad_point() returns a single number of the class lavaan.vector (because it is the output of lavaan::fitMeasures()). This number is the quadratic approximation of the log profile likelihood when the parameter is fixed to theta0.

loglike_range() returns a data frame with these columns:

• theta: The values to which the parameter is fixed to.

• loglike: The log profile likelihood at theta.

• pvalue: The p-values based on the likelihood ratio difference test between the original model and model with the parameter fixed to theta.

Functions

• loglike_compare(): Generates points for log profile likelihood and quadratic approximation, by calling the helper functions loglike_range() and loglike_quad_range().

• loglike_range(): Find the log profile likelihood for a range of values.

• loglike_point(): Find the log likelihood at a value.

• loglike_quad_range(): Find the approximated log likelihood for a range of values.

• loglike_quad_point(): Find the approximated log likelihood at a value.

• loglike_compare_ur(): Generates points for log profile likelihood and quadratic approximation using root finding, by calling the helper functions loglike_range_ur() and loglike_quad_range_ur().

• loglike_range_ur(): Find the log profile likelihood for a range of values using root finding.

• loglike_point_ur(): Find the log likelihood at a value.

• loglike_quad_range_ur(): Find the approximated log likelihood for a range of values using root finding.

• loglike_quad_point_ur(): Find the approximated log likelihood at a value. Support a parameter in the standardized solution.

References

Pawitan, Y. (2013). In all likelihood: Statistical modelling and inference using likelihood. Oxford University Press.

See Also

plot.loglike_compare()

Examples

## loglike_compare

library(lavaan)
data(simple_med)
dat <- simple_med
mod <-
"
m ~ a * x
y ~ b * m
ab := a * b
"
fit <- lavaan::sem(mod, simple_med, fixed.x = FALSE)

# 4 points are used just for illustration
# At least 21 points should be used for a smooth plot
# Remove try_k_more in real applications. It is set
# to zero such that this example does not take too long to run.
# use_pbapply can be removed or set to TRUE to show the progress.
ll_a <- loglike_compare(fit, par_i = "m ~ x", n_points = 4,
                        try_k_more = 0,
                        use_pbapply = FALSE)
plot(ll_a)

# See the vignette "loglike" for an example for the
# indirect effect.


## loglike_range

# Usually not to be used directly.
# Used by loglike_compare().
# 3 points are used just for illustration
ll_1 <- loglike_range(fit, par_i = "y ~ m", n_points = 2)
head(ll_1)


## loglike_point

# Usually not to be used directly.
# Used by loglike_compare().
llp_1 <- loglike_point(theta0 = 0.3, sem_out = fit, par_i = "y ~ m")
llp_1$loglike
llp_1$pvalue
llp_1$lrt



## loglike_quad_range

# Usually not to be used directly.
# Used by loglike_compare().
# 2 points are used just for illustration
lq_1 <- loglike_quad_range(fit, par_i = "y ~ m", n_points = 2)
head(lq_1)



## loglike_quad_point

# Usually not to be used directly.
# Used by loglike_compare().
lqp_1 <- loglike_quad_point(theta0 = 0.3, sem_out = fit, par_i = "y ~ m")
lqp_1

Dataset (SEM, Three Factors, Nine Variables, Mediation)

Description

Generated from a three-factor model with nine variables, n = 150

Usage

mediation_latent

Format

A data frame with 150 rows and nine variables:

x1

x1

x2

x2

x3

x3

x4

x4

x5

x5

x6

x6

x7

x7

x8

x8

x9

x9

Details

This model is used for examples like this one:

mod <-
"
fx =~ x1 + x2 + x3
fm =~ x4 + x5 + x6
fy =~ x7 + x8 + x9
fm ~ a*fx
fy ~ b*fm + cp*fx
ab := a*b
"
fit <- lavaan::sem(mod, mediation_latent)

Examples

print(head(mediation_latent), digits = 3)
nrow(mediation_latent)

Dataset (SEM, Three Factors, Nine Variables, Mediation, Skewed)

Description

Generated from a three-factor model with nine variables, n = 150, with some observed variables positively skewed.

Usage

mediation_latent_skewed

Format

A data frame with 150 rows and nine variables:

x1

x1

x2

x2

x3

x3

x4

x4

x5

x5

x6

x6

x7

x7

x8

x8

x9

x9

Details

This model is used for examples like this one:

mod <-
"
fx =~ x1 + x2 + x3
fm =~ x4 + x5 + x6
fy =~ x7 + x8 + x9
fm ~ a*fx
fy ~ b*fm + cp*fx
ab := a*b
"
fit <- lavaan::sem(mod, mediation_latent)

Examples

print(head(mediation_latent_skewed), digits = 3)
nrow(mediation_latent_skewed)

LBCI Bounds of Nearby Levels of Confidence

Description

Find LBCIs with levels of confidence different from those stored in a semlbci- class object.

Usage

nearby_levels(x, ciperc_levels = c(-0.025, 0.025), ciperc_range = c(0.6, 0.99))

Arguments

Details

It receives a semlbci-class object, gets the original level of confidence, generates one or more levels of confidence different from this level by certain amounts, and repeats the original call to semlbci() with these levels of confidence. The results are returned as a list of class semlbci_list, with the originalsemlbci-class included.

Value

A semlbci_list-class object, which is simply a named list of semlbci-class object, names being the levels of confidence.

Author(s)

Shu Fai Cheung https://orcid.org/0000-0002-9871-9448

See Also

semlbci(), ci_order()

Examples

library(lavaan)
mod <-
"
m ~ x
y ~ m
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
lbci_fit <- semlbci(fit_med)
lbci_fit_nb <- nearby_levels(lbci_fit,
                             ciperc_levels = c(-.050, .050))
names(lbci_fit_nb)
# Check the order of the confidence bounds.
# A confidence interval with a higher level of confidence
# should enclose a confidence interval with
# a lower level of confidence.
ci_order(lbci_fit_nb)

Plot the Output of 'loglike_compare()'

Description

Visualize the log profile likelihood of a parameter fixed to values in a range.

Usage

## S3 method for class 'loglike_compare'
plot(
  x,
  y,
  type = c("ggplot2", "default"),
  size_label = 4,
  size_point = 4,
  nd_theta = 3,
  nd_pvalue = 3,
  size_theta = 4,
  size_pvalue = 4,
  add_pvalues = FALSE,
  ...
)

Arguments

Details

Given the output of loglike_compare(), it plots the log profile likelihood based on quadratic approximation and that based on the original log-likelihood. The log profile likelihood is scaled to have a maximum of zero (at the point estimate) as suggested by Pawitan (2013).

Value

Nothing if type = "default", the generated ggplot2::ggplot() graph if type = "ggplot2".

References

Pawitan, Y. (2013). In all likelihood: Statistical modelling and inference using likelihood. Oxford University Press.

Examples

## loglike_compare

library(lavaan)
data(simple_med)
dat <- simple_med
mod <-
"
m ~ a * x
y ~ b * m
ab := a * b
"
fit <- lavaan::sem(mod, simple_med, fixed.x = FALSE)

# Four points are used just for illustration
# At least 21 points should be used for a smooth plot
# Remove try_k_more in real applications. It is set
# to run such that this example is not too slow.
# use_pbapply can be removed or set to TRUE to show the progress.
ll_a <- loglike_compare(fit, par_i = "m ~ x", n_points = 4,
                        try_k_more = 0,
                        use_pbapply = FALSE)

plot(ll_a)
plot(ll_a, add_pvalues = TRUE)

# See the vignette "loglike" for an example for the
# indirect effect.

Print Method of a 'cibound'-class Object

Description

Print the diagnostic information of a cibound-class object.

Usage

## S3 method for class 'cibound'
print(x, digits = 5, ...)

Arguments

Details

This is the print method for the output of ci_bound_wn_i(), a cibound-class object. It prints the diagnostic information on the bound being found and the search process.

Value

x is returned invisibly. Called for its side effect.

Examples

data(simple_med)
dat <- simple_med

mod <-
"
m ~ x
y ~ m
"

fit_med <- lavaan::sem(mod, simple_med, fixed.x = FALSE)

fn_constr0 <- set_constraint(fit_med)

out1l <- ci_bound_wn_i(i = 1,
                       npar = 5,
                       sem_out = fit_med,
                       f_constr = fn_constr0,
                       which = "lbound")

# Print the output
out1l

Print Method of a 'semlbci' Object

Description

Prints the results of a semlbci object, the output of semlbci().

Usage

## S3 method for class 'semlbci'
print(
  x,
  digits = 3,
  annotation = TRUE,
  time = FALSE,
  verbose = FALSE,
  verbose_if_needed = TRUE,
  drop_no_lbci = TRUE,
  output = c("table", "text", "lavaan"),
  sem_out = NULL,
  lbci_only = drop_no_lbci,
  ratio_digits = 1,
  se = TRUE,
  zstat = TRUE,
  pvalue = TRUE,
  boot.ci.type = "perc",
  ...
)

Arguments

Details

Prints the results of semlbci() as a table.

Value

x is returned invisibly. Called for its side effect.

Author(s)

Shu Fai Cheung https://orcid.org/0000-0002-9871-9448

See Also

semlbci()

Examples

library(lavaan)
mod <-
"
m ~ a*x
y ~ b*m
ab := a * b
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
p_table <- parameterTable(fit_med)
p_table
lbci_med <- semlbci(fit_med,
                    pars = c("ab :="))
lbci_med

print(lbci_med, verbose_if_needed = FALSE)

print(lbci_med, verbose = TRUE)

print(lbci_med, time = TRUE)

print(lbci_med, annotation = FALSE)

print(lbci_med, digits = 4)

# Text output

print(lbci_med, output = "lavaan", sem_out = fit_med)

print(lbci_med, output = "lavaan", sem_out = fit_med, lbci_only = FALSE)

print(lbci_med, output = "lavaan", sem_out = fit_med, lbci_only = FALSE,
      se = FALSE, zstat = FALSE, pvalue = FALSE)

Dataset (Six Variables, One Correlation Close to One)

Description

Generated from a regression model six variables, x4~~x5 correlation close to one.

Usage

reg_cor_near_one

Format

A data frame with 100 rows and six variables:

x1

x1

x2

x2

x3

x3

x4

x4, with correlation with x5 nearly equal to 1

x5

x5, with correlation with x4 nearly equal to 1

y

y, the dependent variable

Details

This model is used for examples like this one:

out <- lm(y ~ x1 + x2 + x3 + x4 + x5, reg_cor_near_one)
summary(out)
cor(reg_cor_near_one[, c("x4", "x5")])

Examples

print(head(reg_cor_near_one), digits = 3)
nrow(reg_cor_near_one)

Likelihood-Based Confidence Interval

Description

Find the likelihood-based confidence intervals (LBCIs) for selected free parameters in an SEM output.

Usage

semlbci(
  sem_out,
  pars = NULL,
  include_user_pars = TRUE,
  remove_variances = TRUE,
  remove_intercepts = TRUE,
  ciperc = 0.95,
  standardized = FALSE,
  method = c("wn", "ur"),
  robust = c("none", "satorra.2000"),
  try_k_more_times = 2,
  semlbci_out = NULL,
  check_fit = TRUE,
  ...,
  parallel = FALSE,
  ncpus = 2,
  use_pbapply = TRUE,
  loadbalancing = TRUE
)

Arguments

Details

semlbci() finds the positions of the selected parameters in the parameter table and then calls ci_i_one() once for each of them. For the technical details, please see ci_i_one() and the functions it calls to find a confidence bound, currently ci_bound_wn_i(). ci_bound_wn_i() uses the approach proposed by Wu and Neale (2012) and illustrated by Pek and Wu (2015).

It supports updating an output of semlbci() by setting semlbci_out. This allows forming LBCIs for some parameters after those for some others have been formed.

If possible, parallel processing should be used (see parallel and ncpus), especially for a model with many parameters.

If the search for some of the confidence bounds failed, with NA for the bounds, try increasing try_k_more_times.

The SEM output will first be checked by check_sem_out() to see whether the model and the estimation method are supported. To skip this test (e.g., for testing or experimenting with some models and estimators), set check_fit to FALSE.

Examples and technical details can be found at Cheung and Pesigan (2023), the website of the semlbci package (https://sfcheung.github.io/semlbci/), and the technical appendices at (https://sfcheung.github.io/semlbci/articles/).

It currently supports lavaan::lavaan outputs only.

Value

A semlbci-class object similar to the parameter table generated by lavaan::parameterEstimates(), with the LBCIs for selected parameters added. Diagnostic information, if requested, will be included in the attributes. See print.semlbci() for options available.

Author(s)

Shu Fai Cheung https://orcid.org/0000-0002-9871-9448

References

Cheung, S. F., & Pesigan, I. J. A. (2023). semlbci: An R package for forming likelihood-based confidence intervals for parameter estimates, correlations, indirect effects, and other derived parameters. Structural Equation Modeling: A Multidisciplinary Journal, 30(6), 985–999. doi:10.1080/10705511.2023.2183860

Falk, C. F. (2018). Are robust standard errors the best approach for interval estimation with nonnormal data in structural equation modeling? Structural Equation Modeling: A Multidisciplinary Journal, 25(2), 244-266. doi:10.1080/10705511.2017.1367254

Pek, J., & Wu, H. (2015). Profile likelihood-based confidence intervals and regions for structural equation models. Psychometrika, 80(4), 1123-1145. doi:10.1007/s11336-015-9461-1

Wu, H., & Neale, M. C. (2012). Adjusted confidence intervals for a bounded parameter. Behavior Genetics, 42(6), 886-898. doi:10.1007/s10519-012-9560-z

Pritikin, J. N., Rappaport, L. M., & Neale, M. C. (2017). Likelihood-based confidence intervals for a parameter with an upper or lower bound. Structural Equation Modeling: A Multidisciplinary Journal, 24(3), 395-401. doi:10.1080/10705511.2016.1275969

See Also

print.semlbci(), confint.semlbci(), ci_i_one(), ci_bound_wn_i()

Examples

library(lavaan)
mod <-
"
m ~ a*x
y ~ b*m
ab := a * b
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
p_table <- parameterTable(fit_med)
p_table
lbci_med <- semlbci(fit_med,
                    pars = c("m ~ x",
                             "y ~ m",
                             "ab :="))
lbci_med

Equality Constraint for Finding the LBCI by Wu-Neale-2012

Description

Create the equality constraint for finding the likelihood-based confidence interval (LBCI) by the Wu-Neale-2012 method.

Usage

set_constraint(sem_out, ciperc = 0.95)

Arguments

Details

Important Notice

This function is not supposed to be used directly by users in typical scenarios. Its interface is user-unfriendly because it should be used through semlbci(). It is exported such that interested users can examine how a confidence bound is found, or use it for experiments or simulations.

Usage

The Wu-Neale-2012 method uses a simple objective function that is optimized with an equality constraint. set_constraint() generates the equality constraint function to be used by ci_bound_wn_i().

It currently supports lavaan::lavaan outputs only.

Value

An equality constraint function to be used by ci_bound_wn_i().

Examples

library(lavaan)
data(simple_med)
dat <- simple_med
mod <-
"
m ~ x
y ~ m
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)

fn_constr0 <- set_constraint(fit_med)
out <- fn_constr0(coef(fit_med), sem_out = fit_med)
out
lavTech(fit_med, "optim")$fx

Dataset (Simple Mediation Model)

Description

Generated from a simple mediation model, n = 200

Usage

simple_med

Format

A data frame with 200 rows and three variables:

x

x, the independent variable

m

m, the mediator

y

y, the dependent variable

Details

This model is used for examples like this one:

library(lavaan)
mod <- "m ~ x
        y ~ m"
fit <- cfa(mod, simple_med)
summary(fit)

Examples

print(head(simple_med), digits = 3)
nrow(simple_med)

Dataset (Simple Mediation Model, Two Groups)

Description

Generated from a simple mediation model, n = 200, two groups, n = 100 each.

Usage

simple_med_mg

Format

A data frame with 500 rows and four variables:

gp

gp, the grouping variable

x

x, the independent variable

m

m, the mediator

y

y, the dependent variable

Details

This model is used for examples like this one:

library(lavaan)
mod <- "m ~ x
        y ~ m"
fit <- sem(mod, simple_med_mg, gp = "group")
summary(fit)

Examples

print(head(simple_med_mg), digits = 3)
nrow(simple_med_mg)
table(simple_med_mg$gp)

Parameter Positions From lavaan Syntax

Description

Converts lavaan syntax to positions in the model parameter table.

Usage

syntax_to_i(syntax, sem_out)

Arguments

Details

syntax_to_i() converts a vector of strings, in lavaan syntax, to the positions in the parameter table of a lavaan::lavaan fit object.

Each element in the vector should have left hand side (lhs), operator (op), and/or right hand side (rhs). For example:all.x

• "m ~ x" denotes the coefficient of the path from x to m.

• "y ~~ x" denotes the covariance between y and x.

For user-defined parameters, only lhs and op will be interpreted. For example:

• To specify the user parameter ab, both "ab := ..." and "ab :=" will do, ... the definition of ab in the model. The right-hand side will be ignored.

To denote a labelled parameters, such as "y ~ a*x", treat it as a user-defined parameters and use ⁠:=⁠, e.g., "a :=" in this example.

For multiple-group models, if a parameter is specified as in a single-group models, then this parameter in all groups will be selected. For example:all.x

• If a model has three groups, "y ~ x" denotes this path parameter in all three groups, and it will be converted to three row numbers.

To select the parameter in a specific group, "multiply" the right-hand-side variable by the group number. For example:

• "y ~ 2*x" denotes the path coefficient from x to y in Group 2.

To denote the parameters in more than one group, multiply the right-hand side variable by a vector of number. For example:all.x

• "f1 =~ c(2,3)*x2" denotes the factor loading of x2 on f1 in Group 2 and Group 3.

Elements that cannot be converted to a parameter in the parameter table will be ignored.

Currently supports lavaan::lavaan outputs only.

Value

A numeric vector of positions (row numbers) in the parameter table.

Examples

library(lavaan)
data(simple_med)
mod <-
"
m ~ a*x
y ~ b*m
ab:= a*b
asq:= a^2
"
fit_med <- sem(mod, simple_med, fixed.x = FALSE)
p_table <- parameterTable(fit_med)

pars <- c("m ~ x",
          "y ~ m",
          "asq := 1",
          "ab  := 2")
out <- syntax_to_i(pars, fit_med)
out
p_table[out, ]