Trajectories of Response to Dorsolateral Prefrontal rTMS in Major Depression: A THREE-D Study

Treatment-resistant major depressive disorder remits for a minority of patients treated with pharmacotherapy or psychotherapy (1). Novel treatments such as repetitive transcranial magnetic stimulation (rTMS) have been shown to be effective for achieving remission in treatment-resistant major depression (2). However, rTMS outcomes are heterogeneous (3, 4), with some studies suggesting subpopulations of fast and slow responders (5).

Previous studies have identified distinct depression response trajectories for pharmacotherapy and psychotherapy (6–8); however, no similar studies have identified response trajectories for rTMS. Understanding the variation in response trajectories with rTMS is important for three reasons. First, it could enable individualized prediction of fast response, slow response, and nonresponse to rTMS using baseline clinical characteristics. Second, it could help determine whether the optimal treatment duration to achieve remission varies among individuals, as does pharmacotherapy (9). Third, it could offer clues regarding the biological heterogeneity of major depression itself, based on pace of response (10).

Using data from a recent rTMS trial, we applied group-based trajectory modeling techniques to identify distinct response trajectories (11). Our primary objective was to characterize the number and pattern of distinct longitudinal response trajectories for adults with treatment-resistant major depression over an acute course of rTMS. Our exploratory objective was to determine baseline clinical characteristics associated with the identified response trajectories.

Methods

Participants

This secondary analysis used data from the THREE-D study (11), a randomized noninferiority trial comparing two rTMS protocols applied to the left dorsolateral prefrontal cortex: conventional high-frequency left (HFL) or intermittent theta burst (iTBS) stimulation. The study was conducted at three Canadian academic hospitals from September 2013 to October 2016: the Centre for Addiction and Mental Health, the University Health Network, and the University of British Columbia. The study was approved by institutional ethics boards at all sites, and participants provided written informed consent.

Participants were outpatients between the ages of 18 and 65 with a diagnosis of unipolar major depressive disorder, confirmed with the Mini International Neuropsychiatric Interview (12). Inclusion criteria were current major depressive episode, with a score ≥18 on the 17-item Hamilton Depression Rating Scale (HAM-D) (13); lack of response to at least one adequate or two inadequate antidepressant trials during the current episode, as assessed by the Antidepressant Treatment History Form (14); and, if taking psychotropic medications, dosages had not been increased for 4 weeks before starting treatment. Exclusion criteria were substance dependence or abuse <3 months before study entry; any unstable medical or neurological illness; acute suicidality; a Mini International Neuropsychiatric Interview diagnosis of bipolar I or II disorder, a primary psychotic disorder, or psychotic symptoms in the current episode; a primary diagnosis of obsessive-compulsive disorder, posttraumatic stress disorder, an anxiety disorder, or a personality disorder; any contraindication to rTMS (history of seizures; intracranial implant); a lifetime history of failure to respond to an adequate course of ECT; previous rTMS treatment; current treatment with any anticonvulsant, or with lorazepam at >2 mg/day; pregnancy; significant laboratory test abnormalities; and failure of more than three adequate antidepressant trials (defined by a score >3 on the Antidepressant Treatment History Form) in the current episode.

Study Design

Participants were randomly assigned in a 1:1 ratio to receive either HFL or iTBS treatment, administered 5 days/week over 4–6 weeks for 20–30 treatments while continuing their psychotropic medications unchanged for the study duration. All participants received an initial course of 20 daily treatments, and participants who achieved a reduction ≥30% from baseline in HAM-D score, but not remission, received an additional 10 treatments over 2 weeks to optimize treatment response and durability (15). Participants who did not achieve a 30% reduction by week 4 exited the study. Participants were withdrawn early if their depression scores (measured by the HAM-D) were >25% higher than at baseline on two consecutive assessments, if they developed significant suicidal ideation, or if they attempted suicide. Participants who missed treatment sessions were rescheduled to achieve the intended course length; however, participants who missed four consecutive treatments were withdrawn. Randomization of participants was stratified by degree of medication resistance (dichotomized as more than one versus one or fewer adequate medication trials without response) and was done using a randomly permuted block method with a random number generator. While the design did not allow the rTMS technician or patient to be blind to treatment allocation, outcome assessors were blind to treatment allocation.

