Coherency of circadian rhythms in the SCN is governed by the interplay of two coupling factors

Experiments were conducted in compliance with the rules and regulations established by the Animal Care and Use Committee of Hokkaido University.

Cry double deficient (Cry1,2^−/−) mice and Vip receptor 2 deficient (Vipr2^−/−) mice were bred with PER2::LUC mice carrying a PER2 luciferase reporter [38]. Wild–type (Cry1,2^+/+/Vipr2^+/+) PER2::LUC transgenic mice on the C57BL/6J background were used as control. Mice were reared in the animal quarters in Hokkaido University, where environmental conditions were controlled (lights–on, 6:00–18:00 h; light intensity, approximately 100 lx at the bottom of cage: humidity, 60 ± 10%).

For the measurement of PER2::LUC bioluminescence from a cultured SCN slice, mice of 8–16 weeks or 2–5 days old, kept under LD condition, were euthanized between 8:00 and 16:00 by cervical dislocation and decapitated. The brain was rapidly removed and a coronal SCN slice of 150 μm or 200 μm was made by a microslicer (D.S.K: DTK–1000; Dosaka EM) or a tissue chopper (McIlwain). The brain slice containing the middle portion of the SCN was selected and trimmed in approximately a 2×2 mm square. The slice was cultured in air at 36.5°C with 1.2 ml Dulbecco’s modified Eagle’s medium (Invitrogen) with 0.2 mM D–luciferin K and 5% supplement solution, the composition of which was described previously [35].

For the measurement of PER2::LUC from dispersed SCN cells, the SCNs from 4–8 neonatal pups (2–5 days old wild-type and Cry double deficient mice) were dissected from hypothalamic slices of 400 μm thick and dissociated using trypsin. Dispersed cells were plated on a 35 mm Petri dish pre–coated with 0.01% Poly–L–ornithine. The cell density was 1100±500 cells/mm². The medium composition was the same as that for the slice culture, except for 5% FBS in dispersed cell culture. In the co–culture experiment, the SCN slices of 150 μm thick were obtained from adult mice carrying the PER2::LUC reporter (recipient). The slice was pre-cultured for 3 or 4 days, and then co–cultured with an SCN slice from mice without the reporter system (donor). The donor SCN slice of 200 μm thick was obtained from WT mice of 7 days old and pre-cultured for one day before the co–culturing. When co–cultured, the graft SCN slice was placed inside out on the surface of recipient SCN slice. Measurement of the bioluminescence was started from the beginning of culturing of the recipient SCN and continued for at least 5 days after the co–culture. AVP receptor antagonists (V1A receptor antagonist: SR49059; TOCRIS, V1B receptor antagonist: SSR149415; Axon Medchem) were dissolved in water (SR49059 and SSR149415: final 2.5 μM). Water (vehicle) or antagonists were applied into the medium 5 to 7 days after co–culturing. The chemicals were either directly added to the culture medium (bath application) or dissolved in the culture medium to exchange with the whole medium in culture.

Bioluminescence at the SCN cell level in cultured slices or in dispersed cells was obtained by DM IRB (Leica), Luminoview 200 (Olympus), or Cellgraph (Atto) equipped with an EMCCD camera cooled at −80°C. The bioluminescence was measured every 60 min with an exposure time of 59 min. The pixel size was 2.3×2.3 μm for DM IRB, 2.0×2.0 μm for Luminoview 200, and 1.6×1.6 μm for Cellgraph. For the measurement of PER2::LUC from dispersed SCN cells, bioluminescence signals were analyzed within a region of interest (ROI). The mean area of a single ROI was about 100 μm², comparable to the size of a single SCN cell. The bioluminescence was expressed with an average intensity of pixels involved in a ROI.

