The bodies of the following functions contain the R implementation of the formulas in V2015. Although many models are covered, there are only four patterns. See the package top page for further references. The code seen here is compacted and lacks the comments. For more redable code, please refer to the Github repo and search for the function names without the preceding regmedint:::. For the type-set mathematical expressions and LaTeX source, please see the supplement.

## mreg linear yreg linear (V2015 p466 Proposition 2.3)

These functions are only used in the setting where both the mediator model and the outcome model are linear regression.

### Point estimates

regmedint:::calc_myreg_mreg_linear_yreg_linear_est
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6)
##     fun_est <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         cde <- (theta1 + theta3 * m_cde + theta5_c) * (a1 - a0)
##         pnde <- (theta1 + theta3 * (beta0 + beta1 * a0 + beta2_c +
##             beta3_c * a0) + theta5_c) * (a1 - a0)
##         tnie <- (theta2 + theta3 * a1 + theta6_c) * (beta1 +
##             beta3_c) * (a1 - a0)
##         tnde <- (theta1 + theta3 * (beta0 + beta1 * a1 + beta2_c +
##             beta3_c * a1) + theta5_c) * (a1 - a0)
##         pnie <- (theta2 + theta3 * a0 + theta6_c) * (beta1 +
##             beta3_c) * (a1 - a0)
##         te <- pnde + tnie
##         pm <- tnie/te
##         c(cde = unname(cde), pnde = unname(pnde), tnie = unname(tnie),
##             tnde = unname(tnde), pnie = unname(pnie), te = unname(te),
##             pm = unname(pm))
##     }
##     return(fun_est)
## }
## <bytecode: 0x7fb889ae7f90>
## <environment: namespace:regmedint>

### Standard error estimates

regmedint:::calc_myreg_mreg_linear_yreg_linear_se
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6, Sigma_beta, Sigma_theta)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6)
##     validate_myreg_vcovs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6,
##         Sigma_beta = Sigma_beta, Sigma_theta = Sigma_theta)
##     Sigma <- Matrix::bdiag(Sigma_beta, Sigma_theta)
##     fun_se <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             pd_cde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_cde_theta5 <- c_cond
##         }
##         Gamma_cde <- matrix(c(0, 0, rep(0, length(beta2)), rep(0,
##             length(beta3)), 0, 1, 0, m_cde, rep(0, length(theta4)),
##             pd_cde_theta5, rep(0, length(theta6))))
##         if (is.null(beta3)) {
##             pd_pnde_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_pnde_beta3 <- theta3 * c_cond
##         }
##         if (is.null(theta5)) {
##             pd_pnde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_pnde_theta5 <- c_cond
##         }
##         Gamma_pnde <- matrix(c(theta3, theta3 * a0, theta3 *
##             c_cond, pd_pnde_beta3, 0, 1, 0, beta0 + beta1 * a0 +
##             beta2_c + beta3_c * a0, rep(0, length(theta4)), pd_pnde_theta5,
##             rep(0, length(theta6))))
##         if (is.null(beta3)) {
##             pd_tnie_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_tnie_beta3 <- c_cond * (theta2 + theta3 * a1 +
##                 theta6_c)
##         }
##         if (is.null(theta6)) {
##             pd_tnie_theta6 <- rep(0, length(theta6))
##         }
##         else {
##             pd_tnie_theta6 <- c_cond * (beta1 + beta3_c)
##         }
##         Gamma_tnie <- matrix(c(0, theta2 + theta3 * a1 + theta6_c,
##             rep(0, length(beta2)), pd_tnie_beta3, 0, 0, beta1 +
##                 beta3_c, a1 * (beta1 + beta3_c), rep(0, length(theta4)),
##             rep(0, length(theta5)), pd_tnie_theta6))
##         if (is.null(beta3)) {
##             pd_tnde_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_tnde_beta3 <- theta3 * a1 * c_cond
##         }
##         if (is.null(theta5)) {
##             pd_tnde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_tnde_theta5 <- c_cond
##         }
##         Gamma_tnde <- matrix(c(theta3, theta3 * a1, theta3 *
##             c_cond, pd_tnde_beta3, 0, 1, 0, beta0 + beta1 * a1 +
##             beta2_c + beta3_c * a1, rep(0, length(theta4)), pd_tnde_theta5,
##             rep(0, length(theta6))))
##         if (is.null(beta3)) {
##             pd_pnie_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_pnie_beta3 <- c_cond * (theta2 + theta3 * a0 +
##                 theta6_c)
##         }
##         if (is.null(theta6)) {
##             pd_pnie_theta6 <- rep(0, length(theta6))
##         }
##         else {
##             pd_pnie_theta6 <- c_cond * (beta1 + beta3_c)
##         }
##         Gamma_pnie <- matrix(c(0, theta2 + theta3 * a0 + theta6_c,
##             rep(0, length(beta2)), pd_pnie_beta3, 0, 0, beta1 +
##                 beta3_c, a0 * (beta1 + beta3_c), rep(0, length(theta4)),
##             rep(0, length(theta5)), pd_pnie_theta6))
##         Gamma_te <- Gamma_pnde + Gamma_tnie
##         pnde <- (theta1 + theta3 * (beta0 + beta1 * a0 + beta2_c +
##             beta3_c * a0) + theta5_c) * (a1 - a0)
##         tnie <- (theta2 + theta3 * a1 + theta6_c) * (beta1 +
##             beta3_c) * (a1 - a0)
##         d_pm <- grad_prop_med_yreg_linear(pnde = unname(pnde),
##             tnie = unname(tnie))
##         Gamma_pm <- (d_pm[["pnde"]] * Gamma_pnde) + (d_pm[["tnie"]] *
##             Gamma_tnie)
##         a1_sub_a0 <- abs(a1 - a0)
##         se_cde <- sqrt(as.numeric(t(Gamma_cde) %*% Sigma %*%
##             Gamma_cde)) * a1_sub_a0
##         se_pnde <- sqrt(as.numeric(t(Gamma_pnde) %*% Sigma %*%
##             Gamma_pnde)) * a1_sub_a0
##         se_tnie <- sqrt(as.numeric(t(Gamma_tnie) %*% Sigma %*%
##             Gamma_tnie)) * a1_sub_a0
##         se_tnde <- sqrt(as.numeric(t(Gamma_tnde) %*% Sigma %*%
##             Gamma_tnde)) * a1_sub_a0
##         se_pnie <- sqrt(as.numeric(t(Gamma_pnie) %*% Sigma %*%
##             Gamma_pnie)) * a1_sub_a0
##         se_te <- sqrt(as.numeric(t(Gamma_te) %*% Sigma %*% Gamma_te)) *
##             a1_sub_a0
##         se_pm <- sqrt(as.numeric(t(Gamma_pm) %*% Sigma %*% Gamma_pm)) *
##             a1_sub_a0
##         c(se_cde = unname(se_cde), se_pnde = unname(se_pnde),
##             se_tnie = unname(se_tnie), se_tnde = unname(se_tnde),
##             se_pnie = unname(se_pnie), se_te = unname(se_te),
##             se_pm = unname(se_pm))
##     }
##     return(fun_se)
## }
## <bytecode: 0x7fb886c1e218>
## <environment: namespace:regmedint>

