Examples for Third-Variable Effect Analysis

Qingzhao Yu & Lin Zhu

2023-08-29

Package installation

R package mma is used for general third-variable effect analysis (Yu and Li 2022). The third-variable effect (TVE) refers to the effect conveyed by intervening variables to an observed relationship between an exposure and a response variable (outcome). In this package, the exposure is also called the predictor, the intervening variables are called mediators or confounders (MC). The TVE include total effect, direct effect, and indirect effects from different third variables. The effects are defined in Yu et al. (2014). The algorithms on how to make inferences on TVEs and how to explain the effects are described in Yu et al. (2014).

To use the R package mma, we first install the package in R (install.packages("mma")) and load it.

library(mma)
#> Loading required package: gbm
#> Loaded gbm 2.1.8.1
#> Loading required package: splines
#> Loading required package: survival
#> Loading required package: car
#> Loading required package: carData
#> Loading required package: gplots
#> 
#> Attaching package: 'gplots'
#> The following object is masked from 'package:stats':
#> 
#>     lowess
#source('C:/Users/qyu/Desktop/Research/mma/version current/R/mma.r')

Data Organization and Identify Potential Mediators/Confounders

We use weight_behavior data as an example to show how to do a preliminary data analysis to identify potential MCs and covariates. The function data.org is used for this purpose. The list of arguments for data.org are defined in Yu and Li (2017). In the following example, the predictor, “pred”, is gender, and the outcome, “y”, is a binary indicator of overweight (1) or not (0). The reference group of the exposure variable is the male as defined by “predref=”M”“. There are 12 variables in”x”, which includes all potential MCs and covariates. As an argument in the function data.org, “mediator=5:12” indicates that the variables from column 5 to column 12 in x should be tested for potential MCs. Other columns in “x” are to be treated as covariates. The “mediator” argument can also be a character vector with names of the potential MCs in x. Two tests are needed to identify a potential MC: first, the variable is significantly related with the predictor adjusting for other covariates. The other covariates here is defined by “cova”. If cova is a data frame by itself, all covariates in cova are for all potential MCs. If cova is a list, the first item will be the covariate data frame and the second item list the names or column numbers of MCs in “x” that should adjust for the covariates. The significance level is set by “alpha2”, whose default value is 0.1. Second, the variable has to be significantly related with the outcome adjusting (testtype=1, by default) or not adjusting (testtype=2) for the predictor and other variables. The significance level is set by “alpha”. In the following example, p-value 1 shows the results for the second test. and p-value 2 are the results for the first test. A variable that pass the second test but not the first test is considered as covariates. Variables do not pass the second test are not adopted in further analysis. Variables in “x” but not in “mediator” are forced in further analysis as covariates.

“jointm” is a list with the first item identify the number of groups of MCs whose joint effects are of interesting, where each group is defined in an item of “jointm”. In the above example, the joint effect of Variables 5, 7, 9 is of interesting. These variables are forced into analysis as potential MCs. A variable can be shown in “jointm” multiple times. The joint effect of each group of variables will be reported. A trick to force in variables as potential MCs is to assign them in a separate group of “jointm”.

One binary/continuous predictor and one binary/continuous outcome

data("weight_behavior")
#binary predictor x
#binary y
 x=weight_behavior[,c(2,4:14)]
 pred=weight_behavior[,"sex"]
 y=weight_behavior[,"overweigh"]
 data.b.b.1<-data.org(x,y,mediator=5:12,jointm=list(n=1,j1=7:9),
                        pred=pred,predref="M", alpha=0.4,alpha2=0.4)
 summary(data.b.b.1)
#> Identified as mediators: 
#> [1] "tvhours"   "cmpthours" "cellhours" "exercises" "sweat"     "sports"   
#> Selected as covariates: 
#> [1] "age"       "race"      "numpeople" "car"      
#> Tests: 
#>             P-Value 1.y P-Value 2.pred
#> age               0.836             NA
#> race              0.850             NA
#> numpeople         0.561             NA
#> car               0.052             NA
#> gotosch           0.710             NA
#> snack             0.724             NA
#> tvhours -         0.830          0.484
#> cmpthours -       0.826          0.760
#> cellhours -       0.067          0.688
#> sports *          0.000          0.003
#> exercises *       0.176          0.203
#> sweat *           0.181          0.046
#> pred              0.178             NA
#> ----
#>  *:mediator,-:joint mediator
#>  P-Value 1:Type-3 tests in the full model (data.org) or estimated coefficients (data.org.big) 
#>  when testtype=1, univariate relationship test with the outcome when testtype=2
#>  P-Value 2:Tests of relationship with the Predictor

