Introduction

Interactions between the hippocampus and the medial prefrontal cortex (mPFC) are essential for cognition. During sleep, hippocampal and cortical activity is continually synchronized by cardinal sleep oscillations such as slow waves (SWs, 0.1-4 Hz), thalamocortical spindles (10-15 Hz), and hippocampal sharp-wave ripples (SWRs, 100-300 Hz) ^1,2. These rhythms have well-known roles in neuronal plasticity ^3,4 and in the consolidation of hippocampus-dependent memories ⁵. In particular, the hierarchical coupling of slow waves, spindles and SWRs is thought to promote the stabilization of short-term memory ⁴. Consistent with this hypothesis, the temporal coupling of hippocampal and mPFC sleep oscillations results in a selective increase in the recall of hippocampus-dependent memory^5,6.

At present, the pathways that coordinate the long-range coupling of slow waves, spindles and SWRs across the prefronto-hippocampal network are not well characterized. Although the hippocampus extends projections to several prefrontal areas^7–11, few regions of the prefrontal cortex reciprocate¹². With respect to medial regions of the prefrontal cortex, no direct projections to the hippocampus have been identified^11,13–16.

The nucleus reuniens of the thalamus (reuniens) is bidirectionally connected to the hippocampus and the mPFC^17,18, making it an ideal hub for coupling hippocampo-cortical interactions during sleep. Reuniens is part of the ventral group of midline thalamic nuclei with widespread connections to cortical, thalamic, hypothalamic and brainstem structures ^18–22. Consistent with its anatomical connections to the hippocampus and the mPFC, reuniens is required for memory function^23,24 and a growing body of evidence points its critical role in processing hippocampus-dependent memory ^23,25. Although evidence for the functional significance of reuniens in cognitive operations is strong, little is known about the nature of reuniens-prefrontal interactions and their interplay with hippocampal events during sleep ²⁶.

Considering the well-established link between sleep and memory consolidation ^4,27,28, we investigated the role of the reuniens in coupling hippocampo-cortical oscillations during non-rapid eye movement (NREM) and REM sleep. Using in vivo intracellular recordings, computational modelling, and spike/local field potential (LFP) recordings, we demonstrate that reuniens is driven by hippocampal oscillations during sleep and that it elicits fast, consistent synaptic responses in mPFC. Hippocampal and reuniens activity precedes the onset of mPFC spindles, suggesting a hippocampo-thalamic trigger for spindle initiation. In turn, the phase of mPFC slow waves determined the amplitude of reuniens spindles and of SWR timing, forming a functional feedback loop through the reuniens.

Results

Reuniens elicits fast, consistent synaptic responses in the mPFC

We obtained intracellular recordings of mPFC neurons in vivo (n=49 neurons) during targeted, electrical stimulation of the ventral hippocampus and reuniens in anesthetized cats (Figure 1A-D). Stimulation of the hippocampus and reuniens elicited synaptic responses in 55% and 51% of recorded mPFC neurons, respectively. Hippocampal stimulation generated variable responses in the same neurons, consisting of complex sequences of either two excitatory postsynaptic potentials (EPSPs) followed by long-lasting hyperpolarization (Figure 1E, top, left), a depolarization-hyperpolarization sequence (Figure 1E, middle, left), or no response (Figure 1E, bottom, left). The polarity of early depolarizing responses elicited by hippocampal stimulation was dependent on membrane potential, reversing with progressive depolarization of the membrane at around -70 mV (Figure 1E, right), suggesting that the initial depolarization was an inhibitory postsynaptic potential (IPSP).

Figure 1

Reuniens stimulation elicited consistent EPSPs followed by long-lasting hyperpolarization (Figure 1C, bottom) at earlier latencies than hippocampal stimulation (Figure 1G). EPSPs elicited by reuniens stimulation had faster rising slopes compared to hippocampal stimulation (mean ± SD) but were of lower amplitude. We observed no significant difference in the duration of the late, hyperpolarizing component in response to reuniens and hippocampal stimulations.

Reuniens stimulation elicited action potentials that were either orthodromic (85%), antidromic (10%) or both (5%), evidenced by jitter in the latency of orthodromic action potentials and consistent latency of antidromic action potentials (Figure 1F). Antidromic responses to reuniens stimulation, confirmed by collision tests, were recorded in cells located 6050-6600 μm ventral to the dorsal aspect of the cortical surface. The mean latency of antidromic spikes was 2.9±0.25 ms (mean ± SD).

To characterize the prefrontal topography of reuniens and hippocampal receptive fields, we descended a 4-shank microelectrode array into the mPFC in 1 mm steps and recorded evoked field potentials in response to reuniens and hippocampal stimulation (Figure 1H). Reuniens stimulation evoked field postsynaptic potential in a wider prefrontal territory compared to hippocampal stimulation which was restricted to the ventromedial wall of the mPFC. In some locations of mPFC, early responses to reuniens stimulation (Figure 1 H, I) corresponded to criteria of field potential antidromic responses²⁹. Areas in the prefrontal cortex that responded to hippocampal stimulation overlapped with areas exhibiting antidromic responses to reuniens stimulation (Figure 1 I, shadowed area). In the following experiments, we targeted this region of the mPFC.

Reuniens and hippocampal activity precede the onset of mPFC spindles

Sleep enhances communication in the hippocampo-cortical pathway by promoting interregional synchrony through highly rhythmic events such as slow waves and spindles. To investigate the role of the reuniens in mediating synchrony in the mPFC-hippocampal network, we performed LFP recordings of the mPFC and hippocampus, as well as targeted spike/LFP recordings of the reuniens during natural sleep (Figure 2A, Figure 2 – figure supplement 1). Recordings were segmented into wake, NREM and REM epochs, based on mPFC delta power, mPFC signal amplitude, and nuchal muscle activity (Figure 2 B-D). During NREM sleep, LFP recordings and corresponding wavelet transforms of the signal suggested correlated activity among mPFC slow waves, spindles, hippocampal SWRs, and reuniens multi-unit activity (Figure 2E). We observed co-occurring reuniens and mPFC spindles (Figure 2F), characterized by significantly earlier onset of reuniens spindles relative to the onset of mPFC spindles (Figure 2G).