## mreg linear yreg non-linear (V2015 p468 Proposition 2.4)

These functions are used in all cases where the mediator model is linear regression and the outcome model is any one of the non-linear models.

### Point estimates

regmedint:::calc_myreg_mreg_linear_yreg_logistic_est
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6, sigma_sq)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6,
##         sigma_sq = sigma_sq)
##     fun_est <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         cde <- (theta1 + theta3 * m_cde + theta5_c) * (a1 - a0)
##         pnde <- ((theta1 + theta5_c + theta3 * (beta0 + beta1 *
##             a0 + beta2_c + beta3_c * a0) + 1/2 * sigma_sq * theta3 *
##             (2 * theta2 + theta3 * a1 + theta3 * a0 + 2 * theta6_c))) *
##             (a1 - a0)
##         tnie <- (theta2 + theta3 * a1 + theta6_c) * (beta1 +
##             beta3_c) * (a1 - a0)
##         tnde <- ((theta1 + theta5_c + theta3 * (beta0 + beta1 *
##             a1 + beta2_c + beta3_c * a1) + 1/2 * sigma_sq * theta3 *
##             (2 * theta2 + theta3 * a1 + theta3 * a0 + 2 * theta6_c))) *
##             (a1 - a0)
##         pnie <- (theta2 + theta3 * a0 + theta6_c) * (beta1 +
##             beta3_c) * (a1 - a0)
##         te <- pnde + tnie
##         pm <- exp(pnde) * (exp(tnie) - 1)/(exp(te) - 1)
##         c(cde = unname(cde), pnde = unname(pnde), tnie = unname(tnie),
##             tnde = unname(tnde), pnie = unname(pnie), te = unname(te),
##             pm = unname(pm))
##     }
##     return(fun_est)
## }
## <bytecode: 0x7fb886d96f38>
## <environment: namespace:regmedint>

### Standard error estimates