The results from “data.org” function is summarized as above.

A potential MC defined in “mediator” argument is identified as categorical or binary if it has only two unique values or if the variable is a vector factors or characters. By default, the reference group is the first level of the variable. If the reference group needs to be changed, the categorical variable should be defined in binmed or catmed, and then the corresponding reference group in binref and catref. Continuous MCs can be defined by contmed. The following example did the same analysis as above, only that it defines that potential continuous mediators are columns 7,8,9,11, and 12 of x, binary mediators are columns 6 and 10, and categorical mediator is column 5 of x with 1 to be the reference group for all categorical and binary mediators.

 data.b.b.2<-data.org(x,y,pred=pred,contmed=c(7:9,11:12),binmed=c(6,10), jointm=list(n=1,j1=7:9),
   binref=c(1,1),catmed=5,catref=1,predref="M",alpha=0.4,alpha2=0.4) 
 summary(data.b.b.2)

Survival outcome

The functions in mma can deal with binary, categorical, continuous, or time-to-event outcomes. If the outcome is time-to-event, it should be defined by the “Surv” function in the survival package. An example of application on breast cancer survival can be found in Yu et al. (2019).

 cgd1<-cgd0
 status<-ifelse(is.na(cgd1$etime1),0,1)
 y=Surv(cgd1$futime,status) 
 x=cgd1[,c(5:7,9:12)]
 pred=cgd1[,4]
 data.b.surv<-data.org(x,y,pred=pred,mediator=1:ncol(x),   
                    alpha=0.4,alpha2=0.4)
 summary(data.b.surv)

Multivariate predictors

In addition, the package can handle multivariate and multicategorical predictors. If the predictor is multicategorical of k levels, “data.org” first transform the predictor to k-1 binary predictors. If a variable significantly relates with any of the k-1 predictors, the variable passes the first test described above. P-value 2 is shown for each predictor. In the following example, the predictor is the race with six levels.

 #multivariate predictor
 x=weight_behavior[,c(2:3,5:14)]
 pred=weight_behavior[,4]
 y=weight_behavior[,15]
 data.mb.b <-data.org(x,y,mediator=5:12,jointm=list(n=1,j1=c(5,7,9)),
                      pred=pred,predref="OTHER", alpha=0.4,alpha2=0.4)
 summary(data.mb.b)
