The purpose of this vignette is to provide users with a step-by-step guide for performing and interpreting the results of a time varying mediation analysis with three treatment (exposure) groups and a continuous outcome. Please note, that this package has been built considering the structure of panel data, where each subject/participant has repeated responses (collected over time) for the outcome and mediator variables. We do not address dynamic treatment regimens in this package. Therefore, we assume the scenario where the treatment (exposure) is time-invariant (i.e., does not change over time). For more details, refer to Cai et al., 2022.

For illustration, we rely on an example dataset, simulated based on
the Wisconsin Smokers’ Health Study 2 data (Baker et. al., 2016) which
includes 1086 individuals assigned to one of three treatment conditions.
One-third of participants received only a `nicotine patch`

;
another one-third received `varenicline`

, and the final third
of participants received a
`combination nicotine replacement therapy (NRT)`

which
included nicotine patch + nicotine mini-lozenge. The outcome of interest
is `cessation fatigue`

; that is, how tired a participant felt
of trying to quit smoking (7-point Likert scale). In addition, mediator
variables were measured by asking participants if they felt a
`negative mood in the last fifteen minutes`

, and whether they
`wanted to smoke in the last 15 minutes`

, also recorded on a
7-point Likert scale. Both the outcome and mediator variables were
assessed two times per day for the first two weeks post-quit day
(rendering 30 time points of response since assessments start on day 0
post target quit day), and every other day (2x per day) for weeks 3 and
4 (rendering 14 time points of response).

A traditional approach to analyzing this type of data would be to use mediation analysis in which the effects are assumed to not vary as a function of time. First, a single (i.e., time-invariant) direct effect would be calculated by regressing the outcome on the treatment condition and mediator. Next, a time-invariant indirect effect would be computed by multiplying the effect of treatment condition on the mediator by the effect of the mediator on the outcome. However, this method potentially misses important information about the dynamic effect that a mediator may have over time. Specifically, we hypothesize that mood changes across and within days and thus, its mediating effect on one’s success of quitting smoking is likely to vary over time. We therefore propose a time varying mediation analysis which estimates the mediation effect as a function that varies over time.

To use the time varying mediation analysis package in R, you must
first install the package and load it. Before that, make sure you have
`R version 4.0.3`

or greater. There are two ways to install
the package from the CRAN (Comprehensive R Archive Network) repository,
by using `install.packages`

or the `devtools`

function.

```
install.packages("tvmediation", dependencies = TRUE)
library(tvmediation)
```

The equivalent code using `devtools`

is:

```
::install_cran("tvmediation", dependencies = TRUE)
devtoolslibrary(tvmediation)
```

If you do not have `devtools`

installed and loaded, you
can do so using the following code:

```
install.packages("devtools", dependencies = TRUE)
library(devtools)
```

Alternatively, if you want to install the package directly from the GitHub repository to access new or revised functions in development, use the following code:

```
::install_github("dcoffman/tvmediation", dependencies = TRUE)
devtoolslibrary(tvmediation)
```

`tvma_3trt`

functionOnce installed, you can type `?tvmediation`

in the console
to view the package documentation, as well as links to the important
functions and data included in the package. The time-varying mediation
analysis for continuous outcomes and three exposure groups, relies on
two user functions `tvma_3trt`

and `LongToWide`

as
well as a number of internal functions of the `tvmediation`

package.

The `tvma_3trt`

function has five required and six
optional arguments.

`T1`

A vector indicating assignment to treatment 1`T2`

A vector indicating assignment to treatment 2`t.seq`

The numeric vector of the time sequence of the measures`mediator`

The matrix of mediator values in wide format`outcome`

The matrix of outcome values in wide format

The optional arguments are:

`t.est`

The time sequence for estimation. This is by default equal to`t.seq`

.`plot`

TRUE or FALSE for plotting the mediation effect. The default value is “FALSE”.`CI`

“none” or “boot” method of deriving confidence intervals (CIs). The default value is “boot”.`replicates`

Number of replicates for bootstrapping CIs. The default value is 1000.`grpname`

Name of the treatment (exposure) groups to be displayed in the results. The default value is “T”.`verbose`

TRUE or FALSE for printing results to screen. The default value is “FALSE”.

