Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

A computational model to characterize the time-course of response to rapid antidepressant therapies

  • Abraham Nunes ,

    Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Project administration, Resources, Software, Validation, Writing – original draft, Writing – review & editing

    nunes@dal.ca

    address: 5909 Veterans’ Memorial Lane, Abbie J. Lane Memorial Building (Room 3083B), QEII Health Sciences Centre, Halifax, Nova Scotia, Canada

    Affiliations Department of Psychiatry, Dalhousie University, Halifax, Nova Scotia, Canada, Faculty of Computer Science, Dalhousie University, Halifax, Nova Scotia, Canada

  • Selena Singh

    Roles Formal analysis, Methodology, Validation, Writing – review & editing

    Affiliation Department of Psychology, Neuroscience & Behaviour, McMaster University, Hamilton, Ontario, Canada

Abstract

Our objective is to propose a method capable of disentangling the magnitude, the speed, and the duration or decay rate of the time course of response to rapid antidepressant therapies. To this end, we introduce a computational model of the time course of response to a single treatment with a rapid antidepressant. Numerical simulation is used to evaluate whether model parameters can be accurately estimated from observed data. Finally, we compare our computational modelling-based approach with linear mixed effects modelling in terms of their ability to detect changes in the magnitude and time-course of response to rapid antidepressant therapies in simulated randomized trials. Simulation experiments show that the parameters of our computational model can be accurately recovered using nonlinear least squares. Parameter estimation accuracy is stable over noise levels reaching as high as 25% of the true antidepressant effect magnitude. Comparison of our approach to mixed effects modelling using simulated randomized controlled trial data demonstrates an inability of linear mixed models to disentangle effect magnitude and time course, while our computational model accurately separates these response components. Our modelling approach may accurately identify the (A) magnitude, (B) speed, and (C) durability or decay rate of response to rapid antidepressant therapies. Future studies should fit this model to data from real clinical trials, and use resulting parameter estimates to uncover predictors and causes of different elements of the temporal course of antidepressant response.

Introduction

Depression is a major public health problem [1], for which there is growing interest in management using rapid antidepressant therapies [2]. These treatments, such as sleep deprivation [3], ketamine [4, 5], and classical psychedelics [6] are notable for their ability to produce robust antidepressant effects after only a single dose. The central problem with these therapies is currently the short duration of response [7].

Understanding the time course of response to these novel antidepressant therapies is central to their development as evidence-based treatments for depression. We must study the predictors and causes of three aspects of treatment response: (1) the magnitude of response, (2) the speed of response (time to peak effect), and (3) the duration or persistence of response.

By predicting which patients are likely to respond significantly, quickly, and for a sustained period of time, we can better personalize treatment through optimized patient selection. By understanding the factors that predict or cause large, rapid, or sustained effects, we will be better able to develop strategies to potentiate or prolong rapid antidepressant action.

To discover factors that will predict the (A) magnitude, (B) speed, and (C) duration of response to rapid antidepressants, we must accurately measure these components of the overall time course, so that they can be used as dependent variables in predictive models. However, given the generally small sample sizes, resource intensity of experimental rapid antidepressant studies, and potential heterogeneity of response time courses across subjects, we must also balance the goal of characterizing temporal response profiles with the need to obtain sufficient statistical power. If we analyze rapid antidepressant response by comparing groups across multiple individual time points, we may not have sufficient statistical power due to multiple comparisons. Using computational modelling may improve statistical power in this setting [8].

Therefore, the present study introduces a nonlinear computational model to facilitate inference of the (A) magnitude, (B) peak response time, and (C) duration of response to rapid antidepressant therapies. We demonstrate that our model’s parameters can be accurately estimated in practice, and outperform linear mixed-effects models in terms of capturing different aspects of the time course of rapid antidepressant response. This latter point is of significance due to the common use of linear mixed-effects models in the analysis of rapid antidepressant response data [5, 914]. Our statistical approach will support the development of real-world studies to uncover biomarkers for each phase of the antidepressant response time-course, so that recovery can be achieved faster, and sustained for longer periods of time.

Materials and methods

Code to replicate all results of this study is included in S1 File.

Computational model