Figure 2

Analysis of reuniens single-unit activity relative to the phase of mPFC and reuniens spindle cycles revealed significant phase-locking of reuniens single-units to mPFC spindle cycles and reuniens spindle cycles (Figure 2H, 2I). The mean phase preference of reuniens single-units was advanced for reuniens spindle cycles relative to mPFC spindle cycles. For pairs of co-occurring reuniens and mPFC spindles, spindle amplitudes were significantly correlated (Figure 2J, 2K).

Lastly, we investigated the change in hippocampal and reuniens signal power during mPFC spindles relative to baseline, which demonstrated increased hippocampal ripple power and increased reuniens spindle power preceding the onset of spindles in the mPFC (Figure 2L, top). The histogram of SWR incidence in the same time window, referenced to the onset of mPFC spindles, indicated increased SWR incidence prior to the onset of mPFC spindles (Figure 2L, middle). Reuniens multi-unit activity increased prior to the onset of mPFC spindles (Figure 2L, bottom).

mPFC slow waves determine spindle amplitude in the reuniens-mPFC network

To examine the modulation of reuniens and hippocampal activity by mPFC slow waves, we analyzed the relationship between the phase of mPFC slow waves with respect to spindles and SWRs. Figure 3A shows a representative example of a coupled slow wave-spindle event from simultaneous recordings of the mPFC and reuniens. We analyzed the average time-frequency representation of mPFC amplitudes locked to the peak of the mPFC slow wave (Figure 3B), as well as the average time-frequency representations of reuniens amplitudes, locked to the peak of the same mPFC slow waves (Figure 3C), which showed significant modulation of spindle amplitudes by the phase of the mPFC slow wave.

Figure 3

The strength of slow wave-spindle coupling, estimated by the absolute value of the synchronization index was higher for mPFC slow wave-reuniens spindle coupling compared to mPFC slow wave-mPFC spindle coupling (Figure 3 D). A comparison of slow wave-spindle coupling strength between mPFC-mPFC pairs and mPFC-reuniens pairs is shown in Figure 3E.

The amplitude of reuniens spindles increased at earlier phases of the mPFC slow wave compared to the amplitude of mPFC spindles, which increased at later phases of the mPFC slow wave (Figure 3B,C,F).

Next, we examined the relative delays in the timing of mPFC and reuniens slow waves. Figure 3G shows a representative example of an overlapping mPFC and reuniens slow wave, with circles indicating slow wave peaks. Single session data in Figure 3H demonstrates the delay between reuniens and mPFC slow waves, showing significant lag in the time of reuniens slow waves relative to mPFC slow waves. The distribution of median delay times for all recording sessions between each reuniens-mPFC slow wave pair are shown in Figure 3I.

Lastly, we analyzed the phase-locking of hippocampal SWRs to mPFC slow waves. Figure 3J shows a representative example of a coupled slow wave-SWR event. The distribution of slow wave phases at the time of SWRs was non-uniform, indicating significant phase-locking of SWRs to mPFC slow waves at a mean phase of 162.7° (Figure 3K). We found no significant phase-locking of SWRs to spindle cycles in the majority of recording sessions (n=3/51 with significant phase for mPFC spindle-SWR and n=10/51 sessions with significant phase for reuniens spindle-SWR coupling, Figure 3 – figure supplement 1 and Table 1).

Table 1

Hippocampal oscillations modulate reuniens activity

Considering the comodulation of reuniens and hippocampal activity prior to mPFC spindles, we asked whether hippocampal oscillations directly influence firing activity in the reuniens. A representative example of LFP recordings of the hippocampus, mPFC, and reuniens during REM sleep is shown in Figure 4A. The wavelet spectrogram of the hippocampal trace demonstrates prominent theta-band activity during REM with periodic locking of reuniens multi-unit activity to hippocampal theta cycles (Figure 4B). Analysis of reuniens firing around the peak of hippocampal theta cycles, shows significant phase-locking of reuniens single-unit activity to the theta rhythm during REM (Figure 4C-F).

Figure 4

Power spectra of mPFC, reuniens, and hippocampal signals during REM are characterized by prominent peaks in the theta range, with coherence between mPFC-hippocampus, reuniens-hippocampus, and mPFC-reuniens recording pairs in the theta range during REM (Figure 4 G, H). mPFC-reuniens recordings were further comodulated by a slow, irregular potential observed in Figure 4A and in the broad low frequency peak (Figure 4H, bottom). Cross-correlations and relative time lags of reuniens-hippocampus and mPFC-hippocampus recording pairs during REM are illustrated in Figure 4I. mPFC signals lagged hippocampal signals by a longer delay than reuniens signals lagged hippocampus (Figure 4J, top) and peak cross-correlation coefficients were significantly higher for reuniens-hippocampus recording pairs than for reuniens-mPFC (Figure 4J, bottom).

Next, we analyzed reuniens firing activity during SWRs in NREM. Figure 4K shows representative LFP recordings of the hippocampus and reuniens multiunit activity during SWR events in NREM sleep. In 20.4 % of reuniens single-units (10/49), firing rate increased significantly following the onset of SWRs in the hippocampus, compared to shuffled events in the same period (Figure 4L, 150 ms post-SWR onset). Firing rates decreased significantly in 32.6% (16/49) of single-units following SWR onset The average thalamic LFP around the onset of the SWR was characterized by a negative, 250 μV deflection coincident with increased reuniens firing.

Hippocampal-thalamocortical model