#> Identified as mediators: 
#> [1] "tvhours"   "cellhours" "sweat"     "sports"    "gotosch"  
#> Selected as covariates: 
#> [1] "age"       "sex"       "numpeople" "car"       "exercises"
#> Tests: 
#>                P-Value 1.y P-Value 2.pred P-Value 2.pred.AFRICAN
#> age                  0.836             NA                     NA
#> sex                  0.178             NA                     NA
#> numpeople            0.561             NA                     NA
#> car                  0.052             NA                     NA
#> gotosch -            0.710          0.112                  0.796
#> snack                0.724             NA                     NA
#> tvhours -            0.830          0.748                  0.535
#> cmpthours            0.826             NA                     NA
#> cellhours -          0.067          0.880                  0.994
#> sports *             0.000          0.587                  0.163
#> exercises            0.176          0.629                  0.731
#> sweat *              0.181          0.647                  0.174
#> pred                 0.543             NA                     NA
#> pred.AFRICAN         0.663             NA                     NA
#> pred.CAUCASIAN       0.242             NA                     NA
#> pred.INDIAN          0.890             NA                     NA
#> pred.MIXED           0.782             NA                     NA
#>                P-Value 2.pred.CAUCASIAN P-Value 2.pred.INDIAN
#> age                                  NA                    NA
#> sex                                  NA                    NA
#> numpeople                            NA                    NA
#> car                                  NA                    NA
#> gotosch -                         0.996                 0.097
#> snack                                NA                    NA
#> tvhours -                         0.334                 0.916
#> cmpthours                            NA                    NA
#> cellhours -                       0.707                 0.555
#> sports *                          0.816                 0.453
#> exercises                         0.723                 0.912
#> sweat *                           0.214                 0.374
#> pred                                 NA                    NA
#> pred.AFRICAN                         NA                    NA
#> pred.CAUCASIAN                       NA                    NA
#> pred.INDIAN                          NA                    NA
#> pred.MIXED                           NA                    NA
#>                P-Value 2.pred.MIXED
#> age                              NA
#> sex                              NA
#> numpeople                        NA
#> car                              NA
#> gotosch -                     0.033
#> snack                            NA
#> tvhours -                     0.092
#> cmpthours                        NA
#> cellhours -                   0.383
#> sports *                      0.091
#> exercises                     0.881
#> sweat *                       0.203
#> pred                             NA
#> pred.AFRICAN                     NA
#> pred.CAUCASIAN                   NA
#> pred.INDIAN                      NA
#> pred.MIXED                       NA
#> ----
#>  *:mediator,-:joint mediator
#>  P-Value 1:Type-3 tests in the full model (data.org) or estimated coefficients (data.org.big) 
#>  when testtype=1, univariate relationship test with the outcome when testtype=2
#>  P-Value 2:Tests of relationship with the Predictor

Multivariate outcomes

Similarly, the package can deal with multivariate outcomes. The following code deals with multivariate predictors and multivariate responses. If the variable is significantly related with any one of the outcomes, it passes the second test described above. The results from “data.org” are summarized for each combination of the exposure-outcome relationship.

 #multivariate responses
 x=weight_behavior[,c(2:3,5:14)]
 pred=weight_behavior[,4]
 y=weight_behavior[,c(1,15)]
 data.mb.mb<-data.org(x,y,mediator=5:12,jointm=list(n=1,j1=c(5,7,9)),
                      pred=pred,predref="OTHER", alpha=0.4,alpha2=0.4)
 summary(data.mb.mb)
#> Identified as mediators: 
#> [1] "tvhours"   "cellhours" "sweat"     "sports"    "gotosch"  
#> Selected as covariates: 
#> [1] "age"       "sex"       "numpeople" "car"       "exercises"
#> Tests: 
#>                P-Value 1.bmi P-Value 1.overweigh P-Value 2.pred
#> age                    0.458               0.836             NA
#> sex                    0.003               0.178             NA
#> numpeople              0.282               0.561             NA
#> car                    0.059               0.052             NA
#> gotosch -              0.527               0.710          0.112
#> snack                  0.830               0.724             NA
#> tvhours -              0.505               0.830          0.748
#> cmpthours              0.676               0.826             NA
#> cellhours -            0.084               0.067          0.880
#> sports *               0.001               0.000          0.587
#> exercises              0.089               0.176          0.629
#> sweat *                0.078               0.181          0.647
#> pred                   0.344               0.543             NA
#> pred.AFRICAN           0.900               0.663             NA
#> pred.CAUCASIAN         0.032               0.242             NA
#> pred.INDIAN            0.446               0.890             NA
#> pred.MIXED             0.833               0.782             NA
#>                P-Value 2.pred.AFRICAN P-Value 2.pred.CAUCASIAN
#> age                                NA                       NA
#> sex                                NA                       NA
#> numpeople                          NA                       NA
#> car                                NA                       NA
#> gotosch -                       0.796                    0.996
#> snack                              NA                       NA
#> tvhours -                       0.535                    0.334
#> cmpthours                          NA                       NA
#> cellhours -                     0.994                    0.707
#> sports *                        0.163                    0.816
#> exercises                       0.731                    0.723
#> sweat *                         0.174                    0.214
#> pred                               NA                       NA
#> pred.AFRICAN                       NA                       NA
#> pred.CAUCASIAN                     NA                       NA
#> pred.INDIAN                        NA                       NA
#> pred.MIXED                         NA                       NA
#>                P-Value 2.pred.INDIAN P-Value 2.pred.MIXED
#> age                               NA                   NA
#> sex                               NA                   NA
#> numpeople                         NA                   NA
#> car                               NA                   NA
#> gotosch -                      0.097                0.033
#> snack                             NA                   NA
#> tvhours -                      0.916                0.092
#> cmpthours                         NA                   NA
#> cellhours -                    0.555                0.383
#> sports *                       0.453                0.091
#> exercises                      0.912                0.881
#> sweat *                        0.374                0.203
#> pred                              NA                   NA
#> pred.AFRICAN                      NA                   NA
#> pred.CAUCASIAN                    NA                   NA
#> pred.INDIAN                       NA                   NA
#> pred.MIXED                        NA                   NA
#> ----
#>  *:mediator,-:joint mediator
#>  P-Value 1:Type-3 tests in the full model (data.org) or estimated coefficients (data.org.big) 
#>  when testtype=1, univariate relationship test with the outcome when testtype=2
#>  P-Value 2:Tests of relationship with the Predictor

