# Implementation of extended formulas when there are effect measure modifiers

#### 2022-04-06

In this document, we demonstrate including effect measure modification (EMM) terms in the mediator or the outcome models. The dataset used in this document is still vv2015.

library(regmedint)
library(tidyverse)
## Prepare dataset
data(vv2015)

# No EMM by covariates

In the first model fit, we do not include any EMM term.

regmedint_obj1 <- regmedint(data = vv2015,
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 3,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
interaction = TRUE,
casecontrol = FALSE)
summary(regmedint_obj1)
## ### Mediator model
##
## Call:
## glm(formula = m ~ x + c, family = binomial(link = "logit"), data = data)
##
## Deviance Residuals:
##     Min       1Q   Median       3Q      Max
## -1.5143  -1.1765   0.9177   1.1133   1.4602
##
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept)  -0.3545     0.3252  -1.090    0.276
## x             0.3842     0.4165   0.922    0.356
## c             0.2694     0.2058   1.309    0.191
##
## (Dispersion parameter for binomial family taken to be 1)
##
##     Null deviance: 138.59  on 99  degrees of freedom
## Residual deviance: 136.08  on 97  degrees of freedom
## AIC: 142.08
##
## Number of Fisher Scoring iterations: 4
##
## ### Outcome model
##
## Call:
## survival::survreg(formula = Surv(y, event) ~ x + m + x:m + c,
##     data = data, dist = "weibull")
##               Value Std. Error     z       p
## (Intercept) -1.0424     0.1903 -5.48 4.3e-08
## x            0.4408     0.3008  1.47    0.14
## m            0.0905     0.2683  0.34    0.74
## c           -0.0669     0.0915 -0.73    0.46
## x:m          0.1003     0.4207  0.24    0.81
## Log(scale)  -0.0347     0.0810 -0.43    0.67
##
## Scale= 0.966
##
## Weibull distribution
## Loglik(model)= -11.4   Loglik(intercept only)= -14.5
##  Chisq= 6.31 on 4 degrees of freedom, p= 0.18
## Number of Newton-Raphson Iterations: 5
## n= 100
##
## ### Mediation analysis
##              est         se         Z          p       lower      upper
## cde  0.541070807 0.29422958 1.8389409 0.06592388 -0.03560858 1.11775019
## pnde 0.505391952 0.21797147 2.3186151 0.02041591  0.07817572 0.93260819
## tnie 0.015988820 0.03171597 0.5041252 0.61417338 -0.04617334 0.07815098
## tnde 0.513662425 0.22946248 2.2385465 0.02518544  0.06392423 0.96340062
## pnie 0.007718348 0.02398457 0.3218047 0.74760066 -0.03929055 0.05472725
## te   0.521380773 0.22427066 2.3247837 0.02008353  0.08181835 0.96094319
## pm   0.039039346 0.07444080 0.5244348 0.59997616 -0.10686194 0.18494063
##
## Evaluated at:
## avar: x
##  a1 (intervened value of avar) = 1
##  a0 (reference value of avar)  = 0
## mvar: m
##  m_cde (intervend value of mvar for cde) = 1
## cvar: c
##  c_cond (covariate vector value) = 3
##
## Note that effect estimates can vary over m_cde and c_cond values when interaction = TRUE.

# EMM by covariates

## There is $$A\times C$$ term in mediator model

Now suppose the covariate $$C$$ modifies the treatment effect on the mediator. We add emm_ac_mreg = c("c") in regmedint(). Although there is only one covariate in our dataset, emm_ac_mreg can take a vector of multiple covariates. Please note that the covariates in emm_ac_mreg should be a subset of the covariates specified in cvar, i.e. if a covariate is an effect measure modifier included in emm_ac_mreg, it must be included in cvar, otherwise an error message will be printed.