To explore the parameters that control interregional synchrony during NREM sleep, we constructed a neural mass firing rate model consisting of two cortical networks (N 1 & 2), representing mPFC layers 5 and 6, respectively, and two hippocampal (CA1-CA3) and three thalamic networks (TRN, MD, and Reu, Figure 5A). Slow waves and sharp waves are spontaneously generated in the isolated CA3 and cortical networks respectively as a result of adaptation and excitatory-to-excitatory recurrent connections. SWRs are reproduced in the CA1 network during the active phase of the sharp wave due to recurrent connections between the CA1 excitatory and inhibitory neurons with the firing rates sharply dependent on the membrane potential in this firing rate model ³⁰. The oscillations in the range of spindle frequency are spontaneously generated in the bursting, recurrently connected MD-TRN network. The beginning and end of the spindles are controlled by cortical inputs in the thalamocortical model with long-range bidirectional projections between cortical and thalamic networks. In the full hippocampo-cortico-thalamic model, the CA1 is directly connected to the cortical network while reuniens, which is bidirectionally connected to both CA1 and cortical networks, mediates the cortical to CA1 influence. With large enough strength of multi-regional projections, the model can spontaneously reproduce the experimentally observed coupling between the slow waves, spindles, and SWRs (Figure 5B). We found that reuniens slow waves occur significantly later than cortical slow waves (Figure 5C). In contrast, thalamic spindles occur significantly earlier than cortical spindles (Figure 5D).

Figure 5

Figure 5E shows the average of cortical slow wave trough-locked time-frequency representations of the membrane potential of mPFC excitatory neurons and reuniens neurons. The slow wave-spindle coupling angles are nonuniformly distributed for both cortical and reuniens spindles indicating strong coupling with slow waves. In quantitative match to the experimental data, reuniens spindles occur mainly around the slow wave peak (active phase) while cortical spindles mainly occur around the slow wave trough. The slow wave-spindle coupling strength is also higher for reuniens than cortical spindles (Figure 5F). Slow waves in the thalamocortical network are also accompanied by large negative deflections and a subsequent peak in CA1 network (Figure 5G). Phase locking analysis revealed that, consistent with the experimental data, SWRs are phase-locked to slow waves, occurring mainly before slow wave peaks (active phase, Figure 5H). The bidirectional projections between reuniens and CA1 coordinate the communication between these two networks so that the membrane potential of the reuniens neurons increases after SWRs (Figure 5I). Reuniens spindles were also preceded by an increase in the SWR rate (Figure 5J, the peak SWR rate occurs 0.32 s before the spindle onset).

Bidirectional projections between the reuniens and CA1 control multiregional interactions

To investigate how bidirectional projections between the reuniens and CA1 control multiregional interactions in the model, we explored hippocampal-thalamocortical coupling during sleep rhythms as a function of the strength of bidirectional projections between reuniens and CA1. The correlation between reuniens and CA1 membrane potentials during SWRs increased by increasing the strength of CA1 to REU projections (Figure 6A). In addition, the SWR activity before spindle onset is an increasing function of bidirectional projections between CA1 and reuniens (Figure 6B). Next, we studied the effect of the strength of reuniens-CA1 bidirectional projections on CA1 activity during cortical slow waves. The peak-to-trough magnitude of CA1 deflections around the cortical slow wave also increases by increasing the reuniens-CA1 bidirectional projections (Figure 6C). The slow wave mean phase at the time of SWR occurrence decreased by increasing the strength of the bidirectional projections between reuniens and CA1 so that SWRs mostly occur before and after slow wave peak for strong and weak projections, respectively (Figure 6D). These results suggest that for strong reuniens-CA1 projection, the depolarizing input from reuniens during down-to-up transition is large enough to induce the SWRs in CA1 before the slow wave peak. Further, slow wave phases at the time of SWR occurrence are more nonuniformly distributed for stronger CA1-reuniens projection, indicating stronger phase locking of SWRs to slow waves.

Figure 6

Discussion

This study examined the functional relationship between mPFC, reuniens, and hippocampal activity in cats using electrophysiological and computational approaches. We found that reuniens activity during natural sleep was strongly modulated by hippocampal oscillations, and that reuniens-mPFC spindles were preceded by increased probability of hippocampal SWRs and ripple power. Electrical stimulation of reuniens generated orthodromic and antidromic spiking in mPFC neurons, indicating bidirectional thalamocortical circuitry between the mPFC and reuniens. Slow waves in reuniens lagged prefrontal slow waves but spindles in reuniens led prefrontal spindles, suggesting a bidirectional, oscillation-dependent dialogue between reuniens and the mPFC. Spindle amplitude and SWR occurrence was modulated by the phase of the prefrontal slow oscillation, consistent with the idea that prefrontal slow waves control the timing of thalamic and hippocampal events. Overall, the study provides physiological confirmation to the hypothesis that reuniens is a core mediator in hippocampo-cortical dialogue.

Hippocampal influence on the activity of prefrontal neurons is a central mechanism in sleep-dependent memory consolidation ^31–33. Previous research has shown that response latency in the hippocampo-prefrontal pathway exhibits substantial variability ³⁴, indicating the involvement of both monosynaptic and polysynaptic mechanisms in prefrontal response to hippocampal input. Our findings support this notion, as we observed inconsistent mPFC responses to hippocampal stimulation, occurring later than those elicited by reuniens stimulation. Notably, the same cell in mPFC could respond or not respond to hippocampal stimuli (Figure 1 C, E) while LFP responses were less variable. This was possibly due to inhibitory mechanisms such as shunting/feedforward inhibition or cholinergic, monoaminergic, or thalamic gating mechanisms ^35–40.