##Third-Variable Effect Analysis Next, we use med function to estimate the TVE using the results from data.org function. The default approach is using a generalized linear model to fit the final full model. If “nonlinear=TRUE”“, multiple Additive Regression Trees (MART) will be used to fit the relationship between outcome(s) and other variables, and smoothing splines used to fit the relationship between the potential MC and exposure variables(s).

Binary predictor and binary outcome

 med.b.b.2<-med(data=data.b.b.2,n=2,nonlinear=TRUE)

Using the linear method and to show the results:

 med.b.b.1<-med(data=data.b.b.2,n=2)
 med.b.b.1
#> 
#> 
#> For the predictor pred :
#>  The estimated total effect:[1] 0.6156
#> 
#>  The estimated indirect effect:
#>       y1.all   y1.tvhours y1.cmpthours y1.cellhours y1.exercises     y1.sweat 
#>       0.2158      -0.0007       0.0207       0.0177      -0.0177       0.0080 
#>    y1.sports        y1.j1 
#>       0.1633       0.0192

We can plot the marginal effect of the MC on the outcome, and the marginal effect of the predictor on the MC. If the predictor is binary, plot.med draws a histogram or boxplot of the marginal density of the MC at each different value of the predictor.

plot(med.b.b.1, data.b.b.2,vari="exercises",xlim=c(0,50))
plot(med.b.b.2, data.b.b.2,vari="sports")

Survival outcome

For survival outcome, the default option is to fit the final full model using Cox proportional hazard model. Note that when use type=“lp”, the linear predictor is the type of predicted value as in predict.coxph, the results is similar to those from MART. Readers are referred to Yu et al. (2019) for more details.

 med.b.surv.1<-med(data=data.b.surv,n=2,type="lp") 
  #close to mart results when use type="lp"
 med.b.surv.2<-med(data=data.b.surv,n=2,nonlinear=TRUE)  
 #results in the linear part unit

Statistical Inference on TVE Analysis

boot.med function is used to make inferences on the TVE. Results return from data.org function is input as the first argument.

Binary predictor and binary outcome

 bootmed.b.b.1<-boot.med(data=data.b.b.2,n=2,n2=50)
#> [1] 1
#> [1] 2
#> [1] 3
#> [1] 4
#> [1] 5
#> [1] 6
#> [1] 7
#> [1] 8
#> [1] 9
#> [1] 10
#> [1] 11
#> [1] 12
#> [1] 13
#> [1] 14
#> [1] 15
#> [1] 16
#> [1] 17
#> [1] 18
#> [1] 19
#> [1] 20
#> [1] 21
#> [1] 22
#> [1] 23
#> [1] 24
#> [1] 25
#> [1] 26
#> [1] 27
#> [1] 28
#> [1] 29
#> [1] 30
#> [1] 31
#> [1] 32
#> [1] 33
#> [1] 34
#> [1] 35
#> [1] 36
#> [1] 37
#> [1] 38
#> [1] 39
#> [1] 40
#> [1] 41
#> [1] 42
#> [1] 43
#> [1] 44
#> [1] 45
#> [1] 46
#> [1] 47
#> [1] 48
#> [1] 49
#> [1] 50
 summary(bootmed.b.b.1)