We model the effect of a single administration of a rapid antidepressant at time t as a normalized difference of exponentials, which is commonly used as a model of electrical conductance across postsynaptic membranes after arrival of discrete presynaptic action potentials [15]. Adopting this model for the present study is appropriate given its qualitative properties and relative simplicity. The model is formally defined as: (1) where g is the magnitude of the peak effect (can be thought of as the magnitude of response) and a > b are time constants. The constant a is the “fall time,” representing the amount of time it takes for the magnitude of Gg,a,b(t) to decline by a factor of e. We adapted the synaptic model in Eq 1 to model change in depressive symptoms, as measured by some rating scale, instead of measuring electrical current. This model was chosen because its three independent parameters (g, a, b) help disentangle the time course of response into its magnitude (g), speed (b), and duration (a).

The time until peak response can be computed as follows (2) unless a = b, in which case tpeak = a = b. The effects of parameter alterations on the value of Gg,a,b(t) and tpeak are shown in Fig 1.

thumbnail
Fig 1. Difference of exponentials model of ketamine response.

Time t = 0 in all cases represents the time of administration, and the units can be set arbitrarily by experimenters. Δ denotes the change in symptoms, conventionally measured using a rating scale. The default parameters for the model in each figure were set to g = −1, a = 1, b = 0.5.

https://doi.org/10.1371/journal.pone.0297708.g001

Using this model, one can show that a rapid and sustained response occurs as a grows and b approaches 0: (3) whereby the response becomes a step function.

Simulation and parameter estimation

We first sought to determine whether, given simulated antidepressant response trajectories, the correct underlying parameters could be recovered statistically. To do this, we first simulated noisy antidepressant response trajectories as follows: (4) where is Gaussian noise with standard deviation σ > 0.

We fit the curve Gg,a,b(t) (Eq 1) to a simulated observed symptom response curve Y = (yt)t=1,2,…,T using nonlinear least squares. Formally, this problem is defined as (5) and can be solved using widely available software packages for nonlinear optimization. In this study, we use the Optim.jl package for the Julia Programming Language (v. 1.8.2) [16]. We fit our model to samples of 100 simulated subjects whose response profiles were sampled according to the following parameters: (6) (7) (8) (9) where Uniform(l, h) denotes a uniform distribution with lower bound l and upper bound h. We repeated this across 25 increasing levels of noise σ ∈ {0.001, 0.011, 0.021, …, 0.241}. Note that the maximal noise level σmax = 0.241 corresponds to an average variation of ≈ 25% of the maximal effect |g| = 1. Parameter ranges were chosen to capture values that may be relevant to real world data. For instance, an effect magnitude of −1 ≤ g ≤ 1 captures the possibility of improvement (g < 0), no change (g ≈ 0), and worsening (g ≥ 0). The ranges of a and b were chosen to allow for a tpeak between 0 and 6 (assuming the relevant time unit is days). We must also maintain the fundamental constraint of a > b.

We quantified the error in parameter estimation across these noise levels using the median absolute deviation (MAD) and interquartile range.

Comparison to linear mixed model

We then sought to evaluate the statistical power of our modelling approach (Eqs 1 and 5) compared to that of linear mixed effects modelling, which is a commonly used analytical strategy in rapid antidepressant research capable of capturing repeated measures data within subjects, while identifying group-level variation [5, 914]. To do this, we simulated a large number of two-group parallel designed randomized trials with 20 subjects per group, in which one group demonstrated an average peak effect of g(0) = −0.5, rise time constant of b(0) = 1 day, and a decay time of a(0) = 5 days. Comparison groups were constructed by systematically changing (g, a, b) by gradually larger amounts. We simulated daily ratings over a 21 day observation period, to model a time horizon and measurement frequency that relevant to the real-world use of rapid antidepressants [5, 9, 17, 18].

When varying one of the parameters (g, a, b), respectively, in the comparison group, all others were held equal to that in the control group. The following perturbation levels were tested: δb ∈ {0.25, 0.5, …, 3.0}, δg ∈ {0.01, 0.03, …, 0.19}, and δa ∈ {0.25, 0.5, …, 3.0}. The comparison group parameter values were subsequently computed as b(1) = b(0) + δb, g(1) = g(0) + δg, and a(1) = a(0) + δa.