The dataset we will use for our illustration is named
`smoker`

and is included in the package.

To load the simulated dataset `smoker.rda`

, type:

`data(smoker)`

The `smoker`

data frame is organized in **long
format** with `SubjectID`

repeating over multiple rows
for each participant. The `tvma_3trt`

function requires that
the data be in **wide format** to estimate the time varying
coefficients. The `tvmediation`

package includes a useful
function `LongToWide`

to help users properly format their
data for analysis.

`LongToWide`

has three required arguments and a fourth
optional argument.

`subject.id`

specifies the column of subject identifiers.`time.sequences`

specifies the column of measurement times.`outcome`

specifies the column for the variable (either the outcome or the mediator) that will be transposed.`verbose`

is an optional argument that if “TRUE” prints the output of`LongToWide`

to the console. The default value is “FALSE”.

The output of `LongToWide`

is a matrix of data in wide
format where columns represent the subjects and rows represent the time
sequence. Thus, each cell contains the j-th subject’s response at the
i-th time point.

The `tvma_3trt`

function requires two matrices, one for
the mediator, and one for the outcome. Thus, we use the
`LongToWide`

function twice as illustrated below:

```
<- LongToWide(smoker$SubjectID, smoker$timeseq,
mediator $NegMoodLst15min)
smoker<- LongToWide(smoker$SubjectID, smoker$timeseq, smoker$cessFatig) outcome
```

`class(mediator)`

`## [1] "matrix" "array"`

`1:16, 1:10] mediator[`

```
## 1 2 3 4 5 6 7 8 9 10
## 0 7 1 1 1 1 NA 4 1 1 3
## 0.5 2 1 1 6 1 NA NA 1 1 1
## 1 1 1 2 NA NA 1 3 1 1 1
## 1.5 1 1 1 1 1 1 2 1 NA 1
## 2 1 1 1 1 NA 1 1 1 1 1
## 2.5 3 1 1 1 1 3 1 1 1 1
## 3 1 6 NA 1 1 1 1 2 NA 1
## 3.5 1 1 2 6 3 7 1 6 1 1
## 4 1 1 1 1 1 1 1 1 2 4
## 4.5 1 1 1 5 1 5 2 NA 1 3
## 5 4 1 2 1 NA 5 2 1 1 1
## 5.5 1 1 2 1 1 2 NA 1 1 2
## 6 2 1 1 4 1 7 1 NA 2 3
## 6.5 1 1 1 1 1 1 1 5 3 1
## 7 1 2 1 NA 1 1 NA 1 1 4
## 7.5 2 2 1 1 1 NA 1 1 1 6
```

`class(outcome)`

`## [1] "matrix" "array"`

`1:16, 1:10] outcome[`

```
## 1 2 3 4 5 6 7 8 9 10
## 0 1 6 1 2 1 NA 1 6 7 3
## 0.5 1 3 1 1 1 NA NA 1 7 1
## 1 1 1 1 NA NA 1 1 7 7 1
## 1.5 6 7 6 1 1 6 6 4 NA 7
## 2 7 7 1 1 NA 1 7 1 1 1
## 2.5 1 2 1 7 1 6 1 2 1 7
## 3 4 1 NA 1 7 2 1 4 NA 5
## 3.5 5 3 7 7 6 1 1 7 5 1
## 4 6 3 5 4 7 2 7 5 7 5
## 4.5 7 7 1 7 7 1 7 NA 4 1
## 5 7 7 1 1 NA 1 1 7 1 1
## 5.5 6 5 1 1 1 3 NA 5 7 1
## 6 1 1 1 7 5 1 7 NA 2 1
## 6.5 5 7 5 7 7 1 1 1 7 7
## 7 1 1 1 NA 7 1 NA 1 7 7
## 7.5 5 1 2 1 1 NA 1 1 1 7
```

If your data are already in wide format, there is no need to use the
`LongToWide`

function and you can simply subset your dataset.
However, `mediator`

and `outcome`

must be of class
`matrix`

; hence make sure you convert the class of the
subsetted `mediator`

and `outcome`

objects to
`matrix`

before proceeding. This can be done using the R
function `as.matrix`

.

The `tvma_3trt`

function requires three more variables
that we have not yet created:

`T1`