rTMS Procedure

Before treatment, all participants underwent anatomical MRI brain scanning. rTMS treatments were delivered using real-time MRI-guided neuronavigation with a Visor2 system (Advanced Neuro Therapeutics, Madison, Wisc.) to optimize coil positioning. The left dorsolateral prefrontal cortex was targeted using the Montreal Neurological Institute’s MNI-152 stereotaxic coordinates (x, y, z: −38, 44, 26) (16). The device used was a MagPro X100/R30 stimulator equipped with a B70 fluid-cooled coil (MagVenture, Farum, Denmark).

The resting motor threshold was determined by visual observation according to published guidelines (17). HFL used treatment settings approved by the U.S. Food and Drug Administration (120% resting motor threshold, 10 Hz, 4 seconds on, 26 seconds off, 3000 pulses/session over 37.5 minutes) (18, 19). iTBS was delivered to the same site with the same intensity but used a different stimulation pattern (triplet 50-Hz bursts, repeated at 5 Hz, 2 seconds on, 8 seconds off, 600 pulses per session over 3 minutes) (20). Further details are provided in the original THREE-D report (11).

Measures

Clinical assessments and prognostic factors.

Clinical assessments were completed by trained research assistants blind to treatment allocation. All characteristics potentially associated with response trajectories were measured at baseline. Validated psychometric scales were used to measure clinical characteristics over time. Depression severity was measured at baseline and then weekly until trial completion with clinician-rated (the HAM-D and the Inventory of Depressive Symptomatology [IDS]) (21) and self-rated (the Quick Inventory of Depressive Symptomatology–Self-Rated [QIDS-SR]) (22) instruments. Response was defined as a reduction ≥50% from baseline in HAM-D score, and remission was defined as a HAM-D score <8, as in the primary THREE-D analysis (11). Supplementary measures were assessed at baseline and at trial completion: anxiety (Brief Symptom Inventory) (23), functional disability (Sheehan Disability Scale) (24), degree of enjoyment and satisfaction experienced in daily functioning (Quality of Life Enjoyment and Satisfaction Questionnaire) (25), and mental well-being (Warwick-Edinburgh Mental Well-Being Scale) (26).

Outcomes.

The primary outcome was to classify study participants into distinct response trajectories, using HAM-D scores, during a course of rTMS. The exploratory outcome was to determine clinical characteristics associated with response trajectories.

Statistical Analysis

All statistical analyses were conducted in SAS, version 9.4 (SAS Institute, Cary, N.C.). For our primary analysis, intended to determine longitudinal response trajectories associated with rTMS treatment, we used a semiparametric group-based trajectory modeling strategy to classify study participants into subgroups based on identifying heterogeneous longitudinal polynomial trajectories. This was implemented via the SAS procedure PROC TRAJ (27). The outcome variable—HAM-D score—is a normally distributed psychometric scale for which we estimated the error structure as a censored normal distribution. We determined the optimal number of response trajectories in the model and optimal polynomial degree in each trajectory, using the Bayesian information criterion (BIC). The BIC measures improvement in model fit gained by inclusion of additional groups or shape parameters, but it also penalizes added complexity. The BIC log Bayes factor approximation, defined as 2×ΔBIC (where ΔBIC is the BIC difference between a more complex and a less complex model), has been shown to be an acceptable approximation to the log Bayes factor criterion (28) and was used to determine the number of response trajectories that best fit the observed data. A log Bayes factor approximation >10 was used as the criterion for favoring the more complex model (27).

We first determined the best-fitting number of response trajectories with the maximum degree of the fitted polynomial fixed at cubic, since previous work has demonstrated that depressive symptoms during treatment generally follow linear, quadratic, or cubic trajectories (8, 29). Next, we determined the polynomial degree in each trajectory by systematically reducing the polynomial degree for each trajectory with the smallest point estimate until all trajectories consisted of linear polynomial degrees. The combination of linear, quadratic, and cubic polynomials that best explained the observed response trajectories (lowest BIC), was considered the best-fitting model.