We implemented linear mixed effects modelling as a baseline comparator method. In R syntax, the mixed effects model was defined as (10)

After inferring the individual level parameters for each subject using nonlinear least squares (Eq 5), denoted , we tested for group differences using the following mixed effects model (11) for . For each perturbation level, we ran 10 simulated experiments. Statistical power was computed for each regression coefficient as the proportion of runs where the p-value of the estimate was under α = 0.05. Inference of all mixed-effects models was implemented using the MixedModels.jl (v4.8.0) package for the Julia Programming Language. We ran the entire experiment on an Apple M1 Max processor (2021; Apple Inc.; Cupertino, CA) in parallel over 16 threads.

Results

Parameter recovery experiment

Results of our parameter recovery experiment are shown in Figs 2 and 3. For low levels of noise, one can observe that the parameter recovery is highly accurate (Fig 2), and increases slowly (sublinearly) for increasing levels of noise (Fig 3). Specifically, at the maximal level of noise (corresponding to 25% of the peak effect magnitude), the median absolute error in estimation of g was under 0.15. Similarly, for parameter estimates of a and b, respectively, the median absolute error levels demonstrate a sub-linear increase with respect to noise levels, suggesting that parameter estimates are relatively stable to noise levels in the underlying data-generating process.

thumbnail
Fig 2. Results of parameter estimation experiment.

The top row demonstrates ground truth (x-axis) vs. estimated (y-axis) parameters. The dashed line along the diagonals indicate the point of perfect parameter recovery. The lower row of plots depict relationships between parameter estimates, to verify absence of correlation in parameter estimates. Note that since the models are constrained such that , the plot of vs. demonstrates what appears to be a positive correlation.

https://doi.org/10.1371/journal.pone.0297708.g002

thumbnail
Fig 3. Median absolute parameter estimation error in relation to response trajectory noise.

Error bars are interquartile ranges.

https://doi.org/10.1371/journal.pone.0297708.g003

Comparison to linear mixed Model

Using our modelling approach to disentangle variation in components of the time-course of rapid antidepressant response generally showed higher statistical power than the typical repeated-measures linear modelling approach (Fig 4).

thumbnail
Fig 4. Statistical power of our model compared to mixed effects modelling of the time course of rapid antidepressant effects.

X-axes denote the magnitude with which control-group parameters were perturbed to establish a comparator (e.g. “Intervention”) group: δg is the difference in peak response magnitude, δa is the difference in decay time, and δb is the difference in rise time. Y-axes denote the proportion of runs in which statistical models detected statistically significant effects (Eqs 10 and 11). Each row depicts results for runs in which a specific parameter was varied (top row is for variation in b, middle row for variation in g, and the bottom row for variation in a). Each column (labeled in the title of each plot) indicates which parameter (using our method) is being tested for group level differences. Blue curves with circular markers represent the statistical power for group level differences in the parameters estimated using our approach (labeled “Ours”). Red curves with diamond markers represent the statistical power for Group effects detected using the classical mixed effects model (labeled “MM-G” in the legend). Green curves with square markers depict the statistical power for Group × Time effects detected using the classical mixed effects model approach (labeled “MM-G×T” in the legend).

https://doi.org/10.1371/journal.pone.0297708.g004

Our modelling approach was better able to detect accurate group level differences across all levels of variation for all parameters compared to the classical modelling approach. For variation in a and b, one should expect that the classical modelling approach should reveal a Group × Time interaction effect, but not a pure Group effect. However, the Group × Time was only statistically significant at extremely large perturbation values, and independent Group effects were detected early. This is a problem, because perturbation of a and b change the time course of response, but not the magnitude of response.

Fig 4 shows that our model-based approach detected group-level differences in variation of a, b, and g, with greater statistical power than the classical approach across virtually all perturbation magnitudes (see plots along leading diagonal of Fig 4, from top left to bottom right). Furthermore, there is very limited risk of false positive results, compared to the mixed modelling method. For instance, group level differences in time course parameters (a, b) are rarely detected when the true underlying variation between groups occurs solely in the magnitude of effect (middle row in Fig 4). One also observes rare detection of statistically significant group differences in g (the magnitude of effect) when the true underlying differences between groups are time course-related (i.e. in a and b).