#> MMA Analysis: Estimated Mediation Effects Using GLM
#> For Predictor/Moderator at pred 
#> $total.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.517     0.660     0.421     1.485    -0.165     1.480    -0.094     1.165 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>    -0.212     1.675    -0.235     1.439    -0.111     0.117     0.120 
#> 
#> $direct.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.369     0.403     0.419     1.224    -0.418     1.120    -0.341     1.120 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>    -0.341     1.538    -0.372     1.154    -0.283     0.336     0.360 
#> 
#> $indirect.effect
#>           y1.all y1.tvhours y1.cmpthours y1.cellhours y1.exercises y1.sweat
#> est        0.148     -0.007        0.051       -0.007       -0.013    0.000
#> mean       0.257     -0.013        0.036        0.027       -0.016    0.011
#> sd         0.149      0.036        0.059        0.053        0.028    0.025
#> upbd       0.549      0.057        0.152        0.130        0.038    0.061
#> lwbd      -0.035     -0.083       -0.080       -0.076       -0.071   -0.039
#> upbd_q     0.528      0.052        0.182        0.139        0.020    0.065
#> lwbd_q    -0.074     -0.085       -0.041       -0.051       -0.086   -0.033
#> upbd_bcbi  0.371      0.078        0.192        0.045        0.010    0.042
#> lwbd_bcbi -0.126     -0.062       -0.036       -0.069       -0.089   -0.039
#> upbd_b     0.529      0.086        0.194        0.145        0.021    0.068
#> lwbd_b    -0.127     -0.087       -0.050       -0.070       -0.089   -0.040
#> upbd_win   0.541      0.052        0.187        0.161        0.013    0.067
#> lwbd_win  -0.074     -0.087       -0.040       -0.036       -0.089   -0.026
#> p_norm     0.085      0.719        0.543        0.610        0.554    0.661
#> p_quan     0.080      0.680        0.520        0.600        0.600    0.800
#>           y1.sports  y1.j1
#> est           0.124  0.011
#> mean          0.169  0.072
#> sd            0.084  0.095
#> upbd          0.334  0.257
#> lwbd          0.003 -0.114
#> upbd_q        0.367  0.246
#> lwbd_q        0.002 -0.094
#> upbd_bcbi     0.172  0.138
#> lwbd_bcbi    -0.031 -0.111
#> upbd_b        0.387  0.308
#> lwbd_b       -0.032 -0.111
#> upbd_win      0.369  0.250
#> lwbd_win      0.022 -0.093
#> p_norm        0.045  0.449
#> p_quan        0.080  0.440

Similarly, the default approach is using a generalized linear model to fit the final full model. If the argument nonlinear=TRUE, multiple Additive Regression Trees (MART) will be used to fit the model in estimating the outcome. Applications for nonlinear Mediation Analysis are described in Yu et al. (2017, 2018).

 bootmed.b.b.2<-boot.med(data=data.b.b.2,n=2,n2=40,nu=0.05,nonlinear=TRUE)
 summary(bootmed.b.b.2)