To analyze spatio–temporal dynamics of the SCN slice movie data, the method of empirical orthogonal functions (EOFs), pioneered by Edward Lorenz in the context of statistical weather prediction [36, 39], was applied. The EOFs extract coherent structures of the spatio–temporal data as empirical eigenfunctions or empirical modes [37, 40, 41]. First, we consider the bioluminescence movie data as T × N matrix A=[a1a2⋯aN],(1) where N and T are the number of pixels in the SCN slice image and the number of time points, respectively. Each column vector a_k = [x_k(1), x_k(2),…, x_k(T)]^T represents time–sequence of the bioluminescence signal at k–th location of the SCN slice image. Interdependence of the dynamics at different locations can be quantified by the covariance matrix R = A^TA, where the (i, j)–element corresponds to covariance of the temporal patterns between locations i and j. The EOFs of the spatio–temporal data A are defined as the eigenvectors e_i of the covariance matrix R, sorted with respect to the size of the eigenvalues Ω_i (in descending order). Time sequence of scalar products between t–th bioluminescence image and i–th eigenvector is called the i–th empirical mode c_i(t). For oscillator network system, spatially coherent patterns are extracted as the major empirical modes, where the normalized eigenvalues, {100×Ωi/∑j=1NΩj[%]:i=1,…,N}, quantify the variance of the corresponding components.

As a model for circadian cells, a generic form of self–sustained oscillators is introduced as follows [42]: dxdt=-λxr(r-α)-ωy+ξx,(2)dydt=-λyr(r-α)+ωx+ξy.(3) The amplitude-phase model is described in Cartesian (x, y)–coordinates with radius r=x2+y2. The system gives rise to a limit cycle attractor with amplitude α and frequency ω, where perturbed dynamics returns to the attractor with a damping ratio of λ. The limit cycle is driven by independent Gaussian noise ξ_x and ξ_y. The single cell model has five unknown parameters {α, ω, λ, D_r, D_φ}, which were estimated for dispersed cell culture data by fitting the autocorrelation function of the model to that of the data [25, 42]. From the estimated parameters, the coefficient of variation CV can be computed, representing the ratio of the standard deviation of the amplitude fluctuations to the oscillator amplitude. The CV provides a criterion to distinguish self–sustained oscillators (CV < 1) from noisy damped oscillators (CV > 1). Detailed procedures of the parameter estimation are described in S1 Text.

By introducing local connections to the single cell models Eqs (2) and (3), which have been fitted to the dispersed data, a cellular network model of the SCN was constructed as dxidt=-λixiri(ri-αi)-ωiyi+∑j∈NiK(xj-xi)+Iavpsin(2π24t)+Ivipsin(2π24(t+ϕ))+ξx,i,(4)dyidt=-λiyiri(ri-αi)+ωixi+ξy,i,(5) where x_i and y_i represent dynamical variables of the i–th cell (i = 1, 2,…, N), ri=xi2+yi2, and N_i stands for neighbors of the i–th cell. The intercellular coupling strength was decomposed into VIP and AVP as K = a_avpK_avp + a_vipK_vip, where K_avp and K_vip stand for default strength of the VIP and AVP couplings. For simulation of the co–culture experiment, external signals from the neonatal wild–type SCN slice (24 h oscillation period) were described by intensities I_avp and I_vip for AVP and VIP signaling, respectively, the inputs of which are phase–delayed by ϕ. The role of Gaussian noise (ξ_x,i, ξ_y,i) is to determine the single cell oscillation property (self-sustained or noisy damped oscillator) and to suppress the network synchrony.

To simulate various types of slices (neonate vs. adult, wild-type vs. knockout), attenuation factors, a_avp, a_vip, were introduced to the AVP and VIP signaling. First, it has been reported that AVP expression in the SCN was significantly reduced in the Cry1 and Cry2 double–knockout mice [11]. Second, VIP expression and release exhibited endogenous circadian rhythms under constant dark condition in the neonatal wild-type SCN, but not in the adult wild-type SCN [28, 30], suggesting that VIP signaling is attenuated in the natural course of development. These findings lead to the following scenario [11]: (1) Through development, the VIP coupling is attenuated in adult; (2) In Cry1,2 double–knockout and Cry1,2 and Vipr2 triple–knockout mice, the AVP coupling is attenuated compared to wild–type; (3) In triple knockout, the VIP coupling is completely inactivated. The actual parameter values were selected based on the synchronization diagrams of S10 Fig, panel a–c, which show dependencies of the network synchrony on the attenuation factors.