Discussion

We have introduced a computational model that can be used to disentangle different components of the time course of rapid antidepressant response: the magnitude, speed, and stability or decay of the effect over time. This model is simple and its parameters are interpretable (magnitude of effect g, decay time a, response time b). Furthermore, these parameters can be accurately recovered from data, and these estimates are robust to noise. Finally, we have demonstrated that an existing standard approach using linear mixed models cannot separate overall antidepressant effects from changes in time course, whereas our model can accurately separate these components. Thus, our approach may be useful for future studies that aim to study predictors, causes, and modifiers of the temporal course of response to rapid antidepressant therapies.

Two major questions about rapid antidepressant therapies concern (1) patient selection and response prediction [19], and (2) prolonging treatment effects [7]. Many investigations are being undertaken to identify predictors of ketamine’s antidepressant effects, given the rapidly widespread use of this treatment modality [19, 20]. Several predictors of response to ketamine have been identified, including intra-infusion dissociation [21], fewer treatment failures [22], family history of alcohol dependence [23], anhedonic or interest-activity symptoms [24], obesity [25], and mixed/anxious features of depression [26]. Electroencephalographic markers have also been identified, most notably increases in gamma power [20]. There is much less evidence available concerning predictors of response to other rapid antidepressant therapies, such as sleep deprivation or psychedelics.

To our knowledge, none of the aforementioned predictors have been specifically associated with either magnitude, speed, or duration of response to ketamine. Yet, it is important to stratify predictors based on which features of the antidepressant response trajectory they predict. It is perhaps most important to find predictors of the duration/durability of rapid antidepressant response (which is encoded in our model using the parameter a). We must also seek to discover adjunctive or maintenance therapies for rapid antidepressants. Using our model-based approach, the effectiveness of these adjunctive or maintenance treatments can be measured directly based on their ability to prolong decay time a. Therefore, future studies should evaluate the degree to which biomarkers are specifically predictive of the magnitude (g), speed (b), and duration (a) of response to rapid antidepressant therapies. By doing so, we may better personalize the application of rapid antidepressant therapies to the needs of individual patients. For example, patients may demonstrate clinical features that predict rapid response (b) to one antidepressant, but sustained response (a) to another antidepressant. Understanding these predictive markers may facilitate devising optimal personalized combinations for such patients, with one agent used for rapid termination of the depressive episode while the other is used for sustaining response.

One limitation of our paper is that all model-fitting experiments used maximum likelihood (point) estimation, which was employed for simplicity and conservatism. As such, the errors in our parameter estimation experiments are likely higher than if we had used maximum a posteriori or Bayesian estimation procedures constraining individual-level parameter estimates using the group-level mean and covariance.

Another limitation of our experiment is the lack of availability of real patient-derived data from rapid-antidepressant therapy trials. However, our group does not yet have access to such data. Furthermore, it is important to conduct internal validation of a modelling approach under controlled circumstances, as we have done here, prior to application with real-world data.

Having demonstrated the ability of our model to disentangle different phases of a rapid antidepressant response trajectory, and having shown its specificity and sensitivity to noise, future studies can apply this model to real-world data. Therefore, future studies should aim to fit our difference of exponentials model to data derived from patients undergoing rapid antidepressant therapies in randomized trials. By regressing the parameters (g, a, b) against various predictors and biomarkers, such studies will be able to discover factors that influence the magnitude, speed, and duration of rapid antidepressant response.

Conclusion

We have developed a computational model of the time course of response to rapid antidepressant therapies. This model can be fit to observed response time-series data in order to disentangle elements responsible for the magnitude, speed, and duration of response. Future studies should apply this model to real world data in order to uncover predictors and causes of these different temporal components of rapid antidepressant response.

Supporting information

S1 File. Code to reproduce the experiments and figures.

https://doi.org/10.1371/journal.pone.0297708.s001

(ZIP)