While the functional impact of the hippocampo-prefrontal pathway has been extensively studied, there are fewer studies on the effects of reuniens on the prefrontal cortex ²⁶. In our study, half of recorded prefrontal cells responded to reuniens stimulation, with a mean response latency of 6 ms and with antidromic spikes with a latency of approximately 3 ms, pointing to bidirectional functional connectivity between the two structures. Similarly, electrical stimulation of the reuniens/rhomboid nucleus in rats ²⁶ elicits prominent evoked potentials in prefrontal cortices that have putatively monosynaptic latencies of 4.5 ms. We found that mPFC responses to reuniens stimulation were qualitatively different from responses elicited by hippocampal stimulation, with significant differences in latency, amplitude, slope of EPSPs, and duration of hyperpolarization. Reuniens-evoked responses occurred earlier and were more consistent than hippocampal-evoked responses, suggesting higher response fidelity in the reuniens-prefrontal pathway than in the hippocampo-cortical pathway. Both pyramidal cells and inhibitory neurons mediate the hippocampo-prefrontal response, with the early component thought to be primarily mediated by AMPA receptors ^34,41 and the early hyperpolarization was likely mediated by GABA_A and the late component by either GABA_B conductances ³⁴ or it was mediated by disfacilitation^42,43. Interneurons also respond to hippocampal stimulation, with a shorter response latency than pyramidal cells, suggesting that the hippocampus can exert feedforward inhibition in the medial prefrontal cortex ⁴⁴. AMPA receptors and GABA inhibition may also mediate the primary and hyperpolarizing responses to reuniens stimulation in medial prefrontal cells as in the mediodorsal-prefrontal pathway ⁴⁵.

Similar to hippocampal afferents to the mPFC, the reuniens nucleus extends glutamatergic afferents to the prefrontal cortex ^46,47 and exerts a strong influence on the prefrontal activity.

Reuniens is believed to play a critical role in connecting the prefrontal cortex and hippocampus due to its bidirectional connections with both structures ²⁶. The prefrontal cortex and hippocampus are both essential for memory processing ³¹, and previous research indicates the reuniens may synchronize the activity of these two structures ⁴⁸. Demonstrations that the nucleus reuniens exerts excitatory actions at the CA1 ^49,50 and also at the prefrontal cortex ²⁶

(intracellularly confirmed and described in our study) provides supporting evidence for this claim. We also found that hippocampal stimulation evoked strong responses in cortical areas coincident with the recording site of antidromic responses to reuniens stimulation (Figure 1). These results support the idea that the reuniens nucleus, hippocampus, and prefrontal cortex form a functional loop, with the prefrontal cortex acting as a mediator of the hippocampo-prefronto-thalamic relay.

In NREM sleep, the coupling of hippocampal SWRs to cortical slow waves and thalamocortical spindles is thought to be the main network mechanisms by which the brain promotes the consolidation of short-term memory into long-term storage ^5,27. However, the mechanism that govern the emergence of this interregional synchrony remain unclear. We found that reuniens single-unit activity exhibited sporadic locking to hippocampal theta cycles and was strongly modulated by hippocampal SWRs. Considering the strong influence of reuniens on prefrontal targets, the hippocampal modulation of reuniens activity points to a thalamic pathway for the coupling of sleep rhythms. Although the minimal requirement to generate spindles is the reticular thalamic nucleus (TRN) ^51,52, our finding of increased SWR probability and ripple power prior to the onset of thalamocortical spindles suggest that hippocampal drive through reuniens may initiate the thalamocortical cascade of events that leads to spindles.

In first-order nuclei of the thalamus, cortical firing can trigger spindles by exciting hyperpolarized neurons of the TRN, eliciting high-frequency spike bursts⁵³. In turn, TRN bursts can drive IPSPs in thalamocortical neurons, eliciting low-threshold calcium spikes ⁵¹ and developing the spindle primarily through recurrent inhibitory-excitatory interactions between the TRN and thalamocortical neurons⁵³. Cortical firing also contributes to spindle synchronization and the maintenance of the spindle through its influence on the TRN ⁵⁴. Thus, cortical activation after the down state can initiate thalamocortical spindling and lock the timing of spindles to the phase of the slow oscillation. However, less is known about spindle generation in higher-order thalamic nuclei, especially in thalamic nuclei with limbic circuitry such as reuniens. In general, although first order thalamocortical cells display inhibitory, spindle-like activities during cortical active state, higher order thalamocortical cells show persistent depolarizing states, excluding their possibility to generate spindle-related rebound spike bursts^55,56. The bidirectional connections of the reuniens to neocortical and archicortical structures make it difficult to discern the possible triggering event that initiates the spindle in this nucleus. We found that hippocampal SWRs trigger thalamic firing and precede the onset of reuniens and mPFC spindles, which points to SWRs as possible candidate events for spindle initiation. In this regard, the strong excitatory drive emerging form the hippocampus during SWRs may represent an ontogenetically different mechanism for the initiation of spindles in the limbic circuitry. However, the generation and maintenance of spindle cycles depend on burst firing from inhibitory structures such as the TRN. To date, it is not yet clear whether the hippocampus influences inhibitory thalamic structures such as the TRN or the zona incerta, although projections from subiculum to the rostral TRN were previously described⁵⁷. We propose that hippocampal SWRs can initiate a cascade of thalamocortical events, beginning with the excitation of reuniens, which in turn excites the mPFC. This mPFC excitation leads to spindle activity generated through cortico-thalamic feedback from the mPFC to other thalamic nuclei, such as the MD nucleus.

Our computational model demonstrates the presence of cross-regional interactions between the prefrontal cortex and its thalamic partners. The model confirmed that thalamic slow waves occur significantly later than cortical slow waves, whereas thalamic spindles occur significantly earlier than cortical spindles, recapturing our experimental findings and pointing to the cortical origin of slow-wave activity and thalamic origin of spindles. The model predicts that, while CA1 to reuniens connection strength plays a critical role in the modulation of reuniens activity by SWRs, bidirectional connections between reuniens and CA1 are critical for SWR-dependent triggering of spindles. Similarly, bidirectional reuniens and CA1 connections are important for transmitting cortical slow wave activity to the CA1 region. With strong bidirectional connections, the SWR occurred primarily before the slow wave active phase peak and with weak connections, the SWR tended to occur after this peak. Given our experimental findings (Figure 3, J-K) the model suggests that the effectiveness of bidirectional connections between reuniens and CA1 is particularly strong during NREM in cat.