Concerning the phase difference ϕ, it determines synergistic or antagonistic interaction between the VIP and AVP signaling. As explained in S1 Text, in-phase (ϕ = 0 h) strengthens the mutual coupling, while out-of-phase (ϕ = 12 h) weakens it. This can be confirmed in the synchronization diagram of S10 Fig, panel d. As a value to realize antagonistic relation between VIP and AVP, their phase difference was empirically determined as ϕ = 11 h. The simulation details are documented in S1 Text.

Synchronized rhythms of SCN neurons are particularly robust in organotypic brain slices from neonatal mice [6, 35]. In Fig 1 (upper graphs), we visualize such rhythms in a preparation from wild–type mice using PER2::LUC bioluminescence recordings. The oscillations appear totally synchronized with a period close to 24 hours and constant amplitudes over a recording time of 6 days. Such spatio–temporal patterns can be analyzed successfully by the EOFs. From the covariance matrix, the dominant spatial modes were extracted and the associated eigenvalues, representing the variance covered by these modes, were computed. Fig 1a shows that about 80% of the variance is represented by the dominant first mode (red color). Interestingly, the second mode (about 10% variance) detects also phase shifted cells in the upper part of the SCN (green). Such an advanced phase of the dorsomedial part of the SCN has been described earlier [6] and seems to be related to shorter period of cells in this area [43, 44]. Higher modes have quite small variances and provide no further information in this case. Neonatal slices typically show more phase coherent patterns than adult slices (S1 Text), pointing to developmental changes of the coupling [25, 35]. Five other slices of the neonate wild–type mice exhibited similar characteristics (sharp peak in period distribution, dominance of first and second modes, and high level of synchrony) as discussed in S1 Text and summarized in S1 Table and S1 Fig, panel d–i.

The locomotor activity of mice without the core clock genes Cry1 and Cry2 appears to be arrhythmic under constant darkness but the rhythmicity can be induced by light–dark cycles [45]. In neonatal slices of Cry1,2 double knockouts, some remaining rhythmicities have been reported [9, 35]. However, amplitudes and periods are quite variable in different slice preparations [11]. In the middle graphs of Fig 1, we analyzed a representative example using EOFs. Here the dominant modes explain about 40% and 12% of the variance. The first mode (red) obeys a period of about 32 hours, whereas the second mode (green) oscillates with a period of about 21 hours. The spatial patterning reveals that such a splitting is induced by a desynchronization of the left and right SCNs. Note that the splitting was confirmed in two slices among eight slices of neonate Cry1,2 double knockouts, where synchronized rhythmicities with fast damping were observed in the other six slices (S1 Text). This finding supports the hypothesis that the coupling between left and right is quite different than the coupling within the nuclei [46, 47]. Antiphase oscillations of the left and right SCNs and bimodal period distributions have been described also in hamsters and mice under constant light conditions [48, 49].

Seven slices of neonate double–knockout mice were further analyzed (see S1 Text, S1 Table, and S1 Fig, panel j–u). As described above, one slice exhibited left–right splitting (S1 Fig, panel j–o). The other six slices showed a single circadian rhythm with a global synchrony in the SCN (S1 Fig, panel p–u), being consistent with the earlier studies [35, 50]. Our interpretation is that the cellular coupling in the neonate double–knockout mice is close to the critical border. Slight difference in the coupling strength may lead either to global synchrony or to multiple clusters in the slice dynamics.

As mentioned in the Introduction section, the neuropeptide VIP is a major coupling agent within the SCN. It was shown that knockouts of the neuropeptide or of its receptor Vipr2 lead to disturbed activity rhythms and broad ranges of single cell rhythms [8]. Thus, neuronal coupling via VIP is essential to establish robust and precise rhythms. In the lower part of Fig 1, we analyzed slice data from a triple knockout, i.e., in addition to the knockouts of Cry1 and Cry2, the gene for the VIP–receptor is lacking. As expected, oscillations and synchrony are largely lost. No eigenvalue exceeds 10% and individual reporter signals appear to be noisy. Still, empirical orthogonal functions can detect weak clusters with periods of about 25 hours and some spatial patterns: the red cells are primarily in the right SCN, whereas most green cells appear on the left side. Two other slices of the neonate triple knockout mice showed also a very noisy behavior (see S1 Text, S1 Table, and S2 Fig).