We assessed model fit by calculating the average posterior probability of group membership (70% minimum for each group), determining the percentage of the total sample within each trajectory (5% minimum for each group), and calculating the odds of correct classification (>5 considered adequate). Because participant data for weeks 5 and 6 were missing not at random and potentially nonignorable, we determined group membership using participant data to week 4. However, because response trajectories may identify protocol-defined early study completers at week 4 (early remission or <30% change in HAM-D score from baseline), we compared response trajectory classifications using week 4 data with classifications using week 6 data. We found that there was an identical number of trajectories and good interrater agreement (kappa=0.74) for assigning response trajectories when week 4 or week 6 data were used. Therefore, we extended participant depressive symptom data to week 6 for descriptive purposes. Furthermore, because this trial used two different rTMS treatment techniques, we repeated the model creation process for each treatment technique separately. We also assigned treatment technique as a covariate in the PROC TRAJ environment to determine whether this was a significant predictor of group membership assignment.

We also completed two secondary analyses: a categorical comparison of HAM-D remission (HAM-D score <8) and response rates (HAM-D change ≥50% from baseline) for weeks 1 to 4 and at trial completion between response trajectories, and a sensitivity analysis determining the response trajectory of participants who received treatment but violated study inclusion criteria.

For the exploratory objective of identifying clinical features associated with response trajectory, we used weighted multinomial multivariate logistic regression to determine associations between potential characteristics (baseline clinical and demographic characteristics and treatment characteristics) and response trajectories. The regressions were weighted by the probability of group membership to account for measurement error introduced by the uncertainty of group membership. The reference group for the regression analysis was chosen a priori to be the response trajectory with the largest membership.

We used the following covariate selection process to determine characteristics associated with response trajectories. First, we removed any clinical, demographic, or treatment variables that measured identical constructs. Second, we assessed covariates for multicollinearity by assessing the variance inflation factor and removing covariates to ensure that the variance inflation factor was <4 for all variables. Third, we assessed covariates for sparse data. Multilevel categorical covariates were collapsed where appropriate. If covariates could not be collapsed, then they were excluded if the expected frequency was <5 within each response trajectory. Fourth, the remaining covariates were selected using a previously described procedure (30): We derived 1,000 bootstrap samples, to which we applied automated backward stepwise selection with variable elimination when p was <0.05. Final covariates were selected if they were included in >75% of the bootstrap samples, to ensure that our selection of covariates would be more stable and robust than from the single apparent sample (31). All bootstrapped regressions were forced to include treatment technique (HFL or iTBS) in the regression, to account for any association with treatment technique. The identified predictors associated with treatment response were then included in a multinomial regression model using the data set to determine their independent association with group membership as measured by odds ratios and corresponding 95% confidence intervals. Model discrimination and fit were assessed in independent logistic models for each response trajectory and are reported using the c-statistic and the Hosmer-Lemeshow test, respectively. Statistical tests were two-tailed with alpha set to 0.05, with the exception of comparing response and remission rates between trajectories, where alpha was set to 0.01 to account for multiple comparisons.

Results

Response Trajectories

A total of 414 participants underwent randomized treatment assignment, with 26 excluded: two before receiving treatment and 24 who received treatment but were subsequently found to have violated study inclusion criteria. Therefore, the analytic cohort for this study consisted of 388 participants who received at least one rTMS treatment. For the primary outcome, we found that four distinct response trajectories adequately fit the observed data (Table 1) and were described by a combination of cubic, quadratic, and linear polynomial components. The longitudinal course of depressive symptoms for each of the four trajectories is depicted in Figure 1. The four response trajectories each had a distinct pattern, which we labeled “nonresponse” (N=43, 11%), with minimal improvement over treatment; “rapid response” (N=73, 19%), with near-maximal improvement by week 2–3, followed by a relative plateau to week 6; “higher baseline symptoms, linear response” (N=118, 30%), with steady linear improvement and no apparent plateau by week 6; and “lower baseline symptoms, linear response” (N=154, 40%), again with steady linear improvement and no apparent plateau by week 6 (Table 2).