A binary numeric vector with treatment1 schedule`T2`

A binary numeric vector with treatment2 schedule`t.seq`

A numeric vector of the time sequence of the measures

When there are three treatment groups, we can create two dummy
variables `NRT1`

and `NRT2`

corresponding to the
treatment (exposure) groups other than the placebo or, in this case, the
`nicotine patch`

only standard of care. `NRT1 = 1`

indicates that the participant has been given `varenicline`

and `NRT1 = 0`

for `no varenicline`

. Similarly,
`NRT2 = 1`

indicates that the participant received
`combination NRT`

and `NRT2 = 0`

indicates that
the participant did not receive the
`combination treatment (comboNRT)`

. Naturally,
`NRT1 = 1`

and `NRT2 = 1`

are mutually exclusive;
that is, a participant cannot have `NRT1 = 1`

and
`NRT2 = 1`

simultaneously. If `NRT1 = 1`

then for
that participant `NRT2`

must be equal to `0`

.
However, `NRT1 = 0`

and `NRT2 = 0`

are not
mutually exclusive. When `NRT1 = 0`

and `NRT2 = 0`

then the participant has received the placebo or the standard of care,
in this case `nicotine patch`

only. Thus if there are three
treatment groups, two dummy variables need to be created indicating
`NRT1`

or `treatment1`

and `NRT2`

or
`treatment2`

. In the `smoker.rda`

dataset, three
columns indicating the assignment of the three treatments are present,
such that we can use any two of the columns and the omitted column
becomes the reference group.

We create two dummy variables to indicate whether a participant was
given `varenicline`

or `combination NRT`

. Because
`treatment`

is time-invariant, we need only one instance of
each subject’s response for the treatment group
(e.g. `varenicline`

or `combination NRT`

). We then
convert it to a numeric value and subtract 1 to yield a vector of zeros
and ones, as shown below.

```
# Step 1: Since each subject has multiple rows of data, extract the unique response of each subject to receiving varenicline. The data is still in dataframe format.
<- unique(smoker[ , c("SubjectID","varenicline")])
trtv 1:10,] trtv[
```

```
## SubjectID varenicline
## 1 1 No
## 2 2 Yes
## 3 3 Yes
## 4 4 No
## 5 5 Yes
## 6 6 Yes
## 7 7 Yes
## 8 8 Yes
## 9 9 Yes
## 10 10 Yes
```

```
# Step 2: `2` to those subjects who received varenicline and `1` to the rest. The data is now in vector format.
<- as.numeric(trtv[ , 2])
trtv2 1:10] trtv2[
```

`## [1] 2 1 1 2 1 1 1 1 1 1`

```
# Step 3: subtract 1 from these numeric responses and procure a vector of zeros and ones
<- trtv2 -1
NRT1 1:10] NRT1[
```

`## [1] 1 0 0 1 0 0 0 0 0 0`

```
# Step 1: Since each subject has multiple rows of data, extract the unique response of each subject to receiving comboNRT. The data is still in dataframe format.
<- unique(smoker[ , c("SubjectID","comboNRT")])
trtc 1:10,] trtc[
```

```
## SubjectID comboNRT
## 1 1 Yes
## 2 2 No
## 3 3 Yes
## 4 4 Yes
## 5 5 No
## 6 6 Yes
## 7 7 Yes
## 8 8 Yes
## 9 9 Yes
## 10 10 Yes
```

```
# Step 2: `2` to those subjects who received comboNRT and `1` to the rest. The data is now in vector format.
<- as.numeric(trtc[ , 2])
trtc2 1:10] trtc2[
```

`## [1] 1 2 1 1 2 1 1 1 1 1`

```
# Step 3: subtract 1 from these numeric responses and procure a vector of zeros and ones
<- trtc2 -1
NRT2 1:10] NRT2[
```

`## [1] 0 1 0 0 1 0 0 0 0 0`

This steps can be alternatively collated into a single step and written as follows:

```
<- as.numeric(unique(smoker[ ,c("SubjectID","varenicline")])[,2])-1
NRT1 <- as.numeric(unique(smoker[ ,c("SubjectID","comboNRT")])[,2])-1
NRT2 1:10] NRT1[
```

`## [1] 1 0 0 1 0 0 0 0 0 0`