In the same manner as the neonate slices, the cultured SCN slice data from adult mice were analyzed (six slices of wild–type mice, four slices of Cry1 and Cry2 double–knockout mice, and four slices of Cry1, Cry2, and Vipr2 triple–knockout mice). The results of the slice analyses are shown in detail in S1 Text and summarized in S2 Table. Representative graphs are also shown in S3 Fig.

Briefly, the adult wild–type slices showed clear circadian rhythms with global phase coherence (, panel a–f). Concentration of the phase–advanced pixels around innermost part of the dorsomedial SCN, observed in the neonate wild-type slice, was not recognized in the adult slice (, panel b), due to its slightly different configuration of the phase waves (see acrophase mapping of, panel c). Since phase waves and tides in the SCN are rather variable in different experimental settings, the EOF cannot be always expected to extract the same pattern of phase waves. S3 Fig S3 Fig S1 Fig

Adult double–knockout mice exhibited noisy and desynchronized rhythms (S3 Fig, panel g–l). As reported in [35], qualitative dynamics of the double–knockout mice changed significantly through development from neonate to adult. Adult triple–knockout mice showed even noisier behavior (S3 Fig, panel m–r). This is expected, because the VIP coupling was further diminished in the knockout slice.

To examine the four quantities (average and standard deviation of cellular periods, sum of principle eigenvalues, and synchronization index) that characterized the thirty slices from neonate and adult SCN, one-way analysis of variance (ANOVA) was carried out with respect to six groups (neonate wild-type, neonate double–knockout, neonate triple–knockout, adult wild-type, adult double–knockout, and adult triple–knockout). Statistically significant effect (p < 0.01) was detected for all the four quantities. According to post hoc comparisons using Fisher’s least significant difference, pairs of groups, whose means differ significantly (p < 0.01), were extracted. Although the results were similar among the four quantities, different pairs were also detected from one quantity to the other (see S1 Text). This indicates that the cellular periods, EOFs, and synchronization index capture similar but somewhat different features of the slice. These quantities should be utilized in a complementary fashion to detect the group differences.

The SCN can be regarded as a network of coupled oscillators and has been modeled extensively [24, 46, 51]. In most network models, the individual oscillator is based on transcriptional/translational feedback loops of the core clock genes [52–54]. When dealing with phenotypes displaying complex slice behaviors, such detailed biomechanical modeling approach may face difficulties, since many models might reproduce such experiments [55, 56]. Moreover, tedious optimization procedure of biochemical parameters is needed for the gene regulatory networks [56]. Our amplitude–phase model [25, 42], on the other hand, does not rely on complex gene networks. It simply connects dynamical properties of individual cells, which are quantified from dispersed cells, via inter-cellular coupling. Our former study [25] showed that such network of amplitude-phase oscillators can produce essentially the same results as those of complex gene regulatory network models. Although the amplitude–phase models do not provide a straightforward interpretation of specific gene mechanism, it has a generality of being independent of the choice of single cell models. As explained in detail in the Methods section, parameters of our single cell model in Fig 2a–2f were estimated from dispersed cells of neonate wild–type SCNs. For simulations of knockouts (Fig 2g–2r), we fitted our single cell models to SCN slices from Cry1 and Cry2 double–knockout mice. Our simulated cells are locally coupled via VIP and AVP terms (in Eq (4)). The corresponding coupling terms may exhibit different phases reflecting complex rhythmicities of VIP, AVP, and their corresponding receptors [8, 28, 50, 57]. Our simulations of the interplay of two different coupling terms reproduced observed counter-intuitive effects as discussed below.

Fig 2 shows that the observed spatio–temporal patterns described in Fig 1 can be simulated using the data–based stochastic single cell oscillators, local coupling, and imposed period differences. In Fig 2a–2f, we implemented the observation of Noguchi et al. [43] that the dorsomedial cells exhibit shorter periods. Even though all periods are locked, the second mode (green) indicates a different phase as found experimentally (compare Fig 1a–1f). Different periods of the left and right SCNs allow the simulation of splitting in Fig 2g–2l, that is comparable to the experimental data in Fig 1g–1l. Finally, we simulated triple knockouts in Fig 2m–2r by a reduced VIP–coupling and found largely random periods with small clusters that resemble the corresponding EOF analysis in Fig 1m–1r.

