A new phase model of the spatiotemporal relationships between three circadian oscillators in the brainstem

Analysis of ex vivo Per2 bioluminescent rhythm previously recorded in the mouse dorsal vagal complex reveals a characteristic phase relationship between three distinct circadian oscillators. These signals represent core clock gene expression in the area postrema (AP), the nucleus of the solitary tract (NTS) and the ependymal cells surrounding the 4th ventricle (4Vep). Initially, the data suggests a consistent phasing in which the AP peaks first, followed shortly by the NTS, with the 4Vep peaking 8–9 h later. Wavelet analysis reveals that this pattern is not consistently maintained throughout a recording, however, the phase dynamics strongly imply that oscillator interactions are present. A simple phase model of the three oscillators is developed and it suggests that realistic phase dynamics occur between three model oscillators with coupling close to a synchronisation transition. The coupling topology suggests that the AP bidirectionally communicates phase information to the NTS and the 4Vep to synchronise the three structures. A comparison of the model with previous experimental manipulations demonstrates its feasibility to explain DVC circadian phasing. Finally, we show that simulating steadily decaying coupling improves the model’s ability to capture experimental phase dynamics.


AP-4Vep PD analysis and estimation of coupling parameters
Supplementary Figure S3. AP-4Vep PD analysis and estimation of coupling parameters. (A) The AV-PD dynamics of the five experiments that were selected to fit the AP-4Vep phase oscillator model. Each data trajectory (grey) has a period of constant PD between the vertical dashed lines, defined such that |dθ av /dt| < 0.01. After estimating the parameters of an AP-4Vep coupled oscillator model, the simulated trajectory (orange) matches well with the data. During an epoch of constant PD, the two oscillators are approximately at the same frequency, indicated by nearly identical period dynamics in (B). The average period over this epoch for each oscillator in all five experiments is plotted in (C). The difference between the oscillator's average period within this epoch is plotted in (D). The small differences suggest that the oscillators are approximately at the same period, hence their collective period is calculated as the mean period between the oscillators throughout this epoch. The AN-PD dynamics of the five experiments that were selected to fit the AP-NTS phase oscillator model. Each data trajectory (grey) has a period of constant PD between the vertical dashed lines, defined such that |dθ an /dt| < 0.01. After estimating the parameters of an AP-NTS coupled oscillator model, the simulated trajectory (blue) matches well with the data. To estimate all the parameters of the AP-NTS model, the linear decay rate of the PD data was estimated by fitting an exponential curve (red dots) to the PD trajectory when it is sufficiently close to its constant PD state. During an epoch of constant PD, the two oscillators are approximately at the same frequency, indicated by nearly identical period dynamics in (B). The average period over this epoch for each oscillator in all five experiments is plotted in (C). The difference between the oscillator's average period within this epoch is plotted in (D). The small differences suggest that the oscillators are approximately at the same period, hence their collective period is calculated as the mean period between the oscillators throughout this epoch.

4/7
Analysis of the stability of phase differences Supplementary Figure S5. The prominence of a peak in the histogram of phase differences is an indication of the stability of the PD time series. The prominence is measured using the expression in (A), which is 1 for a histogram with all its weight within one bin (the d-distribution) and close to zero (1/N, where N is the number of bins) when the distribution is uniform. This metric provides a measure of the stability of a PD trajectory. Simulations of the model system (1) with realistic periods and initial conditions (the same used in the simulations of Figure 4C and D) and zero coupling (s = 0) are shown in (B). PD trajectories of uncoupled oscillators are linear and display no tendency to remain at a particular PD value, which is indicated in the histograms in (C) (averaged over all simulations in (B)). Such distributions are nearly uniform and have a low y value, indicated as 'S0' in (Da) and (Ea). In (Da) and (Ea), the peak prominence is displayed for all 11 experiments, and example time series and histograms (15min-binned) are shown in parts b and c respectively. Three AN-PD traces are multistable (red 'MS'), as shown by the trajectories in Supplementary Figure S1 and their PD histograms displayed bimodality. Despite their high score for stability, these traces could not be used to fit the model since it was ambiguous which near-constant PD should be considered the as the steady-state. Furthermore, similar occurrences of multiple constant PD regions arise in simulations with slowly decaying coupling ( Figure 4C; black curve), hence this could be an effect of reducing the oscillators coupling.