For the secondary outcomes, we observed significant differences between trajectories in response rates by week 1 and remission rates by week 3, and these differences were maintained until trial completion (Table 3). Compared with the included participants (N=388), the excluded participants (N=24) had a higher proportion in the nonresponding trajectory (42% compared with 10%), but formal statistical testing was not possible because of the small number of excluded participants (see Table S1 in the online supplement).

For the sensitivity analysis of rTMS treatment technique (HFL versus iTBS), we found that the optimal number of trajectory groups was three—likely because of smaller sample size (29)—but the models demonstrated a qualitatively similar pattern (a nonresponding group, a rapid-responding group achieving most improvement by week 2–3, and an intermediate group showing steady gains and no plateau by treatment end) (see Figure S1A,B in the online supplement). Notably, including treatment condition in the PROC TRAJ framework did not reveal any significant association of the treatment protocol (HFL or iTBS) with response trajectories.

Characteristics Associated With Each Symptom Trajectory

The results of the multinomial logistic regression model for group membership are presented in Table 4. After adjusting for treatment technique, the following characteristics were significantly associated with response trajectory group: clinician-rated baseline depression severity (HAM-D score), self-rated baseline depression severity (QIDS-SR score), age, and benzodiazepine use (present or absent). Characteristics that were associated with membership in the rapid response group were older age (odds ratio=1.04, 95% CI=1.01, 1.07), lower baseline QIDS-SR score (odds ratio=0.79, 95% CI=0.71, 0.87), and absence of benzodiazepine use (odds ratio=0.40, 95% CI=0.18, 0.90). In contrast, characteristics that were significantly associated with membership in the nonresponse group were higher baseline HAM-D score (odds ratio=1.31, 95% CI=1.17, 1.47) and QIDS-SR score (odds ratio=1.20, 95% CI=1.05, 1.38). The association of benzodiazepine use with the nonresponse trajectory fell short of statistical significance (odds ratio=2.25, 95% CI=0.99, 5.11). Sensitivity analyses in which benzodiazepine use was modeled as a continuous covariate (total daily dose) yielded a similar pattern of findings, with greater total daily benzodiazepine doses associated with reduced odds of membership in the rapid response group and greater odds of membership in the nonresponse group, although the findings were not significant (see Table S2 in the online supplement). Individual logistic regression models for each response trajectory indicated adequate model discrimination (c-statistic, 0.67–0.84) and fit (Hosmer-Lemeshow test p>0.05).

Discussion

This is, to our knowledge, the first study to describe depression response trajectories with 4 to 6 weeks of rTMS. We identified four distinct and clinically relevant trajectories. Differences in remission and response rates between trajectories were evident by week 1 and persisted for the duration of treatment. The rapid response group demonstrated dramatic improvement by week 2, whereas the nonresponse group demonstrated no improvement throughout treatment. Two intermediate groups demonstrated slower, linear improvement with no apparent plateau of improvement but with only a minority remitting by the end of treatment. At treatment completion, there were significantly different outcomes, with response rates ranging from 0% in the nonresponse group to nearly 90% in the rapid response group. Our exploratory analysis also identified four clinical characteristics independently associated with response trajectories: baseline HAM-D score, QIDS-SR score, age, and benzodiazepine use.

Our finding of distinct fast and slow response trajectories with rTMS agrees with previous trajectory-based analyses of pharmacotherapy in depression in older adults (8) and younger adults (6) and combined pharmacotherapy and psychotherapy in younger adults (7). However, unlike in trajectory-based analyses of pharmacotherapy (8, 9), we did not identify trajectories with delayed response to rTMS. This may be due to the fact that the pharmacotherapy trajectory-based analyses examined treatment outcomes up to 12 weeks, compared with a maximum of 6 weeks in the present study. Previous studies have suggested that some patients may benefit from longer rTMS treatment courses (5, 15, 32), and this could have been observed had our study allowed for prolonged treatment. However, delayed responders represent a minority of patients in pharmacotherapy trajectory analyses (5%−15%) (8, 9), and given the costs associated with rTMS (33), the cost-benefit analysis of prolonged treatment courses will necessitate careful consideration on an individual basis. Furthermore, our results raise the possibility that the identified response trajectories represent distinct neurophysiological phenotypes of major depression with preferential response to neurostimulation (34).