#> MMA Analysis: Estimated Mediation Effects Using MART
#> For Predictor/Moderator at pred 
#> $total.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.164     0.194     0.084     0.359     0.030     0.287     0.123     0.287 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>     0.123     0.291     0.126     0.287     0.123     0.020     0.000 
#> 
#> $direct.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.070     0.021     0.074     0.167    -0.125     0.122    -0.020     0.132 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>    -0.018     0.132    -0.017     0.122    -0.018     0.776     0.500 
#> 
#> $indirect.effect
#>           y1.all y1.tvhours y1.cmpthours y1.cellhours y1.exercises y1.sweat
#> est        0.094      0.002        0.003        0.010       -0.003    0.006
#> mean       0.173      0.015        0.028        0.008       -0.016    0.000
#> sd         0.050      0.025        0.014        0.007        0.021    0.006
#> upbd       0.272      0.063        0.056        0.022        0.024    0.011
#> lwbd       0.075     -0.034        0.000       -0.006       -0.057   -0.011
#> upbd_q     0.241      0.035        0.041        0.017        0.005    0.006
#> lwbd_q     0.143     -0.014        0.011        0.002       -0.042   -0.005
#> upbd_bcbi  0.143      0.035        0.010        0.018        0.006    0.007
#> lwbd_bcbi  0.143     -0.015        0.010        0.003       -0.022   -0.005
#> upbd_b     0.248      0.036        0.038        0.018        0.006    0.007
#> lwbd_b     0.143      0.003        0.010        0.004       -0.044   -0.005
#> upbd_win   0.241      0.035        0.041        0.017        0.005    0.006
#> lwbd_win   0.143     -0.014        0.011        0.003       -0.042   -0.005
#> p_norm     0.001      0.556        0.050        0.276        0.430    0.992
#> p_quan     0.000      0.500        0.000        0.000        0.500    1.000
#>           y1.sports  y1.j1
#> est           0.078  0.015
#> mean          0.147  0.029
#> sd            0.055  0.029
#> upbd          0.255  0.086
#> lwbd          0.038 -0.028
#> upbd_q        0.194  0.061
#> lwbd_q        0.075 -0.003
#> upbd_bcbi     0.175  0.045
#> lwbd_bcbi     0.068 -0.004
#> upbd_b        0.196  0.063
#> lwbd_b        0.068  0.017
#> upbd_win      0.194  0.061
#> lwbd_win      0.068 -0.003
#> p_norm        0.008  0.321
#> p_quan        0.000  0.500

Plots of the fitted mma object from boot.med

plot.mma() plots the marginal effect of the selected variable on the outcome, and the marginal effect of the predictor on the selected variable. Interpretations for similar plots can be found in Yu et al. (2017).

 plot(bootmed.b.b.1,vari="exercises",xlim=c(0,50))

 plot(bootmed.b.b.1,vari="sports")

Combined function for multiple TVE analysis

The mma function is a combined function that automatically identify potential MCs, based on which to make statistical inference on the mediation effects.

 x=weight_behavior[,c(2,4:14)]
 pred=weight_behavior[,3]
 y=weight_behavior[,15]
 mma.b.b.glm<-mma(x,y,pred=pred,contmed=c(7:9,11:12),binmed=c(6,10),binref=c(1,1), catmed=5,catref=1,predref="M",alpha=0.4,alpha2=0.4,n=2,n2=2)
#> [1] 1
#> [1] 2
 summary(mma.b.b.glm)
#> MMA Analysis: Estimated Mediation Effects Using GLM
#> For Predictor/Moderator at pred 
#> $total.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.572     0.956     0.304     1.551     0.361     1.160     0.752     0.742 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>     0.742     1.171     1.171     1.160     0.752     0.002     0.000 
#> 
#> $direct.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.417     0.785     0.304     1.381     0.189     0.989     0.581     0.570 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>     0.570     1.000     1.000     0.989     0.581     0.010     0.000 
#> 
#> $indirect.effect
#>           y1.all y1.exercises y1.sweat y1.sports
#> est        0.155        0.014    0.000     0.091
#> mean       0.171       -0.024    0.008     0.174
#> sd         0.001        0.011    0.019     0.003
#> upbd       0.173       -0.002    0.045     0.181
#> lwbd       0.170       -0.046   -0.029     0.168
#> upbd_q     0.172       -0.016    0.020     0.177
#> lwbd_q     0.171       -0.031   -0.005     0.172
#> upbd_bcbi  0.171       -0.016    0.020     0.172
#> lwbd_bcbi  0.171       -0.016   -0.005     0.172
#> upbd_b     0.172        0.172    0.172     0.172
#> lwbd_b    -0.016       -0.016   -0.016    -0.016
#> upbd_win   0.172       -0.016    0.020     0.177
#> lwbd_win   0.171       -0.031   -0.005     0.172
#> p_norm     0.000        0.034    0.677     0.000
#> p_quan     0.000        0.000    1.000     0.000