regmedint:::calc_myreg_mreg_linear_yreg_logistic_se
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6, sigma_sq, Sigma_beta, Sigma_theta,
##     Sigma_sigma_sq)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6,
##         sigma_sq = sigma_sq)
##     validate_myreg_vcovs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6,
##         sigma_sq = sigma_sq, Sigma_beta = Sigma_beta, Sigma_theta = Sigma_theta,
##         Sigma_sigma_sq = Sigma_sigma_sq)
##     Sigma <- Matrix::bdiag(Sigma_beta, Sigma_theta, Sigma_sigma_sq)
##     fun_se <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             pd_cde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_cde_theta5 <- c_cond
##         }
##         Gamma_cde <- matrix(c(0, 0, rep(0, length(beta2)), rep(0,
##             length(beta3)), 0, 1, 0, m_cde, rep(0, length(theta4)),
##             pd_cde_theta5, rep(0, length(theta6)), 0))
##         if (is.null(beta3)) {
##             pd_pnde_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_pnde_beta3 <- theta3 * a1 * c_cond
##         }
##         if (is.null(theta5)) {
##             pd_pnde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_pnde_theta5 <- c_cond
##         }
##         if (is.null(theta6)) {
##             pd_pnde_theta6 <- rep(0, length(theta6))
##         }
##         else {
##             pd_pnde_theta6 <- theta3 * sigma_sq * c_cond
##         }
##         Gamma_pnde <- matrix(c(theta3, theta3 * a0, theta3 *
##             c_cond, pd_pnde_beta3, 0, 1, theta3 * sigma_sq, a0 *
##             (beta1 + beta3_c) + beta0 + beta2_c + sigma_sq *
##             ((theta2 + theta6_c) + theta3 * (a0 + a1)), rep(0,
##             length(theta4)), pd_pnde_theta5, pd_pnde_theta6,
##             1/2 * theta3 * (2 * (theta2 + theta6_c) + theta3 *
##                 (a0 + a1))))
##         if (is.null(beta3)) {
##             pd_tnie_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_tnie_beta3 <- c_cond * (theta2 + theta3 * a1 +
##                 theta6_c)
##         }
##         if (is.null(theta6)) {
##             pd_tnie_theta6 <- rep(0, length(theta6))
##         }
##         else {
##             pd_tnie_theta6 <- c_cond * (beta1 + beta3_c)
##         }
##         Gamma_tnie <- matrix(c(0, theta2 + theta3 * a1 + theta6_c,
##             rep(0, length(beta2)), pd_tnie_beta3, 0, 0, beta1 +
##                 beta3_c, a1 * (beta1 + beta3_c), rep(0, length(theta4)),
##             rep(0, length(theta5)), pd_tnie_theta6, 0))
##         if (is.null(beta3)) {
##             pd_tnde_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_tnde_beta3 <- theta3 * a1 * c_cond
##         }
##         if (is.null(theta5)) {
##             pd_tnde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_tnde_theta5 <- c_cond
##         }
##         if (is.null(theta6)) {
##             pd_tnde_theta6 <- rep(0, length(theta6))
##         }
##         else {
##             pd_tnde_theta6 <- theta3 * sigma_sq * c_cond
##         }
##         Gamma_tnde <- matrix(c(theta3, theta3 * a1, theta3 *
##             c_cond, pd_tnde_beta3, 0, 1, theta3 * sigma_sq, a1 *
##             (beta1 + beta3_c) + beta0 + beta2_c + sigma_sq *
##             ((theta2 + theta6_c) + theta3 * (a0 + a1)), rep(0,
##             length(theta4)), pd_tnde_theta5, pd_tnde_theta6,
##             1/2 * theta3 * (2 * (theta2 + theta6_c) + theta3 *
##                 (a0 + a1))))
##         if (is.null(beta3)) {
##             pd_pnie_beta3 <- rep(0, length(beta3))
##         }
##         else {
##             pd_pnie_beta3 <- c_cond * (theta2 + theta3 * a0 +
##                 theta6_c)
##         }
##         if (is.null(theta6)) {
##             pd_pnie_theta6 <- rep(0, length(theta6))
##         }
##         else {
##             pd_pnie_theta6 <- c_cond * (beta1 + beta3_c)
##         }
##         Gamma_pnie <- matrix(c(0, theta2 + theta3 * a0 + theta6_c,
##             rep(0, length(beta2)), pd_pnie_beta3, 0, 0, beta1 +
##                 beta3_c, a0 * (beta1 + beta3_c), rep(0, length(theta4)),
##             rep(0, length(theta5)), pd_pnie_theta6, 0))
##         Gamma_te <- Gamma_pnde + Gamma_tnie
##         pnde <- ((theta1 + theta5_c + theta3 * (beta0 + beta1 *
##             a0 + beta2_c + beta3_c * a0) + 1/2 * sigma_sq * theta3 *
##             (2 * theta2 + theta3 * a1 + theta3 * a0 + 2 * theta6_c))) *
##             (a1 - a0)
##         tnie <- (theta2 + theta3 * a1 + theta6_c) * (beta1 +
##             beta3_c) * (a1 - a0)
##         d_pm <- grad_prop_med_yreg_logistic(pnde = unname(pnde),
##             tnie = unname(tnie))
##         Gamma_pm <- d_pm[["pnde"]] * Gamma_pnde + d_pm[["tnie"]] *
##             Gamma_tnie
##         a1_sub_a0 <- abs(a1 - a0)
##         se_cde <- sqrt(as.numeric(t(Gamma_cde) %*% Sigma %*%
##             Gamma_cde)) * a1_sub_a0
##         se_pnde <- sqrt(as.numeric(t(Gamma_pnde) %*% Sigma %*%
##             Gamma_pnde)) * a1_sub_a0
##         se_tnie <- sqrt(as.numeric(t(Gamma_tnie) %*% Sigma %*%
##             Gamma_tnie)) * a1_sub_a0
##         se_tnde <- sqrt(as.numeric(t(Gamma_tnde) %*% Sigma %*%
##             Gamma_tnde)) * a1_sub_a0
##         se_pnie <- sqrt(as.numeric(t(Gamma_pnie) %*% Sigma %*%
##             Gamma_pnie)) * a1_sub_a0
##         se_te <- sqrt(as.numeric(t(Gamma_te) %*% Sigma %*% Gamma_te)) *
##             a1_sub_a0
##         se_pm <- sqrt(as.numeric(t(Gamma_pm) %*% Sigma %*% Gamma_pm)) *
##             a1_sub_a0
##         c(se_cde = unname(se_cde), se_pnde = unname(se_pnde),
##             se_tnie = unname(se_tnie), se_tnde = unname(se_tnde),
##             se_pnie = unname(se_pnie), se_te = unname(se_te),
##             se_pm = unname(se_pm))
##     }
##     return(fun_se)
## }
## <bytecode: 0x7fb886ef6e90>
## <environment: namespace:regmedint>