We also identified four demographic and clinical characteristics associated with membership in each trajectory: baseline HAM-D score, QIDS-SR score, age, and benzodiazepine use. On self- and clinician-rated scales, low baseline depression severity was associated with the rapid response trajectory and high baseline severity with the nonresponse trajectory. Similar observations were made between baseline severity and response trajectories in two previous pharmacotherapy trajectory analyses of depression (8, 35). This suggests that patients with higher baseline severity are less likely to respond quickly and may require longer durations of rTMS treatment. Interestingly, older age was associated with increased odds of rapid rTMS response trajectory, although our sample only included adults under age 65. This finding is consistent with previous work finding that rTMS was more effective for older adults under age 65 and for late-life depression (36, 37) when rTMS coils at higher stimulus intensities were used, as in the present study (i.e., at 120% resting motor threshold). This suggests that the hypothesis of older age predicting poor response to rTMS may be due to early clinical trials that used insufficient stimulus intensities (38). Further studies using modern rTMS technologies to treat individuals across the entire lifespan will be required to clarify the association between age and rTMS outcomes.

Of particular clinical relevance is the association between poorer response trajectories and benzodiazepine use, which was the only modifiable characteristic identified. Although high-dose benzodiazepine use was excluded, low-dose benzodiazepine use (32% of study participants) was associated with 60% lower odds of membership in the rapid response trajectory and more than twice the odds of membership in the nonresponse trajectory. This finding is likely independent of comorbid anxiety, which has been associated with worse outcomes (40), because clinician-rated anxiety (using the Brief Symptom Inventory) was not identified during our variable selection process. This is an important finding because benzodiazepines are positive allosteric modulators at GABA_A receptors and have been demonstrated to interfere with cortical excitability (41). This suggests that medications that affect cortical excitability, such as benzodiazepines or antiepileptic medications, may have a negative impact on rTMS treatment outcomes, as has been suggested in the literature on electroconvulsive therapy (42–44). From a clinical perspective, our work suggests that clinicians should consider discontinuing even low-dose benzodiazepines, if possible, before pursuing rTMS treatment.

Limitations

Some important limitations of this study should be noted. First, the selection criteria prevent the generalization of our results to individuals over age 65, individuals with bipolar depression, and individuals with significant psychiatric comorbidity. Second, two different treatments were used in this study (HFL rTMS and iTBS rTMS), which was a noninferiority trial, and although we performed multiple sensitivity analyses, and despite the fact that the original study found nearly identical longitudinal response trajectories (11), residual trajectory differences between treatment techniques remain possible. Given that these two treatments likely result in distinct neurophysiological changes, further work will be required to determine whether individual clinical or biological characteristics may be useful in the choice of treatment with iTBS or HFL rTMS. However, until biomarkers of response to either form of treatment have been identified, the increased treatment capacity offered by using iTBS remains a compelling advantage over standard HFL (3 minutes compared with 37.5 minutes). Third, our analysis of the characteristics associated with response trajectories was exploratory and data driven. While we applied conservative significance thresholds to mitigate this issue, the characteristics identified could still be a result of model overfitting and therefore require replication in independent cohorts in order to develop clinically useful predictive models. Fourth, this analysis only considered clinical characteristics associated with response trajectories and did not use any biological markers such as baseline anatomical or functional MRI scans or neurophysiological markers. These characteristics will be assessed in future studies for any potential further predictive value.

Conclusions

In this study of depression response trajectories with rTMS treatment, to our knowledge the first of its kind, we identified four distinct trajectories: nonresponse, rapid response, and two intermediate response trajectories. Greater age, absence of benzodiazepine use, and lower baseline depression severity were associated with the rapid response trajectory, and greater baseline depression severity was associated with the nonresponse trajectory. While this study provides clinically relevant information, further work will be needed to replicate these findings as well as to determine the utility of prolonged treatment courses and to explore the possibility that these trajectories represent biologically distinct major depression subtypes.