`1:10] NRT2[`

`## [1] 0 1 0 0 1 0 0 0 0 0`

To generate `t.seq`

we found the unique instance of each
time point and then sorted from smallest to largest. There are 44 unique
time points in the dataset where `0`

after the decimal
indicates the morning measurement and `5`

after the decimal
indicates the evening measurement recorded for that day.

```
<- sort(unique(smoker$timeseq))
t.seq t.seq
```

```
## [1] 0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0 5.5 6.0 6.5 7.0
## [16] 7.5 8.0 8.5 9.0 9.5 10.0 10.5 11.0 11.5 12.0 12.5 13.0 13.5 14.0 14.5
## [31] 16.0 16.5 18.0 18.5 20.0 20.5 22.0 22.5 24.0 24.5 26.0 26.5 28.0 28.5
```

We are now ready to perform the time varying mediation analysis.

`tvma_3trt`

functionAs discussed earlier, the `tvma_3trt`

function has five
required and six optional arguments.

`T1`

A vector indicating assignment to treatment 1`T2`

A vector indicating assignment to treatment 2`t.seq`

The numeric vector of the time sequence of the measures`mediator`

The matrix of mediator values in wide format`outcome`

The matrix of outcome values in wide format

The optional arguments are:

`t.est`

The time sequence for estimation. This is by default equal to`t.seq`

.`plot`

TRUE or FALSE for plotting the mediation effect. The default value is “FALSE”.`CI`

“none” or “boot” method of deriving CIs. The default value is “boot”.`replicates`

Number of replicates for bootstrapping CIs. The default value is 1000.`grpname`

Name of the treatment (exposure) groups to be displayed in the results. The default value is “T”.`verbose`

TRUE or FALSE for printing results to screen. The default value is “FALSE”.

The argument `grpname`

is important if the user wants the
results from the bootstrapped CIs for the mediation effects to be
displayed with a particular name in the `Estimates`

dataframe
and plots. For example, if the user defines
`grpname = "exposure"`

, instead of the default value
`T`

, the effect of this choice will be reflected in the
column names; `CI.upper.T1`

will be
`CI.upper.exposure1`

and so forth. The plot titles will also
change accordingly. For our illustration we use `NRT`

since
the example is drawn from a Nicotine Replacement Therapy study.

We will call the function with additional optional arguments
`plot=TRUE`

and `replicates = 250`

. We decreased
the number of bootstrap replicates so that this vignette compiles faster
but we suggest using at least 500 replicates in an actual analysis. The
remaining optional arguments are left to their respective default
values.

```
<- tvma_3trt(NRT1, NRT2, t.seq, mediator, outcome, plot = TRUE,
results_tvma_3trt replicates = 250, grpname = "NRT")
```

The `tvma_3trt`

function returns a list of results that
include:

`hat.alpha1`

estimated time-varying effect of treatment 1 on the mediator`CI.lower.alpha1`

CI lower limit for`hat.alpha1`

`CI.upper.alpha1`

CI upper limit for`hat.alpha1`

`hat.alpha2`

estimated time-varying effect of treatment 2 on the mediator`CI.lower.alpha2`

CI lower limit for`hat.alpha2`

`CI.upper.alpha2`

CI upper limit for`hat.alpha2`

`hat.gamma1`

estimated time-varying direct effect of treatment 1 on the outcome`CI.lower.gamma1`

CI lower limit for`hat.gamma1`

`CI.upper.gamma1`

CI upper limit for`hat.gamma1`

`hat.gamma2`

estimated time-varying direct effect of treatment 2 on the outcome`CI.lower.gamma2`

CI lower limit for`hat.gamma2`

`CI.upper.gamma2`

CI upper limit for`hat.gamma2`

`hat.tau1`

estimated time-varying total effect of treatment 1 on the outcome`CI.lower.tau1`

CI lower limit for`hat.tau1`

`CI.upper.tau1`

CI upper limit for`hat.tau1`

`hat.tau2`

estimated time-varying total effect of treatment 2 on the outcome`CI.lower.tau2`

CI lower limit for`hat.tau2`

`CI.upper.tau2`

CI upper limit for`hat.tau2`

`hat.beta`

estimated time-varying effect of the mediator on the outcome`CI.lower.beta`