Materials and Methods

Experimental model and subject details

Experiments were carried out in accordance with the guidelines of the Canadian Council on Animal Care and were approved by the Committee for Animal Care of Université Laval and are consistent with ARRIVE guidelines. In 20 adult cats (16 anesthetized during recording sessions and 4 chronically implanted, not anesthetized during recording sessions), aged 1-1.5 years, electrophysiological recordings were obtained from the mPFC, midline thalamus, hippocampus and/or various cortical areas including somatosensory, motor, auditory, associative, and visual cortex. Prefrontal intracellular recordings (n=20) and LFP-array recordings (n=9) were obtained from anesthetized animals and spike/local field potential (n=4) recordings were obtained from non-anesthetized animals during sleep and wake cycles.

Surgeries

Anesthetized Animals

For acute intracellular recordings, cats were anesthetized with a mixture of ketamine and xylazine (10-15 and 2-3 mg/kg i.m., respectively) and maintained in anesthetic state during the experiment by supplemental doses of ketamine (5 mg/kg). Body temperature was maintained at 37°C via a feedback-controlled heating pad and the heart rate was continuously monitored (90-110 bpm). All pressure points and incision sites were infiltrated with lidocaine (0.5%).

Animals were paralyzed with gallamine triethiodide and artificially ventilated to an end-tidal CO₂ concentration of 3.5-3.8%. Recording stability was maintained by cannulation of the cisterna magna, bilateral pneumothorax, hip suspension, and by sealing craniotomies with agar solution (4% in 0.9% saline).

Non-anesthetized animals

For chronic spike/LFP recordings, cats were implanted with indwelling electrodes and a head-restraint system according to surgical procedures described previously in Timofeev et al. ⁵⁹ and Chauvette et al. ⁶⁰. Briefly, cats were pre-anesthetized with an intramuscular injection of ketamine (3 mg/kg), buprenorphine (0.02 mg/kg), and dexmedetomidine (5 µg/kg). The incision site was shaved, and the cats were intubated for gaseous anesthesia. Lidocaine (0.5%) and bupivacaine (0.25%) was injected at the site of the incision and in all pressure points where the head contacted the stereotaxic frame. Four stainless-steel electrodes were implanted under stereotaxic guidance in the mPFC (from the interaural line: AP 23 to 24 mm, ML 1 mm, DV 5 mm, θ=10°) according to the Reinoso-Suarez stereotaxic atlas of the cat brain ⁶¹. Additionally, single stainless-steel electrodes were implanted in deep cortical layers (1.2 mm from cortical surface) in the ipsi-or contralateral hemisphere in the marginal gyrus [visual cortex (areas 17 and 18)], suprasylvian gyrus [associative cortex (areas 5, 7, and 21)], ectosylvian gyrus [auditory cortex (areas 22 and 50)], postcruciate gyrus [somatosensory cortex (area 3)], precruciate gyrus [motor cortex (areas 4 and 6)]. Tetrodes, constructed out of four twisted microwires (platinum-iridium, 12.5 µm, gold-plated to ∼200 kΩ), were implanted into the midline thalamus (from the interaural line: AP 11.5, ML 0.5, DV 1.5, θ=10°) and into the CA1 region of the hippocampus (from the interaural line: AP 7, ML 10.5, DV -3.5, θ=10°).

For electro-oculography (EOG), a silver electrode was implanted through the orbital bone and for electromyography (EMG), two teflon-insulated stainless-steel electrodes were implanted into the neck muscle. All recordings were referenced to a silver electrode fixed to the skull above the cerebellum. To allow for subsequent head-restrained recordings, head-fixation holders were also attached and encased into the dental cement. Throughout the surgery, the body temperature was maintained at 37°C using a water-circulating thermo-regulated blanket. Heartbeat, oxygen saturation and blood pressure were continuously monitored using a pulse oximeter (Rad-8; MatVet), and the level of anesthesia was adjusted to maintain a heartbeat at 130-150 per minute. Lactate Ringer’s solution was given intravenously (5 ml/ kg/h, i.v.) during the surgery. After the surgery, cats received meloxicam (0.05 mg/kg) once a day for 3 days.

Electrophysiological recordings

Acute recordings in anesthetized animals

Intracellular recordings of the mPFC were obtained in vivo using glass micropipettes pulled on a vertical puller (Narishige PP-830, Tokyo, Japan) from borosilicate capillary tubes (WPI; P-97, Sutter Instrument) and filled with 3 M potassium acetate (DC resistances of 30-60 MΩ). The prefrontal region around the cruciate gyrus was targeted according to the stereotaxic cat atlas of Reinoso-Suárez ⁶¹ (from the interaural line: AP 23 to 24 mm, ML 1 mm, DV 5-14.5 mm). The micropipette was advanced dorsoventrally until a cell was penetrated, confirmed in situ by electrophysiological observations and classification. A high-impedance amplifier (Neurodata IR-283 amplifiers; Cygnus Technology) with active bridge circuitry was used to record the membrane potential (sampled at 20 kHz) and to inject current into neurons.

Stimuli (0.2 ms, 0.3-1.5 mA) were delivered in single pulses at 1 Hz through bipolar (coaxial) electrodes implanted in the ventral midline thalamus (from the interaural line: AP 12, ML 0.0, DV 1.0) and hippocampus (from the interaural line: AP 7.2, ML 10.5, DV -4.5). Bandpassed (0.1 Hz-10 kHz) field potentials were recorded through the same bipolar electrodes at 20 kHz sampling.

