This article is a brief illustration of how
to use indirect_effect()
from the package
manymome
(Cheung & Cheung, 2023)
to estimate the indirect effects
when the model parameters are estimated by
ordinary least squares (OLS) multiple regression
using lm()
.
This is the sample dataset used for illustration:
library(manymome) dat <- data_med_complicated print(round(head(dat), 2)) #> x1 x2 m11 m12 m2 y1 y2 c1 c2 #> 1 10.16 4.00 5.84 6.78 7.11 6.06 10.38 5.01 8.20 #> 2 10.89 4.79 4.95 5.81 7.92 5.15 8.52 3.92 9.90 #> 3 9.99 5.79 4.95 4.47 8.21 3.25 7.81 5.91 11.36 #> 4 12.36 4.80 5.56 6.21 8.88 6.27 9.41 6.06 10.49 #> 5 10.85 6.39 6.19 5.39 8.76 5.36 9.84 4.28 9.81 #> 6 10.28 4.58 4.88 4.28 8.50 4.57 10.42 4.71 11.38
This dataset has 9 variables: 2
predictors (x1
and x2
),
three mediators (m11
, m12
, and m2
),
two outcome variables (y1
and y2
),
and two control variables (c1
and c2
).
Suppose this is the model to be fitted:
Despite the apparent complexity, the path parameters can be estimated by five multiple regression models:
lm_m11 <- lm(m11 ~ x1 + x2 + c1 + c2, dat) lm_m12 <- lm(m12 ~ m11 + x1 + x2 + c1 + c2, dat) lm_m2 <- lm(m2 ~ x1 + x2 + c1 + c2, dat) lm_y1 <- lm(y1 ~ m12 + m2 + m11 + x1 + x2 + c1 + c2, dat) lm_y2 <- lm(y2 ~ m12 + m2 + m11 + x1 + x2 + c1 + c2, dat)
These are the regression coefficient estimates of the paths (those of control variables omitted):
#> m11 m12 m2 y1 y2 #> x1 0.352 -0.212 0.022 -0.078 0.115 #> x2 -0.045 -0.072 0.289 0.003 0.062 #> m11 0.454 0.147 0.024 #> m12 0.234 0.135 #> m2 -0.433 -0.436
Although not mandatory, it is recommended to combine these
five models into one object (a system of regression models)
using lm2list()
:
fit_lm <- lm2list(lm_m11, lm_m12, lm_m2, lm_y1, lm_y2) fit_lm #> #> The model(s): #> m11 ~ x1 + x2 + c1 + c2 #> m12 ~ m11 + x1 + x2 + c1 + c2 #> m2 ~ x1 + x2 + c1 + c2 #> y1 ~ m12 + m2 + m11 + x1 + x2 + c1 + c2 #> y2 ~ m12 + m2 + m11 + x1 + x2 + c1 + c2
Simply use the lm()
outputs as arguments. Order does not
matter. To ensure that
the regression outputs can be validly combined,
lm2list()
will also check:
a. whether the same sample is used in all regression analysis (not just same sample size, but the same set of cases), and
b. whether the models are "connected", to ensure that the regression outputs can be validly combined.
To form nonparametric bootstrap confidence interval for
indirect effects to be computed, do_boot()
can be used
to generate bootstrap estimates for all regression
coefficients first. These estimates can be reused for
any indirect effects to be estimated.
boot_out_lm <- do_boot(fit_lm, R = 100, seed = 54532, ncores = 1)
Please see vignette("do_boot")
or
the help page of do_boot()
on how
to use this function. In real research,
R
, the number of bootstrap samples,
should be set to 2000 or even 5000.
The argument ncores
can usually be omitted
unless users want to manually control
the number of CPU cores used in
parallel processing.
We can now use indirect_effect()
to
estimate the indirect effect and form its
bootstrap confidence interval for any path
in the model. By reusing the generated bootstrap
estimates, there is no need to repeat the
resampling.
Suppose we want to estimate the indirect
effect from x1
to y1
through m11
and m12
:
(Refer to vignette("manymome")
and the help page
of indirect_effect()
on the arguments.)
out_x1m11m12y1 <- indirect_effect(x = "x1", y = "y1", m = c("m11", "m12"), fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm) out_x1m11m12y1 #> #> == Indirect Effect == #> #> Path: x1 -> m11 -> m12 -> y1 #> Indirect Effect: 0.037 #> 95.0% Bootstrap CI: [0.003 to 0.077] #> #> Computation Formula: #> (b.m11~x1)*(b.m12~m11)*(b.y1~m12) #> #> Computation: #> (0.35204)*(0.45408)*(0.23402) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m11~x1 0.352 #> m12~m11 0.454 #> y1~m12 0.234
The indirect effect is 0.037, with 95% confidence interval [0.003, 0.077].
Similarly, we can estimate the indirect
effect from x2
to y2
through m2
:
out_x2m2y2 <- indirect_effect(x = "x2", y = "y2", m = "m2", fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm) out_x2m2y2 #> #> == Indirect Effect == #> #> Path: x2 -> m2 -> y2 #> Indirect Effect: -0.126 #> 95.0% Bootstrap CI: [-0.233 to -0.043] #> #> Computation Formula: #> (b.m2~x2)*(b.y2~m2) #> #> Computation: #> (0.28901)*(-0.43598) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m2~x2 0.289 #> y2~m2 -0.436
The indirect effect is -0.126, with 95% confidence interval [-0.233, -0.043].
Note that any indirect path in the model can be
estimated this way. Suppose, after doing the regression analysis,
we want to estimate the indirect effect from x2
to m12
through m11
, we just call indirect_effect()
:
out_x2m11m12 <- indirect_effect(x = "x2", y = "m12", m = "m11", fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm) out_x2m11m12 #> #> == Indirect Effect == #> #> Path: x2 -> m11 -> m12 #> Indirect Effect: -0.020 #> 95.0% Bootstrap CI: [-0.139 to 0.110] #> #> Computation Formula: #> (b.m11~x2)*(b.m12~m11) #> #> Computation: #> (-0.04471)*(0.45408) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m11~x2 -0.0447 #> m12~m11 0.4541
The indirect effect is -0.020, with 95% confidence interval [-0.139, 0.110].
There is no limit on the path to be estimated, as long
as all required path coefficients are in the model.
indirect_effect()
will also check whether a path is valid.
Therefore, estimating the effect from x1
to m2
through
m11
will result in an error because this path does not
exist in the model defined by the regression outputs.
The standardized indirect
effect from x1
to y1
through m11
and m12
can be estimated by setting
standardized_x
and standardized_y
to `TRUE:
std_x1m11m12y1 <- indirect_effect(x = "x1", y = "y1", m = c("m11", "m12"), fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm, standardized_x = TRUE, standardized_y = TRUE) std_x1m11m12y1 #> #> == Indirect Effect (Both 'x1' and 'y1' Standardized) == #> #> Path: x1 -> m11 -> m12 -> y1 #> Indirect Effect: 0.039 #> 95.0% Bootstrap CI: [0.004 to 0.085] #> #> Computation Formula: #> (b.m11~x1)*(b.m12~m11)*(b.y1~m12)*sd_x1/sd_y1 #> #> Computation: #> (0.35204)*(0.45408)*(0.23402)*(1.11605)/(1.06579) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m11~x1 0.352 #> m12~m11 0.454 #> y1~m12 0.234 #> #> NOTE: #> - The effects of the component paths are from the model, not #> standardized.
The standardized indirect effect is 0.039, with 95% confidence interval [0.004,0.085].
Similarly, we can estimate the standardized indirect
effect from x1
to y1
through m2
:
std_x1m2y1 <- indirect_effect(x = "x1", y = "y1", m = "m2", fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm, standardized_x = TRUE, standardized_y = TRUE) std_x1m2y1 #> #> == Indirect Effect (Both 'x1' and 'y1' Standardized) == #> #> Path: x1 -> m2 -> y1 #> Indirect Effect: -0.010 #> 95.0% Bootstrap CI: [-0.069 to 0.067] #> #> Computation Formula: #> (b.m2~x1)*(b.y1~m2)*sd_x1/sd_y1 #> #> Computation: #> (0.02233)*(-0.43300)*(1.11605)/(1.06579) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m2~x1 0.0223 #> y1~m2 -0.4330 #> #> NOTE: #> - The effects of the component paths are from the model, not #> standardized.
The standardized indirect effect is -0.010, with 95% confidence interval [-0.069, 0.067].
Note that the results of indirect_effect()
can be added using +
.
For example, to find the total indirect effect of
x1
on y1
, we need to compute
the indirect effects along the following paths:
a. x1
to m11
to m12
to y1
b. x1
to m11
to y1
c. x1
to m12
to y1
d. x1
to m2
to y1
The indirect effects along Path a has already been computed. We compute the indirect effects along Paths b, c, and d below:
out_x1m11y1 <- indirect_effect(x = "x1", y = "y1", m = "m11", fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm) out_x1m11y1 #> #> == Indirect Effect == #> #> Path: x1 -> m11 -> y1 #> Indirect Effect: 0.052 #> 95.0% Bootstrap CI: [-0.036 to 0.103] #> #> Computation Formula: #> (b.m11~x1)*(b.y1~m11) #> #> Computation: #> (0.35204)*(0.14694) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m11~x1 0.352 #> y1~m11 0.147
out_x1m12y1 <- indirect_effect(x = "x1", y = "y1", m = "m12", fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm) out_x1m12y1 #> #> == Indirect Effect == #> #> Path: x1 -> m12 -> y1 #> Indirect Effect: -0.050 #> 95.0% Bootstrap CI: [-0.125 to -0.000] #> #> Computation Formula: #> (b.m12~x1)*(b.y1~m12) #> #> Computation: #> (-0.21182)*(0.23402) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m12~x1 -0.212 #> y1~m12 0.234
out_x1m2y1 <- indirect_effect(x = "x1", y = "y1", m = "m2", fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm) out_x1m2y1 #> #> == Indirect Effect == #> #> Path: x1 -> m2 -> y1 #> Indirect Effect: -0.010 #> 95.0% Bootstrap CI: [-0.077 to 0.064] #> #> Computation Formula: #> (b.m2~x1)*(b.y1~m2) #> #> Computation: #> (0.02233)*(-0.43300) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m2~x1 0.0223 #> y1~m2 -0.4330
We can now compute the total indirect effect:
out_x1y1_total <- out_x1m11m12y1 + out_x1m11y1 + out_x1m12y1 + out_x1m2y1 out_x1y1_total #> #> == Indirect Effect == #> #> Path: x1 -> m11 -> m12 -> y1 #> Path: x1 -> m11 -> y1 #> Path: x1 -> m12 -> y1 #> Path: x1 -> m2 -> y1 #> Function of Effects: 0.030 #> 95.0% Bootstrap CI: [-0.088 to 0.132] #> #> Computation of the Function of Effects: #> (((x1->m11->m12->y1) #> +(x1->m11->y1)) #> +(x1->m12->y1)) #> +(x1->m2->y1) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples.
The total effect of f1
on f4
is
0.030,
with 95% confidence interval
[-0.088, 0.132].
See help("math_indirect")
for further details on addition
for indirect_effect()
outputs.
Subtraction can also be conducted using -
. For
example, we can compute the difference between
the indirect effect of x1
on y1
through m11
and m12
and the indirect effect of x1
on y1
through m2
:
out_x1_diff <- out_x1m11m12y1 - out_x1m2y1 out_x1_diff #> #> == Indirect Effect == #> #> Path: x1 -> m11 -> m12 -> y1 #> Path: x1 -> m2 -> y1 #> Function of Effects: 0.047 #> 95.0% Bootstrap CI: [-0.038 to 0.118] #> #> Computation of the Function of Effects: #> (x1->m11->m12->y1) #> -(x1->m2->y1) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples.
The difference in effects is 0.047, with 95% confidence interval [-0.038, 0.118].
See help("math_indirect")
for further details on subtraction
for indirect_effect()
outputs.
If there are several indirect paths in
a model, the function all_indirect_paths()
can
be used to automatically identify all
indirect paths (a path with at least
one mediator) in a model:
all_paths <- all_indirect_paths(fit = fit_lm) all_paths #> Call: #> all_indirect_paths(fit = fit_lm) #> Path(s): #> path #> 1 m11 -> m12 -> y1 #> 2 m11 -> m12 -> y2 #> 3 x1 -> m11 -> m12 #> 4 x1 -> m11 -> m12 -> y1 #> 5 x1 -> m11 -> y1 #> 6 x1 -> m12 -> y1 #> 7 x1 -> m2 -> y1 #> 8 x1 -> m11 -> m12 -> y2 #> 9 x1 -> m11 -> y2 #> 10 x1 -> m12 -> y2 #> 11 x1 -> m2 -> y2 #> 12 x2 -> m11 -> m12 #> 13 x2 -> m11 -> m12 -> y1 #> 14 x2 -> m11 -> y1 #> 15 x2 -> m12 -> y1 #> 16 x2 -> m2 -> y1 #> 17 x2 -> m11 -> m12 -> y2 #> 18 x2 -> m11 -> y2 #> 19 x2 -> m12 -> y2 #> 20 x2 -> m2 -> y2 #> 21 c1 -> m11 -> m12 #> 22 c1 -> m11 -> m12 -> y1 #> 23 c1 -> m11 -> y1 #> 24 c1 -> m12 -> y1 #> 25 c1 -> m2 -> y1 #> 26 c1 -> m11 -> m12 -> y2 #> 27 c1 -> m11 -> y2 #> 28 c1 -> m12 -> y2 #> 29 c1 -> m2 -> y2 #> 30 c2 -> m11 -> m12 #> 31 c2 -> m11 -> m12 -> y1 #> 32 c2 -> m11 -> y1 #> 33 c2 -> m12 -> y1 #> 34 c2 -> m2 -> y1 #> 35 c2 -> m11 -> m12 -> y2 #> 36 c2 -> m11 -> y2 #> 37 c2 -> m12 -> y2 #> 38 c2 -> m2 -> y2
The initial list is very long because
control variables (c1
and c2
) are
included in the search. Moreover, paths
that start from a mediator or end at
a mediator are also included.
Users can customize the search:
all_paths <- all_indirect_paths(fit = fit_lm, x = c("x1", "x2"), y = c("y1", "y2"), exclude = c("c1", "c2")) all_paths #> Call: #> all_indirect_paths(fit = fit_lm, exclude = c("c1", "c2"), x = c("x1", #> "x2"), y = c("y1", "y2")) #> Path(s): #> path #> 1 x1 -> m11 -> m12 -> y1 #> 2 x1 -> m11 -> y1 #> 3 x1 -> m12 -> y1 #> 4 x1 -> m2 -> y1 #> 5 x1 -> m11 -> m12 -> y2 #> 6 x1 -> m11 -> y2 #> 7 x1 -> m12 -> y2 #> 8 x1 -> m2 -> y2 #> 9 x2 -> m11 -> m12 -> y1 #> 10 x2 -> m11 -> y1 #> 11 x2 -> m12 -> y1 #> 12 x2 -> m2 -> y1 #> 13 x2 -> m11 -> m12 -> y2 #> 14 x2 -> m11 -> y2 #> 15 x2 -> m12 -> y2 #> 16 x2 -> m2 -> y2
x
is a vector of names. Only paths
start from these variables will be
included.
y
is a vector of names. Only paths
end at these variables will be
included.
exclude
is a vector of names. Paths
that involve these variables will be
excluded.
The output is a all_paths
-class object.
It can be used in many_indirect_effects()
out_all <- many_indirect_effects(paths = all_paths, fit = fit_lm, boot_ci = TRUE, boot_out = boot_out_lm)
The first argument, paths
, is the output
of all_indirect_paths()
. The other
arguments
will be passed to indirect_effect()
.
The output is an indirect_list
-class
object, which is a list of the outputs
of indirect_effects()
. If printed,
a summary of the indirect effects will
be printed:
out_all #> #> == Indirect Effect(s) == #> ind CI.lo CI.hi Sig #> x1 -> m11 -> m12 -> y1 0.037 0.003 0.077 Sig #> x1 -> m11 -> y1 0.052 -0.036 0.103 #> x1 -> m12 -> y1 -0.050 -0.125 -0.000 Sig #> x1 -> m2 -> y1 -0.010 -0.077 0.064 #> x1 -> m11 -> m12 -> y2 0.022 -0.015 0.065 #> x1 -> m11 -> y2 0.009 -0.086 0.093 #> x1 -> m12 -> y2 -0.029 -0.089 0.014 #> x1 -> m2 -> y2 -0.010 -0.074 0.068 #> x2 -> m11 -> m12 -> y1 -0.005 -0.044 0.017 #> x2 -> m11 -> y1 -0.007 -0.055 0.034 #> x2 -> m12 -> y1 -0.017 -0.054 0.049 #> x2 -> m2 -> y1 -0.125 -0.258 -0.031 Sig #> x2 -> m11 -> m12 -> y2 -0.003 -0.021 0.022 #> x2 -> m11 -> y2 -0.001 -0.021 0.052 #> x2 -> m12 -> y2 -0.010 -0.056 0.045 #> x2 -> m2 -> y2 -0.126 -0.233 -0.043 Sig #> #> - [CI.lo to CI.hi] are 95.0% percentile confidence intervals by #> nonparametric bootstrapping with 100 samples. #> - The 'ind' column shows the indirect effects. #>
The output of many_indirect_effects()
is a named list, names being the path
name as appeared in the output.
Individual indirect
effects can be extracted using either
the indices or the path names
An example using index:
out1 <- out_all[[1]] out1 #> #> == Indirect Effect == #> #> Path: x1 -> m11 -> m12 -> y1 #> Indirect Effect: 0.037 #> 95.0% Bootstrap CI: [0.003 to 0.077] #> #> Computation Formula: #> (b.m11~x1)*(b.m12~m11)*(b.y1~m12) #> #> Computation: #> (0.35204)*(0.45408)*(0.23402) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m11~x1 0.352 #> m12~m11 0.454 #> y1~m12 0.234
An example using path name (though not recommended because the name is usually long):
out2 <- out_all[["x2 -> m2 -> y2"]] out2 #> #> == Indirect Effect == #> #> Path: x2 -> m2 -> y2 #> Indirect Effect: -0.126 #> 95.0% Bootstrap CI: [-0.233 to -0.043] #> #> Computation Formula: #> (b.m2~x2)*(b.y2~m2) #> #> Computation: #> (0.28901)*(-0.43598) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples. #> #> Coefficients of Component Paths: #> Path Coefficient #> m2~x2 0.289 #> y2~m2 -0.436
The extracted element can be used just
like the outputs of indirect_effect()
in previous section.
See the help page of all_indirect_paths()
and many_indirect_effects()
for other
arguments available.
The total indirect effect between two
variables in a list of paths can be
computed by total_indirect_effect()
.
total_x1_y1 <- total_indirect_effect(out_all, x = "x1", y = "y1") total_x1_y1 #> #> == Indirect Effect == #> #> Path: x1 -> m11 -> m12 -> y1 #> Path: x1 -> m11 -> y1 #> Path: x1 -> m12 -> y1 #> Path: x1 -> m2 -> y1 #> Function of Effects: 0.030 #> 95.0% Bootstrap CI: [-0.088 to 0.132] #> #> Computation of the Function of Effects: #> (((x1->m11->m12->y1) #> +(x1->m11->y1)) #> +(x1->m12->y1)) #> +(x1->m2->y1) #> #> #> Percentile confidence interval formed by nonparametric bootstrapping #> with 100 bootstrap samples.
The first argument is the output
of many_indirect_effects()
or a list
of indirect
-class object. x
is the
name of the variable that starts the paths.
y
is the name of the variable that
ends the paths.
For further information on do_boot()
and indirect_effect()
,
please refer to their help pages.
Cheung, S. F., & Cheung, S.-H. (2023). manymome: An R package for computing the indirect effects, conditional effects, and conditional indirect effects, standardized or unstandardized, and their bootstrap confidence intervals, in many (though not all) models. Behavior Research Methods. https://doi.org/10.3758/s13428-023-02224-z
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.