Our simulations in Fig 2 illustrate that rather few assumptions are required to reproduce quite complex spatio–temporal patterns in the SCN. Noisy single cell oscillators close to the Hopf bifurcation can be synchronized efficiently [44] and imposed period differences lead to phase and frequency clusters as observed experimentally.

Slices from adult Cry1 and Cry2 double knockout mice lose synchrony [25, 35]. Along the lines of Maywood et al. [9], SCN slices were co–cultured with neonatal wild–type SCN slices, which do not carry a bioluminescence reporter [11]. From a dynamical systems point of view, this protocol corresponds to a periodic forcing via paracrine signaling. Consequently, we extended our model by adding periodic forcing terms that represent external VIP and AVP signaling (see Eq (4)). Fig 3c–3e and S5 Fig, panel a–c,g–i, display a representative double knockout slice with co–culture. We find a partial rescue with a wide range of periods ranging from 15 to 37 hours with a broad peak around 24 hours. The dominant mode represents a variance of about 14%. In order to study the interplay of VIP and AVP coupling, a cocktail of AVP receptor antagonists has been added [11]. Unexpectedly, the cocktail enhanced significantly the amplitudes and the synchrony of adult Cry1 and Cry2 double knockout SCNs (see Fig 3f–3h and S5 Fig, panel d–f,j–l). The period distribution is much narrower and the dominant mode has an increased variance of 25%. The same feature was observed in two other slices (see S1 Text, S5 and S6 Figs, S3 Table). According to paired t-test applied to n = 3 slices, significant difference (p = 0.001) between control and AVP antagonists was detected using average period as the statistical quantity. Moreover, the 24 h-period component was strengthened by the antagonists treatment with a significant difference (p < 0.01) between antagonists and vehicle (Fig 3b). This is a counter-intuitive observation, since the weakening of coupling via AVP receptor antagonists improved synchrony.

Our modeling provides insight into the combinatorial effects of multiple coupling agents. Since the VIP and AVP coupling terms in our model exhibit different phases, their effect can by synergistic or competitive depending on their phase relationship (ϕ in Eq (4)). Fig 3l–3n demonstrates that the inhibition of one coupling agent can indeed improve synchrony. Thus, experimental data and simulations indicate that, in the preparations from adult double knockouts, VIP and AVP couplings act antagonistically. This explains why inhibition of AVP coupling can improve rhythmicity.

If there is indeed an antagonistic relationship between VIP and AVP couplings, perturbations of VIP signaling alone might also improve rhythmicity in periodically forced SCN slices. In order to test such situations, triple knockouts of Cry1, Cry2, and the VIP receptor Vipr2 were studied [11]. Surprisingly, the slices from triple knockout mice exhibited indeed improved rescue behavior compared to those from the double knockout mice [11]. Fig 4c–4e shows an example of such a rescued rhythmicity. The periods center around 24 hours and the first mode has a variance of more than 20%. Simulations confirmed that very weak single cell oscillators (compare Fig 1m–1r) can be synchronized efficiently with external forcing (Fig 4i–4k).

Finally, we studied the combined perturbation of both coupling agents. In Fig 4f–4h, the triple knockouts were further inhibited by the AVP antagonist. This implies that both major coupling factors were no longer acting, because the AVP signaling from the co–culture was inhibited. At the end, the synchrony was lost, a wide range of periods were observed, and all the EOFs have their variances below 5%. The reduced level of synchrony was observed also in two other slices (see supplementary S1 Text, S7 and S8 Figs, and S3 Table). The difference between control and AVP antagonists was significant (p = 0.02) using average period as the statistical quantity. Furthermore, the 24 h-period component was weakened by the antagonists with a significant difference (p < 0.01) from vehicle control (Fig 4b). Such a loss of synchrony is also visible in the associated simulations in Fig 4l–4n (also in S12 Fig, panel f,h).