CI lower limit for`hat.beta`

`CI.upper.beta`

CI upper limit for`hat.beta`

`hat.mediation1`

time varying mediation effect of treatment 1 on the outcome`hat.mediation2`

time varying mediation effect of treatment 2 on the outcome

Optional returns based on argument `CI = "boot"`

include:

`CI.lower.T1`

CI lower limit for the time varying mediation effect of treatment 1`CI.upper.T1`

CI upper limit for the time varying mediation effect of treatment 1`CI.lower.T2`

CI lower limit for the time varying mediation effect of treatment 2`CI.upper.T2`

CI upper limit for the time varying mediation effect of treatment 2`SE_MedEff1`

estimated standard error of the time varying mediation effect for treatment 1`SE_MedEff2`

estimated standard error of the time varying mediation effect for treatment 2

The above estimates are compiled in a single dataframe which can be
accessed using `nameOfStoredResultsObj$Estimates`

. The
following line of code displays only the estimates at the first 6
time-points.

`head(results_tvma_3trt$Estimates)`

```
## timeseq hat.alpha1 CI.lower.alpha1 CI.upper.alpha1 hat.alpha2 CI.lower.alpha2
## 1 0.0 -0.7157486 -0.3763980 -1.0326121 -0.8871672 -0.3725386
## 2 0.5 -0.7000508 -0.4451674 -0.9045140 -0.7883056 -0.4075918
## 3 1.0 -0.6924259 -0.5103652 -0.8447188 -0.7115165 -0.4556211
## 4 1.5 -0.6913309 -0.5643743 -0.8142764 -0.6621873 -0.5055364
## 5 2.0 -0.6949645 -0.5909751 -0.7884795 -0.6423338 -0.5381122
## 6 2.5 -0.7014576 -0.6088315 -0.7766434 -0.6489869 -0.5757931
## CI.upper.alpha2 hat.gamma1 CI.lower.gamma1 CI.upper.gamma1 hat.gamma2
## 1 -1.1991911 0.5192259 0.8313453 0.1153684 0.6184312
## 2 -1.0054830 0.5979069 0.8345743 0.3041455 0.6632065
## 3 -0.8975101 0.6577964 0.8263311 0.4609218 0.6978976
## 4 -0.8014465 0.6981700 0.8404256 0.5706859 0.7219649
## 5 -0.7671274 0.7199237 0.8494571 0.6042949 0.7358780
## 6 -0.7534555 0.7253630 0.8449502 0.6123135 0.7408749
## CI.lower.gamma2 CI.upper.gamma2 hat.tau1 CI.lower.tau1 CI.upper.tau1
## 1 0.9241760 0.2589904 -0.271779695 0.2223223 -1.2098842
## 2 0.8697178 0.3985754 -0.174314576 0.1440487 -0.8789909
## 3 0.8746704 0.4986232 -0.097974873 0.1029291 -0.6203830
## 4 0.8855653 0.5739183 -0.042454107 0.1212315 -0.3597596
## 5 0.8836782 0.6142476 -0.006592278 0.1560464 -0.2033455
## 6 0.8607226 0.6242374 0.011546099 0.1811074 -0.1225190
## hat.tau2 CI.lower.tau2 CI.upper.tau2 hat.beta CI.lower.beta
## 1 -0.379043226 0.2985753 -0.9370762 0.1854940 0.2674796
## 2 -0.258942738 0.1713527 -0.6369531 0.1919473 0.2560043
## 3 -0.149034271 0.1376944 -0.4517728 0.1949410 0.2456874
## 4 -0.057886622 0.1704225 -0.2904095 0.1965561 0.2413538
## 5 0.008665267 0.2123508 -0.1711534 0.1983277 0.2411077
## 6 0.048551637 0.2301981 -0.1074584 0.2012113 0.2435668
## CI.upper.beta hat.mediation1 SE_MedEff1 CI.lower.NRT1 CI.upper.NRT1
## 1 0.06438462 -0.1327670 0.04097862 -0.2055110 -0.05669306
## 2 0.12760761 -0.1343729 0.02827195 -0.1823950 -0.08504908
## 3 0.14725170 -0.1349822 0.02087980 -0.1737506 -0.09842275
## 4 0.16119942 -0.1358853 0.01710684 -0.1687136 -0.10387303
## 5 0.16529845 -0.1378307 0.01549510 -0.1686552 -0.10710785
## 6 0.16929354 -0.1411412 0.01484733 -0.1712198 -0.11110177
## hat.mediation2 SE_MedEff2 CI.lower.NRT2 CI.upper.NRT2
## 1 -0.1645642 0.04824542 -0.2434645 -0.06742130
## 2 -0.1513131 0.03214720 -0.2047666 -0.08408166
## 3 -0.1387037 0.02249151 -0.1831488 -0.09283067
## 4 -0.1301569 0.01758947 -0.1634173 -0.09545963
## 5 -0.1273926 0.01564932 -0.1586542 -0.09759681
## 6 -0.1305835 0.01473999 -0.1598271 -0.10352024
```