regmedint_obj2 <- regmedint(data = vv2015,
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
emm_ac_mreg = c("c"),
emm_ac_yreg = NULL,
emm_mc_yreg = NULL,
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 3,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
interaction = TRUE,
casecontrol = FALSE)
summary(regmedint_obj2)
## ### Mediator model
##
## Call:
## glm(formula = m ~ x + c + x:c, family = binomial(link = "logit"),
##     data = data)
##
## Deviance Residuals:
##     Min       1Q   Median       3Q      Max
## -1.5689  -1.1585   0.8925   1.1242   1.4342
##
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.32727    0.34979  -0.936    0.349
## x            0.30431    0.56789   0.536    0.592
## c            0.24085    0.24688   0.976    0.329
## x:c          0.09216    0.44624   0.207    0.836
##
## (Dispersion parameter for binomial family taken to be 1)
##
##     Null deviance: 138.59  on 99  degrees of freedom
## Residual deviance: 136.04  on 96  degrees of freedom
## AIC: 144.04
##
## Number of Fisher Scoring iterations: 4
##
## ### Outcome model
##
## Call:
## survival::survreg(formula = Surv(y, event) ~ x + m + x:m + c,
##     data = data, dist = "weibull")
##               Value Std. Error     z       p
## (Intercept) -1.0424     0.1903 -5.48 4.3e-08
## x            0.4408     0.3008  1.47    0.14
## m            0.0905     0.2683  0.34    0.74
## c           -0.0669     0.0915 -0.73    0.46
## x:m          0.1003     0.4207  0.24    0.81
## Log(scale)  -0.0347     0.0810 -0.43    0.67
##
## Scale= 0.966
##
## Weibull distribution
## Loglik(model)= -11.4   Loglik(intercept only)= -14.5
##  Chisq= 6.31 on 4 degrees of freedom, p= 0.18
## Number of Newton-Raphson Iterations: 5
## n= 100
##
## ### Mediation analysis
##             est         se         Z          p       lower      upper
## cde  0.54107081 0.29422958 1.8389409 0.06592388 -0.03560858 1.11775019
## pnde 0.50404000 0.21666437 2.3263632 0.01999919  0.07938564 0.92869435
## tnie 0.02377050 0.05679639 0.4185213 0.67556602 -0.08754838 0.13508937
## tnde 0.51632801 0.23444392 2.2023519 0.02764046  0.05682637 0.97582965
## pnie 0.01148248 0.03882957 0.2957149 0.76744784 -0.06462208 0.08758704
## te   0.52781049 0.22811645 2.3137765 0.02067998  0.08071046 0.97491053
## pm   0.05727853 0.13042523 0.4391676 0.66054011 -0.19835021 0.31290728
##
## Evaluated at:
## avar: x
##  a1 (intervened value of avar) = 1
##  a0 (reference value of avar)  = 0
## mvar: m
##  m_cde (intervend value of mvar for cde) = 1
## cvar: c
##  c_cond (covariate vector value) = 3
##
## Note that effect estimates can vary over m_cde and c_cond values when interaction = TRUE.

## There is $$A\times C$$ term in both mediator and outcome models

Now suppose in addition to the EMM on mediator, the covariate $$C$$ also modifies the treatment effect on the outcome We add emm_ac_yreg = c("c") in regmedint(). Please note that the covariates in emm_ac_yreg should be a subset of the covariates specified in cvar, i.e. if a covariate is an effect measure modifier included in emm_ac_yreg, it must be included in cvar, otherwise an error message will be printed.

regmedint_obj3 <- regmedint(data = vv2015,
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
emm_ac_mreg = c("c"),
emm_ac_yreg = c("c"),
emm_mc_yreg = NULL,
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 3,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
interaction = TRUE,
casecontrol = FALSE)
summary(regmedint_obj3)
## ### Mediator model
##
## Call:
## glm(formula = m ~ x + c + x:c, family = binomial(link = "logit"),
##     data = data)
##
## Deviance Residuals:
##     Min       1Q   Median       3Q      Max
## -1.5689  -1.1585   0.8925   1.1242   1.4342
##
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.32727    0.34979  -0.936    0.349
## x            0.30431    0.56789   0.536    0.592
## c            0.24085    0.24688   0.976    0.329
## x:c          0.09216    0.44624   0.207    0.836
##
## (Dispersion parameter for binomial family taken to be 1)
##
##     Null deviance: 138.59  on 99  degrees of freedom
## Residual deviance: 136.04  on 96  degrees of freedom
## AIC: 144.04
##
## Number of Fisher Scoring iterations: 4
##
## ### Outcome model
##
## Call:
## survival::survreg(formula = Surv(y, event) ~ x + m + x:m + c +
##     x:c, data = data, dist = "weibull")
##                Value Std. Error     z       p
## (Intercept) -1.04148    0.19261 -5.41 6.4e-08
## x            0.43626    0.33092  1.32    0.19
## m            0.09138    0.26954  0.34    0.73
## c           -0.06844    0.10315 -0.66    0.51
## x:m          0.09681    0.43437  0.22    0.82
## x:c          0.00725    0.22300  0.03    0.97
## Log(scale)  -0.03473    0.08104 -0.43    0.67
##
## Scale= 0.966
##
## Weibull distribution
## Loglik(model)= -11.4   Loglik(intercept only)= -14.5
##  Chisq= 6.31 on 5 degrees of freedom, p= 0.28
## Number of Newton-Raphson Iterations: 5
## n= 100
##
## ### Mediation analysis
##             est         se         Z         p       lower      upper
## cde  0.55481498 0.51519657 1.0768996 0.2815251 -0.45495174 1.56458171
## pnde 0.51905802 0.51048271 1.0167984 0.3092493 -0.48146970 1.51958574
## tnie 0.02345019 0.05730239 0.4092359 0.6823666 -0.08886042 0.13576081
## tnde 0.53092081 0.50659327 1.0480218 0.2946285 -0.46198375 1.52382537
## pnie 0.01158740 0.03904582 0.2967641 0.7666466 -0.06494101 0.08811581
## te   0.54250821 0.50656000 1.0709654 0.2841850 -0.45033114 1.53534756
## pm   0.05535403 0.13968838 0.3962680 0.6919074 -0.21843016 0.32913822
##
## Evaluated at:
## avar: x
##  a1 (intervened value of avar) = 1
##  a0 (reference value of avar)  = 0
## mvar: m
##  m_cde (intervend value of mvar for cde) = 1
## cvar: c
##  c_cond (covariate vector value) = 3
##
## Note that effect estimates can vary over m_cde and c_cond values when interaction = TRUE.

