In any association study, it is important to distinguish direct and indirect effects in order to build truly functional models. For this purpose, we consider a directed acyclic graph (DAG) setting with an exposure variable, primary and intermediate outcome variables, and confounding factors. In order to make valid statistical inference on the direct effect of the exposure on the primary outcome, it is necessary to consider all potential effects in the graph, and we propose to use the estimating equation method with robust Huber-White sandwich standard errors. Then, a large-sample Wald-type test statistic is computed for testing the absence of the direct effect. In this package, the proposed causal inference method based on estimating equations (CIEE) is implemented for both the analysis of continuous and time-to-event primary outcomes subject to censoring for the model in Figure 1. Additionally, standard multiple regression, regression of residuals, and the structural equation approach are implemented for fitting the same model.
Results from simulation studies (Konigorski et al., 2018) showed that CIEE successfully removes the effect of intermediate outcomes from the primary outcome and is robust against measured and unmeasured confounding of the indirect effect through observed factors. Also, an application in a genetic association study in the same study showed that CIEE can identify genetic variants that would be missed by traditional regression methods. Both multiple regression methods and the structural equation method fail in some scenarios where their corresponding test statistics lead to inflated type I errors. An alternative approach for the analysis of continuous traits is the sequential G-estimation method (Vansteelandt et al., 2009).
In this package, CIEE is implemented for the model described in the DAG in Figure 1, which includes the direct effect \(\alpha_{XY}\) of an exposure X on the primary outcome Y and an indirect effect of X on Y through a secondary outcome K. The model further includes measured and unmeasured factors L and U, respectively, which potentially confound the effect of K on Y. CIEE can also be applied to different models with different error distributions. The goal is to estimate and test the direct effect \(\alpha_{XY}\), while removing the indirect effect of X on Y through K, and with robustness against effects of L and U. Without restriction of generality, it is assumed that there aren’t any factors affecting X and that any such factors are included as covariates in the analysis or have been dealt with using other approaches. Also, we generally assume that either \(\alpha_{LY}=0\) (L is a factor influencing K) or \(\alpha_{XL}=0\) (L is a measured confounder of \(K \to Y\)). Otherwise, the effect of L as intermediate outcome could be removed from Y in the analysis analogously to K.
Two traditional methods for the aim to estimate and test \(\alpha_{XY}\) are (i) to include the intermediate outcomes and factors as covariates in a multiple regression (MR) model of the primary outcome on the exposure, or (ii) to first regress the primary outcome on the intermediate outcome and factors, and then regress the extracted residuals on the exposure (regression of residuals, RR). In more detail, estimates of \(\alpha_{XY}\) are obtained from fitting the following models in the quantitative outcome setting for a normally-distributed Y (GLM setting):
MR: Obtain the least squares (LS) estimate of \(\alpha_{XY}\) by fitting \[Y_i = \alpha_0 + \alpha_{XY} x_i + \alpha_1 k_i + \alpha_2 l_i + \varepsilon_i, \ \varepsilon_i \sim N(0,\sigma_1^2)\].
RR: First, obtain residuals \(\hat{\epsilon}_{1i} = y_i - \hat{\alpha}_0 -
\hat{\alpha}_1 k_i - \hat{\alpha}_2 l_i\) by fitting
\[Y_i = \alpha_0 + \alpha_1 k_i + \alpha_2
l_i + \varepsilon_{1i}, \ \varepsilon_{1i} \sim N(0,\sigma_1^2
)\] using the LS estimation. Second, obtain the LS estimate of
\(\alpha_{XY}\) by fitting \[\hat{\varepsilon}_{1i} = \alpha_3 + \alpha_{XY}
x_i + \varepsilon_{2i}, \ \varepsilon_{2i} \sim N(0,
\sigma_2^2).\]
Then, \(H_0: \alpha_{XY} = 0\)
versus \(H_A: \alpha_{XY} \neq 0\) is
tested using the default t-test in the lm() function in R.
For the analysis of a censored time-to-event primary trait Y
(accelerated failure time setting; AFT), only the MR approach is
implemented. Here, the equivalent censored log-linear regression model
is fitted using the survreg() function in the
survival R package to obtain the maximum likelihood
estimate of \(\alpha_{XY}\), and a
Wald-type test is performed for testing the null hypothesis \(H_0: \alpha_{XY} = 0\). Both approaches are
implemented in the functions mult_reg() and
res_reg() and can be used as follows. For this
illustration, data is first generated using the
generate_data() function, which generates data for the
quantitative outcomes Y and K, a genetic marker X (single nucleotide
polymorphism, SNP, taking values 0, 1, 2) as exposure, and observed as
well as unobserved confounders L, U.
dat <- generate_data(setting="GLM", n = 1000, maf = 0.2, cens = 0.3, a = NULL, b = NULL,
aUL = 0, aXL = 0, aXK = 0.2, aLK = 0, aUY = 0, aKY = 0.3, aXY = 0.1,
aLY = 0, mu_U = 0, sd_U = 1, X_orth_U = TRUE, mu_X = NULL,
sd_X = NULL, mu_L = 0, sd_L = 1, mu_K = 0, sd_K = 1, mu_Y = 0,
sd_Y = 1)
head(dat)## Y K X L U
## 1 -0.5393958 -1.4012617 1 1.5752288 -0.45602230
## 2 0.9072012 -0.7548686 0 -0.7398913 0.03192698
## 3 -0.1980821 1.1968833 1 -0.4357265 1.47753546
## 4 -1.2654489 -2.4168304 0 -0.6935100 0.95051614
## 5 0.6307854 0.9748694 1 -0.4385381 -0.24502667
## 6 -0.1440354 0.6684717 0 0.1108722 -0.50603496## $point_estimates
## alpha_0 alpha_1 alpha_XY alpha_2
## 0.001251528 0.280530841 0.128234963 -0.003962981
##
## $SE_estimates
## alpha_0 alpha_1 alpha_XY alpha_2
## 0.03822808 0.03097598 0.05722875 0.03044872
##
## $pvalues
## alpha_0 alpha_1 alpha_XY alpha_2
## 9.738897e-01 6.954293e-19 2.526244e-02 8.964720e-01## $point_estimates
## alpha_0 alpha_1 alpha_2 alpha_3 alpha_XY
## 0.050745720 0.284617731 -0.001622808 -0.049648634 0.127631450
##
## $SE_estimates
## alpha_0 alpha_1 alpha_2 alpha_3 alpha_XY
## 0.03126371 0.03098453 0.03049211 0.03813408 0.05703737
##
## $pvalues
## alpha_0 alpha_1 alpha_2 alpha_3 alpha_XY
## 1.048733e-01 2.306295e-19 9.575668e-01 1.932341e-01 2.546185e-02As another approach for modeling DAGs, the structural equation
modeling method (SEM; Bollen, 1989) can be used. Among others, it is
implemented in the sem() function of the
lavaan package (Rosseel, 2012). For a comparison of the
results, the function sem_appl() applies the SEM method to
the DAG in Figure 1 based on the following model equations:
\[L_i = \alpha_0 + \alpha_1 x_i + \varepsilon_{1i}, \ \varepsilon_{1i} \sim N(0,\sigma_1^2 ) \] \[K_i = \alpha_2 + \alpha_3 x_i + \alpha_2 l_i + \varepsilon_{2i}, \ \varepsilon_{2i} \sim N(0,\sigma_2^2 ) \] \[Y_i = \alpha_5 + \alpha_6 k_i + \alpha_{XY} x_i + \varepsilon_{3i}, \ \varepsilon_{3i} \sim N(0,\sigma_3^2 ) \]
## $point_estimates
## alpha_1 alpha_3 alpha_4 alpha_6 alpha_XY
## 0.06617423 0.10878398 0.01503394 0.28046918 0.12797949
##
## $SE_estimates
## alpha_1 alpha_3 alpha_4 alpha_6 alpha_XY
## 0.05930231 0.05832235 0.03108087 0.03091062 0.05708106
##
## $pvalues
## alpha_1 alpha_3 alpha_4 alpha_6 alpha_XY
## 2.644738e-01 6.215067e-02 6.285961e-01 1.151959e-19 2.495713e-02Further proposed approaches for fitting the model in Figure 1 include the two-stage sequential G-estimation method (Vansteelandt et al., 2009). It first removes the effect of K from the primary outcome Y, and then tests the association of X with the adjusted primary outcome.
In more detail, for the analyis of a quantitative outcome Y with n independent observations, in the first stage, the effect of K on Y, \(\alpha_1\), is estimated and \(\hat{\alpha}_1\) is obtained using the LS estimation method under the model
\[\begin{equation} Y_i = \alpha_0 + \alpha_1 k_i + \alpha_2 x_i + \alpha_3 l_i + \varepsilon_i, \ \varepsilon_i \sim N(0,\sigma_1^2), \quad i=1,…,n \end{equation}\]
Then, to block all indirect paths from X to Y, the adjusted outcome \(\tilde{Y}\) is obtained by removing the effect of K on Y with
\[\begin{equation} \tilde{y}_i = y_i - \bar{y} - \hat{\alpha}_1 (k_i - \bar{k}) \end{equation}\]
where \(\bar{y} = \sum_{i=1}^n y_i\) and \(\bar{k} = \sum_{i=1}^n k_i\). In the second stage, the significance of the direct effect of X on Y, \(\alpha_{XY}\), is tested under the model
\[\begin{equation} \tilde{Y}_i = \alpha_4 + \alpha_{XY} x_i + \varepsilon_i, \ \varepsilon_i \sim N(0,\sigma_2^2 ) \end{equation}\]
using the proposed test statistic in Vansteelandt and colleagues
(2009). The approach is implemented in the CGene package,
which can be obtained from cran.r-project.org/src/contrib/Archive/CGene/.
CIEE follows the general idea of the two-stage sequential G-estimation method with the major difference that the approach is one-stage and obtains coefficient estimates of all parameters simultaneously by solving estimating equations. This also allows building on existing asymptotic properties of the estimator and obtaining robust sandwich standard error estimates considering the additional variability of the estimates from the outcome adjustment.
In more detail, for the analyis of a quantitative outcome Y, we formulate unbiased estimating equations \(U(\theta)=0\) for a consistent estimation of the unknown parameter vector \(\theta = (\alpha_0, \alpha_1, \alpha_2, \alpha_3, \sigma_1^2, \alpha_4, \alpha_{XY}, \sigma_2^2)^T\) where
\[ U(\theta) = \left( \frac{\partial l_1(\theta)}{\partial \alpha_0}, \frac{\partial l_1(\theta)}{\partial \alpha_1}, \frac{\partial l_1(\theta)}{\partial \alpha_2}, \frac{\partial l_1(\theta)}{\partial \alpha_3}, \frac{\partial l_1(\theta)}{\partial \sigma_1^2}, \frac{\partial l_1(\theta)}{\partial \alpha_4}, \frac{\partial l_1(\theta)}{\partial \alpha_{XY}}, \frac{\partial l_1(\theta)}{\partial \sigma_2^2} \right)^T \] with
\[ l_1(\theta) = \sum_{i=1}^n \left[ -log(\sigma_1) + log\left( \phi\left( \frac{y_i - \alpha_0 - \alpha_1 k_i - \alpha_2 x_i - \alpha_3 l_i}{\sigma_1} \right) \right) \right]\]
and
\[ l_2(\theta) = \sum_{i=1}^n \left[ -log(\sigma_2) + log\left( \phi\left( \frac{y_i - \bar{y} - \alpha_1 (k_i - \bar{k}) - \alpha_4 - \alpha_{XY} x_i}{\sigma_2} \right) \right) \right],\]
where \(\phi\) is the probability density function of the standard normal distribution.
By solving the first five estimating equations based on \(l_1(\theta)\), we are hence obtaining estimates of \(\alpha_0, \alpha_1, \alpha_2, \alpha_3, \sigma_1^2\). Analogously, solving the last three estimating equations based on \(l_2(\theta)\) yields estimates of \(\alpha_4, \alpha_{XY}, \sigma_2^2\). To give an intuition on how these estimating equations are obtained, \(l_1(\theta)\) is the log-likelihood function under the model in (1) and \(l_2(\theta)\) is the log-likelihood function under the model in (3) given that \(\alpha_1\) is known. All parameters in \(\theta\) are estimated simultaneously and the additional variability obtained in the outcome adjustment in (2) is considered by using the robust Huber-White sandwich estimator of the standard error of \(\hat{\theta}\). Then, the large sample Wald-type test statistic \(\hat{\theta}/\widehat{SE}(\hat{\theta})\), which has an asymptotic standard normal distribution, is computed to test the absence of the exposure effect \(\alpha_{XY}\).
For the analysis of a censored time-to-event primary outcome Y, the estimating equations can be constructed as described above for a quantitative primary phenotype (for the same parameters except for \(\sigma_1\) instead of \(\sigma_1^2\)), but in order to remove the effect of K from Y, the true underlying log survival times \(Y_{est}\) need to be estimated for censored survival times. For uncensored survival times, \(Y_{est}\) equals the observed log-survival time Y. Then, the estimating equations are constructed accordingly, the robust Huber-White sandwich estimator of the standard error is obtained and the large-sample Wald-type tests are computed. For more statistical details, see Konigorski et al. (2018).
In the implementation of CIEE in this package, the
est_funct_expr() function contains the estimating equations
as an expression.
## $logL1
## expression(log((1/sqrt(sigma1sq)) * dnorm((y_i - alpha0 - alpha1 *
## k_i - alpha2 * x_i - alpha3 * l_i)/sqrt(sigma1sq), mean = 0,
## sd = 1)))
##
## $logL2
## expression(log((1/sqrt(sigma2sq)) * dnorm((y_i - y_bar - alpha1 *
## (k_i - k_bar) - alpha4 - alphaXY * x_i)/sqrt(sigma2sq), mean = 0,
## sd = 1)))## $logL1
## expression(-c_i * log(sigma1) + c_i * log(dnorm((y_i - alpha0 -
## alpha1 * k_i - alpha2 * x_i - alpha3 * l_i)/sigma1, mean = 0,
## sd = 1)) + (1 - c_i) * log(1 - pnorm((y_i - alpha0 - alpha1 *
## k_i - alpha2 * x_i - alpha3 * l_i)/sigma1, mean = 0, sd = 1)))
##
## $logL2
## expression(log((1/sqrt(sigma2sq)) * dnorm(((c_i * y_i + (1 -
## c_i) * ((alpha0 + alpha1 * k_i + alpha2 * x_i + alpha3 *
## l_i) + (sigma1 * dnorm((y_i - alpha0 - alpha1 * k_i - alpha2 *
## x_i - alpha3 * l_i)/sigma1, mean = 0, sd = 1)/(1 - pnorm((y_i -
## alpha0 - alpha1 * k_i - alpha2 * x_i - alpha3 * l_i)/sigma1,
## mean = 0, sd = 1))))) - y_adj_bar - alpha1 * (k_i - k_bar) -
## alpha4 - alphaXY * x_i)/sqrt(sigma2sq), mean = 0, sd = 1)))The function get_estimates() obtains estimates of the
parameters in the models (1)-(3) by using the lm() and
survreg() functions for computational purposes. The
estimates are identical to estimates obtained by solving the estimating
equations.
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4
## 0.001251528 0.280530841 0.128234963 -0.003962981 0.959974401 -0.049781387
## alpha_XY sigma_2_sq
## 0.127972716 0.959990732In order to compute the robust Huber-White sandwich estimator of the
parameters, in a first step, the deriv_obj() function
computes the expression of all first and second derivatives for the 8
parameters \(\alpha_0, \alpha_1, \alpha_2,
\alpha_3, \sigma_1^2, \alpha_4, \alpha_{XY}, \sigma_2^2\) by
using the expressions from the est_funct_expr() function as
input. Then, the numerical values of all first and second derivatives
are obtained for the observed data and parameter point estimates for all
observed individuals, using the scores() and
hessian() functions.
derivobj <- deriv_obj(setting = "GLM", logL1 = estfunct$logL1, logL2 = estfunct$logL2,
Y = dat$Y, X = dat$X, K = dat$K, L = dat$L, estimates = estimates)
names(derivobj)## [1] "logL1_deriv" "logL2_deriv"## alpha0 alpha1 alpha2 alpha3 sigma1sq alpha4 alphaXY
## [1,] -0.2807810 0.3934477 -0.2807810 -0.44229432 -0.4814282 0 0
## [2,] 1.1612616 -0.8766000 0.0000000 -0.85920733 0.1534170 0 0
## [3,] -0.6927873 -0.8291856 -0.6927873 0.30186575 -0.2808701 0 0
## [4,] -0.6161136 1.4890421 0.0000000 0.42728094 -0.3310492 0 0
## [5,] 0.2355063 0.2295879 0.2355063 -0.10327850 -0.4931156 0 0
## [6,] -0.3462327 -0.2314468 0.0000000 -0.03838758 -0.4609087 0 0
## sigma2sq
## [1,] 0
## [2,] 0
## [3,] 0
## [4,] 0
## [5,] 0
## [6,] 0## alpha0 alpha1 alpha2 alpha3 sigma1sq alpha4
## [1,] -0.2807810 0.3934477 -0.2807810 -0.44229432 -0.4814282 -0.2871468
## [2,] 1.1612616 -0.8766000 0.0000000 -0.85920733 0.1534170 1.1641553
## [3,] -0.6927873 -0.8291856 -0.6927873 0.30186575 -0.2808701 -0.6908446
## [4,] -0.6161136 1.4890421 0.0000000 0.42728094 -0.3310492 -0.6133812
## [5,] 0.2355063 0.2295879 0.2355063 -0.10327850 -0.4931156 0.2374449
## [6,] -0.3462327 -0.2314468 0.0000000 -0.03838758 -0.4609087 -0.3468255
## alphaXY sigma2sq
## [1,] -0.2871468 -0.4796117
## [2,] 0.0000000 0.1567904
## [3,] -0.6908446 -0.2822053
## [4,] 0.0000000 -0.3327201
## [5,] 0.2374449 -0.4926483
## [6,] 0.0000000 -0.4606944## num [1:1000, 1:8, 1:8] -1.04 -1.04 -1.04 -1.04 -1.04 ...
## - attr(*, "dimnames")=List of 3
## ..$ : NULL
## ..$ : chr [1:8] "alpha0" "alpha1" "alpha2" "alpha3" ...
## ..$ : chr [1:8] "alpha0" "alpha1" "alpha2" "alpha3" ...The robust Huber-White sandwich estimator of the standard error can
then be obtained using the sandwich_se() function:
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4 alpha_XY
## 0.03839344 0.03162343 0.05723497 0.03101816 0.04340919 0.03856425 0.05729095
## sigma_2_sq
## 0.04342725Alternatively, bootstrap standard error estimates can be computed
using the bootstrap_se() function. Also, for comparison,
the function naive_se() computes naive standard error
estimates of the parameter estimates of \(\alpha_0, \alpha_1, \alpha_2, \alpha_3, \alpha_4,
\alpha_{XY}\) without accounting for the additional variability
due to the two stages in the model in (1)-(3):
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4 alpha_XY
## 0.03920897 0.03148622 0.05681736 0.03036809 0.04496922 0.02232444 0.05669928
## sigma_2_sq
## 0.04500299## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4 alpha_XY
## 0.03822808 0.03097598 0.05722875 0.03044872 NA 0.03813395 0.05703718
## sigma_2_sq
## NAFinally, the functions ciee() and
ciee_loop() allow an easy integrated use of all above
functions and a simultaneous computation of the estimating equations
approach using either standard error computation, the traditional
regression-based approaches, and the SEM method. ciee()
fits the model in equations (1)-(3) (e.g. the model in Figure 1) and
yields parameter estimates, standard error estimates, and p-values for
all parameters. ciee_loop() provides an extension of
ciee() and allows the input of multiple exposure variables
(e.g. multiple SNPs) which are tested sequentially. In the output of
ciee_loop(), only the coefficient estimates, standard error
estimates, and p-values with respect to the direct effect \(\alpha_{XY}\) are provided.
results_ciee <- ciee(setting = "GLM", Y = dat$Y, X = dat$X, K = dat$K, L = dat$L,
estimates = c("ee", "mult_reg", "res_reg", "sem"),
ee_se = "sandwich")
results_ciee## $results_ee
## $results_ee$point_estimates
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4
## 0.001251528 0.280530841 0.128234963 -0.003962981 0.959974401 -0.049781387
## alpha_XY sigma_2_sq
## 0.127972716 0.959990732
##
## $results_ee$SE_estimates
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4 alpha_XY
## 0.03839344 0.03162343 0.05723497 0.03101816 0.04340919 0.03856425 0.05729095
## sigma_2_sq
## 0.04342725
##
## $results_ee$wald_test_stat
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq alpha_4
## 0.03259746 8.87098157 2.24050016 -0.12776326 22.11454423 -1.29086878
## alpha_XY sigma_2_sq
## 2.23373342 22.10572002
##
## $results_ee$pvalues
## alpha_0 alpha_1 alpha_2 alpha_3 sigma_1_sq
## 9.739956e-01 7.250396e-19 2.505847e-02 8.983363e-01 2.290100e-108
## alpha_4 alpha_XY sigma_2_sq
## 1.967492e-01 2.550062e-02 2.784581e-108
##
##
## $results_mult_reg
## $results_mult_reg$point_estimates
## alpha_0 alpha_1 alpha_XY alpha_2
## 0.001251528 0.280530841 0.128234963 -0.003962981
##
## $results_mult_reg$SE_estimates
## alpha_0 alpha_1 alpha_XY alpha_2
## 0.03822808 0.03097598 0.05722875 0.03044872
##
## $results_mult_reg$pvalues
## alpha_0 alpha_1 alpha_XY alpha_2
## 9.738897e-01 6.954293e-19 2.526244e-02 8.964720e-01
##
##
## $results_res_reg
## $results_res_reg$point_estimates
## alpha_0 alpha_1 alpha_2 alpha_3 alpha_XY
## 0.050745720 0.284617731 -0.001622808 -0.049648634 0.127631450
##
## $results_res_reg$SE_estimates
## alpha_0 alpha_1 alpha_2 alpha_3 alpha_XY
## 0.03126371 0.03098453 0.03049211 0.03813408 0.05703737
##
## $results_res_reg$pvalues
## alpha_0 alpha_1 alpha_2 alpha_3 alpha_XY
## 1.048733e-01 2.306295e-19 9.575668e-01 1.932341e-01 2.546185e-02
##
##
## $results_sem
## $results_sem$point_estimates
## alpha_1 alpha_3 alpha_4 alpha_6 alpha_XY
## 0.06617423 0.10878398 0.01503394 0.28046918 0.12797949
##
## $results_sem$SE_estimates
## alpha_1 alpha_3 alpha_4 alpha_6 alpha_XY
## 0.05930231 0.05832235 0.03108087 0.03091062 0.05708106
##
## $results_sem$pvalues
## alpha_1 alpha_3 alpha_4 alpha_6 alpha_XY
## 2.644738e-01 6.215067e-02 6.285961e-01 1.151959e-19 2.495713e-02
##
##
## attr(,"class")
## [1] "ciee"maf <- 0.2
n <- 1000
dat <- generate_data(n = n, maf = maf)
datX <- data.frame(X = dat$X)
names(datX)[1] <- "X1"
for(i in 2:10){
X <- rbinom(n, size = 2, prob = maf)
datX$X <- X
names(datX)[i] <- paste("X", i, sep="")
}
results_ciee_loop <- ciee_loop(setting = "GLM", Y = dat$Y, X = datX, K = dat$K, L = dat$L)
results_ciee_loop## $results_ee
## $results_ee$point_estimates
## X1 X2 X3 X4 X5 X6
## 0.063077824 0.017219872 0.063068442 0.103458279 -0.030938837 0.020921676
## X7 X8 X9 X10
## 0.002527108 0.072445972 0.164248588 -0.003469517
##
## $results_ee$SE_estimates
## X1 X2 X3 X4 X5 X6 X7
## 0.05173170 0.05134489 0.05678258 0.05214891 0.05194370 0.05632630 0.05182672
## X8 X9 X10
## 0.05431966 0.05265740 0.05328431
##
## $results_ee$wald_test_stat
## X1 X2 X3 X4 X5 X6
## 1.21932637 0.33537658 1.11070051 1.98390095 -0.59562253 0.37143707
## X7 X8 X9 X10
## 0.04876071 1.33369714 3.11919290 -0.06511329
##
## $results_ee$pvalues
## X1 X2 X3 X4 X5 X6
## 0.222720345 0.737341029 0.266697283 0.047266883 0.551427425 0.710312021
## X7 X8 X9 X10
## 0.961109991 0.182303132 0.001813472 0.948083798
##
##
## $results_mult_reg
## $results_mult_reg$point_estimates
## X1 X2 X3 X4 X5 X6
## 0.062628442 0.017644859 0.060999591 0.104372791 -0.031949076 0.021189455
## X7 X8 X9 X10
## 0.003105817 0.072015147 0.163694259 -0.003600446
##
## $results_mult_reg$SE_estimates
## X1 X2 X3 X4 X5 X6 X7
## 0.05320910 0.05236999 0.05579270 0.05231128 0.05466858 0.05439856 0.05196940
## X8 X9 X10
## 0.05347421 0.04996759 0.05269134
##
## $results_mult_reg$pvalues
## X1 X2 X3 X4 X5 X6
## 0.239466722 0.736242990 0.274515196 0.046290384 0.559074221 0.696973076
## X7 X8 X9 X10
## 0.952356843 0.178374662 0.001089184 0.945535953
##
##
## $results_res_reg
## $results_res_reg$point_estimates
## X1 X2 X3 X4 X5 X6
## 0.061436896 0.017640202 0.060380185 0.104154953 -0.031897823 0.021183820
## X7 X8 X9 X10
## 0.003103701 0.071981218 0.162710874 -0.003597357
##
## $results_res_reg$SE_estimates
## X1 X2 X3 X4 X5 X6 X7
## 0.05264836 0.05231059 0.05545340 0.05220449 0.05456996 0.05433680 0.05189961
## X8 X9 X10
## 0.05340804 0.04976895 0.05261593
##
## $results_res_reg$pvalues
## X1 X2 X3 X4 X5 X6
## 0.243517852 0.736021598 0.276484888 0.046299820 0.558995364 0.696722257
## X7 X8 X9 X10
## 0.952325297 0.178041484 0.001114884 0.945504713
##
##
## $results_sem
## $results_sem$point_estimates
## X1 X2 X3 X4 X5 X6
## 0.063229060 0.017216199 0.063028121 0.103417368 -0.030916939 0.020936983
## X7 X8 X9 X10
## 0.002543383 0.072460179 0.164185505 -0.003438170
##
## $results_sem$SE_estimates
## X1 X2 X3 X4 X5 X6 X7
## 0.05309423 0.05226550 0.05546768 0.05218873 0.05453222 0.05429461 0.05186050
## X8 X9 X10
## 0.05336572 0.04985986 0.05259176
##
## $results_sem$pvalues
## X1 X2 X3 X4 X5 X6
## 0.2336991512 0.7418541850 0.2558295692 0.0475236488 0.5707494289 0.6997795567
## X7 X8 X9 X10
## 0.9608852001 0.1745259155 0.0009914583 0.9478756754
##
##
## attr(,"class")
## [1] "ciee"Both ciee() and ciee_loop() return
ciee objects as output, and the implemented
summary.ciee() function can be used through the generic
summary() to provide a reader-friendly formatted output of
the results.
## [1] "Results based on estimating equations."
## point_estimates SE_estimates wald_test_stat pvalues
## CIEE_alpha_0 0.001251528 0.03839344 0.03259746 9.739956e-01
## CIEE_alpha_1 0.280530841 0.03162343 8.87098157 7.250396e-19
## CIEE_alpha_2 0.128234963 0.05723497 2.24050016 2.505847e-02
## CIEE_alpha_3 -0.003962981 0.03101816 -0.12776326 8.983363e-01
## CIEE_sigma_1_sq 0.959974401 0.04340919 22.11454423 2.290100e-108
## CIEE_alpha_4 -0.049781387 0.03856425 -1.29086878 1.967492e-01
## CIEE_alpha_XY 0.127972716 0.05729095 2.23373342 2.550062e-02
## CIEE_sigma_2_sq 0.959990732 0.04342725 22.10572002 2.784581e-108
## [1] "Results based on traditional multiple regression."
## point_estimates SE_estimates pvalues
## MR_alpha_0 0.001251528 0.03822808 9.738897e-01
## MR_alpha_1 0.280530841 0.03097598 6.954293e-19
## MR_alpha_XY 0.128234963 0.05722875 2.526244e-02
## MR_alpha_2 -0.003962981 0.03044872 8.964720e-01
## [1] "Results based on traditional regression of residuals."
## point_estimates SE_estimates pvalues
## RR_alpha_0 0.050745720 0.03126371 1.048733e-01
## RR_alpha_1 0.284617731 0.03098453 2.306295e-19
## RR_alpha_2 -0.001622808 0.03049211 9.575668e-01
## RR_alpha_3 -0.049648634 0.03813408 1.932341e-01
## RR_alpha_XY 0.127631450 0.05703737 2.546185e-02
## [1] "Results based on structural equation modeling."
## point_estimates SE_estimates pvalues
## SEM_alpha_1 0.06617423 0.05930231 2.644738e-01
## SEM_alpha_3 0.10878398 0.05832235 6.215067e-02
## SEM_alpha_4 0.01503394 0.03108087 6.285961e-01
## SEM_alpha_6 0.28046918 0.03091062 1.151959e-19
## SEM_alpha_XY 0.12797949 0.05708106 2.495713e-02## [1] "Results based on estimating equations."
## point_estimates SE_estimates wald_test_stat pvalues
## CIEE_X1 0.063077824 0.05173170 1.21932637 0.222720345
## CIEE_X2 0.017219872 0.05134489 0.33537658 0.737341029
## CIEE_X3 0.063068442 0.05678258 1.11070051 0.266697283
## CIEE_X4 0.103458279 0.05214891 1.98390095 0.047266883
## CIEE_X5 -0.030938837 0.05194370 -0.59562253 0.551427425
## CIEE_X6 0.020921676 0.05632630 0.37143707 0.710312021
## CIEE_X7 0.002527108 0.05182672 0.04876071 0.961109991
## CIEE_X8 0.072445972 0.05431966 1.33369714 0.182303132
## CIEE_X9 0.164248588 0.05265740 3.11919290 0.001813472
## CIEE_X10 -0.003469517 0.05328431 -0.06511329 0.948083798
## [1] "Results based on traditional multiple regression."
## point_estimates SE_estimates pvalues
## MR_X1 0.062628442 0.05320910 0.239466722
## MR_X2 0.017644859 0.05236999 0.736242990
## MR_X3 0.060999591 0.05579270 0.274515196
## MR_X4 0.104372791 0.05231128 0.046290384
## MR_X5 -0.031949076 0.05466858 0.559074221
## MR_X6 0.021189455 0.05439856 0.696973076
## MR_X7 0.003105817 0.05196940 0.952356843
## MR_X8 0.072015147 0.05347421 0.178374662
## MR_X9 0.163694259 0.04996759 0.001089184
## MR_X10 -0.003600446 0.05269134 0.945535953
## [1] "Results based on traditional regression of residuals."
## point_estimates SE_estimates pvalues
## RR_X1 0.061436896 0.05264836 0.243517852
## RR_X2 0.017640202 0.05231059 0.736021598
## RR_X3 0.060380185 0.05545340 0.276484888
## RR_X4 0.104154953 0.05220449 0.046299820
## RR_X5 -0.031897823 0.05456996 0.558995364
## RR_X6 0.021183820 0.05433680 0.696722257
## RR_X7 0.003103701 0.05189961 0.952325297
## RR_X8 0.071981218 0.05340804 0.178041484
## RR_X9 0.162710874 0.04976895 0.001114884
## RR_X10 -0.003597357 0.05261593 0.945504713
## [1] "Results based on structural equation modeling."
## point_estimates SE_estimates pvalues
## SEM_X1 0.063229060 0.05309423 0.2336991512
## SEM_X2 0.017216199 0.05226550 0.7418541850
## SEM_X3 0.063028121 0.05546768 0.2558295692
## SEM_X4 0.103417368 0.05218873 0.0475236488
## SEM_X5 -0.030916939 0.05453222 0.5707494289
## SEM_X6 0.020936983 0.05429461 0.6997795567
## SEM_X7 0.002543383 0.05186050 0.9608852001
## SEM_X8 0.072460179 0.05336572 0.1745259155
## SEM_X9 0.164185505 0.04985986 0.0009914583
## SEM_X10 -0.003438170 0.05259176 0.9478756754Bollen KA (1989). Structural equations with latent variables. New York: John Wiley & Sons.
Konigorski S, Wang Y, Cigsar C, Yilmaz YE (2018). Estimating and testing direct genetic effects in directed acyclic graphs using estimating equations. Genetic Epidemiology, 42: 174-186.
Rosseel Y (2012). lavaan: an R package for structural equation modeling. Journal of Statistical Software, 48(2), 1–36.
Vansteelandt S, Goetgeluk S, Lutz S, et al. (2009). On the adjustment for covariates in genetic association analysis: a novel, simple principle to infer direct causal effects. Genetic Epidemiology, 33, 394-405.