At each time point of interest `timeseq = t.est`

, which in
this case is equal to `t.seq`

, the effects of the treatment
on the mediator, the treatment on the outcome (adjusted and not adjusted
for mediator), and the mediator on the outcome are estimated along with
the respective 95% CIs. The CIs are computed via a non-parametric
bootstrap method (Efron and Tibshirani, 1986), drawing samples of size
1086 from the original sample with replacement, estimating the sample
mean and then applying the percentile method to compute the 95% CIs.
Note that the CIs for the `alpha`

, `gamma`

,
`beta`

and `tau`

coefficients
`(hat.alpha.1, hat.alpha.2, hat.gamma1, hat.gamma2, hat.tau1, hat.tau2, hat.beta)`

are computed regardless of the value of `CI`

argument in the
function. `hat.mediation1`

and `hat.mediation2`

are the estimated mediation effects of `NRT1 (varenicline)`

and `NRT2 (combination NRT)`

compared to the
`placebo (nicotine patch only)`

that varies over
`t.est`

. For `CI = "boot"`

(which is the default
option unless the user chooses otherwise) the standard errors of the
estimated mediation effects and 95% CIs are estimated via a similar
bootstrapping technique described earlier for the coefficients.

If `plot = TRUE`

, the results will also include the
following figures:

`plot1_a1`

plot for`hat.alpha1`

with 95% CIs across`timeseq`

`plot2_a2`

plot for`hat.alpha2`

with 95% CIs across`timeseq`

`plot3_g1`

plot for`hat.gamma1`

with 95% CIs across`timeseq`

`plot4_g2`

plot for`hat.gamma2`

with 95% CIs across`timeseq`

`plot5_t1`

plot for`hat.tau1`

with 95% CIs across`timeseq`

`plot6_t2`

plot for`hat.tau2`

with 95% CIs across`timeseq`

`plot7_b`

plot for`hat.beta`

with 95% CIs across`timeseq`

`MedEff_T1`

plot for`hat.mediation1`

across`timeseq`

`MedEff_T2`

plot for`hat.mediation2`

across`timeseq`

`MedEff_CI_T1`

plot for`hat.mediation1`

with 95% CIs across`timeseq`

`MedEff_CI_T2`

plot for`hat.mediation2`

with 95% CIs across`timeseq`

We recommend using the plots to interpret your findings as it may be
difficult to derive meaning from the numerical values alone. To display
the plots, use `nameOfStoredResultsObj$`

followed by the name
of the plot to access the required plot accordingly. For example:

`$plot1_a1 results_tvma_3trt`

In the above plot, the effect of `varenicline`

on
subjects’ feeling of
`negative mood in the last fifteen minutes`

varies somewhat
over time. The effect is negative compared to the effect of
`nicotine patch`

only. That is, those assigned to
`varenicline`

have a less negative mood in the last 15 min.
compared to those assigned to the `nicotine patch`

only
group.

`$plot2_a2 results_tvma_3trt`

In the above plot, the effect of `combination NRT`

on
subjects’ feeling of
`negative mood in the last fifteen minutes`

varies over time.
The effect is negative compared to the effect of
`nicotine patch`

only. That is, those assigned to
`combination NRT`

have a less negative mood in the last 15
min. compared to those assigned to the `nicotine patch`

only
group.

`$plot3_g1 results_tvma_3trt`

The above plot shows the time-varying direct effect of
`varenicline`

