R/ci.R
In quantPsych/RMediation: Mediation Analysis Confidence Intervals
Defines functions
ci
Documented in
ci
#' CI for a nonlinear function of coefficients estimates
#'
#' This function returns a (\eqn{1-\alpha})% confidence interval (CI) for a
#' well--defined nonlinear function of the coefficients in single--level and
#' multilevel structural equation models. The \code{ci} function uses the Monte
#' Carlo (\code{type="MC"}) and the asymptotic normal theory
#' (\code{type="asymp"}) with the multivariate delta standard error
#' (Asymptotic--Delta) method (Sobel, 1982) to compute a CI. In addition, for
#' each of the methods, when a user specifies \code{plot=TRUE} and
#' \code{plotCI=TRUE}, a plot of the sampling distribution of the quantity of
#' interest in the \code{quant} argument and an overlaid plot of the CI will be
#' produced. When \code{type="all"} and \code{plot=TRUE}, two overlaid plots of
#' the sampling distributions corresponding to each method will be produced;
#' when \code{plotCI=TRUE}, then the overlaid plots of the CIs for both methods
#' will be displayed as well.
#'
#' @param mu (1) a \link{vector} of means (e.g., coefficient estimates) for the
#' normal random variables. A user can assign a name to each mean value, e.g.,
#' \code{mu=c(b1=.1,b2=3)}; otherwise, the coefficient names are assigned
#' automatically as follows: \code{b1,b2,...}. Or, (2) a \link{lavaan} object.
#' @param Sigma either a covariance matrix or a \link{vector} that stacks all
#' the columns of the lower triangle variance--covariance matrix one
#' underneath the other.
#' @param quant quantity of interest, which is a nonlinear/linear function of
#' the model parameters. Argument \code{quant} is a \link{formula} that
#' \strong{must} start with the symbol "tilde" (\code{~}): e.g.,
#' \code{~b1*b2*b3*b4}. The names of coefficients must conform to the names
#' provided in the argument \code{mu} or to the default names, i.e.,
#' \code{b1,b2,...}.
#' @param alpha significance level for the CI. The default value is .05.
#' @param type method used to compute a CI. It takes on the values \code{"MC"}
#' (default) for Monte Carlo, \code{"asymp"} for Asymptotic--Delta, or
#' \code{"all"} that produces CIs using both methods.
#' @param plot when \code{TRUE}, plot the approximate sampling distribution of
#' the quantity of interest using the specified method(s) in the argument
#' \code{type}. The default value is \code{FALSE}. When \code{type="all"},
#' superimposed density plots generated by both methods are displayed.
#' @param plotCI when \code{TRUE}, overlays a CI plot with error bars on the
#' density plot of the sampling distribution of \code{quant}. When
#' \code{type="all"}, the superimposed CI plots generated by both methods are
#' added to the density plots. Note that to obtain a CI plot, one must also
#' specify \code{plot="TRUE"}. The default value is \code{FALSE}.
#' @param n.mc Monte Carlo sample size. The default sample size is 1e+6.
#' @param H0 False. If \code{TRUE}, it will estimate the sampling distribution
#' of \eqn{H_{0}:f(\bm b)=0}. See the arguments \code{mu0} and \code{Sigma0}.
#' @param mu0 a \link{vector} of means (e.g., coefficient estimates) for the
#' normal random variables that satisfy the null hypothesis \eqn{H_{0}:f(\bm
#' b)=0}. If it is not provided, smallest z value of \code{mu} is zet to zero.
#' @param Sigma0 either a covariance matrix or a \link{vector} that stacks all
#' the columns of the lower triangle variance--covariance matrix one
#' underneath the other. If it is not provided, then \code{Sigma} is used
#' instead.
#' @param ... additional arguments.
#' @return When \code{type} is \code{"MC"} or \code{"asymp"}, \code{ci} returns
#' a \link{list} that contains: \item{(\eqn{1-\alpha})% CI}{a vector of lower
#' and upper confidence limits,} \item{Estimate}{a point estimate of the
#' quantity of interest,} \item{SE}{standard error of the quantity of
#' interest,} \item{MC Error}{When \code{type="MC"}, error of the Monte Carlo
#' estimate.} When \code{type="all"}, \code{ci} returns a \link{list} of two
#' objects, each of which a \link{list} that contains the results produced by
#' each method as described above.
#' @keywords regression distribution
#'
#' @examples
#' ci(
#' mu = c(b1 = 1, b2 = .7, b3 = .6, b4 = .45),
#' Sigma = c(.05, 0, 0, 0, .05, 0, 0, .03, 0, .03),
#' quant = ~ b1 * b2 * b3 * b4, type = "all", plot = TRUE, plotCI = TRUE
#' )
#' # An Example of Conservative Null Sampling Distribution
#' ci(c(b1 = .3, b2 = .4, b3 = .3), c(.01, 0, 0, .01, 0, .02),
#' quant = ~ b1 * b2 * b3, type = "mc", plot = TRUE, plotCI = TRUE,
#' H0 = TRUE, mu0 = c(b1 = .3, b2 = .4, b3 = 0)
#' )
#' # An Example of Less Conservative Null Sampling Distribution
#' ci(c(b1 = .3, b2 = .4, b3 = .3), c(.01, 0, 0, .01, 0, .02),
#' quant = ~ b1 * b2 * b3, type = "mc", plot = TRUE, plotCI = TRUE,
#' H0 = TRUE, mu0 = c(b1 = 0, b2 = .4, b3 = 0.1)
#' )
#' @author Davood Tofighi \email{dtofighi@@gmail.com}
#' @references Tofighi, D. and MacKinnon, D. P. (2011). RMediation: An R
#' package for mediation analysis confidence intervals. \emph{Behavior
#' Research Methods}, \bold{43}, 692--700. \doi{doi:10.3758/s13428-011-0076-x}
#' @seealso \code{\link{medci}} \code{\link{RMediation-package}}
#' @importFrom lavaan lav_matrix_vech_reverse
#' @importFrom MASS mvrnorm
#' @importFrom e1071 skewness kurtosis
#' @note A shiny web application for Monte Carlo method of this function is available at \url{https://amplab.shinyapps.io/MEDMC/}
#' @export
ci <- function(mu, Sigma, quant, alpha = 0.05, type = "MC", plot = FALSE, plotCI = FALSE, n.mc = 1e+06, H0 = FALSE, mu0 = NULL, Sigma0 = NULL, ...) {
if (missing(mu) | is.null(mu)) stop(paste("argument", sQuote("mu"), "must be specified"))
if (missing(quant) | is.null(quant)) stop(paste("argument", sQuote("quant"), "must be specified"))
if (!type %in% c("MC", "mc", "asymp", "Asymp", "all")) stop(paste("Please enter a ccorrect", sQuote("type"), "argument"))
if (!isa(mu, "lavaan")) {
if (is.null(Sigma) | missing(Sigma)) stop(paste("argument", sQuote("Sigma"), "cannot be a NULL value"))
if (!is.matrix(Sigma)) {
if (length(mu) != (sqrt(1 + 8 * length(Sigma)) - 1) / 2) stop(paste("Please check the length of", sQuote("Sigma"), "and", sQuote("mu"), ". If the length(dimension) of the", sQuote("mu"), "vector (", length(mu), ") is correct, the stacked lower triangle matrix", sQuote("Sigma"), "must have ", ((2 * length(mu) + 1)^2 - 1) / 8, "elements, instead of", length(Sigma)))
Sigma <- lavaan::lav_matrix_vech_reverse(Sigma) # converts to a symmetric matrix
}
if (is.null(names(mu))) names(mu) <- paste("b", 1:length(mu), sep = "") # if mu names is NULL
if (!all(all.vars(quant) %in% names(mu))) stop(paste("The parameters names in formula", sQuote("quant"), "must match the parameters names provided in", sQuote("mu"), "."))
} else {
fm1 <- mu
# if(!fm1@Options[['do.fit']]) fm1 <- update(mu, do.fit=TRUE) # If it's not fitted, refit.
pEstM1 <- coef(fm1) # parameter estimate
name1 <- all.vars(quant)
if (!all(name1 %in% names(pEstM1))) stop(paste("The parameters names in formula", sQuote("quant"), "must match the parameters names provided in lavaan object."))
mu <- pEstM1[name1]
## Cov
covM1 <- (vcov(fm1)) # covariance of the coef estimates
Sigma <- covM1[name1, name1]
}
if (length(mu) * n.mc > .Machine$integer.max) {
n.mc <- 1e6
warning(paste("n.mc is too large. It is reset to", n.mc))
}
# This line sets the plot margins and other plot parameters
if (plot) {
op <- par(mar = c(7, 4, 7, 4) + 0.1, xpd = TRUE, ask = FALSE)
}
if (type %in% c("mc", "MC")) res <- confintMC(mu = mu, Sigma = Sigma, quant = quant, alpha = alpha, plot = plot, plotCI = plotCI, n.mc = n.mc, H0 = H0, mu0 = mu0, Sigma0 = Sigma0, ...)
if (type %in% c("asymp", "Asymp")) res <- confintAsymp(mu = mu, Sigma = Sigma, quant = quant, alpha = alpha, plot = plot, plotCI = plotCI)
if (type %in% "all") {
res1 <- confintMC(mu = mu, Sigma = Sigma, quant = quant, type = type, alpha = alpha, plot = plot, plotCI = plotCI, n.mc = n.mc)
res2 <- confintAsymp(mu = mu, Sigma = Sigma, quant = quant, type = type, alpha = alpha, plot = plot, plotCI = plotCI)
res <- list(MC = res1, Asymptotic = res2)
}
if (plot) {
# par(mar = c(5, 4, 4, 2) + 0.1, xpd=FALSE)
par(op)
}
return(res)
}
quantPsych/RMediation documentation built on March 4, 2024, 6 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions ci
Documented in ci
#' CI for a nonlinear function of coefficients estimates
#'
#' This function returns a (\eqn{1-\alpha})% confidence interval (CI) for a
#' well--defined nonlinear function of the coefficients in single--level and
#' multilevel structural equation models. The \code{ci} function uses the Monte
#' Carlo (\code{type="MC"}) and the asymptotic normal theory
#' (\code{type="asymp"}) with the multivariate delta standard error
#' (Asymptotic--Delta) method (Sobel, 1982) to compute a CI. In addition, for
#' each of the methods, when a user specifies \code{plot=TRUE} and
#' \code{plotCI=TRUE}, a plot of the sampling distribution of the quantity of
#' interest in the \code{quant} argument and an overlaid plot of the CI will be
#' produced. When \code{type="all"} and \code{plot=TRUE}, two overlaid plots of
#' the sampling distributions corresponding to each method will be produced;
#' when \code{plotCI=TRUE}, then the overlaid plots of the CIs for both methods
#' will be displayed as well.
#'
#' @param mu (1) a \link{vector} of means (e.g., coefficient estimates) for the
#' normal random variables. A user can assign a name to each mean value, e.g.,
#' \code{mu=c(b1=.1,b2=3)}; otherwise, the coefficient names are assigned
#' automatically as follows: \code{b1,b2,...}. Or, (2) a \link{lavaan} object.
#' @param Sigma either a covariance matrix or a \link{vector} that stacks all
#' the columns of the lower triangle variance--covariance matrix one
#' underneath the other.
#' @param quant quantity of interest, which is a nonlinear/linear function of
#' the model parameters. Argument \code{quant} is a \link{formula} that
#' \strong{must} start with the symbol "tilde" (\code{~}): e.g.,
#' \code{~b1*b2*b3*b4}. The names of coefficients must conform to the names
#' provided in the argument \code{mu} or to the default names, i.e.,
#' \code{b1,b2,...}.
#' @param alpha significance level for the CI. The default value is .05.
#' @param type method used to compute a CI. It takes on the values \code{"MC"}
#' (default) for Monte Carlo, \code{"asymp"} for Asymptotic--Delta, or
#' \code{"all"} that produces CIs using both methods.
#' @param plot when \code{TRUE}, plot the approximate sampling distribution of
#' the quantity of interest using the specified method(s) in the argument
#' \code{type}. The default value is \code{FALSE}. When \code{type="all"},
#' superimposed density plots generated by both methods are displayed.
#' @param plotCI when \code{TRUE}, overlays a CI plot with error bars on the
#' density plot of the sampling distribution of \code{quant}. When
#' \code{type="all"}, the superimposed CI plots generated by both methods are
#' added to the density plots. Note that to obtain a CI plot, one must also
#' specify \code{plot="TRUE"}. The default value is \code{FALSE}.
#' @param n.mc Monte Carlo sample size. The default sample size is 1e+6.
#' @param H0 False. If \code{TRUE}, it will estimate the sampling distribution
#' of \eqn{H_{0}:f(\bm b)=0}. See the arguments \code{mu0} and \code{Sigma0}.
#' @param mu0 a \link{vector} of means (e.g., coefficient estimates) for the
#' normal random variables that satisfy the null hypothesis \eqn{H_{0}:f(\bm
#' b)=0}. If it is not provided, smallest z value of \code{mu} is zet to zero.
#' @param Sigma0 either a covariance matrix or a \link{vector} that stacks all
#' the columns of the lower triangle variance--covariance matrix one
#' underneath the other. If it is not provided, then \code{Sigma} is used
#' instead.
#' @param ... additional arguments.
#' @return When \code{type} is \code{"MC"} or \code{"asymp"}, \code{ci} returns
#' a \link{list} that contains: \item{(\eqn{1-\alpha})% CI}{a vector of lower
#' and upper confidence limits,} \item{Estimate}{a point estimate of the
#' quantity of interest,} \item{SE}{standard error of the quantity of
#' interest,} \item{MC Error}{When \code{type="MC"}, error of the Monte Carlo
#' estimate.} When \code{type="all"}, \code{ci} returns a \link{list} of two
#' objects, each of which a \link{list} that contains the results produced by
#' each method as described above.
#' @keywords regression distribution
#'
#' @examples
#' ci(
#' mu = c(b1 = 1, b2 = .7, b3 = .6, b4 = .45),
#' Sigma = c(.05, 0, 0, 0, .05, 0, 0, .03, 0, .03),
#' quant = ~ b1 * b2 * b3 * b4, type = "all", plot = TRUE, plotCI = TRUE
#' )
#' # An Example of Conservative Null Sampling Distribution
#' ci(c(b1 = .3, b2 = .4, b3 = .3), c(.01, 0, 0, .01, 0, .02),
#' quant = ~ b1 * b2 * b3, type = "mc", plot = TRUE, plotCI = TRUE,
#' H0 = TRUE, mu0 = c(b1 = .3, b2 = .4, b3 = 0)
#' )
#' # An Example of Less Conservative Null Sampling Distribution
#' ci(c(b1 = .3, b2 = .4, b3 = .3), c(.01, 0, 0, .01, 0, .02),
#' quant = ~ b1 * b2 * b3, type = "mc", plot = TRUE, plotCI = TRUE,
#' H0 = TRUE, mu0 = c(b1 = 0, b2 = .4, b3 = 0.1)
#' )
#' @author Davood Tofighi \email{dtofighi@@gmail.com}
#' @references Tofighi, D. and MacKinnon, D. P. (2011). RMediation: An R
#' package for mediation analysis confidence intervals. \emph{Behavior
#' Research Methods}, \bold{43}, 692--700. \doi{doi:10.3758/s13428-011-0076-x}
#' @seealso \code{\link{medci}} \code{\link{RMediation-package}}
#' @importFrom lavaan lav_matrix_vech_reverse
#' @importFrom MASS mvrnorm
#' @importFrom e1071 skewness kurtosis
#' @note A shiny web application for Monte Carlo method of this function is available at \url{https://amplab.shinyapps.io/MEDMC/}
#' @export
ci <- function(mu, Sigma, quant, alpha = 0.05, type = "MC", plot = FALSE, plotCI = FALSE, n.mc = 1e+06, H0 = FALSE, mu0 = NULL, Sigma0 = NULL, ...) {
if (missing(mu) | is.null(mu)) stop(paste("argument", sQuote("mu"), "must be specified"))
if (missing(quant) | is.null(quant)) stop(paste("argument", sQuote("quant"), "must be specified"))
if (!type %in% c("MC", "mc", "asymp", "Asymp", "all")) stop(paste("Please enter a ccorrect", sQuote("type"), "argument"))
if (!isa(mu, "lavaan")) {
if (is.null(Sigma) | missing(Sigma)) stop(paste("argument", sQuote("Sigma"), "cannot be a NULL value"))
if (!is.matrix(Sigma)) {
if (length(mu) != (sqrt(1 + 8 * length(Sigma)) - 1) / 2) stop(paste("Please check the length of", sQuote("Sigma"), "and", sQuote("mu"), ". If the length(dimension) of the", sQuote("mu"), "vector (", length(mu), ") is correct, the stacked lower triangle matrix", sQuote("Sigma"), "must have ", ((2 * length(mu) + 1)^2 - 1) / 8, "elements, instead of", length(Sigma)))
Sigma <- lavaan::lav_matrix_vech_reverse(Sigma) # converts to a symmetric matrix
}
if (is.null(names(mu))) names(mu) <- paste("b", 1:length(mu), sep = "") # if mu names is NULL
if (!all(all.vars(quant) %in% names(mu))) stop(paste("The parameters names in formula", sQuote("quant"), "must match the parameters names provided in", sQuote("mu"), "."))
} else {
fm1 <- mu
# if(!fm1@Options[['do.fit']]) fm1 <- update(mu, do.fit=TRUE) # If it's not fitted, refit.
pEstM1 <- coef(fm1) # parameter estimate
name1 <- all.vars(quant)
if (!all(name1 %in% names(pEstM1))) stop(paste("The parameters names in formula", sQuote("quant"), "must match the parameters names provided in lavaan object."))
mu <- pEstM1[name1]
## Cov
covM1 <- (vcov(fm1)) # covariance of the coef estimates
Sigma <- covM1[name1, name1]
}
if (length(mu) * n.mc > .Machine$integer.max) {
n.mc <- 1e6
warning(paste("n.mc is too large. It is reset to", n.mc))
}
# This line sets the plot margins and other plot parameters
if (plot) {
op <- par(mar = c(7, 4, 7, 4) + 0.1, xpd = TRUE, ask = FALSE)
}
if (type %in% c("mc", "MC")) res <- confintMC(mu = mu, Sigma = Sigma, quant = quant, alpha = alpha, plot = plot, plotCI = plotCI, n.mc = n.mc, H0 = H0, mu0 = mu0, Sigma0 = Sigma0, ...)
if (type %in% c("asymp", "Asymp")) res <- confintAsymp(mu = mu, Sigma = Sigma, quant = quant, alpha = alpha, plot = plot, plotCI = plotCI)
if (type %in% "all") {
res1 <- confintMC(mu = mu, Sigma = Sigma, quant = quant, type = type, alpha = alpha, plot = plot, plotCI = plotCI, n.mc = n.mc)
res2 <- confintAsymp(mu = mu, Sigma = Sigma, quant = quant, type = type, alpha = alpha, plot = plot, plotCI = plotCI)
res <- list(MC = res1, Asymptotic = res2)
}
if (plot) {
# par(mar = c(5, 4, 4, 2) + 0.1, xpd=FALSE)
par(op)
}
return(res)
}
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.