To study the topography of prefrontal responses to thalamic and/or hippocampal stimulation, an array of 4 tungsten microelectrodes (9-12 MΩ) was positioned in the mPFC in 9 animals. The array was positioned close to the midline in the pericruciate gyrus sigmoidus and advanced ventrally in 1 mm steps. At each step, thalamic or hippocampal stimuli (0.2 ms, 0.3-1.5 mA) were delivered at 1 Hz or in 5-pulse trains every 2 seconds at 5Hz, 10Hz, 20 Hz and 100 Hz through bipolar (coaxial) electrodes implanted in the ventral midline thalamus (from the interaural line: AP 12, ML 0.0, DV 1.0) and hippocampal formation (from the interaural line: AP 7.2, ML 10.5, DV -4.5). Antidromic events were detected according to methods described previously in Lipski ⁶².

To label stimulation sites, electrolytic lesions were made through thalamic and hippocampal electrodes via 0.75 mA, 1 s current pulses, delivered every 3 s for 5 minutes.

Chronic recordings in non-anesthetized animals

For non-anesthetized cats, a post-operative recovery time of 1 week was maintained prior to the first recording session. Cats were trained over 2-3 days to remain in head-restrained position for 2– 4 h and cycle through periods of quiet wakefulness, NREM, and REM sleep.

Recordings were obtained over two weeks following the post-operative recovery period. All recordings were conducted within a Faraday chamber using AM 3000 amplifiers (A-M Systems), bandpassed 0.1 Hz to 10 kHz with a 1k gain. All signals were sampled at 20 kHz and digitized with PowerLab (ADInstruments - Data Acquisition Systems for Life Science, RRID:SCR_001620).

Histology

At the conclusion of experiments, animals were deeply anesthetized and perfused intracardially with 4°C saline (0.9%) followed by 4% paraformaldehyde (PFA) in 0.1 M phosphate buffer in saline (PBS). Brains were extracted and cryoprotected by immersion in several sucrose-PFA solutions of increasing sucrose concentration (10, 20 and 30% sucrose in PFA over 1 week).

Brains were then sectioned in a freezing microtome (80 or 100 µm thickness), mounted on gelatinized glass slides, dried, and subsequently labelled with DAPI (300 nM, Invitrogen D3571) or stained with cresyl violet by standard Nissl-staining procedures. The coverslipped sections were visualized in brightfield (Nissl) or fluorescence (DAPI) microscopy and sections with electrode tracks were photographed in 10 and 20 x magnification.

Analysis

Signals were analyzed using Igor Pro 4.0 (Wavemetrics, Inc.), built-in and custom routines in Spike2 (Cambridge Electronic Design, Ltd.) and custom-written and open-source code in MATLAB (Mathworks Inc.)

Intracellular data analysis

Stimuli were detected manually, and minimal binomial smoothing was performed to diminish noise on occasion. The latency of orthodromic responses was calculated as the time when the response reached 10% of the maximal amplitude of the sigmoidal fit of the rise of the EPSP. The amplitude of the EPSP was calculated as the maximum of this sigmoidal fit, and the slope of the EPSP calculated as the maximum of the differentiation of this fit. A horizontal line was drawn across the response from the base of the early EPSP and the interception of this line by the hyperpolarizing potential was used to calculate the duration of this latter component.

Responses to sweeps of five pulse stimulation trains were averaged with an IgorPro algorithm (n=30 per stimulation protocol). The EPSP response magnitude elicited by the stimulation of the reuniens nucleus was computed by computing the area below the curve. Average response areas evoked by pulses 2 through 5 of the train were subsequently normalized by dividing by the response evoked by the first pulse.

Wake, NREM and REM detection

Detection of wake, NREM and REM sleep was based on automated analysis of the mPFC and EMG signals, followed by manual REM scoring using mPFC, EMG and EOG data. For automated detection, the Spike2 script RatSleepAuto, based on rodent recordings ⁶³, was adapted to cat mPFC recordings. Briefly, the mPFC signal was downsampled to 100 Hz and bandpassed in the delta range (1-4 Hz). The mPFC, EMG and the delta-filtered mPFC signals were then squared and the mean of each signal was calculated in 5 second bins. Each 5-second bin was classified as NREM if the mPFC and delta signals exceed the mean + 1 standard deviation (SD) of the signal and if the EMG fell below mean - 1 s.d. All bins other than NREM were initially classified as wake. In the subsequent manual correction, “wake” bins with EMG below mean – 1 s.d. and EOG above mean + 1 s.d. were classified as REM. In the final stage, the 5-second bins were combined into groups of 4 to create 20 second epochs of wake, REM or NREM according to the majority of 5-second bins that comprised the 20 second epoch. Epochs without a majority of either NREM, REM or wake were labelled as “Transition”.

Slow wave, spindle and SWR detection

Slow wave, spindle and SWR detection was based on methodology described previously in Alizadeh et al. ⁶⁴ and respectively in Zugaro et al. ⁶⁵, Klinzing et al. ⁶⁶ and Levenstein et al. ⁶⁷. The detection protocol for all three events was restricted to periods of NREM sleep.

The Freely Moving Animal (FMA) toolbox⁶⁵ was used for slow wave and spindle detection. For slow waves, mPFC and reuniens LFP signals were bandpass filtered between 0.5 to 4 Hz using a FIR filter with a filter order corresponding to 3 cycles of the low frequency cut-off. Tentative slow wave events whose positive peaks were greater than 2.5 standard deviations of the signal and whose difference between positive and negative peaks was greater than 4 standard deviations were detected. Slow waves were finally defined as periods when consecutive positive-to-negative zero crossings of the signal filtered in 0.5-2 Hz occurred between 0.5 and 2 seconds of each other.