## mreg logistic yreg linear (V2015 p471 Proposition 2.5)

These functions are only used in the setting where the mediator model is logistic regression and the outcome model is non-linear regression.

### Point estimates

regmedint:::calc_myreg_mreg_logistic_yreg_linear_est
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6)
##     fun_est <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         expit <- function(x) {
##             exp(x)/(1 + exp(x))
##         }
##         cde <- (theta1 + theta3 * m_cde + theta5_c) * (a1 - a0)
##         pnde <- (theta1 + theta3 * expit(beta0 + beta1 * a0 +
##             beta2_c + beta3_c * a0) + theta5_c) * (a1 - a0)
##         tnie <- (theta2 + theta3 * a1 + theta6_c) * (expit(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1) - expit(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0))
##         tnde <- (theta1 + theta3 * expit(beta0 + beta1 * a1 +
##             beta2_c + beta3_c * a1) + theta5_c) * (a1 - a0)
##         pnie <- (theta2 + theta3 * a0 + theta6_c) * (expit(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1) - expit(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0))
##         te <- pnde + tnie
##         pm <- tnie/te
##         c(cde = unname(cde), pnde = unname(pnde), tnie = unname(tnie),
##             tnde = unname(tnde), pnie = unname(pnie), te = unname(te),
##             pm = unname(pm))
##     }
##     return(fun_est)
## }
## <bytecode: 0x7fb88964bbf8>
## <environment: namespace:regmedint>

### Standard error estimates

