title: Estimation and Inference for the Mediation Effect in a Time-varying Mediation Model
authors: Xizhen Cai Donna L. Coffman Megan E. Piper Runze Li
text: 
  
  
  Traditional mediation analysis typically examines the relations among an intervention, a time-invariant mediator, and a time-invariant outcome variable. Although there may be a direct effect of the intervention on the outcome, there is a need to understand the process by which the intervention affects the outcome (i.e. the indirect effect through the mediator). This
  indirect effect is frequently assumed to be time-invariant. With improvements in data collection technology, it is possible to obtain repeated assessments over time resulting in intensive longitudinal data. This calls for an extension of traditional mediation analysis to incorporate time-varying variables as well as time-varying effects. In this paper, we focus on estimation and inference for the time-varying mediation model, which allows mediation effects to vary as a function of time. We propose a two-step approach to estimate the time-varying mediation effect. Moreover, we use a simulation based approach to derive the corresponding point-wise confidence band for the time-varying mediation effect. Simulation studies show that the proposed procedures perform well when comparing the confidence band and the true underlying model. We further apply the proposed model and the statistical inference procedure to real-world data collected from a smoking cessation study. 
  
  
  
  
  
  SmartSON:A Smart contract driven incentive management framework for Self-Organizing Networks
      Abdullah Yousafzai, and Choong Seon Hong, Senior Member, IEEE
         
          
  Abdullah Yousafzai is with the Networking Lab, Kyunghee University, Republic of Korea.
  
  
  
  E-mail: see dr.abdullah@khu.ac.kr
  Choong Seon Hong is Professor at the School of Electronics and Information, Kyung Hee University, Republic of Korea.
  Manuscript received August 26, 2020.
      December 30, 2023
  ========================================================================================================================================================================================================================================================================================================================================================================
  
  
  
  
  
  § INTRODUCTION
  
  Developments in mobile and wearable device technology have enabled the collection of intensive longitudinal data <cit.>, such as ecological momentary assessment (EMA) <cit.>, in which data on variables such as craving, withdrawal symptoms, or stress, are measured in real-time, real-world contexts. EMA is particularly useful in health behavior change studies, for example, smoking cessation studies. Often, the variables that are collected during EMA are variables that are targets of a behavior change intervention and are also thought to affect the health outcomes of interest. In other words, they are mediators, variables that lie on the pathway between the intervention and the outcome.
  
  
  
  
  As the collection of data using EMA has grown, so have methods for analyzing and making the most of the temporal density of measurements, such as the mixed-effects location scale model <cit.> and the time-varying effect model <cit.>. However, there have been very few proposed methods for assessing mediation using this type of intensively measured data. This paper aims to propose an approach to mediation in which both the mediator and outcome are collected using EMA, and therefore can be examined as they change over time. Specifically, we first proposed a two-step approach to estimate the time-varying mediation effect. 
  To make statistical inference on the time-varying mediation effect, we develop a simulation-based approach to derive the corresponding point-wise 
  confidence band for the time-varying mediation effect.
  
  The rest of this paper is organized as follows. In section <ref>, we present relevant background material on varying-coefficient models and the proposed model for time-varying mediation, including estimation and bootstrap inference. In section <ref>, we present simulation studies to examine the performance of the bootstrap confidence intervals. In section <ref>, we apply the proposed methods to data from a smoking cessation intervention study. In section <ref>, we discuss limitations, future directions, and conclusions.
  
  
  
  § VARYING-COEFFICIENT MODELS AND PROPOSED EXTENSION TO MEDIATION
   
  Time-varying coefficient models <cit.> have been used to model time-varying effects of an independent variable on a dependent variable <cit.>. These are essentially varying-coefficient models <cit.> applied to intensive longitudinal data. For each individual, i, the independent variable and the outcome variable are measured at multiple time points {t_ij, j=1,2,…, T_i}. The data collected are 
  
      {t_ij, X_i(t_ij), Y_i(t_ij)},      i=1,2,…,n,    j=1,2,…,T_i
  
  and the model can be written as 
  
      Y_i(t_ij)=β_0(t_ij)+X_i(t_ij)β_1(t_ij)+ ϵ_i(t_ij),
  
  where β_0(t) and β_1(t) are smoothing functions of time, and therefore are called the time-varying coefficient functions, and ϵ(t) is a zero-mean stochastic process with covariance function, γ(s,t). Not only are the effects (i.e., coefficients) of the predictor variables time-varying, but the values of the variables themselves can also change over time. Different types of estimation procedures have been well summarized <cit.>. There are essentially two estimation approaches for time-varying effect models: splines and local smoothing methods. In this paper, we focus on 
  local smoothing methods, which locally approximate coefficient functions by linear or polynomial functions <cit.>.
  
  <cit.> proposed a powerful two-step procedure that uses kernel methods to estimate the time-varying coefficients and their corresponding standard errors. Both simulations and real data applications showed the efficiency of their method over other previous proposals. Later,  <cit.> extended the model and its estimation to the case where the outcome variable depends not only on the current but also past values of the independent variables, 
  
      Y_ij(t_ij)=β_0(t_ij)+∑_r=1^pβ_r(t_ij)X_i,j-q-(r-1)(t_i,j-q-(r-1))+ϵ_ij(t_ij),
   where p is the number of past time points that are assumed to affect the current response and q incorporates the possible time lag. The estimation is done through a variant of the two-step procedure and has been shown to perform well. This two-step procedure <cit.> provides an important foundation for our proposed estimation procedure for time-varying mediation effects which combines the traditional linear mediation model estimation procedure and local polynomial smoothing.
  
  Although time-varying coefficient models are relatively common for examining the time-varying effect of an independent variable on a dependent variable, relatively little work has examined time-varying effects for mediation.  <cit.> first introduced functional (or time-varying) mediation effects in which the independent and dependent variables were measured at a single point in time but the mediator was measured intensively over time using fMRI. More recently, <cit.> proposed a mediation g-formula, which allows time-varying treatments, time-varying mediators, and an end-of-study point outcome. They mention the possibility of time-varying effects, but their study does not directly address them. In our application to a smoking cessation study, the mediator and outcome are both measured repeatedly over time (i.e., time-varying mediator and time-varying outcome) and the independent variable is random assignment to the intervention (not a time-varying treatment). Thus, neither of these previous models apply directly to our smoking cessation study.
  
  
  
  
  Traditional methods of assessing mediation, shown in Figure <ref>, generally specify the direct effect (i.e., that does not go through the mediator) as γ, and the indirect or mediated effect as the product of paths α (effect of intervention on mediator) and β (effect of mediator on outcome). The standard error of this product term, αβ, is obtained either asymptotically <cit.> or via bootstrap procedures in order to test the statistical significance of the mediated effect. Several simulation studies have shown bootstrap standard errors to be superior, especially in smaller samples because αβ may not be normally distributed <cit.>.
  
  As in traditional mediation analysis with time-invariant effects, the time-varying mediation effect is the product of two effects, but in this case, both effects are time-varying. That is, the two effects are no longer single numbers such as α and β; rather, they are functions of time, and the product term is also a function of time. Figure <ref> is extended in Figure <ref> to include time-varying effects.
  
  
  
  
  
  
  
  
  
  
  
  In this paper, we propose to estimate the time-varying mediation model by extending the two-step approach <cit.>, followed by bootstrapping to obtain confidence intervals for the indirect effect (i.e., the effect of the intervention on the outcome through the mediator). Mediation is inherently about causal pathways - the intervention changes the mediator, which in turn has an effect on the outcome. To infer causality, we will need to assume that there are no unmeasured confounders, additivity (no interactions or non-linearities), and no time-varying confounders. In addition, we assume temporal ordering such that the intervention occurs before the mediator which occurs before the outcome. These are the standard assumptions needed for a linear structural equation model to estimate a “causal" effect. These are strong assumptions, to which we return when discussing future directions.
  
  
  
  
  
   §.§ The proposed model
  
  Extending the mediation framework in Figure <ref> to take advantage of the temporal density of intensive longitudinal data allows us to estimate time-varying effects as shown in the dynamic mediation diagram in Figure <ref>. In this model, we consider the measurement timing of the variables consistent with modeling mediation as a process that unfolds over time (i.e., intervention must precede change on the mediator, and mediator must precede change on the outcome). The intervention or independent variable, denoted X, is time-invariant, and assigned at time t_0. Across time, the effect of X on the value of the mediator M at any time t>t_0 (i.e., M(t)) is denoted by α(t). The value of the outcome variable Y at time t (i.e., Y(t)) is affected by the value of the mediator at a small window before time t, i.e., M(t-Δ t). Here Δ t is a small constant which represents the time-lag of the effect of the mediator on the outcome. More discussion on the choice of Δ t will be presented shortly.
  
  The diagram in Figure <ref> leads to the following time-varying mediation model: 
  
      M(t)    =    α_0 (t)+I (t ≥ t_0) α (t)   X+ ϵ^M(t),
      
      Y(t)     =   β_0 (t)+I (t ≥ t_0) {γ(t)  X + β(t) M(t-Δ t) }+ ϵ^Y(t),
  
  where ϵ^M (t) and ϵ^Y (t) are both zero-mean stochastic processes. The time-varying mediation effect of interest is α(t-Δ t)β(t). Suppose there are repeated measurements of N subjects at multiple time points {t_ij}, then the observed data are 
  
      {X_i, (t_ij, M_i(t_ij), Y_i(t_ij))},    i=1,2,…, N,     j=1,2,…,T_i,
  
  and the model is
  
      M_i(t_ij)    =    α_0 (t_ij)+I (t_ij≥ t_0) α (t_ij)   X_i+ ϵ^M_i (t_ij)
      
      Y_i(t_ij)     =   β_0 (t_ij)+I (t_ij≥ t_0) {γ (t_ij)  X_i + β(t_ij) M_i(t_ij -Δ t) }+ ϵ^Y_i (t_ij).
  
  Note that all effects of the intervention, X_i, are controlled by a post-intervention indicator, I (t_ij≥ t_0), because the intervention is assigned at time t_0. Since the mediation (i.e., indirect) effect is of primary interest, we focus on the time points after the intervention is assigned, thus, the indicator term can be dropped and the final model is, 
  
      M_i(t_ij)    =    α_0 (t_ij)+ α (t_ij)   X_i+ ϵ^M_i (t_ij) 
      
      Y_i(t_ij)     =   β_0 (t_ij)+ γ (t_ij)  X_i + β (t_ij) M_i(t_ij-Δ t)+  ϵ^Y_i (t_ij).
  
  
  To use all information in the collected data effectively, the model in equation (<ref>) indicates that the window size Δ t should be smaller than the distance between the current time point t_ij and the previous time point t_i,j-1. If the observation times are not equally spaced, the window size can be chosen to be a value that is smaller than the smallest distance. In our simulation and applied analysis, we assume that the value of the mediator at time t_ij-Δ t can be substituted by its value at the previous time point M_i(t_i,j-1). Additionally, models (<ref>)
  and (<ref>) assume that there are only two intervention groups (e.g., treatment versus control) such that X_i is a binary indicator of the treatment condition, and the time-varying effect α(·) is the effect of the treatment as compared to the control group. For more than two intervention groups, the proposed model can be easily extended by adding more indicator variables (see the smoking cessation study in Section <ref> as an example). Without loss of generality, we present the following proposed estimation procedure and bootstrap inference for the models in equations (<ref>) and (<ref>). 
  
  
  
   §.§ Estimation of the time-varying mediation effect
  
  
  We could estimate the time-varying effects in model equations (<ref>)  and  (<ref>)  separately using the two-step estimation procedure <cit.>. Here we propose a variant of that approach to estimate them simultaneously. 
  
  Let {t_1, t_2, …, t_T} be the distinct time points when the data are measured.  For any fixed time point t_j∈{t_2, …, t_T}, we observe complete data from N_j subjects (N_j does not necessarily equal N). Then for any individual i at this fixed time point t_j, the observed data are 
  
      (X_i, M_ij, Y_ij),    i=1,2,…, N_j,
  
  where M_ij=M_i(t_ij) and Y_ij=Y_i(t_ij). Similar to the first step of the two-step procedure <cit.>, at any fixed time t_j, model equations (<ref>)  and  (<ref>)  become the traditional linear mediation model. So we can estimate the value of the varying coefficient functions α(t_j), β(t_j), and γ(t_j), which are treated as three parameters rather than three functions, by the least squares method, namely, by solving the following two optimization problems,
  
      min_α∑_i=1^N_j (M_ij-α(t_j) X_i)^2 min_β, γ∑_i=1^N_j (Y_ij-γ(t_j) X_i - β(t_j) M_i,j-1)^2.
  
  Suppose the outcome and independent variables are mean centered or standardized at the fixed time point so that we can safely drop the intercept terms. To derive a joint distribution of the estimated coefficients, we propose to combine the two least squares problems together to create a new least squares problem given as
  
      min_α, β, γ{∑_i=1^N_j (M_ij-α (t_j) X_i)^2 + ∑_i=1^N_j (Y_ij-γ (t_j) X_i - β (t_j) M_i,j-1)^2}
         ⇔   min_δ∑_i=1^2N_j (Y^*_ij-δ^⊤(t_j)X^*_ij)^2
  
  where δ(t_j)=(α(t_j), γ(t_j), β(t_j))^⊤, and Y^*_ij, X^* _ij in matrix forms are, 
  
      𝐘^* _j = ([     M_1j;     M_2j;        ⋮;  M_N_j,j;     Y_1j;     Y_2j;        ⋮; Y_ N_j,j ])_2N_j× 1,
      𝐗^* _j = ([       X_1         0         0;       X_2         0         0;         ⋮         ⋮         ⋮;     X_N_j         0         0;         0       X_1   M_1,j-1;         0       X_2   M_2,j-1;         ⋮         ⋮         ⋮;         0     X_N_j M_N_j,j-1 ])_2N_j× 3
  
  Denote the solution to the least squares problem in (<ref>) as d (t_j) = (a(t_j), c(t_j),  b(t_j))^⊤,  i.e. d(t_j) is an estimate of δ(t_j) and is a 3× (T-1) dimensional vector, which includes values of the estimated time-varying coefficient functions at all time points, 
  
      d=(a(t_2),c(t_2),b(t_2), a(t_3),c(t_3),b(t_3), ⋯, a(t_T),c(t_T),b(t_T))^⊤
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  Similar to the second step of the two-step procedure <cit.>, the smoothed coefficient functions α̂ and β̂ in model equations (<ref>) and  (<ref>) are further calculated by local polynomial regression using (a(t_j), b(t_j)) as, 
  
      α̂ (t- Δ t)    =   ∑ _l=2^T w(t_l, t- Δ t) a (t_l)
      β̂ (t)    =   ∑ _l=2^T w(t_l, t) b (t_l)
  
  where w(t_j, t) can be weights from any linear smoothing techniques. Here, we use local polynomial weights. Then the desired mediation effect is 
  
      α̂ (t- Δ t)β̂ (t)  = {∑ _l=2^T w(t_l, t- Δ t) a (t_l) }{∑ _l=2^T w(t_l, t) b (t_l)},
  
  and each part can be rewritten as linear combinations of d: 
  
      α̂ (t- Δ t)β̂ (t)  = (𝐰_a^T𝐝)(𝐰_b^T𝐝),
  
  where
  
  
      𝐰_𝐚=
      ([ w(t_2, t- Δ t);              0;              0; w(t_3, t- Δ t);              0;              0;              ⋮; w(t_T, t- Δ t);              0;              0;                ]), 𝐰_𝐛=([         0;         0; w(t_2, t);         0;         0; w(t_3, t);         ⋮;         0;         0; w(t_T, t);           ]).
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
   §.§ Estimating point-wise confidence interval of mediation effect through bootstrap
  
  
  
  
  To identify a statistically significant mediation effect in the time-varying setting, we consider the following hypothesis:
  
  
      H_0:α (t- Δ t)β (t) = 0t
         vs.    H_A: .
  
  Since the distribution of the mediation effect is not necessarily normal, we use a bootstrap approach to construct confidence intervals. Specifically, the lower and upper bounds of the 1-α% confidence interval are taken to be the corresponding lower (α/2) and upper (1-α/2) percentiles of the distribution of the estimated mediation effect from a large number of bootstrapped samples <cit.>.
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  For any fixed time t, the above bootstrap percentile method creates a point-wise confidence interval for the mediation effect at that t. Connecting all confidence intervals yields a confidence band, but this is different from a simultaneous confidence band throughout the entire time interval, since  the nominal confidence level is only satisfied at each fixed time point t. We return to this point in Section <ref>. The estimation procedure and bootstrapped confidence intervals are implemented in an R package, tvmediation, that is available on GitHub.
  
  
  
  
  
  
  
  
  
  
  
  § SIMULATION STUDIES
  
  To examine the performance of the proposed point-wise confidence interval, we consider the following two simulation models,
  
  
    * α_1(t)=10+12t^3, γ(t)=-20-18t,
  β(t)=50+150t^2, γ(s,t)=15exp(-0.3|s-t|)
  
  
    * α(t)=15+8.7sin(0.5π t),
  γ(t)=4-17(t-1/2)^2,
  β(t)=1+2t^2+11.3(1-t)^3,
   γ(s,t)=15exp(-0.3|s-t|)
  
  
  The first model includes polynomial functions of different orders, and the second model incorporates a sin function to increase the complexity of the mediation effect. The two models are similar to those in <cit.> and <cit.>
  
  Without loss of generality, observation times are generated as 50 equally spaced time points between 0 and 1. And the time lag Δ t is chosen to be half of the length between any two consecutive time points. To generate the simulated data, we first randomly assign intervention and control group, each with probability of 0.5. The error term is generated from multivariate normal distribution with mean zero and covariance γ(s, t). The value of the mediator and the outcome variables are generated according to equations (<ref>) and (<ref>). In the second step of the estimation procedure, local linear regression is used and the bandwidth is chosen by the rule of thumb formula in  section 4.2 of <cit.>. 
  
  We consider three simulation settings with N=100, 200, and 500, separately. To verify the nominal level for 95% confidence intervals, we calculated the coverage rate of the proposed point-wise confidence interval at t=0.2, 0.4, 0.6, and 0.8, separately. Table <ref> summarizes the results based on 500 simulation replications. 
  
  
  
  
  
  
  
  
  Except for a few settings, the coverage rates are all near a 95% confidence level at these points. We also evaluated the performance of the estimated time-varying mediation effect by the mean absolute deviation error (MADE) and weighted average squared error (WASE) <cit.>, defined as follows,
  
  
      MADE =  (4T)^-1∑_j=1^T|η(t_j)-η̂(t_j)|/(η),   
      WASE = (4T)^-1∑_j=1^T{η(t_j)-η̂(t_j)}^2/^2(η)
  
  
  where η(t)=α (t- Δ t)β (t) is the time-varying mediation effect. 
  Figure <ref> presents boxplots of these measurements for the two models. Not surprisingly, both of them show similar patterns with different sample sizes and models.  Specifically, both MADE and WASE decrease as sample size increases for a particular model, and the error for model ii is slightly higher than that of model i.
  
  
  
  To present a typical fit of the proposed procedure, we selected the simulation sample with MADE closest to the median value among all 500 replications.  The estimated time-varying mediation effect and the  corresponding confidence intervals are plotted in Figure <ref>, as compared to the true effect. 
  
  
  The red solid line is the true time-varying mediation effect, and the blue solid line is the estimated effect. For both models, the estimated effect is close to the true underlying effect. The blue dashed lines are the limits of the point-wise confidence band estimated by the proposed method. For both models, the width of the confidence band is not constant throughout the whole time range, but at each time, the true effect is fully contained in the confidence band. As the sample size increases, the confidence band becomes narrower. 
  
  
  
  
  § APPLICATION: THE WISCONSIN SMOKERS' HEALTH STUDY 2
   
  
  We applied the proposed method to conduct an empirical analysis of data collected from a smoking cessation study, the Wisconsin Smoker's Health Study 2 <cit.>, which used EMA to assess negative affect and cessation fatigue during a smoking cessation trial. The study was a randomized comparative efficacy trial <cit.> directly comparing the two most effective smoking cessation therapies (varenicline vs. combination nicotine replacement therapy - nicotine patch + nicotine mini-lozenge) with one another and with an active comparator treatment (nicotine patch only). In total, 1086 smokers recruited from Madison and Milwaukee, WI were randomly assigned to one of the three 12-week pharmacotherapies. Participants completed one morning EMA prompt, and one evening prompt every day for one week prior to the quit day and for two weeks after the quit day and then every other day for the remaining two weeks of the EMA period (i.e., total of one week pre-quit and four weeks post-quit). Thus, there are 14 EMAs prior to the quit day and 42 after the quit day.  The goal of our empirical analysis is to examine whether the intervention has an effect on cessation fatigue that is mediated by negative affect (see Figure <ref> for the mediation diagram).  
  
  
  
  Cessation fatigue, defined as tiredness of trying to quit smoking <cit.>, and negative affect, measured by asking participants if they were in a negative mood in the last 15 minutes, were both measured on 7-point Likert scales. Previous studies have found that negative affect and cessation fatigue are positively related and related to cessation failure <cit.>.
  
  
  
  
  We use data from the 42 EMAs after the quit day. Both the time-varying outcome, cessation fatigue, and time-varying mediator, negative affect, are assessed at each EMA prompt. Unlike in the simulation studies, and as is common in most empirical studies, especially with wearable and mobile devices, there are intermittent missing values in the data. Excluding individuals with no data at all, we have 1047 individuals in total, and the observed data are
  
      {_i, _i,(t_ij, _ij, _ij)},   i=1,2,…,1047,    j=1,2,…,42.
  
  There are two indicator variables for the intervention: Varen_i indicates assignment to the varenicline group and cNRT_i indicates assignment to the combination nicotine replacement therapy group. The nicotine patch only condition is the reference group as it is considered the standard of care. Additionally, the observation times are not equally spaced (i.e. everyday for the first two weeks and every other day for the remaining two weeks). The previous proposed model can be modified to incorporate the additional intervention condition as follows:
  
      NA_i(t_ij)    =    α_0 (t_ij)+ α_1(t_ij)   Varen_i+ α_2(t_ij)   cNRT_i + ϵ^NA _i (t_ij) 
      
      CF_i(t_ij)     =   β_0 (t_ij)+ γ_1 (t_ij) Varen_i + γ_2(t_ij) cNRT_i + β(t_ij) NA_i(t_ij-Δ t)+  ϵ^CF_i (t_ij) .
  
  Using the proposed method, the estimated mediation effects, α_1(t-Δ t)β(t) and α_2(t-Δ t)β(t), with the corresponding confidence bands are presented in Figure <ref>. 
  
  
  
  As compared to the nicotine patch only, the effect of varenicline on cessation fatigue that is mediated by negative affect becomes more negative shortly after quitting. More specifically, the magnitude of this negative mediation effect increases quickly in the first week after quit day, has a slower decrease in the second week, becomes stable in week 3, and slightly increases in week 4. The pattern for the effect of combination nicotine replacement therapy (vs. the nicotine patch only) was similar although the initial increase in the magnitude was not as pronounced but the mediation effect was statistically significant throughout the four weeks post-quit. Examination of the two time-varying effects that make up the mediation effect are also informative. Examination of the bottom row of Figure <ref> shows that varenicline (vs. the nicotine patch only) has a negative effect that becomes stronger over the course of the first week. This effect then begins to diminish during the following three weeks. In contrast, the strong negative effect of combination nicotine replacement therapy (vs. the nicotine patch only) on negative affect is apparent at the beginning of week one but then, similar to the varenicline group, diminishes during the following three weeks. Examination of the time-varying effect of negative affect on cessation fatigue reveals that there is a strong positive relationship (i.e., more negative affect results in more cessation fatigue) initially during the first week that diminishes over the following three weeks.
  
  Additionally, the mediation effects for both interventions (compared to the nicotine patch only group) are time-varying. Compared to using nicotine patch only, the effects of varenicline and combination nicotine replacement therapy on cessation fatigue, as mediated by negative affect, are not only negative, but also, time-varying for the four weeks post quit day. Both mediation effects have a narrow confidence band and thus, we can rule out a constant mediation effect over time because we would not be able to fit a flat line over time within the confidence interval.
  
  
  
  § DISCUSSION
  
  In this paper, we have described a model for assessing mediation in the context of intensive longitudinal data in which both the mediator and outcome variables are time-varying. This model allows for estimation of time-varying mediation effects. Intensive longitudinal data often arise from the collection of EMA data but may also arise from the collection of data from mobile devices, such as wrist-worn or hip-worn accelerometers. The temporal density of these data allow for more nuanced research questions that cannot be addressed by, for example, averaging over the EMA data and/or assuming that the mediated effect does not vary as a function of time. By allowing mediated effects to vary as a function of time, research questions such as the timing of important mediation effects can be assessed. Thus, our approach may prove useful to other researchers who wish to conduct mediation analysis in the context of intensive longitudinal data. 
  
  The simulation study showed that the proposed bootstrap pointwise confidence intervals contained the true time-varying mediation effect and that the estimated time-varying mediation effect was close to the true time-varying mediation effect. We then applied our approach to examine the mediation effect of three smoking cessation treatments (i.e., varenicline, combination nicotine replacement therapy, and nicotine patch only) on cessation fatigue via negative affect. The results indicated that the mediated effect 1) did indeed vary as a function of time, 2) was statistically different from zero throughout the four weeks post-quit day, and 3) that the effect was strongest in the first week post-quit for the varenicline group (vs. nicotine patch only). That is, the varenicline group experienced decreased negative affect during the first week, leading to decreased cessation fatigue. Interestingly, the effect was also strongest in the first week for the patch only group and the effect was immediate whereas the varenicline effect improved over the first half of the first week. The mediated effect for both treatments, compared to the patch alone, appeared to dissipate over the course of the first four weeks of the quit attempt. This information may lead to modifications and/or adaptations of the intervention to, for example, implement a behavioral component to address negative affect, with a specific focus on reducing negative affect in the first week of the quit attempt. This information would not have been evident had we assumed that the mediation effect was invariant across the four week post-quit period.
  
  There are several limitations of the current approach. First, the proposed method only constructs a point-wise confidence interval. For inference at a fixed time point, a point-wise confidence interval is useful. However, a simultaneous confidence band is needed to make inferences over the entire time span. Thus, an obvious future direction is developing and estimating a simultaneous confidence band. Second, although the current approach does not require observations from all participants at all time points <cit.>, the algorithm will not work if  (𝐗^*_j) <3 (or <d in the general case). In such cases, one can implement the four methods discussed in <cit.> (see their Remark 1).
  
  A third limitation is that as with all mediation analyses, we need to make strong assumptions regarding no unmeasured confounding, temporal order, and additivity as mentioned previously. In our particular application, individuals were randomly assigned to the smoking cessation treatments; however, they are not randomly assigned to the mediator and therefore, there may be confounders of the mediator and the outcome. In addition, due to the intensive longitudinal nature of the study, we cannot rule out the possibility of time-varying confounding. Of particular concern is the possibility of time-varying confounders of the mediator and outcome that have themselves been affected by the smoking cessation treatments. In future work, we will propose sensitivity analyses to address potential violations of these assumptions. Additional future work will also address binary outcomes, for example, daily smoking, and count outcomes, for example, daily number of cigarettes smoked.
  
  
  
  In conclusion, we have presented a model for estimating time-varying mediation effects which builds on previous work <cit.> to allow a time-varying outcome as well as a time-varying mediator. We also presented a method for obtaining point-wise confidence intervals for the product of two time-varying coefficient functions (i.e., a time-varying mediation effect), evaluated its feasibility in a small simulation study, and applied the method to evaluate the time-varying mediation effects of three pharmacotherapy smoking cessation interventions. We have implemented the estimation and bootstrap procedure in a user-friendly R package, tvmediation, available on GitHub. Although we cannot share the actual data, the R package contains data simulated to mimic the real data along with tutorials on how to use the functions to fit the model described above. We believe that this approach will be useful to those collecting frequent data using mobile devices for self-reported EMA and who wish to examine mediation effects.
  
  *
  rss
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  