Major Depressive Episodes and Random Mood

Context  Mathematical models describing changes in mood in affective disorders may assist in the identification of underlying pathologic and neurobiologic mechanisms and in differentiating between alternative interpretations of psychiatric data.

Objective  Using time-to-event data from a large epidemiologic survey on recovery from major depression, to model the survival probability, in terms of an underlying process, with parameters which might be recognized and influenced in clinical practice.

Design  We present a sequential-phase model for survival analysis, which describes depression as a state with or without an additional incubation phase. Recovery is seen as the transition to a nondepressive state. We show that this sequential-phase model finds a microscopic realization in a dynamic description, the random-mood model, which depicts mood as governed by an Ornstein-Uhlenbeck type of stochastic process, driven by intermittent gaussian noise.

Results  For reversible depression (80%), the fractional probability of recovery is remarkably independent of the history of the depression. Analysis with the sequential-phase model suggests single exponential decay in this group, possibly with a short incubation phase. Within the random-mood model, the data for this reversibly depressed cohort are compatible with an intermittent noise pattern of stimuli with average spacing of 4 months and incompatible with nonintermittent noise.

Conclusions  Time-to-event data from psychiatric epidemiologic studies can be conceptualized through modeling as intrasubject processes. The proposed random-mood model reproduces the time-to-event data and explains the incubation phase as an artifact due to the inclusion criterion of 14 days in most current psychiatric diagnostic systems. Depression is found to result more often from pileup of negative stimuli than from single life events. Time sequences, generated using the random-mood model, produce power plots, phase-space trajectories, and pair-correlation sums, similar to recent results for individual patients. This suggests possible clinical relevance along with the model's use as a tool in survival analysis.

Mathematical models describing changes in mood in affective disorders may assist in the identification of underlying pathologic and neurobiologic mechanisms and in distinguishing between alternative interpretations of psychiatric data. Current mood models are based on data collected in long-term observational studies. Mood in individual persons is subject to change, which may occur within a relatively short time, and a recorded time sequence of mood may look haphazard to an outside observer. Stochastic behavior and, on the other hand, chaotic and deterministic mechanisms have been suggested.^1-7 Mood models, based on nonlinear differential equations, have been proposed, and the addition of a noise component has been studied.^1,8-10

The models presented herein are derived from and applied to epidemiologic data. We modeled time-to-event data from the Netherlands Mental Health Survey and Incidence Study (NEMESIS),^11,12 a prospective psychiatric epidemiologic survey of nonhospitalized subjects in the Dutch population. We investigate whether the underlying hazard rate can be understood in terms of sequential phases or states. Such an approach has also been proposed by Aalen and Gjessing,^13,14 who describe a model based on diffusion-type transitions between states. Economizing on the complexity of such a model, we specifically address the question whether a 2-state description, depressed vs nondepressed, is a sufficient basis (1-step model) or whether the data provide evidence for 2 distinct states of depression, the first of which would act as an incubation phase preceding the second phase, from which recovery occurs (2-step model).

The analysis with this sequential-phase model gives values for the fraction of subjects who remain depressed and the fraction of those who do eventually recover. For the latter group, the decay times (mean durations) of the phases are determined. These are the global parameters that characterize the average time course of depression. One would like to get an understanding of the mechanisms that drive the time sequences of mood in individual subjects or, one step less ambitious, the average time course of depression in a group (ensemble) of subjects. We propose here that this purpose is served by a random-mood model, which is based on a linear, noise-driven, equation describing an Ornstein-Uhlenbeck process.^15-17 To simulate an ensemble of subjects, a very long time sequence of mood is generated and characteristic parameters, such as average time between noise stimuli and the relaxation time by which mood is restored to normal, are regularly randomized. The time-to-event pseudodata derived from this simulation are then compared with the NEMESIS data and the sequential-phase analysis of these data.

We further investigated whether this random-mood model holds any promise for application in individual subjects. Power plots and phase-space trajectories are presented and compared with available data, such as those of Gottschalk et al^2,3 and Heiby et al.⁵ As these data have been interpreted as possible evidence for low-dimensional chaotic and deterministic behavior, we investigated whether simulation data of the (largely nondeterministic) random-mood model, when analyzed with the method of Grassberger and Procaccia,^18,19 suggest a higher dimensionality than the 1-dimensional space on which it has been defined.