5/7
Additional mathematical details The NTS is modelled as a single oscillator For demonstrating the assumptions underlying the description of the whole NTS PER2::LUC rhythm by a single phase variable, we shall consider an AP-NTS system. Consider an AP (θ a ) oscillator connected to two NTS oscillators (θ n 1 , θ n 2 ). The phase dynamics are given bẏ θ a = ω a + K an 1 sin(θ n 1 − θ a ) + K an 2 sin(θ n 2 − θ a ) θ n 1 = ω n 1 + K n 1 a sin(θ a − θ n 1 ) + K n 1 n 2 sin(θ n 2 − θ n 1 ) θ n 2 = ω n 2 + K n 2 a sin(θ a − θ n 2 ) + K n 2 n 1 sin(θ n 1 − θ n 2 ), First, we assume that the oscillatory signals from both NTS are identical, such that ω n 1 = ω n 2 = ω n and θ n 1 = θ n 2 = θ n . Next, we assume that AP-NTS interactions are symmetrical across hemispheres, such that K n 1 a = K n 2 a = K na and K an 1 = K an 2 1 = K an . This results iṅ θ a = ω a + 2K an sin(θ n − θ a ) θ n = ω n + K na sin(θ a − θ n ) θ n = ω n + K na sin(θ a − θ n ). (2) Clearly, these assumptions result in either NTS being described by the same equation, hence it makes sense only to study one variable to represent the NTS oscillations.
This mathematical consideration also motivates the protocol we use for our simulated NTS disconnect experiment. Disconnecting one NTS results in removing AP-NTS interactions from one hemisphere, K an 2 = K n 2 a = 0, and removing any possible communication between either side of the NTS, K n 1 n 2 = K n 2 n 1 = 0. Applying this disconnection to equation 1, whilst assuming that the intrinsic properties are unchanged giveṡ θ a = ω a + K an 1 sin(θ n 1 − θ a ) θ n 1 = ω n + K n 1 a sin(θ a − θ n 1 ) θ n 2 = ω n . (3) Under the assumption of bilaterally symmetric AP-NTS coupling, K an 1 = K an 2 = K an , we see that removing one NTS results in a reduction of the NTS to AP coupling (K an = K an 1 ) by one half (compare the first equations of systems 2 and 3).

6/7
A phase lag, γ, is required in the AP-NTS coupling term To motivate the inclusion of the phase lag, γ, in the AP-NTS coupling term consider an isolated AP-NTS system, the phase difference of which is given bẏ θ an = ω an −K an sin θ an , whereK an = K an + K na , and ω an = ω a − ω n < 0. The steady-state solution, θ * an , and the eigenvalue at this solution, λ , is given by sin θ * an = ω an /K an and λ = −K an cos θ * an , respectively. In our experiments, we observe small phase differences around zero between the AP and NTS, and so we want a stable (biologically observable) solution to be small, which we define as θ * an ∈ [6, 6] hrs ([−π/2, π/2] rad). For this solution to be stable we require λ < 0, and since cos θ * an > 0 this implies thatK an > 0. Hence for small phase differences to be stable, we require that the sum of the coupling is positive. However, this immediately leads to a problem since forK an > 0, equation 4 only allows positive solutions; the minimum phase difference possible, in the limit of infinite coupling, is zero. This does not agree with our observations of constant phase differences below zero (see Supplementary Figure S1, experiments 9, 12 and 13). Negative, stable solutions would be possible if there were an overlap in the intrinsic frequency distributions of the AP and the NTS, but our data reveals that this is not the case. We include a phase lag into the AP-NTS interaction terms to allow phase difference states around zero to be stable.