When the argument nonlinear=TRUE, multiple Additive Regression Trees (MART) will be used to fit the final full model in estimating the outcome.

 mma.b.b.mart<-mma(x,y,pred=pred,contmed=c(7:9,11:12),binmed=c(6,10),binref=c(1,1), catmed=5,catref=1,predref="M",alpha=0.4,alpha2=0.4,nonlinear=TRUE,n=2,n2=5)
 summary(mma.b.b.mart)
 summary(mma.b.b.mart)
#> MMA Analysis: Estimated Mediation Effects Using MART
#> For Predictor/Moderator at pred 
#> $total.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.123     0.362     0.075     0.509     0.214     0.448     0.257     0.246 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>     0.246     0.387     0.246     0.453     0.249     0.000     0.000 
#> 
#> $direct.effect
#>       est      mean        sd      upbd      lwbd    upbd_q    lwbd_q upbd_bcbi 
#>     0.067     0.262     0.104     0.467     0.058     0.322     0.098     0.078 
#> lwbd_bcbi    upbd_b    lwbd_b  upbd_win  lwbd_win    p_norm    p_quan 
#>     0.078     0.322     0.078     0.322     0.078     0.012     0.000 
#> 
#> $indirect.effect
#>           y1.all y1.exercises y1.sweat y1.sports
#> est        0.056       -0.001    0.003     0.065
#> mean       0.099        0.002    0.009     0.110
#> sd         0.067        0.015    0.018     0.046
#> upbd       0.231        0.032    0.044     0.201
#> lwbd      -0.032       -0.028   -0.026     0.020
#> upbd_q     0.175        0.021    0.033     0.156
#> lwbd_q     0.036       -0.016   -0.010     0.045
#> upbd_bcbi  0.173        0.021    0.030     0.125
#> lwbd_bcbi  0.035       -0.013   -0.011     0.038
#> upbd_b     0.168        0.014    0.035     0.157
#> lwbd_b     0.034       -0.017   -0.011     0.038
#> upbd_win   0.176        0.021    0.035     0.157
#> lwbd_win   0.040       -0.013   -0.008     0.040
#> p_norm     0.139        0.873    0.602     0.017
#> p_quan     0.000        0.800    0.400     0.000

plot(mma.b.b.mart,vari="exercises")

plot(mma.b.b.glm,vari="sweat")

References

Yu, Qingzhao, Ying Fan, and Xiaocheng Wu. 2014. “General Multiple Mediation Analysis with an Application to Explore Racial Disparities in Breast Cancer Survival.” Journal of Biometrics & Biostatistics 5 (2): 1–9. https://doi.org/10.4172/2155-6180.1000189.
Yu, Qingzhao, and Bin Li. 2017. “Mma: An r Package for Mediation Analysis with Multiple Mediators.” Journal of Open Research Software 5 (1): 11. https://doi.org/10.5334/jors.160.
———. 2022. Statistical Methods for Mediation, Confounding and Moderation Analysis Using r and SAS. Chapman & Hall/CRC. https://www.routledge.com/Statistical-Methods-for-Mediation-Confounding-and-Moderation-Analysis/Yu-Li/p/book/9780367365479.
Yu, Qingzhao, Kaelen L. Medeiros, Xiaocheng Wu, and Roxanne E. Jensen. 2018. “Nonlinear Predictive Models for Multiple Mediation Analysis: With an Application to Explore Ethnic Disparities in Anxiety and Depression Among Cancer Survivors.” Psychometrika 83 (4): 991–1006. https://doi.org/10.1007/s11336-018-9612-2.
Yu, Qingzhao, Richard A. Scribner, Claudia Leonardi, Lu Zhang, Chi Park, Liwei Chen, and Neal R. Simonsen. 2017. “Exploring Racial Disparity in Obesity: A Mediation Analysis Considering Geo-Coded Environmental Factors.” Spatial and Spatio-Temporal Epidemiology 21: 13–23. https://doi.org/10.1016/j.sste.2017.02.001.
Yu, Qingzhao, Xiaocheng Wu, Bin Li, and Richard A. Scribner. 2019. “Multiple Mediation Analysis with Survival Outcomes: With an Application to Explore Racial Disparity in Breast Cancer Survival.” Statistics in Medicine 38 (3): 398–412. https://doi.org/10.1002/sim.7977.