Data used are from NEMESIS.^11,12 This prospective psychiatric epidemiologic survey of the Dutch population aged 18 to 64 years was conducted in 3 waves, in 1996, 1997, and 1999. From each selected household, 1 respondent was randomly chosen. In the first wave, data were collected for 7076 subjects; 1458 of these were lost to attrition in the second wave, and 822 were lost in the third wave. The remaining 4796 subjects were interviewed during all 3 waves. The respondents were interviewed using the Composite International Diagnostic Interview version 1.1.²⁰ Using the DSM-III-R,²¹ a variety of diagnoses was identified, including bipolar disorder, major depression, and dysthymia. Further details are described elsewhere.^11,12

The focus of the present study is on major depression. To include only new (first or recurrent) cases, respondents with a diagnosis of major depression in the period between T1 (1997) and T2 (1999) but no diagnosis of major depression during the month before T1 (1997), were identified. Subjects with bipolar disorder and primary psychotic disorder were excluded. The state of the major depressive episodes was assessed at 3-month intervals using the Life Chart Interview²² and was retrospectively discretized into 6-week intervals from recall.

Ten respondents were classified as having had a major depressive episode of 0.5-month duration. The duration of major depressive episodes in 250 respondents with depression was determined for the first depressive episode recorded on the Life Chart Interview. Data from NEMESIS are summarized in Table 1, which also gives the nonparametric estimates of the survival or Kaplan-Meier²³ curve and the conditional recovery probability estimates per time interval. We analyzed all respondents as a single group. The role of risk factors such as severity of depression, comorbid anxiety, comorbid dysthymia, somatic disease, recurrence vs nonrecurrence, sex, and treatment status are discussed in a forthcoming article (S.Y.W., K.I.K., P.J., J.S., R.G., J.K., unpublished data).

Data are analyzed in terms of a sequential model, which depicts recovery from depression as a transition from one state, A (depressed), to another state, B (nondepressed). In the general population, the prevalence of depression is rather constant, giving rise to a quasi-equilibrium, and the rates at which subjects become depressed (B→ A) or recover from depression (A→ B) are constant.

The assumption of exponential decay as a starting point may be justified by the following simple argument: just as in a chemical reaction, the equilibrium condition is d[A]/dt = −k[A]+k′[B] = 0, where k and k′ are the forward and backward rate constants and [A] and [B] are the numbers of subjects in states A and B. Considering only the process A→ B, that is, selectively following the time course of depression to recovery, without taking into account its starting point in time, the process is governed by d[A]/dt=−k[A], which describes exponential decay.

For the model to represent the general population, compatible with NEMESIS, it needs a few refinements. Depression appears irreversible in about 20% of the cases, as is already evident from inspection of Table 1. We, therefore, dichotomize the subjects retrospectively into a nonrecoverable fraction, S(∞), and a recoverable fraction, [1 − S(∞)], where S(∞) is the value to which the survival curve converges for large times. Exponential decay is assumed for the recoverable cases. In this 1-step model, the time-to-event probability reads

 

where τ is a relaxation time or decay time and the decay rate constant is k = 1/τ.

The possible existence of an incubation time is suggested by the low initial probability for recovery, as seen in Table 1. This can be modeled by splitting the depressed state, A, into 2 substates, A₁ and A₂. When initially all subjects are in state A₁ and recovery can only proceed from A₂, then A₁ acts as an incubation phase. The whole process reads as A₁→ A₂→ B, with decay times τ₁ and τ₂, respectively. For this 2-step model, the time-to-event probability reads

 

The hazard rate, h, is the partial probability or conditional probability for recovery. The conditional survival probability is (1 − h). These may easily be obtained from equations 1 and 2.

Maximum likelihood fits have been conducted using the 1-step model, with adjustable parameters τ and S(∞), and the 2-step model, with adjustable parameters τ₁, τ₂, and S(∞). The fitted data are the time-to-event probabilities, Ŝ_i, and the conditional or fractional decay probabilities, ĥ_i. The latter are statistically independent. The time-to-event probabilities, by their product form, are strongly correlated. The quality of the fit is judged, therefore, from the partial fit to the fractional decay probabilities only.