regmedint:::calc_myreg_mreg_logistic_yreg_linear_se
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6, Sigma_beta, Sigma_theta)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6)
##     validate_myreg_vcovs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6,
##         Sigma_beta = Sigma_beta, Sigma_theta = Sigma_theta)
##     Sigma <- Matrix::bdiag(Sigma_beta, Sigma_theta)
##     fun_se <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         expit <- function(x) {
##             exp(x)/(1 + exp(x))
##         }
##         if (is.null(theta5)) {
##             pd_cde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_cde_theta5 <- c_cond
##         }
##         Gamma_cde <- matrix(c(0, 0, rep(0, length(beta2)), rep(0,
##             length(beta3)), 0, 1, 0, m_cde, rep(0, length(theta4)),
##             pd_cde_theta5, rep(0, length(theta6))))
##         .e3 <- a0 * (beta1 + beta3_c) + beta0 + beta2_c
##         .e4 <- exp(.e3)
##         .e5 <- 1 + .e4
##         .e6 <- 1 - .e4/.e5
##         pnde_d1 <- theta3 * .e6 * .e4/.e5
##         pnde_d2 <- a0 * theta3 * .e6 * .e4/.e5
##         pnde_d3 <- c_cond * theta3 * .e6 * .e4/.e5
##         if (is.null(beta3)) {
##             pnde_d4 <- rep(0, length(beta3))
##         }
##         else {
##             pnde_d4 <- a0 * c_cond * theta3 * .e6 * .e4/.e5
##         }
##         pnde_d5 <- 0
##         pnde_d6 <- 1
##         pnde_d7 <- 0
##         pnde_d8 <- expit(.e3)
##         pnde_d9 <- rep(0, length(theta4))
##         if (is.null(theta5)) {
##             pnde_d10 <- rep(0, length(theta5))
##         }
##         else {
##             pnde_d10 <- c_cond
##         }
##         pnde_d11 <- rep(0, length(theta6))
##         Gamma_pnde <- matrix(c(pnde_d1, pnde_d2, pnde_d3, pnde_d4,
##             pnde_d5, pnde_d6, pnde_d7, pnde_d8, pnde_d9, pnde_d10,
##             pnde_d11))
##         .e1 <- beta1 + beta3_c
##         .e2 <- beta2_c
##         .e5 <- a0 * .e1 + beta0 + .e2
##         .e8 <- a1 * .e1 + beta0 + .e2
##         .e9 <- exp(.e5)
##         .e10 <- exp(.e8)
##         .e11 <- 1 + .e9
##         .e12 <- 1 + .e10
##         .e13 <- 1 - .e9/.e11
##         .e14 <- 1 - .e10/.e12
##         .e17 <- a1 * theta3 + theta6_c + theta2
##         .e20 <- expit(.e8) - expit(.e5)
##         .e25 <- .e14 * .e10/.e12 - .e13 * .e9/.e11
##         .e32 <- a1 * .e14 * .e10/.e12 - a0 * .e13 * .e9/.e11
##         tnie_d1 <- .e25 * .e17
##         tnie_d2 <- .e32 * .e17
##         tnie_d3 <- c_cond * .e25 * .e17
##         if (is.null(beta3)) {
##             tnie_d4 <- rep(0, length(beta3))
##         }
##         else {
##             tnie_d4 <- c_cond * .e32 * .e17
##         }
##         tnie_d5 <- 0
##         tnie_d6 <- 0
##         tnie_d7 <- .e20
##         tnie_d8 <- a1 * .e20
##         tnie_d9 <- rep(0, length(theta4))
##         tnie_d10 <- rep(0, length(theta5))
##         if (is.null(theta6)) {
##             tnie_d11 <- rep(0, length(theta6))
##         }
##         else {
##             tnie_d11 <- c_cond * .e20
##         }
##         Gamma_tnie <- matrix(c(tnie_d1, tnie_d2, tnie_d3, tnie_d4,
##             tnie_d5, tnie_d6, tnie_d7, tnie_d8, tnie_d9, tnie_d10,
##             tnie_d11))
##         .e3 <- a1 * (beta1 + beta3_c) + beta0 + beta2_c
##         .e4 <- exp(.e3)
##         .e5 <- 1 + .e4
##         .e6 <- 1 - .e4/.e5
##         tnde_d1 <- theta3 * .e6 * .e4/.e5
##         tnde_d2 <- a1 * theta3 * .e6 * .e4/.e5
##         tnde_d3 <- c_cond * theta3 * .e6 * .e4/.e5
##         if (is.null(beta3)) {
##             tnde_d4 <- rep(0, length(beta3))
##         }
##         else {
##             tnde_d4 <- a1 * c_cond * theta3 * .e6 * .e4/.e5
##         }
##         tnde_d5 <- 0
##         tnde_d6 <- 1
##         tnde_d7 <- 0
##         tnde_d8 <- expit(.e3)
##         tnde_d9 <- rep(0, length(theta4))
##         if (is.null(theta5)) {
##             tnde_d10 <- rep(0, length(theta5))
##         }
##         else {
##             tnde_d10 <- c_cond
##         }
##         tnde_d11 <- rep(0, length(theta6))
##         Gamma_tnde <- matrix(c(tnde_d1, tnde_d2, tnde_d3, tnde_d4,
##             tnde_d5, tnde_d6, tnde_d7, tnde_d8, tnde_d9, tnde_d10,
##             tnde_d11))
##         .e1 <- beta1 + beta3_c
##         .e2 <- beta2_c
##         .e5 <- a0 * .e1 + beta0 + .e2
##         .e8 <- a1 * .e1 + beta0 + .e2
##         .e9 <- exp(.e5)
##         .e10 <- exp(.e8)
##         .e11 <- 1 + .e9
##         .e12 <- 1 + .e10
##         .e13 <- 1 - .e9/.e11
##         .e14 <- 1 - .e10/.e12
##         .e17 <- a0 * theta3 + theta6_c + theta2
##         .e20 <- expit(.e8) - expit(.e5)
##         .e25 <- .e14 * .e10/.e12 - .e13 * .e9/.e11
##         .e32 <- a1 * .e14 * .e10/.e12 - a0 * .e13 * .e9/.e11
##         pnie_d1 <- .e25 * .e17
##         pnie_d2 <- .e32 * .e17
##         pnie_d3 <- c_cond * .e25 * .e17
##         if (is.null(beta3)) {
##             pnie_d4 <- rep(0, length(beta3))
##         }
##         else {
##             pnie_d4 <- c_cond * .e32 * .e17
##         }
##         pnie_d5 <- 0
##         pnie_d6 <- 0
##         pnie_d7 <- .e20
##         pnie_d8 <- a0 * .e20
##         pnie_d9 <- rep(0, length(theta4))
##         pnie_d10 <- rep(0, length(theta5))
##         if (is.null(theta6)) {
##             pnie_d11 <- rep(0, length(theta6))
##         }
##         else {
##             pnie_d11 <- c_cond * .e20
##         }
##         Gamma_pnie <- matrix(c(pnie_d1, pnie_d2, pnie_d3, pnie_d4,
##             pnie_d5, pnie_d6, pnie_d7, pnie_d8, pnie_d9, pnie_d10,
##             pnie_d11))
##         Gamma_te <- ((a1 - a0) * Gamma_pnde) + Gamma_tnie
##         pnde <- (theta1 + theta3 * expit(beta0 + beta1 * a0 +
##             beta2_c + beta3_c * a0) + theta5_c) * (a1 - a0)
##         tnie <- (theta2 + theta3 * a1 + theta6_c) * (expit(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1) - expit(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0))
##         d_pm <- grad_prop_med_yreg_linear(pnde = unname(pnde),
##             tnie = unname(tnie))
##         Gamma_pm <- (d_pm[["pnde"]] * (a1 - a0) * Gamma_pnde) +
##             (d_pm[["tnie"]] * Gamma_tnie)
##         a1_sub_a0 <- abs(a1 - a0)
##         se_cde <- sqrt(as.numeric(t(Gamma_cde) %*% Sigma %*%
##             Gamma_cde)) * a1_sub_a0
##         se_pnde <- sqrt(as.numeric(t(Gamma_pnde) %*% Sigma %*%
##             Gamma_pnde)) * a1_sub_a0
##         se_tnie <- sqrt(as.numeric(t(Gamma_tnie) %*% Sigma %*%
##             Gamma_tnie))
##         se_tnde <- sqrt(as.numeric(t(Gamma_tnde) %*% Sigma %*%
##             Gamma_tnde)) * a1_sub_a0
##         se_pnie <- sqrt(as.numeric(t(Gamma_pnie) %*% Sigma %*%
##             Gamma_pnie))
##         se_te <- sqrt(as.numeric(t(Gamma_te) %*% Sigma %*% Gamma_te))
##         se_pm <- sqrt(as.numeric(t(Gamma_pm) %*% Sigma %*% Gamma_pm))
##         c(se_cde = unname(se_cde), se_pnde = unname(se_pnde),
##             se_tnie = unname(se_tnie), se_tnde = unname(se_tnde),
##             se_pnie = unname(se_pnie), se_te = unname(se_te),
##             se_pm = unname(se_pm))
##     }
##     return(fun_se)
## }
## <bytecode: 0x7fb8897ae250>
## <environment: namespace:regmedint>

