Nothing
#' @title Create a Wrapper To Be Used
#' in 'lavaan' Models
#'
#' @description Make a function on
#' `lavaan` object usable in a `lavaan`
#' model syntax.
#'
#' @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.
#'
#' @param func A function that receives a
#' `lavaan`-object and returns a scalar.
#' See Details on the restriction on this
#' function.
#'
#' @param sem_out A `lavaan`-class object
#' to be modified.
#'
#' @return
#'
#' ## `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.
#'
#' @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)
#'
#'
#' @rdname gen_userp
#'
#' @export
gen_userp <- function(func,
sem_out) {
stopifnot(is.function(func))
stopifnot(inherits(sem_out, "lavaan"))
out <- function(.x., ...) {
if (missing(.x.)) {
.x. <- get(".x.", parent.frame())
}
# Just in case ...
force(sem_out)
sem_out@Model <- lavaan::lav_model_set_parameters(sem_out@Model, .x.)
sem_out@implied <- lavaan::lav_model_implied(sem_out@Model)
f0 <- func(sem_out)
if (is.na(f0) || is.nan(f0)) {
f0 <- Inf
}
f0
}
out
}
# @title Generate a Function to Fix a
# User-Defined Parameter
#
# @description Add a function to a
# `lavaan` object as a user-defined
# parameter and return a function which
# can fix this user-defined parameter
# to a value.
#
#' @details
#'
#' ## `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.
#'
#' @param userp A function that is
#' generated by [gen_userp()].
#'
#' @param userp_name The name of the
#' function `userp` to be used in the
#' `lavaan` model. It does not have to
#' be the name of the function in
#' `userp`. Should be changed only if it
#' conflicts with another object in the
#' parent environment, which should not
#' happen if the model is always fitted
#' in a clean R session.
#'
#' @param fix If `TRUE`, the default,
#' the function generated is used to fix
#' the value of `userp` to a target
#' value using an equality constraint.
#' If `FALSE`, then the function simply
#' fits the model to the data.
#'
#' @param control_args To be passed to
#' the argument of the same name in
#' [lavaan::lavaan()]. Default is
#' `list()`. Can be used to set the
#' default values of this argument
#' in the generated function.
#'
#' @param iter.max The maximum number of
#' iteration when the generated function
#' fit the model. Default is 10000.
#'
#' @param max_attempts If the initial
#' fit with the equality constraint
#' fails, how many more attempts
#' will be made by the generated
#' function. Default is 5.
#'
#' @return
#'
#' ## `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.
#'
#' @rdname gen_userp
#'
#' @export
gen_sem_out_userp <- function(userp,
sem_out,
userp_name = "semlbciuserp1234",
fix = TRUE,
control_args = list(),
iter.max = 10000,
max_attempts = 5) {
# Generate a function that:
# - refits the model with the user-parameter fixed to a target value.
# CHECKED: Identical to the experimental version
# For add_fun() using callr.
iter.max <- max(lavaan::lavInspect(sem_out, "optim")$iterations,
iter.max)
ptable <- lavaan::parameterTable(sem_out)
# Create a unique label: userp_label
labels <- unique(ptable$label)
labels <- labels[nchar(labels) > 0]
userp_label <- make.unique(c(labels, "user"))
userp_label <- userp_label[length(userp_label)]
# Generate the user parameter row using userp()
user_pt1 <- data.frame(lhs = userp_label[length(userp_label)],
op = ":=",
rhs = paste0(userp_name, "()"),
user = 1,
block = 0,
group = 0,
free = 0,
ustart = NA,
label = userp_label,
exo = 0)
user_pt1 <- lavaan::lav_partable_complete(partable = user_pt1,
start = TRUE)
# Add the user parameter to the model
ptable1 <- lavaan::lav_partable_add(partable = ptable,
add = user_pt1)
ptable1$se[which(ptable1$label == userp_label)] <- NA
# Extract slots to be used
slot_opt <- sem_out@Options
slot_dat <- sem_out@Data
if (fix) {
# Create the row of equality constraint
user_pt2 <- data.frame(lhs = userp_label,
op = "==",
rhs = NA,
user = 1,
block = 0,
group = 0,
free = 0,
ustart = NA,
exo = 0)
user_pt2 <- lavaan::lav_partable_complete(partable = user_pt2,
start = TRUE)
# Add the equality constraint to the parameter tahle
ptable2 <- lavaan::lav_partable_add(partable = ptable1,
add = user_pt2)
# Store the position of the constraint
i_eq <- which((ptable2$lhs == userp_label) &
(ptable2$op == "=="))
# Fix the ids of the rows
ptable2$id <- seq_along(ptable2$id)
# Set the options
slot_opt$do.fit <- TRUE
# The "mplus" version does not take into account the user-parameter
slot_opt$start <- "simple"
# Heywood cases are allowed
slot_opt$check.start <- FALSE
slot_opt$check.post <- FALSE
# Disable some options for efficiency
slot_opt$se <- "none"
# slot_opt$h1 <- FALSE
# slot_opt$baseline <- FALSE
out <- function(target,
verbose = FALSE,
control = list(),
seed = NULL) {
if (!is.null(seed)) set.seed(seed)
# Just in case ...
force(control_args)
force(slot_opt)
force(slot_dat)
force(ptable2)
force(userp_label)
force(i_eq)
ptable3 <- ptable2
ptable3$rhs[i_eq] <- target
control_args <- utils::modifyList(control,
list(iter.max = iter.max))
slot_opt$control <- control_args
# Fit once to update the model slot
attempted <- 0
fit_new_converged <- FALSE
if (verbose) {
cat("Target:", target, "\n")
}
if (verbose) {
cat("First attempt ...\n")
}
# Fail when
# eq.constraints.K is 0x0
# Succeed when
# fit_new@Model@eq.constraints, or
# length(fit_new@Model@ceq.nonlinear.idx) > 0
slot_opt_tmp <- slot_opt
slot_opt_tmp$do.fit <- FALSE
eq_not_ok <- TRUE
while (eq_not_ok) {
fit_new <- lavaan::lavaan(model = ptable3,
slotOptions = slot_opt_tmp,
slotData = slot_dat)
eq_not_ok <- !fit_new@Model@eq.constraints &&
!(length(fit_new@Model@ceq.nonlinear.idx) > 0)
}
# Extract the updated slots
slot_pat_new <- fit_new@ParTable
slot_model_new <- fit_new@Model
slot_smp_new <- fit_new@SampleStats
# Fit the model with the user parameter constrained
# to the target value
fit_new <- tryCatch(lavaan::lavaan(slotParTable = slot_pat_new,
slotModel = slot_model_new,
slotSampleStats = slot_smp_new,
slotOptions = slot_opt,
slotData = slot_dat),
error = function(e) e)
if (inherits(fit_new, "error")) {
stop("Something's with the function.",
"Error occurred when fitting the modified model.")
}
fit_new_converged <- lavaan::lavInspect(fit_new, "converged")
# Try max_attempts more times
if (!fit_new_converged) {
if (verbose) {
cat("More attempts ...\n")
}
while (!((attempted > max_attempts) ||
fit_new_converged)) {
if (verbose) {cat("Attempt", attempted, "\n")}
ptable_new <- fit_new@ParTable
i_free <- ptable_new$free > 0
ptable_new$est[i_free] <- jitter(ptable_new$est[i_free],
factor = 50)
slot_opt$start <- as.data.frame(ptable_new)
fit_old <- fit_new
fit_new <- tryCatch(lavaan::lavaan(slotParTable = ptable_new,
slotModel = slot_model_new,
slotSampleStats = slot_smp_new,
slotOptions = slot_opt,
slotData = slot_dat),
error = function(e) e)
if (inherits(fit_new, "error")) {
fit_new <- fit_old
fit_new_converged <- FALSE
} else {
fit_new_converged <- lavaan::lavInspect(fit_new, "converged")
}
attempted <- attempted + 1
}
}
if (!fit_new_converged) {
if (verbose) {
cat("The modified model failed to converge.")
}
}
if (verbose) {
cat("Done!\n")
}
fit_new
}
} else {
# The arguments are not used in this version
out <- function(target,
verbose = FALSE,
control = list(),
seed = NULL) {
# Generate a model with the user parameter added.
# For debugging
# Just in case ...
force(slot_opt)
force(slot_dat)
force(ptable1)
slot_opt$do.fit <- TRUE
fit_new <- lavaan::lavaan(model = ptable1,
slotOptions = slot_opt,
slotData = slot_dat)
fit_new
}
}
out
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.