We attempt a description of how mood, y, might develop in time. Our approach is to investigate here a simple linear, noise-driven, equation:

 

The left-hand side of equation 3, by itself, produces exponential decay with relaxation time, T. When left alone, mood will therefore tend to zero, and it will stay there when nothing drives it away. On the right-hand side, g_w(D) is gaussian white noise of unit variance, occurring as sharp stimuli randomly distributed in time with average spacing, D. Z is a susceptibility or frailty coefficient,²⁴ which determines the extent to which the system follows the imposed stimuli: Z  = 0 would be immunity and Z  = 1 maximal sensitivity. The relaxation time, T, determines how fast the effect of a stimulus is digested or overcome.

The left-hand side of equation 3 is linear. Since different sources of gaussian white noise add up to gaussian white noise again, equation 3 is linear as a whole. Its autocorrelation function is just C(t) = exp(−t/T). This model is similar to the Ornstein-Uhlenbeck process.^15-17 The suggestion that it might be used in survival analysis has recently been put forward by Aalen and Gjessing.¹⁷

Simulation data are produced with the above model, discretized in Δt  = 1 (day) steps

 

Intermittent gaussian white noise is generated using a generalization of the Box-Muller algorithm²⁵

 

Here, a, b, and d are random numbers between 0 and 1. This intermittent noise has variance σ²  = 1/D.

Intermittent gaussian white noise (in time space) has the important property that Fourier transforms into nonintermittent gaussian white noise (in frequency space) with a variance reduced by a factor 1/D. Suppose mood is sampled in steps of 1 day, from m  = 0 to (N − 1). The mood-frequency amplitudes, c_n, can be obtained in closed form. Their absolute square (power plot) is given by

 

where the corresponding frequencies are f_n = n/N. The mean of this power plot is obtained by replacing |g_w(1)|²→ 1.

The discretized Fourier transform is mirror-symmetric around f  = 0.5, the Nyquist limit. This is not true for the continuous-time sampling, whose power plot reads

 

and to which equation 6 converges for low frequencies.

The continuous-sampling power plot exhibits 1/f² noise towards high frequencies. The discrete-sampling power plot follows this behavior only in the midfrequency region, typically for f < 0.2, and levels off towards f  = 0.5 because of its symmetry around this value. The spectrum represents exponentially correlated colored noise for T > 0 and becomes completely random (white noise) in the limit T→ 0.

Figure 1A, B, D, and E show the NEMESIS data as they have been published, with the first 1.5-month interval subdivided into 0.5 and 1.0 month. The time-to-event curve is not displayed in its more common form, as a stepwise decreasing function, but by data points and their standard deviations. The data points of the fractional decay probability (Figure 1D-F) are statistically independent. Those of the survival curve (Figure 1A-C) are not, and their interpretation is different: they rather define a 68% confidence band (shading in Figure 1A-C) within which the fit as a whole should be confined.

Details about the fitting results, using the 1-step and 2-step models, are summarized in Table 2. An irreversible fraction of 20% is found for both models. The relaxation time of the depressed phase is τ = 3.80 ± 0.31 months in the 1-step model. The 2-step model splits this over an incubation phase, τ₁ =0.29 ± 0.09 months, and a second phase from which recovery proceeds, with τ₂ =3.49 ± 0.28 months.

The conspicuously low values of the first data points of the conditional decay rate, that is, low recovery rate, are suggestive for the existence of an incubation phase. This conclusion hinges strongly, however, on the subdivision of the first 1.5 months into intervals of 0.5 and 1.0 month. It is noted that 0.5 month is close to the inclusion criterion of 2 weeks. We, therefore, recombined the data for the first 2 subintervals into 1 interval of 1.5 months (Figure 1C and F). Leaving τ₁ as an adjustable parameter, a value of zero is preferred, thus reducing the 2-step model to the 1-step version. The best fit gives τ = 3.71 ± 0.30 months and a nonrecovering fraction of 0.20 ± 0.03.

The analysis has dichotomized the subjects, retrospectively, into 2 groups: those who eventually recover and those who remain depressed. The partial survival and decay curves of the subgroup of those who eventually recover are shown in Figure 1C and F as dashed curves. For the nonrecovering group, the conditional decay probability remains zero.