## mreg logistic yreg non-linear (V2015 p473 Proposition 2.6)

These functions are used in all cases where the mediator model is logistic regression and the outcome model is any one of the non-linear models.

### Point estimates

regmedint:::calc_myreg_mreg_logistic_yreg_logistic_est
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6)
##     fun_est <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         expit <- function(x) {
##             exp(x)/(1 + exp(x))
##         }
##         cde <- (theta1 + theta3 * m_cde + theta5_c) * (a1 - a0)
##         pnde <- ((theta1 + theta5_c) * (a1 - a0)) + log(1 + exp(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0 + theta2 + theta3 *
##             a1 + theta6_c)) - log(1 + exp(beta0 + beta1 * a0 +
##             beta2_c + beta3_c * a0 + theta2 + theta3 * a0 + theta6_c))
##         tnie <- log(1 + exp(beta0 + beta1 * a1 + beta2_c + beta3_c *
##             a1 + theta2 + theta3 * a1 + theta6_c)) - log(1 +
##             exp(beta0 + beta1 * a0 + beta2_c + beta3_c * a0 +
##                 theta2 + theta3 * a1 + theta6_c)) + log(1 + exp(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0)) - log(1 + exp(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1))
##         tnde <- ((theta1 + theta5_c) * (a1 - a0)) + log(1 + exp(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1 + theta2 + theta3 *
##             a1 + theta6_c)) - log(1 + exp(beta0 + beta1 * a1 +
##             beta2_c + beta3_c * a1 + theta2 + theta3 * a0 + theta6_c))
##         pnie <- log(1 + exp(beta0 + beta1 * a1 + beta2_c + beta3_c *
##             a1 + theta2 + theta3 * a0 + theta6_c)) - log(1 +
##             exp(beta0 + beta1 * a0 + beta2_c + beta3_c * a0 +
##                 theta2 + theta3 * a0 + theta6_c)) + log(1 + exp(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0)) - log(1 + exp(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1))
##         te <- pnde + tnie
##         pm <- (exp(pnde) * (exp(tnie) - 1))/(exp(te) - 1)
##         c(cde = unname(cde), pnde = unname(pnde), tnie = unname(tnie),
##             tnde = unname(tnde), pnie = unname(pnie), te = unname(te),
##             pm = unname(pm))
##     }
##     return(fun_est)
## }
## <bytecode: 0x7fb889b38b00>
## <environment: namespace:regmedint>

### Standard error estimates