For spindles, mPFC and reuniens LFP signals were bandpassed in the spindle frequency band (7-15 Hz) using finite-impulse-response (FIR) filters from the EEGLAB toolbox ⁶⁸. The instantaneous amplitude of the filtered signal was then computed from its Hilbert transform and the results were smoothed using a 300 ms Gaussian window. Potential spindles were detected as periods when this signal crossed 3 standard deviations (SD) above its mean for 0.5-3 s A threshold of 2.5 SD above the mean was set to detect spindle onset and offset. Events that occurred within 0.5 seconds of each other were merged into a single event. Events longer than 0.5 s and shorter than 3 s were defined as spindles. The largest peak in spindle amplitude was defined as the time of the spindle event.

For SWRs, hippocampal LFP signals were bandpassed in the ripple frequency band (150-200 Hz, FIR filter from the EEGLAB toolbox). The root mean square (RMS) of the filtered signal was then computed and the results were smoothed using a 50 ms Gaussian windows. SWRs were defined as periods when the smoothed RMS crossed 4 standard deviations above the median and stayed crossed for at least 30 ms. The largest peak in SWR amplitude was defined as the time of the SWR event.

Phase-amplitude coupling

Computations of phase-amplitude coupling (PAC) were based on methodology published previously in Alizadeh et al.⁶⁴, Dehnavi et al.⁶⁹ and Staresina et al.⁷⁰ Time-frequency representations (TFRs) of every spindle/slow wave event were obtained using the mtmconvol function in the FieldTrip toolbox⁷¹. Spindle power around the slow wave was estimated using sliding Hanning tapered windows (10 ms steps) with a variable length that included five cycles. To compute coupling, TFRs were then normalized as percentage change from pre-event baseline (−3.5 to −2.5 second before event) for all events. To calculate the modulation of spindle power by slow wave phase, TFR bins around slow waves were averaged across the spindle frequency range and then filtered in the slow wave frequency range. Phase values of this signal were calculated from their Hilbert transform.

Synchronization index and coupling strength

We defined a synchronization index (SI) for each event:

where m is the number of time points, 𝜃_{slow wave}(𝑗) represents slow wave phase value at time j and 𝜃_spindle(𝑗) represents phase values of spindle power fluctuations at time j. Preferred phase was estimated in the interval -0.25 to 0.25s around slow wave peak. The strength of phase-amplitude coupling was estimated from the absolute value of the SI with 1 representing the most consistent phase shift between spindle amplitude and slow wave event. The phase of SI represents the amount of phase shift between the spindle amplitude and slow wave event.

Spindle-SWR coupling was computed using a similar procedure. Briefly, the normalized TFRs in the SWR frequency range around each spindle event (baseline defined from -2.5 to -1.5 second before spindle event) were filtered in the spindle frequency range and the SI was defined as the vector mean of unit vectors, each showing the phase difference between SWR amplitude and slow wave event at different time points around the spindle peak.

Average wavelet transforms

Wavelet transforms of the LFP were calculated using the continuous one-dimensional wavelet transform, estimated by Morse wavelets defined in the cwt function, available in the MATLAB wavelet toolbox. Average wavelet transforms around spindle onsets (Figure 2L) were computed using custom-written MATLAB code. First, 5-second LFP epochs, centered on spindle onsets were extracted from the continuous recordings and downsampled to 250 Hz. A filter bank, based on the known epoch length (5 s) and sampling frequency (250 Hz), was created with voices per octave set to 20 and frequency limits set to 0.1 to 50 Hz. Using the filter bank, wavelet transforms of each epoch around a spindle onset were calculated and stored in a 3-dimensional array corresponding to time, frequency, and spindle event. To obtain the magnitude of the signal, the absolute value of the wavelet transform was then computed. To obtain the average wavelet transform of the LFP signals around the spindle event onsets, the mean of the magnitude matrices along the third dimension (spindle events) was calculated.

Spike detection

For spike detection, reuniens LFP signals were bandpass filtered in the 300-8000 Hz range using finite impulse response (FIR) filters from the EEGLAB toolbox. Epochs of ± 250 ms around times when the signal crossed ± 16 SD were considered as artifacts and removed from the analysis. Spikes were detected using UltraMegaSort2000 ⁷² with a negative threshold of -4SD. During a shadow time of 0.75 ms after the threshold crossing, no new event was detected to avoid multiple threshold crossings of one event. A window of 1.5 ms was considered around each detected spike from 0.6 ms before the threshold crossing time to 0.9 ms after this time.

Spike sorting was based on a method described previously in Souza et al.⁷³. We first extracted spike features for each channel. The features consisted of peak-to-trough value, trough value and also an index which was set to one for the channel with maximum trough amplitude and zero for the other channels. We used principal component analysis (PCA) to combine 12 features (3 features times 4 channels) and selected the two first components. Gaussian distributions were fitted to PCA components with gaussian mixture model (GMM) algorithm. Each distribution was determined as a potential spike cluster. We further limited the distribution by removing spikes with large Mahalanobis distance from the center of the distribution. We removed, combined, and limited the distributions manually by visual inspection of the plots showing PCA components versus each other and the plots showing the features for different channel pairs (Figure 2 – figure supplement 1). For multi-unit analysis, spikes were detected according to threshold crossing and template matching algorithms in Spike2⁷⁴. Briefly, fast signals (∼1-3 ms) that exceeded the threshold were first detected as putative spikes and timestamped. Thresholds were manually determined according to the signal-to-noise ratio of spikes in each recording. Next, a temporary spike template was formed according to the shape of the putative spike and a “template width” was estimated (twice the mean difference between the template and the spike that created it). Each new spike was then compared against the template and added to it if the spike’s sample points fell within the template width. The template was modified with the addition of each new spike up to a maximum of eight spikes after which the template was checked against previous templates. If a match was found, the temporary template was merged with the existing template. Otherwise, a new confirmed template was generated.