The analysis demonstrates that an incubation phase, if it exists, is short, at least shorter than the inclusion criterion of 2 weeks. Hence, time to recovery from reversible major depressive episodes seems to be governed by a single decay constant. This implies that the fractional decay probability, that is, the probability for recovering from depression, is constant. More specifically, the probability of recovery seems to be independent of the time that the depressed condition has existed. This suggests that at any time the remaining group of subjects with depression is always just a scaled-down replica of the initial cohort.

Independence of the fractional decay probability of the history, found in the reversibly depressed cohort, is a natural property of the random-mood model. The model (equations 3-5) follows mood on a day-by-day basis over a stretch of time equivalent to many typical lifetimes of a single person. The relaxation parameter, T, and the susceptibility, Z, are individual properties that must be randomized to make the model represent an ensemble. T is randomly chosen from a Poisson distribution with mean <T> and then kept constant for 1000 integration steps (≈3 years). For the moment, we will adopt equal susceptibility, Z = 1, for all subjects and postpone the discussion of a multivariate distribution until later.

We adopt a depression line at a value of −2. This means that only 5% of the negative stimuli is strong enough to initiate a depression by themselves, in which case they are life events. When mood drops below the −2 limit, this is taken as the start of an episode of depression and the duration of the depression is recorded. An inclusion criterion is applied, that the depressed state should last at least 14 days, in line with the current DSM-III-R guidelines. Time-to-event curves, derived from these simulation data, are then compared with the NEMESIS data, and the parameters D and <T> may be adjusted to achieve a best fit.

Figure 2 shows a mood pattern, modeled after the NEMESIS data, with parameters explained in the following subsection “Time-to-Event Simulation Data From the Random-Mood Model.” Figure 2A shows a simulated mood pattern extending over 80 years, the typical lifetime of one person and illustrating that depressive episodes come in various lengths. Figure 2B follows mood over 1000 years. On this extended scale, mood exhibits a global pattern that closely resembles noise. The simulation data are based on 10 million days, equivalent to following a cohort of more than 9000 subjects for 3 years.

The parameters of the random mood model are D, the average spacing between stimuli, and <T>, the mean of the Poisson distribution from which the relaxation time is taken. These are rather well determined by requiring that the simulation data should best reproduce the NEMESIS time-to-event curve, for those subjects who will recover only. The result of a (nonautomated) search is D = 120 days, <T> = 365 days. The stimuli come, on average, once in 4 months and, when strong enough to initiate a depression, they are life events that require a long time to cope with.

It has been verified that the analysis and its results are changed in no way when a low-amplitude nonintermittent noise term is added in equations 3 and 4. In other words: the addition of small daily hassles, with SD ≤ .03, does not affect the results when at the same time these small irritations and pleasures are quickly forgotten (<T′> ≤ 1 day).

The simulation data are shown in Figure 3, together with the fits that were obtained for the recoverable cohort in the NEMESIS data. The agreement is excellent. The exponential decay persists even longer than the NEMESIS data can be followed. Three more results of this model deserve attention. First, we find that in our model depression results more often (>80% of cases) from pileup of negative stimuli than from single life events. Second, the fractional decay probability tends to zero for short duration. Under the inclusion criterion of 14 days (DSM-III-R), the simulation data reproduce this drop. Changing the inclusion criterion to only 1 day, however, the first points of the fractional decay probability go up and the drop at short duration vanishes altogether. The effect is only significant for the first 2 points (asterisks in Figure 3D). This suggests that the evidence, for a possible incubation phase, found above, is most likely an artifact, caused by the DSM-III-R inclusion criterion. Third, the fraction of time that mood falls below the adopted depression line is 4% to 5%, which is 3 to 4 times higher than the prevalence of depression in the general Dutch population. This is easily remedied by the ad hoc introduction of a bivariate distribution for the susceptibility factor: Z = 1 for approximately 30% of the cohort and Z = 0 (immunity) for the rest. This is known as a cure model.²⁴ Alternatively, and probably more realistically, one might take Z as a random number between 0 and 1, in which case <T> must be raised to 500 days to get the same survival curve.