regmedint:::calc_myreg_mreg_logistic_yreg_logistic_se
## function (beta0, beta1, beta2, beta3, theta0, theta1, theta2,
##     theta3, theta4, theta5, theta6, Sigma_beta, Sigma_theta)
## {
##     validate_myreg_coefs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6)
##     validate_myreg_vcovs(beta0 = beta0, beta1 = beta1, beta2 = beta2,
##         beta3 = beta3, theta0 = theta0, theta1 = theta1, theta2 = theta2,
##         theta3 = theta3, theta4 = theta4, theta5 = theta5, theta6 = theta6,
##         Sigma_beta = Sigma_beta, Sigma_theta = Sigma_theta)
##     Sigma <- Matrix::bdiag(Sigma_beta, Sigma_theta)
##     fun_se <- function(a0, a1, m_cde, c_cond) {
##         if (is.null(beta2)) {
##             assertthat::assert_that(is.null(c_cond))
##             beta2_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(beta2))
##             beta2_c <- sum(t(matrix(beta2)) %*% matrix(c_cond))
##         }
##         if (is.null(beta3)) {
##             beta3_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(beta3))
##             beta3_c <- sum(t(matrix(beta3)) %*% matrix(c_cond))
##         }
##         if (is.null(theta4)) {
##             assertthat::assert_that(is.null(c_cond))
##             theta4_c <- 0
##         }
##         else {
##             assertthat::assert_that(!is.null(c_cond))
##             assertthat::assert_that(length(c_cond) == length(theta4))
##             theta4_c <- sum(t(matrix(theta4)) %*% matrix(c_cond))
##         }
##         if (is.null(theta5)) {
##             theta5_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta5))
##             theta5_c <- sum(t(matrix(theta5)) %*% matrix(c_cond))
##         }
##         if (is.null(theta6)) {
##             theta6_c <- 0
##         }
##         else {
##             assertthat::assert_that(length(c_cond) == length(theta6))
##             theta6_c <- sum(t(matrix(theta6)) %*% matrix(c_cond))
##         }
##         expit <- function(x) {
##             exp(x)/(1 + exp(x))
##         }
##         if (is.null(theta5)) {
##             pd_cde_theta5 <- rep(0, length(theta5))
##         }
##         else {
##             pd_cde_theta5 <- c_cond
##         }
##         Gamma_cde <- matrix(c(0, 0, rep(0, length(beta2)), rep(0,
##             length(beta3)), 0, (a1 - a0), 0, (a1 - a0) * m_cde,
##             rep(0, length(theta4)), pd_cde_theta5, rep(0, length(theta6))))
##         pnde_expit_a1 <- expit(a0 * (beta1 + beta3_c) + a1 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         pnde_expit_a0 <- expit(a0 * (beta1 + beta3_c) + a0 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         pnde_d1 <- pnde_expit_a1 - pnde_expit_a0
##         pnde_d2 <- a0 * pnde_d1
##         pnde_d3 <- c_cond * pnde_d1
##         if (is.null(beta3)) {
##             pnde_d4 <- rep(0, length(beta3))
##         }
##         else {
##             pnde_d4 <- c_cond * a0 * pnde_d1
##         }
##         pnde_d5 <- 0
##         pnde_d6 <- a1 - a0
##         pnde_d7 <- pnde_d1
##         pnde_d8 <- a1 * pnde_expit_a1 - a0 * pnde_expit_a0
##         pnde_d9 <- rep(0, length(theta4))
##         if (is.null(theta5)) {
##             pnde_d10 <- rep(0, length(theta5))
##         }
##         else {
##             pnde_d10 <- c_cond * (a1 - a0)
##         }
##         if (is.null(theta6)) {
##             pnde_d11 <- rep(0, length(theta6))
##         }
##         else {
##             pnde_d11 <- c_cond * pnde_d1
##         }
##         Gamma_pnde <- matrix(c(pnde_d1, pnde_d2, pnde_d3, pnde_d4,
##             pnde_d5, pnde_d6, pnde_d7, pnde_d8, pnde_d9, pnde_d10,
##             pnde_d11))
##         tnie_expit_q1 <- expit(a0 * (beta1 + beta3_c) + beta0 +
##             beta2_c)
##         tnie_expit_q2 <- expit(a0 * (beta1 + beta3_c) + a1 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         tnie_expit_q3 <- expit(a1 * (beta1 + beta3_c) + beta0 +
##             beta2_c)
##         tnie_expit_q4 <- expit(a1 * (beta1 + beta3_c) + a1 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         tnie_d1 <- tnie_expit_q1 - tnie_expit_q2 - tnie_expit_q3 +
##             tnie_expit_q4
##         tnie_d2 <- a0 * (tnie_expit_q1 - tnie_expit_q2) + a1 *
##             (-tnie_expit_q3 + tnie_expit_q4)
##         tnie_d3 <- c_cond * tnie_d1
##         if (is.null(beta3)) {
##             tnie_d4 <- rep(0, length(beta3))
##         }
##         else {
##             tnie_d4 <- c_cond * tnie_d2
##         }
##         tnie_d5 <- 0
##         tnie_d6 <- 0
##         tnie_d7 <- tnie_expit_q4 - tnie_expit_q2
##         tnie_d8 <- a1 * tnie_d7
##         tnie_d9 <- rep(0, length(theta4))
##         tnie_d10 <- rep(0, length(theta5))
##         if (is.null(theta6)) {
##             tnie_d11 <- rep(0, length(theta6))
##         }
##         else {
##             tnie_d11 <- c_cond * tnie_d7
##         }
##         Gamma_tnie <- matrix(c(tnie_d1, tnie_d2, tnie_d3, tnie_d4,
##             tnie_d5, tnie_d6, tnie_d7, tnie_d8, tnie_d9, tnie_d10,
##             tnie_d11))
##         tnde_expit_a1 <- expit(a1 * (beta1 + beta3_c) + a1 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         tnde_expit_a0 <- expit(a1 * (beta1 + beta3_c) + a0 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         tnde_d1 <- tnde_expit_a1 - tnde_expit_a0
##         tnde_d2 <- a1 * tnde_d1
##         tnde_d3 <- c_cond * tnde_d1
##         if (is.null(beta3)) {
##             tnde_d4 <- rep(0, length(beta3))
##         }
##         else {
##             tnde_d4 <- c_cond * a1 * tnde_d1
##         }
##         tnde_d5 <- 0
##         tnde_d6 <- a1 - a0
##         tnde_d7 <- tnde_d1
##         tnde_d8 <- a1 * tnde_expit_a1 - a0 * tnde_expit_a0
##         tnde_d9 <- rep(0, length(theta4))
##         if (is.null(theta5)) {
##             tnde_d10 <- rep(0, length(theta5))
##         }
##         else {
##             tnde_d10 <- c_cond * (a1 - a0)
##         }
##         if (is.null(theta6)) {
##             tnde_d11 <- rep(0, length(theta6))
##         }
##         else {
##             tnde_d11 <- tnde_d3
##         }
##         Gamma_tnde <- matrix(c(tnde_d1, tnde_d2, tnde_d3, tnde_d4,
##             tnde_d5, tnde_d6, tnde_d7, tnde_d8, tnde_d9, tnde_d10,
##             tnde_d11))
##         pnie_expit_q1 <- expit(a0 * (beta1 + beta3_c) + beta0 +
##             beta2_c)
##         pnie_expit_q2 <- expit(a0 * (beta1 + beta3_c) + a0 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         pnie_expit_q3 <- expit(a1 * (beta1 + beta3_c) + beta0 +
##             beta2_c)
##         pnie_expit_q4 <- expit(a1 * (beta1 + beta3_c) + a0 *
##             theta3 + beta0 + beta2_c + theta6_c + theta2)
##         pnie_d1 <- pnie_expit_q1 - pnie_expit_q2 - pnie_expit_q3 +
##             pnie_expit_q4
##         pnie_d2 <- a0 * (pnie_expit_q1 - pnie_expit_q2) + a1 *
##             (-pnie_expit_q3 + pnie_expit_q4)
##         pnie_d3 <- c_cond * pnie_d1
##         if (is.null(beta3)) {
##             pnie_d4 <- rep(0, length(beta3))
##         }
##         else {
##             pnie_d4 <- c_cond * pnie_d2
##         }
##         pnie_d5 <- 0
##         pnie_d6 <- 0
##         pnie_d7 <- pnie_expit_q4 - pnie_expit_q2
##         pnie_d8 <- a0 * pnie_d7
##         pnie_d9 <- rep(0, length(theta4))
##         pnie_d10 <- rep(0, length(theta5))
##         if (is.null(theta6)) {
##             pnie_d11 <- rep(0, length(theta6))
##         }
##         else {
##             pnie_d11 <- c_cond * pnie_d7
##         }
##         Gamma_pnie <- matrix(c(pnie_d1, pnie_d2, pnie_d3, pnie_d4,
##             pnie_d5, pnie_d6, pnie_d7, pnie_d8, pnie_d9, pnie_d10,
##             pnie_d11))
##         Gamma_te <- Gamma_pnde + Gamma_tnie
##         pnde <- ((theta1 + theta5_c) * (a1 - a0)) + log(1 + exp(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0 + theta2 + theta3 *
##             a1 + theta6_c)) - log(1 + exp(beta0 + beta1 * a0 +
##             beta2_c + beta3_c * a0 + theta2 + theta3 * a0 + theta6_c))
##         tnie <- log(1 + exp(beta0 + beta1 * a1 + beta2_c + beta3_c *
##             a1 + theta2 + theta3 * a1 + theta6_c)) - log(1 +
##             exp(beta0 + beta1 * a0 + beta2_c + beta3_c * a0 +
##                 theta2 + theta3 * a1 + theta6_c)) + log(1 + exp(beta0 +
##             beta1 * a0 + beta2_c + beta3_c * a0)) - log(1 + exp(beta0 +
##             beta1 * a1 + beta2_c + beta3_c * a1))
##         d_pm <- grad_prop_med_yreg_logistic(pnde = unname(pnde),
##             tnie = unname(tnie))
##         Gamma_pm <- (d_pm[["pnde"]] * Gamma_pnde) + (d_pm[["tnie"]] *
##             Gamma_tnie)
##         se_cde <- sqrt(as.numeric(t(Gamma_cde) %*% Sigma %*%
##             Gamma_cde))
##         se_pnde <- sqrt(as.numeric(t(Gamma_pnde) %*% Sigma %*%
##             Gamma_pnde))
##         se_tnie <- sqrt(as.numeric(t(Gamma_tnie) %*% Sigma %*%
##             Gamma_tnie))
##         se_tnde <- sqrt(as.numeric(t(Gamma_tnde) %*% Sigma %*%
##             Gamma_tnde))
##         se_pnie <- sqrt(as.numeric(t(Gamma_pnie) %*% Sigma %*%
##             Gamma_pnie))
##         se_te <- sqrt(as.numeric(t(Gamma_te) %*% Sigma %*% Gamma_te))
##         se_pm <- sqrt(as.numeric(t(Gamma_pm) %*% Sigma %*% Gamma_pm))
##         c(se_cde = unname(se_cde), se_pnde = unname(se_pnde),
##             se_tnie = unname(se_tnie), se_tnde = unname(se_tnde),
##             se_pnie = unname(se_pnie), se_te = unname(se_te),
##             se_pm = unname(se_pm))
##     }
##     return(fun_se)
## }
## <bytecode: 0x7fb888d8ffd8>
## <environment: namespace:regmedint>

## Bibliography

• V2015: VanderWeele (2015) Explanation in Causal Inference