## There are $$A\times C$$ term in both mediator and outcome models, and $$M\times C$$ term in outcome model

Now suppose in addition to the EMM of treatment effect, the covariate $$C$$ also modifies the mediator effect on the outcome. We add emm_mc_yreg = c("c") in regmedint(). Please note that the covariates in emm_mc_yreg should be a subset of the covariates specified in cvar, i.e. if a covariate is an effect measure modifier included in emm_mc_yreg, it must be included in cvar, otherwise an error message will be printed.

regmedint_obj4 <- regmedint(data = vv2015,
## Variables
yvar = "y",
avar = "x",
mvar = "m",
cvar = c("c"),
emm_ac_mreg = c("c"),
emm_ac_yreg = c("c"),
emm_mc_yreg = c("c"),
eventvar = "event",
## Values at which effects are evaluated
a0 = 0,
a1 = 1,
m_cde = 1,
c_cond = 3,
## Model types
mreg = "logistic",
yreg = "survAFT_weibull",
interaction = TRUE,
casecontrol = FALSE)
summary(regmedint_obj4)
## ### Mediator model
##
## Call:
## glm(formula = m ~ x + c + x:c, family = binomial(link = "logit"),
##     data = data)
##
## Deviance Residuals:
##     Min       1Q   Median       3Q      Max
## -1.5689  -1.1585   0.8925   1.1242   1.4342
##
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.32727    0.34979  -0.936    0.349
## x            0.30431    0.56789   0.536    0.592
## c            0.24085    0.24688   0.976    0.329
## x:c          0.09216    0.44624   0.207    0.836
##
## (Dispersion parameter for binomial family taken to be 1)
##
##     Null deviance: 138.59  on 99  degrees of freedom
## Residual deviance: 136.04  on 96  degrees of freedom
## AIC: 144.04
##
## Number of Fisher Scoring iterations: 4
##
## ### Outcome model
##
## Call:
## survival::survreg(formula = Surv(y, event) ~ x + m + x:m + c +
##     x:c + m:c, data = data, dist = "weibull")
##               Value Std. Error     z       p
## (Intercept) -0.9959     0.2071 -4.81 1.5e-06
## x            0.4185     0.3354  1.25    0.21
## m           -0.0216     0.3112 -0.07    0.94
## c           -0.1339     0.1405 -0.95    0.34
## x:m          0.0905     0.4265  0.21    0.83
## x:c          0.0327     0.2242  0.15    0.88
## m:c          0.1275     0.1861  0.69    0.49
## Log(scale)  -0.0406     0.0814 -0.50    0.62
##
## Scale= 0.96
##
## Weibull distribution
## Loglik(model)= -11.1   Loglik(intercept only)= -14.5
##  Chisq= 6.78 on 6 degrees of freedom, p= 0.34
## Number of Newton-Raphson Iterations: 5
## n= 100
##
## ### Mediation analysis
##             est         se         Z         p      lower     upper
## cde  0.60705735 0.52594922 1.1542128 0.2484129 -0.4237842 1.6378989
## pnde 0.57902523 0.51447701 1.1254638 0.2603926 -0.4293312 1.5873816
## tnie 0.05333600 0.10591830 0.5035579 0.6145721 -0.1542601 0.2609321
## tnde 0.58889505 0.51488644 1.1437377 0.2527324 -0.4202638 1.5980539
## pnie 0.04346618 0.09107534 0.4772552 0.6331804 -0.1350382 0.2219706
## te   0.63236123 0.52776615 1.1981845 0.2308452 -0.4020414 1.6667639
## pm   0.11082259 0.20960355 0.5287248 0.5969964 -0.2999928 0.5216380
##
## Evaluated at:
## avar: x
##  a1 (intervened value of avar) = 1
##  a0 (reference value of avar)  = 0
## mvar: m
##  m_cde (intervend value of mvar for cde) = 1
## cvar: c
##  c_cond (covariate vector value) = 3
##
## Note that effect estimates can vary over m_cde and c_cond values when interaction = TRUE.