The frequency pattern of mood may be studied from its Fourier transform. Figure 4A shows the power plot of a time sequence for the simulation data, obtained with the random-mood model. Also shown in Figure 4 is the continuous-time sampling (equation 7), which towards higher frequencies keeps descending as 1/f², in contrast to the discrete sampling data, which are symmetric around the Nyquist frequency.

Figure 4B shows the power plots and their means for different single relaxation times. It shows that the flat region at low frequencies grows with decreasing relaxation time, leaving an increasingly shorter 1/f region.

The analysis has thus far been concerned with an ensemble of subjects. We extend the model to individuals by considering T, Z, and, in part, also D, as parameters characteristic for a particular individual. This extension allows us to compare our model with results from longitudinal observations.^2,3,5 The curves in Figure 4B might then be interpreted as individual power plots. Suppose one is presented with an experimental time sequence and wishes to analyze it in the context of the random-mood model. The power plot may give the relaxation time, T, but it contains no information on whether the stimulus pattern is intermittent: equations 6 and 7 tell us that daily occurring stimuli of variance 1/D give the same power plot as intermittent noise of variance 1 and average time spacing, D. A related problem is estimating the absolute scale of mood, especially if the sequence has been presented on a self-rated scale. One person may, in Goethe's words, swing from “Himmelhoch jauchzend bis zum Tode betrübt,”²⁶ while another may only be driven by daily ups and downs. Yet, when an absolute scale is missing, their power plots may be identical as long as their relaxation times are the same.

Unlike power plots, phase-space trajectories are sensitive to intermittency. Figure 5 gives such trajectories, derived with the random-mood model, driven by intermittent noise of average time spacing (D  = 120 days) and variance 1 (Figure 5A, C, and E) and by nonintermittent noise (D = 1) with variance 1/120 (Figure 5B, D, and F). Phase-space has been defined here as the 2-dimensional space, spanned by mood and its derivative, displayed as the mood difference between 2 consecutive days. The scale is the same in Figure 5A through F.

Because its variance is reduced by a factor 1/D, the phase-space plots for nonintermittent noise (Figure 5B, D, and F) exhibit significantly lower spread in their day-to-day fluctuations compared with the intermittent noise model (Figure 5A, C, and E). Nonintermittent stimuli occur, however, D times more often, and mood itself may “random-walk away” about equally far in both cases. Thus, for <T> = 1 year, mood is found outside the depression-manic limits, adopted as −2 and +2, with about the same probability, roughly 4% to 5%.

A number of studies has reported analyses of mood in terms of chaos theory.^2-7 The idea is rather simple: one may suspect that mood changes along with some other, unknown quantities that are all interconnected in a deterministic scheme, for example, a set of equations, and which together span a multidimensional space. Solutions of these equations may then be found on a subspace, the attractor, the dimension of which is called the Hausdorff dimension. Estimating this Hausdorff dimension is, however, rather intricate. The method that is probably used most often is that of Grassberger and Procaccia,^18,19 in which one considers smaller time sequences (y_i, y_i + Δ,  . . . , y_{i+ (f-1)Δ}), all taken a time lag Δ apart. One then determines how much 2 different such mini–time sequences may differ from each other by making up their so-called pair-correlation sum, C(r). This sum should, for small distances r, behave like r^ν. When the value of ν, read from a log-log plot, shows convergence on increasing the embedding space (f, the length of the mini–time sequences), this might be an indication for chaotic, deterministic behavior, and ν would be a good estimator of the Hausdorff dimension.

Figure 6 shows simulation data for intermittent noise with D = 16 and T  = 64. A small component of nonintermittent noise with variance 0.01 and T = 1 has been added to simulate small daily ups and downs. The embedding dimension was taken up to f  = 16 and the time lag Δ = 5, as in the study by Gottschalk et al.² For the smallest distances, the apparent correlation dimension is about equal to the embedding dimension itself; there is no convergence, and the process is practically random. On the other hand, the large-amplitude intermittent noise produces a maximum slope just before the correlation sum saturates. This region shows near-convergence towards a value of 4, close to that found by Gottschalk et al^2,3 and by Heiby et al,⁵ who evaluated their correlation dimension in this local linear scaling region or in the flat area just before that bump.