on the outcome `cessation fatigue`

is positive compared to the direct effect of `nicotine patch`

only. That is, those assigned to `varenicline`

have greater
cessation fatigue (due to factors other than negative mood) compared to
those in the `nicotine patch`

only group.

`$plot4_g2 results_tvma_3trt`

The above plot shows the time-varying direct effect of
`combination NRT`

on the outcome
`cessation fatigue`

is positive compared to the direct effect
of `nicotine patch`

only. That is, those assigned to
`combination NRT`

have greater cessation fatigue (due to
factors other than negative mood) compared to those in the
`nicotine patch`

only group.

`$plot5_t1 results_tvma_3trt`

The above plot shows the time-varying total effect of
`varenicline`

on the outcome `cessation fatigue`

,
which is not statistically significant over much of the time. The
estimated 95% CI covers 0 (no effect), except during days 6-9.

`$plot6_t2 results_tvma_3trt`

The above plot shows the time-varying total effect of
`combination NRT`

on the outcome
`cessation fatigue`

, which is not statistically significant
over time. The estimated 95% CI covers 0 (no effect).

`$plot7_b results_tvma_3trt`

In the above plot, the time-varying effect of subjects’
`negative mood in the last fifteen minutes`

on the outcome
`cessation fatigue`

increases slightly until around day 6-7
(end of week 1) and then decreases beginning around day 20-21 (end of
week 3).

`$MedEff_T1 results_tvma_3trt`

`$MedEff_CI_T1 results_tvma_3trt`

In the above plots, the time-varying effect of
`varenicline`

on `cessation fatigue`

that is
mediated by the `negative mood in the last fifteen minutes`

is negative indicating that in comparison to the
`nicotine patch`

only group, `varenicline`

reduces
`cessation fatigue`

due to a decrease in negative mood. This
effect decreases during the first week, increases slightly from day
8-10, and then continues to weaken further until week 3. The effect
increases after week 3 with a slight decrease around days 24-26. Given
that the CIs do not cross zero, we can conclude that this effect is
statistically significant.

`$MedEff_T2 results_tvma_3trt`

`$MedEff_CI_T2 results_tvma_3trt`

In the above plots, the effect of `combination NRT`

on
`cessation fatigue`

that is mediated by the
`negative mood in the last fifteen minutes`

is weaker than
the effect of using the `nicotine patch`

only. This effect is
sinusoidal over time, with occasional spikes at days 2, 11, 24 and 28.
The effect takes a decreasing trend between week 2 and week 3. The
effect increases after week 3 and continues its upward trend with
occasional decrease. Given that the confidence intervals do not cross
zero, we can conclude that this effect is statistically significant.

The `tvma_3trt`

function computes bootstrap confidence
intervals by default. Therefore, if the user decides to not bootstrap
CIs for the mediation effect by specifying `CI = "none"`

, but
by mistake also specifies `replicates = 500`

, the function
will not display an error, but will simply execute without computing the
CIs for mediation effect. Note that the CIs for the effects of the
treatment on the mediator and the mediator on the outcome, and for the
direct and total effects are computed even if the user passes the
argument `CI = "none"`

.

The `tvmediation`

package provides a set of functions for
estimating mediation effects that vary over time for both binary and
continuous time-varying outcomes. Currently, the package only allows for
a time-invariant treatment. The mediator and outcome are assumed to be
time-varying, such as intensive longitudinal measurements obtained via
Ecological Momentary Assessment or via wearable and mobile devices. The
development of this tool has widespread application for use in human
behavior research, clinical trials, addiction research, and others by
allowing specification of mediation effects that vary as a function of
time.

Cai X, Coffman DL, Piper ME, Li R. Estimation and inference for the mediation effect in a time-varying mediation model. BMC Med Res Methodol. 2022;22(1):1-12.

Baker TB, Piper ME, Stein JH, et al. Effects of Nicotine Patch vs Varenicline vs Combination Nicotine Replacement Therapy on Smoking Cessation at 26 Weeks: A Randomized Clinical Trial. JAMA. 2016;315(4):371.

B. Efron, R. Tibshirani. Bootstrap Methods for Standard Errors, Confidence Intervals, and Other Measures of Statistical Accuracy. Statistical Science. 1986;1(1):54-75.