Peri-event histograms

To analyze reuniens single-unit activity around the onset of SWRs, peri-event spike histograms (PETHs) were calculated for a 1-second period around each event onset and the results were smoothed with a 10 ms Gaussian window. For comparison, surrogate PETHs were obtained by randomly shuffling SWR events detected within 1 second windows. This procedure was repeated 100 times. True-even PETHs were z-score normalized to the firing rate of surrogate PETHs. Single-units with activity greater than 2 SDs and smaller than -2 SDs in the 150 ms window following SWR onset were considered activated and suppressed, respectively.

Phase-locking analysis

For phase-locking analysis, we calculated spindle phase at the time of single-unit spikes occurring within ± 0.25 s around spindle peak (Figure 2H), slow wave phase at SWR time (Figure 3K) and theta phase at single-unit spike time occurring within ± 0.5 s around theta peak (Figure 4C). First, mPFC and reuniens LFPs were filtered in the frequency band of interest (spindle, slow wave, or theta) and the phase of the filtered signal was extracted from its Hilbert transform. These phase values were used to obtain phase histograms and calculate circular means.

Computational model

Our model represents a minimal neural mass model that describes the network architecture required for generating three key oscillatory activities during NREM: the slow oscillation, spindles, and SWRs. The hippocampal-thalamocortical model represents CA1-CA3 networks, two cortical networks (CX, CX’), and three thalamic networks consisting of mediodorsal nucleus of the thalamus (MD), thalamic reticular nucleus (TRN) and nucleus reuniens (REU). The details of the CA1-CA3, MD-TRN and cortical networks were described previously ^30,75,76. In brief, in this neural mass firing rate model, each of cortical and hippocampal networks consists of one identical group of excitatory and one identical group of inhibitory neurons. V_k^m and N_k^m indicate the membrane potential and the number of type k (either e, excitatory or i, inhibitory) neurons of network m, respectively. The firing rates of hippocampal and cortical neurons as well as the REU neurons are sigmoid functions of the membrane potential of the neurons with the sharpness and threshold of g^m, V^*m respectively:

where, r₁= 70 Hz (r₀ = 0.1 Hz) is the maximal (minimum) firing rate. T-type calcium currents is modeled only for MD and TRN neurons by defining a bursting variable, u_m so that these neurons show burst mode in addition to the tonic mode:

Here show the sharpness of the firing rate dependence on the membrane potential, the maximum firing rate, and the threshold potential during tonic (burst) mode. The sharpness of switch between the two modes is controlled by L_k (a positive parameter).

The probability and the strength of connections from neuron k of network m to neuron j of network n are shown by respectively. m and n can be CA1, CA3, CX, MD, TRN or REU. k and j can take e or i (representing excitatory and inhibitory neurons, respectively). For the short-range connections since m and n are equal, only one of them is shown. E-E connections of cortical and hippocampal excitatory neurons are subject to dendritic spike frequency adaptation so that is a sigmoid function of an adaptation variable, c^m with sharpness and threshold of g_c and c*respectively:

Here is the maximal synaptic strength. The full model is described by 17 rate equations:

The nonlinear equations of the model were solved using the fourth-order Runge–Kutta method with a time step of 1 ms by MATLAB R2017b software.

Quantification and statistical analysis

Statistical analyses were performed using the MATLAB Statistics and Machine Learning Toolbox (Mathworks Inc., Natick, MA, USA). Phase analysis was performed using the circular statistics tool box developed by Berens ⁷⁷. The Rayleigh test was performed to determine the uniformity in the phase preference of events. The Watson-Williams test was used to test the inequality of phases. Polar plots showing phases, and mean phases were plotted using MATLAB.

Intracellular data was analyzed using JMP 5.01 software. Normality was assessed by normal quantile plots and by the Shapiro-Wilk test. Homogeneity of variance was assessed by Levene’s test. Confidence intervals were set at 95% (p<0.05). ANOVA was used when data was normally distributed and variance homogeneous, Welch ANOVA when variance was non-homogeneous. To determine group differences in normally distributed data, Tukey’s range test was used.

Otherwise, differences were assessed by non-parametric testing (van der Waerden test). Standard deviations are reported as SD and the mean absolute deviation is reported as MAD.

Data and code availability

Further information and requests for resources and reagents should be directed to and will be fulfilled by the Lead Contact, Igor Timofeev (igor.timofeev@fmed.ulaval.ca). Original electrophysiological recordings have been deposited at Mendeley Data and are publicly available as of the date of publication.

The following data sets were generated

Basha, Diellor; Timofeev, Igor (2022), “Natural sleep electrophysiology: mPFC, thalamus, hippocampus”, Mendeley Data, V2, doi: 10.17632/4×27jjp9mv.258. In the public domain at Mendeley Data, https://data.mendeley.com/datasets/4x27jjp9mv/2.

MATLAB code generated for this study is available at https://github.com/DiellorBasha/MATLAB_for_Electrophysiology and https://github.com/azarmehri/brain-network-interactions.

Tetrodes and stainless-steel electrodes manufactured in-house for this study are available from the lead contact with a completed materials transfer agreement. This study did not generate new unique reagents. Any additional information required to reanalyze the data reported in this paper is available from the lead contact upon request.

Key Resources Table

Contribution

Conceptualization, Methodology, Investigation, Validation, Writing - Review & Editing, Supervision, Project administration, Funding Acquisition.

Competing interests

No competing interests declared.

Centre de recherche CERVO, 2601 de la Canardiere Quebec, QC, Canada

Funding

National Sciences and Engineering Research Council of Canada (grant RGPIN-2018-06291) Igor Timofeev Canadian Institutes of Health Research (grants PTJ - 183862) Igor Timofeev

Acknowledgements

We thank Serge Ftomov for technical assistance.