We want to stress that the random-mood model might be extended to individual subjects. We, therefore, need to make a comparison of random-mood simulation data and individual time sequences as available in the literature, for which we chose those of Gottschalk et al^2,3 and Heiby et al.⁵ Raw time sequences cannot be compared directly, but power plots, phase-space plots, and correlation sums derived from them can. These are found largely similar to our simulation data. This proves our point.

What one learns from this in addition is a caveat: that unfortunately, for relatively short time sequences of N ≈ 1000, the proper r→ 0 limit remains elusive by lack of statistics and is otherwise obscured by low-amplitude noise, which cannot be avoided; mood, determined by rating on some scale, has an inherent uncertainty caused by rounding-off errors. There is bound to be a noise component. Estimating the Hausdorff dimension from a higher r region, as Gottschalk et al^2,3 and Heiby et al⁵ do, comes down to stretching a mathematical theorem to a region where its application was not intended.

Analysis of time-to-event data on major depressive episodes in the general population shows that about 75% of subjects recover within a year and about 50% within 3 months. For about 20% of the subjects, however, recovery takes longer than 2 years. When analyzing the subjects who eventually recover as a group by itself, one observes a remarkable constancy in hazard rate, suggesting that the probablility of recovery is largely independent of the preceding history.

The random-mood model, designed to reproduce the average time course of depression among an ensemble of subjects, seems to provide a natural microscopic description of the survival data. It also explains the evidence for an initial incubation phase as an artifact, introduced by the inclusion criterion of 14 days, which is dictated by most current diagnostic systems (eg, DSM-III-R, International Classification of Diseases).

On the level of individual subjects, the random-mood model holds promise: it produces the essentials of power plots and phase-space trajectories as reported in the literature. We have alluded to the possibility that mood might be split into a long-term component driven by a large-amplitude, intermittent stimulus pattern and with long relaxation time and, on the other hand, a short-term component arising from nonintermittent, low-amplitude stimuli and characterized by a short relaxation time. The short-term component would describe small day-to-day fluctuations, whereas the long-term component corresponds to more serious events and may give rise to major depression.

When analyzing the pair-correlation sum in terms of chaos theory, the r→0 limit is obscured by a small noise component. Associating the correlation dimension of the intermittent large noise component with a locally straight section on a log-log plot at higher r values suggests near-convergence to values that have been reported in the literature.^2,3,5 Yet the random-mood model is linear and 1-dimensional by its construction. This illustrates that convergence of the apparent correlation dimension does not automatically imply a deterministic nonlinear mechanism, especially not if the analysis is done in a range that does not represent the r→0 limit.

We propose that the random-mood model might serve as a basis to analyze mood sequences. Its simple interpretation is that mood is driven by an intermittent stimulus pattern of gaussian noise with average rate 1/D. In clinical terms, the susceptibility, Z, determines how strongly mood follows the imposed stimulus pattern. The relaxation time, T, determines how long afterwards the effect of a stimulus is felt. We studied the influence of both of these parameters on the probability for getting depressed and on the probability for recovering from depression.

By analogy of equation 7 with an electric circuit, the quantity

 

may be interpreted as an impedance or resistance against imposed mood swings of frequency, f. Both the parameters Z and T can be influenced by medication and by care. Thus, the random-mood model may offer a framework for testing the effects of therapeutic interventions or of intrinsic factors (including somatic comorbidity) on the time structure of major depressive disorders. For example, antidepressant drugs may lessen the consequences of acute negative experiences by lowering the susceptibility and shortening the relaxation time, hence increasing the impedance. An unfortunate aspect is that medication mostly lessens the effects of positive experiences as well. An example of differences in susceptibility would be that lowering of brain serotonin concentration with a diet low in tryptophan induces depressive symptoms in genetically vulnerable individuals but not in others.^27,28

Correspondence: Siebren Y. van der Werf, PhD, Kernfysisch Versneller Instituut of the University of Groningen, Zernikelaan 25, 9747AA Groningen, the Netherlands (vdwerf@kvi.nl or s.y.van.der.werf@rug.nl).

Submitted for Publication: July 11, 2005; final revision received September 15, 2005; accepted September 21, 2005.