References

  1. 1. Palay J, Taillieu TL, Afifi TO, Turner S, Bolton JM, Enns MW, et al. Prevalence of Mental Disorders and Suicidality in Canadian Provinces. Canadian Journal of Psychiatry Revue Canadienne de Psychiatrie. 2019;64(11):761–769. pmid:31619055
  2. 2. Swainson J, McGirr A, Blier P, Brietzke E, Richard-Devantoy S, Ravindran N, et al. The Canadian Network for Mood and Anxiety Treatments (CANMAT) Task Force Recommendations for the Use of Racemic Ketamine in Adults with Major Depressive Disorder: Recommandations Du Groupe De Travail Du Réseau Canadien Pour Les Traitements De L’humeur Et De L’anxiété (Canmat) Concernant L’utilisation De La Kétamine Racémique Chez Les Adultes Souffrant De Trouble Dépressif Majeur. The Canadian Journal of Psychiatry. 2021;66(2):113–125. pmid:33174760
  3. 3. Ioannou M, Wartenberg C, Greenbrook JTV, Larson T, Magnusson K, Schmitz L, et al. Sleep deprivation as treatment for depression: Systematic review and meta-analysis. Acta Psychiatrica Scandinavica. 2021;143(1):22–35. pmid:33145770
  4. 4. Berman RM, Cappiello A, Anand A, Oren DA, Heninger GR, Charney DS, et al. Antidepressant effects of ketamine in depressed patients. Biological Psychiatry. 2000;47(4):351–354. pmid:10686270
  5. 5. Phillips JL, Norris S, Talbot J, Birmingham M, Hatchard T, Ortiz A, et al. Single, Repeated, and Maintenance Ketamine Infusions for Treatment-Resistant Depression: A Randomized Controlled Trial. American Journal of Psychiatry. 2019;176(5):401–409. pmid:30922101
  6. 6. Goodwin GM, Aaronson ST, Alvarez O, Arden PC, Baker A, Bennett JC, et al. Single-Dose Psilocybin for a Treatment-Resistant Episode of Major Depression. New England Journal of Medicine. 2022;387(18):1637–1648. pmid:36322843
  7. 7. McMullen EP, Lee Y, Lipsitz O, Lui LMW, Vinberg M, Ho R, et al. Strategies to Prolong Ketamine’s Efficacy in Adults with Treatment-Resistant Depression. Advances in Therapy. 2021;38(6):2795–2820. pmid:33929660
  8. 8. Huys QJM, Maia TV, Frank MJ. Computational psychiatry as a bridge from neuroscience to clinical applications. Nature Neuroscience. 2016;19(3):404–413. pmid:26906507
  9. 9. Zarate CA, Brutsche NE, Ibrahim L, Franco-Chaves J, Diazgranados N, Cravchik A, et al. Replication of Ketamine’s Antidepressant Efficacy in Bipolar Depression: A Randomized Controlled Add-On Trial. Biological Psychiatry. 2012;71(11):939–946. pmid:22297150
  10. 10. Ionescu DF, Luckenbaugh DA, Niciu MJ, Richards EM, Zarate CA Jr. A single infusion of ketamine improves depression scores in patients with anxious bipolar depression. Bipolar Disorders. 2015;17(4):438–443. pmid:25400146
  11. 11. Saligan LN, Luckenbaugh DA, Slonena EE, Machado-Vieira R, Zarate CA. An assessment of the anti-fatigue effects of ketamine from a double-blind, placebo-controlled, crossover study in bipolar disorder. Journal of Affective Disorders. 2016;194:115–119. pmid:26807672
  12. 12. Wilkinson ST, Ballard ED, Bloch MH, Mathew SJ, Murrough JW, Feder A, et al. The Effect of a Single Dose of Intravenous Ketamine on Suicidal Ideation: A Systematic Review and Individual Participant Data Meta-Analysis. American Journal of Psychiatry. 2018;175(2):150–158. pmid:28969441
  13. 13. Farmer CA, Gilbert JR, Moaddel R, George J, Adeojo L, Lovett J, et al. Ketamine metabolites, clinical response, and gamma power in a randomized, placebo-controlled, crossover trial for treatment-resistant major depression. Neuropsychopharmacology. 2020;45(8):1398–1404. pmid:32252062
  14. 14. Gilbert JR, Ballard ED, Galiano CS, Nugent AC, Zarate CA. Magnetoencephalographic Correlates of Suicidal Ideation in Major Depression. Biological Psychiatry: Cognitive Neuroscience and Neuroimaging. 2020;5(3):354–363. pmid:31928949
  15. 15. Rothman JS. Modeling Synapses. In: Jaeger D, Jung R, editors. Encyclopedia of Computational Neuroscience. New York, NY: Springer New York; 2014. p. 1–15. Available from: http://link.springer.com/10.1007/978-1-4614-7320-6_240-1.
  16. 16. Mogensen PK, Riseth AN. Optim: A mathematical optimization package for Julia. Journal of Open Source Software. 2018;3(24):615.
  17. 17. Nugent AC, Ballard ED, Gould TD, Park LT, Moaddel R, Brutsche NE, et al. Ketamine has distinct electrophysiological and behavioral effects in depressed and healthy subjects. Molecular Psychiatry. 2019;24(7):1040–1052. pmid:29487402
  18. 18. Diazgranados N, Ibrahim L, Brutsche NE, Newberg A, Kronstein P, Khalife S, et al. A Randomized Add-on Trial of an N-methyl-D-aspartate Antagonist in Treatment-Resistant Bipolar Depression. Archives of General Psychiatry. 2010;67(8):793–802. pmid:20679587
  19. 19. Rong C, Park C, Rosenblat JD, Subramaniapillai M, Zuckerman H, Fus D, et al. Predictors of Response to Ketamine in Treatment Resistant Major Depressive Disorder and Bipolar Disorder. International Journal of Environmental Research and Public Health. 2018;15(4):771. pmid:29673146
  20. 20. Gilbert JR, Zarate CA. Electrophysiological biomarkers of antidepressant response to ketamine in treatment-resistant depression: Gamma power and long-term potentiation. Pharmacology Biochemistry and Behavior. 2020;189:172856. pmid:31958471
  21. 21. Luckenbaugh DA, Niciu MJ, Ionescu DF, Nolan NM, Richards EM, Brutsche NE, et al. Do the dissociative side effects of ketamine mediate its antidepressant effects? Journal of Affective Disorders. 2014;159:56–61. pmid:24679390
  22. 22. Jesus-Nunes AP, Leal GC, Correia-Melo FS, Vieira F, Mello RP, Caliman-Fontes AT, et al. Clinical predictors of depressive symptom remission and response after racemic ketamine and esketamine infusion in treatment-resistant depression. Human Psychopharmacology: Clinical and Experimental. 2022;37(4):e2836. pmid:35179810
  23. 23. Luckenbaugh DA, Ibrahim L, Brutsche N, Franco-Chaves J, Mathews D, Marquardt CA, et al. Family history of alcohol dependence and antidepressant response to an N-methyl-D-aspartate antagonist in bipolar depression. Bipolar Disorders. 2012;14(8):880–887. pmid:22978511
  24. 24. Lally N, Nugent AC, Luckenbaugh DA, Ameli R, Roiser JP, Zarate CA. Anti-anhedonic effect of ketamine and its neural correlates in treatment-resistant bipolar depression. Translational Psychiatry. 2014;4(10):e469–e469. pmid:25313512
  25. 25. Lipsitz O, McIntyre RS, Rodrigues NB, Lee Y, Gill H, Subramaniapillai M, et al. Does body mass index predict response to intravenous ketamine treatment in adults with major depressive and bipolar disorder? Results from the Canadian Rapid Treatment Center of Excellence. CNS Spectrums. 2022;27(3):322–330. pmid:33267928
  26. 26. McIntyre RS, Rodrigues NB, Lee Y, Lipsitz O, Subramaniapillai M, Gill H, et al. The effectiveness of repeated intravenous ketamine on depressive symptoms, suicidal ideation and functional disability in adults with major depressive disorder and bipolar disorder: Results from the Canadian Rapid Treatment Center of Excellence. Journal of Affective Disorders. 2020;274:903–910. pmid:32664031