A diversity of interneurons and Hebbian plasticity facilitate rapid compressible learning in the hippocampus

The hippocampus is able to rapidly learn incoming information, even if that information is only observed once. Furthermore, this information can be replayed in a compressed format in either forward or reverse modes during sharp wave–ripples (SPW–Rs). We leveraged state-of-the-art techniques in training recurrent spiking networks to demonstrate how primarily interneuron networks can achieve the following: (1) generate internal theta sequences to bind externally elicited spikes in the presence of inhibition from the medial septum; (2) compress learned spike sequences in the form of a SPW–R when septal inhibition is removed; (3) generate and refine high-frequency assemblies during SPW–R-mediated compression; and (4) regulate the inter-SPW interval timing between SPW–Rs in ripple clusters. From the fast timescale of neurons to the slow timescale of behaviors, interneuron networks serve as the scaffolding for one-shot learning by replaying, reversing, refining, and regulating spike sequences. How does the hippocampus learn with only a single presentation of a stimulus? Nicola and Clopath show how interneuron networks can dilate pre-existing spike sequences during sharp waves into theta sequences for one-shot learning.


Introduction
Shortly after watching a movie, we can recall many scenes in detail and remember long stretches of dialogue. However, as time passes, our memories may only preserve some aspects of the plot and a few memorable lines. This example illustrates how memories operate on separate timescales. On the short timescale, we remember recent events clearly with only a single viewing, while we lose finer details on longer timescales. This observation, among others, led to the two-stage model of memory formation [1,2]. In the initial stage, a labile form of the memory is imprinted onto the hippocampus. This initial acquisition is accompanied by a theta oscillation in the hippocampal Local Field Potential (LFP). The second stage of memory formation occurs during consummatory behaviours such as sleeping or eating. In this stage, strongly correlated activity is initiated in the hippocampus as a Sharp-Wave-Ripple Complex (SPW-R). The spike sequence elicited by the stimulus during waking is replayed in compressed time during SPW-Rs [3]. The SPW-R propagates out of the hippocampus, presumably for long-term memory consolidation in cortex.
This two-stage model requires some mechanism for the initial memory formation. This initial mechanism should handle the fact that we often only see a stimulus once [4]. Furthermore, this hypothesized fast-learning mechanism should allow for compressible spiking, as observed in SPW-Rs. Finally, this mechanism should allow for the reversal of spike sequences. Spike sequences can be compressed in reverse order, even for mice navigating novel mazes in the initial lap [5].
What known hippocampal dynamics allows for learning subject to these constraints? Narrowing our focus to episodic memory formation yields one possible candidate, the Internally Generated Theta Sequence (IGTS) [6,7]. First reported in CA1, these cells partition time during episodic memory tasks [6]. The IGTSs are dependent on the Medial Septum (MS) and the hippocampal theta oscillation for normal functioning [7]. Here, we explore the hypothesis that IGTSs may serve as a compressible backbone to bind information to during single-trial-learning.
Leveraging recent advances in training spiking neural networks [8], we construct networks of phase precessing IGTSs using a dual oscillator model [9,10]. Unlike previous dual oscillator models, we consider a novel implementation: one oscillator corresponds to inhibition from the medial septum, while another oscillator corresponds to interneuron control of spike-times during SPW-Rs. In effect, the interference between these oscillators serves to dilate spike sequences during SPW-Rs to create IGTSs. The IGTS binds new, stimulusevoked sequences in a single-trial to form memory engrams that can be compressed and reversed in time. Our results support a single interneuron mechanism for gating, controlling, transforming and refining spike sequences for rapid learning.

Internally Generated Theta Assemblies with Network-Based Dual Oscillators
The network model consists of three layers which we construct and analyze incrementally: the inputs, the Septal-Hippocampal-Oscillatory-Theta or SHOT network, and the readout layer. The inputs serve as one of the oscillators. The next layer, the SHOT network, is the critical nucleus of the circuit that creates IGTS and SPW sequences. The SHOT network contains the second oscillator. The final layer is the readout layer where sequences must be learned in a single-trial. These to-be-learned-sequences are elicited by an externally applied supervisor.
Thus, given the architecture of this network model, we interpret these components to correspond to the Medial Septum (INP-MS, input layer), hippocampal area CA3 (SHOT-CA3), and CA1 (RO-CA1). Finally, we interpret the supervisor which elicits spiking in the RO-CA1 layer to be the entorhinal cortex (SUP-EC, supervisor). The interference pattern or "beat" between the SHOT-CA3 network oscillator and the INP-MS oscillator creates IGTSs, while removal of INP-MS creates compressed sequences (SPWs). This is our interpretation of an inherently abstract model which we implement in a spiking network. However, this implementation is not unique. For example, the sources of these oscillations could entirely be intra-hippocampal (see [11,12]). Additionally, given that only inhibitory connections are required for sequence generation in the SHOT network, it may be contained in CA1.
To model the oscillatory and inhibitory nature of MS, the SHOT-CA3 interneurons (SHOT-CA3I) receive the INP-MS signal: I M S (t) = I GABA (κ + cos(2πθ M S t)) . (1) The parameter θ M S is the theta oscillation frequency of INP-MS, κ ≥ 1 controls the baseline inhibition levels, and I GABA < 0 controls the amplitude of this oscillation. Next, the IGTSs observed in [6] generate sequential firing-fields on a behavioural timescale, while neurons oscillate at a theta frequency inside the field [6]. The spiking rate within a field is slightly faster than the LFP theta frequency, an effect termed phase precession [6,10]. Unlike previous works [11][12][13], we take the perspective that firing-fields are created by septal-hippocampal interference of oscillators. This is motivated by recent work that supports intra-hippocampal theta oscillations [14][15][16].
Thus, we hypothesized that INP-MS could interact with a recurrent interneuron network that oscillates at a frequency of θ IN T (INT for Interneuron) to generate a robust interference pattern. To avoid exhaustive parameter searches, we used supervised training methods ([8]) to couple the interneurons. We successfully

Theta-Phase and SPW-R Compression in the SHOT-CA3 Network
When INP-MS is on, the SHOT-CA3E fire IGTSs while the SHOT-CA3I fire bursts ( Figure 1D). When INP-MS is off, the SHOT-CA3E also fire sequential bursts when provided with a sufficiently large (but constant) background current ( Figure 1E). The burst sequence was identical to the IGTS, but compressed in time. Thus, we naturally interpret these two network states as the theta and SPW state, respectively.
Sequence compression is a consequence of the dual oscillator model. For example, when INP-MS is off, the non-dimensionalized (Materials and Methods) current is where ψ denotes the phase preference of a neuron with respect to the θ IN T oscillator. The phase preference dictates a neurons firing order in the compressed SPW sequence. When INP-MS is on, the current is . (4) The sum becomes a product of a carrier (C(t)) and envelope (E(t)) function. Note, again, that the currents are dimensionless and lack amplitude terms (e.g. I GABA ). The envelope oscillates on a long timescale ( θ IN T −θ M S 2 ) and controls the IGTS firing-fields. The envelope inherits an identical order of firing from the θ IN T oscillator through its phase, ψ 2 when θ IN T > θ M S . This allows for SPW-R-like compression of the IGTS by removing INP-MS (Supplementary Figure 2).
In order to determine how the SHOT-CA3I control spike sequences, we analyzed the currents that the neurons receive ( Figure 1F). The total current to the SHOT-CA3E displayed characteristics of an interference pattern, as expected. However, the SHOT-CA3I did not. Their current was dominated by a single oscillation, θ IN T . This was surprising as the SHOT-CA3I directly receive INP-MS and generate the θ IN T oscillation, unlike the excitatory neurons.
We investigated this further by decomposing the total current into all sources: the background current, the SHOT-CA3I-synaptic current, and INP-MS ( Figure 1G). For both populations, the inhibitory synaptic currents display an interference pattern. One oscillator is θ IN T . The second oscillator, however, does not correspond to INP-MS directly, but instead to a π phase-shifted or vertically "flipped" version of INP-MS. As the SHOT-CA3I also receive INP-MS directly, INP-MS destructively interferes with its phase-shifted counterpart, allowing SHOT-CA3I to oscillate at θ IN T .
Next, we investigated the structure of the FORCE trained weight matrices. At first, the weights appeared to be random (Supplementary Figure 3A-B). However after sorting neurons according to their phase preferences in θ IN T , a spatial structure emerged (Supplementary Figure 3A-B). We observed recurring light and dark bands in the weight matrix, in addition to random, phase-independent inhibitory connections. The random connections sample from the entire SHOT-CA3I. By receiving these connections, a neuron can detect the population level oscillation in θ M S for destructive interference (Materials and Methods). The banding structure, however, implies interneurons preferentially inhibit interneurons that already fired in θ IN T , thus advancing the θ IN T oscillation forward. This band reappeared in multiple trials of FORCE training (Supplementary Figure 3C-F).
Finally, by using either the INP-MS oscillation (not shown) or the population activity ( Figure 1H,I) as a reference theta oscillation, both excitatory neurons and interneurons precess in phase. This is a classical result of dual oscillator models as the carrier, C(t), controls the spike timing in the IGTS [10]. As the carrier oscillates faster than INP-MS ( θ IN T +θ M S 2 > θ M S ), the neurons precess in phase relative to the population activity (average rate precession ≈ − π 2 /s, SD of 0.09π/s n = 100, SHOT-CA3E, and ≈ −0.91π/s, SD of 0.4π/s, n = 100, SHOT-CA3I).

Firing Field Generation is Robust
For the IGTSs to serve as a compressible backbone for learning, the interference mechanism should be robust. To that end, we retrained the network under different conditions such as with 1) smaller integration time steps (Supplementary Figure 4), 2) parameter heterogeneity (Supplementary Figure 5) 3) larger networks (Supplementary Figure 6C), and 4) different SHOT-CA3E/SHOT-CA3I ratios (Supplementary Figure 6A-B,D-E) and 5) synaptic failure (Supplementary Figure 7). In all cases, the network generated ordered IGTSs that compressed when INP-MS was removed and a weight matrix banding structure (Supplementary Figure  5,6).
Next, we considered how the trained network (in Figure 1) would extrapolate to varied θ M S frequencies. If SHOT-CA3 was operating as a linear oscillator, we should expect that the firing-field periodicity scales like |θ M S − θ IN T | −1 near θ IN T , while the population activity has a period of θ −1 M S . Furthermore, our model predicts that IGTS order reverses when θ M S > θ IN T (Materials and Methods). This regime should also induce phase recession. Finally, we also predict firing fields to emerge near harmonics of θ IN T . To test this hypothesis, we slowly swept θ M S (2 Hz to 12 Hz, Supplementary Figure 8). We found that SHOT-CA3 behaves like a linear oscillator for θ M S tested. This surprising ability to extrapolate to new inputs is due to the SHOT-CA3I ability to create destructive interference of any input to maintain the θ IN T oscillator. The firing fields also reverse near the harmonics of θ IN T . Thus, θ M S directly determines the frequency of the population activity, as predicted from a dual oscillator model. Finally, ramping the background current to the interneurons has no effect on the frequency of the population activity and a small yet statistically significant effect on the period of the interference pattern (Supplementary Figure 9). However, a stronger reduction in the background current destroys θ IN T and causes SHOT-CA3I to lock onto INP-MS (not shown). This mirrors recent results where optogenetic suppression of interneurons reduces phase precession rates but increases intracellular theta amplitude [17]. Ramping the INP-MS oscillation amplitude (I GABA ) has a large and statistically significant effect on the periodicity of the IGTSs and a small yet significant effect on the frequency of the population activity (Supplementary Figure 9). Thus, the network is sensitive to large modulations of the INP-MS amplitude but robust to θ M S variations.

Sharp Waves Can Be Initiated Stochastically
So far, the SPWs induced in our model (when INP-MS is deactivated) are periodic and must be activated by an externally applied current as INP-MS dis-inhibits the SHOT-CA3E. We investigated if SPWs could be initiated stochastically when INP-MS is off by adding recurrent excitation to the SHOT-CA3E. Once again, we will do this incrementally to mechanistically understand how ordered SPWs are generated in the SHOT-CA3 network.
First, in lieu of constant background current, we injected SHOT-CA3I with white noise (Figure 2A-D). All inhibitory connections are identical as before (as in Figure 1). A current pulse to SHOT-CA3I can transiently activate the θ IN T oscillator. In the absence of this current, the θ IN T oscillator is destroyed by the noise as the oscillator is no longer detected in FORCE oscillator components ( Figure 2B), voltage traces ( Figure 2C), and voltage auto-correlations ( Figure 2D). However, the spikes are still biased towards a particular order as noise percolates through the inhibitory weight structure for both populations ( Figure 2B).
Next, we investigated potential mechanisms that could initiate and terminate a SPW in lieu of externally applied currents (in Figure 2B). To initiate a SPW, we densely coupled a subset of 50 SHOT-CA3E (the initiators). The initiators projected to the rest of the excitatory population. To terminate SPWs, all excitatory neurons had spike frequency adaptation currents. As noise induces spiking throughout the network, a random subset of initiators are stimulated to spike ( Figure 2E-G). This triggers a positive feedback loop through recurrent excitation which causes a synchronized burst throughout SHOT-CA3E. However, this triggers the accumulation of adaptation currents which terminate the burst. This induces a refractory period where SHOT-CA3E cannot spike until the adaptation current has decayed. As the θ IN T oscillator has not been activated, the spikes in the bursts have no sequential ordering and are pathologically synchronized.
Finally, we add SHOT-CA3E to SHOT-CA3I connections ( Figure 2H). The SHOT-CA3E self-activate through recurrent excitation and now also activate the θ IN T oscillator ( Figure 2I-J). The θ IN T oscillator protracts the recruitment SHOT-CA3E into the burst, thereby creating an ordered sequence during the SPW rather than synchronized bursts. The adaptation variable then terminates the burst. The resulting inter-SPWinterval distributions can be uni-modal with varying amounts of skew, or multi-modal (Supplementary Figure  10). Thus, SPWs can be activated stochastically with interneurons protracting the recruitment of excitatory neurons into the ordered burst.

Reverse Compressed Replay
Given that the IGTS is compressible, we investigated if reverse compression was also possible, as in other models [3,18]. We discovered three potential mechanisms for reverse compression (Figure 3). The first two mechanisms (harmonic reversion, and mirror reversion) are based on manipulations to the INP-MS frequency to reverse the order during the IGTS. The third mechanism utilizes a secondary population of interneurons (interneuron reversion) to reverse the order of spiking in SHOT-CA3E during the SPW.
In mirror reversion, the frequencies satisfy θ M S > θ IN T which reverses IGTS order ( Figure 3A-E) as E(t) now inherits − ψ 2 from the θ IN T oscillation. Removing of INP-MS reverses burst order during compressed replay ( Figure 3C). However, this is not a plausible mechanism as this results in phase recession and reversed theta sequences ( Figure 3D), seldom observed experimentally [19,20]. Harmonic reversion is possible when θ M S is near harmonics of θ IN T , where IGTS order also reverses. For harmonics smaller than θ IN T , this does not induce phase recession. However, memories encoded onto the IGTS under mirror/harmonic Reversion only display reverse replays as the SPW sequence order is fixed, while both forward and reverse sequences are observed in SPWs [3,5].
However, a more plausible mechanism for reverse compression utilizes another population of interneurons. During SPWs, excitatory neurons receive the non-dimensionalized currents z F OR (t) from the θ IN T oscillator: The phase preference, ψ, orders the neurons in a sequence. Thus, a reverse replay is equivalent to replacing ψ with −ψ.
To perform this transform, we utilize another interneuron current (the Reversion Current) where the total current during reverse replay is now: The reversion current we add to Equation (5) is a synchronized pulse of the θ IN T oscillator. However, each neuron receives more or less of this pulse through the amplitude term 2 sin(ψ). Thus, we implement the reversion current as another population of interneurons that receives SHOT-CA3I ( Figure 3F-J). This reverses the excitatory SHOT-CA3 sequences when INP-MS is removed ( Figure 3H-I, see Materials and Methods). The pulse amplitude 2 sin(ψ) need not be exact for (Supplementary Figure 11). The reversion interneurons may also receive INP-MS to prevent their firing during theta states (see Supplementary Figure 11). Further, phase recession is not observed in Interneuron Reversion ( Figure 3J). Thus, interneurons with selective connectivity can reverse the timing of spike-timing during SPWs in a plausible way.

Online and Immediate Learning Using a Theta Backbone
Compressible theta sequences are one component of fast learning in this model. The other component is a mechanism to write new information onto these sequences. Reading out information falls to a layer of readout (RO-CA1) neurons, which we interpret to be CA1 pyramidal neurons ( Figure 4A).
We derived a local learning rule that allows us to store new sequences in RO-CA1 onto the IGTSs in SHOT-CA3 with only one presentation of the to-be-learned RO-CA1 sequences. The learning rule governing the excitatory weight (ω CA1E,CA3E pre,post ) from a SHOT-CA3E to a RO-CA1 neuron is given by: This learning rule is Hebbian (λ · r CA3E pre (t)r CA1E post (t)), with a forgetting term ( λ · r CA3E pre (t)r CA1E post (t − τ )), where λ dictates the learning rate and r CA3E/CA1E pre/post (t) are the synaptically filtered spikes. The learned weights are interpreted to be Schaffer-Collateral connections. We call this the local Fourier rule based on its origins in function approximation theory (see Supplementary Figure 12, Supplementary Section S1). Like other rules, this rule is phenomenological [21,22].
To test the local Fourier rule, we triggered to-be-learned sequences in RO-CA1 while simultaneously presenting INP-MS ( Figure 4B,C) to SHOT-CA3, (Equation (7)). We interpret the external activity which causes the to-be-learned sequences to come from the Entorhinal Cortex (SUP-EC, SUP for supervisor). Critically, the ω CA3E,CA3I pre,post weights are too weak to trigger spikes in the presence of INP-MS and do not interfere with SUP-EC triggered spikes. However, if INP-MS is removed, RO-CA1 sequences are replayed in compressed sequential order after a single-trial of learning ( Figure 4D).
Surprisingly, this compression does not capture all spikes in RO-CA1. In fact, the local Fourier rule combined with IGTSs creates high-frequency assemblies in RO-CA1 during compressed replay ( Figure 4D,F). This is due to the Hebbian learning rule as only RO-CA1 neurons that fired in phase with SHOT-CA3E are remembered. The frequency of assembly firing is Γ = θ M S + Our result implies that theta oscillations and Hebbian plasticity can induce the formation of high-frequency assemblies, which appear during SPWs. Unlike other mechanisms of high-frequency oscillations [23-25], the assemblies here are learned through Hebbian plasticity operating on phase precessing IGTSs rather than recurrence, or different populations serving as separate gamma and theta oscillators [26]. Stimulating the reversion interneurons triggers a reverse replay of these assemblies ( Figure 4E). More complicated spiking patterns can also be learned ( Figure 4G).
Finally, we show that the local Fourier rule is robust to synchronized supervisors (Supplementary Figure  13  The RO-CA1E were (as in Figure 4) stimulated to spike in a to-be-learned sequence by SUP-EC ( Figure  5C). The local Fourier rule constructed a memory engram in a single trial. The INP-MS was removed and triggered a compressed replay of the now-learned sequence ( Figure 5B-D). The population activity ( Figure 5F) in RO-CA1 displayed ripple-like high-frequency oscillations during compression. The activity in SHOT-CA3, however, demonstrated a sharp increase when INP-MS was removed ( Figure 5G). This is a combined effect of removing INP-MS from SHOT-CA3I and applying a super-threshold current to SHOT-CA3E. Compression resulted in high-frequency assembly formation in RO-CA1E ( Figure 5C,D) while RO-CA1I fired high-frequency synchronous volleys ( Figure 5E). The number of RO-CA1E assemblies was less than the number of bursts triggered in RO-CA1I (8 vs. 10).
The majority of the activity is due to RO-CA1I which are locked to the population activity ripples ( Figure  5H). If we regard the population activity as a proxy (albeit a poor one) for the LFP, then our results mirror PVbasket cells locking to LFP ripples [27,28]. The RO-CA1E assemblies also display some locking to the ripples ( Figure 5H). However, the assemblies precede the ripple peaks as RO-CA1E synapses onto RO-CA1I trigger the synchronized ripples. This is similar to the Pyramidal Interneuron Network Gamma (PING) mechanism [23, 24].
As SHOT-CA3 and Reversion interneurons coordinated spike timing, we investigated if RO-CA1I performed a similar function. To test this hypothesis, we disabled RO-CA1I post-learning while triggering a compressed replay ( Figure 5I) by removing INP-MS. The assemblies in RO-CA1E were no longer segregated to discrete moments in time and overlapped extensively. In fact, the first assembly overlapped with the last assembly during replay. This implies that RO-CA1I may be refining the duration of learned pyramidal assemblies and segregating them to discrete intervals of time ( Figure 5J) with strong, synchronized pulses of inhibition. Indeed, in vivo, it appears that pulses of inhibition dominate during ripples in CA1 over excitation [29]. This may help to transmit information reliably out of the hippocampus. It is possible to manually separate these assemblies, for example (as in Figure 4D). However, adding interneurons in RO-CA1 in these cases will still refine assembly segregation further (as in Figure 5D,I). As with SHOT-CA3 and Reversion interneurons, our results support the idea that interneurons in CA1 also regulate spike-timing by refining the width of learned assemblies in CA1 via a PING-like mechanism [23, 24].

High-Frequency Assemblies Nested on a Theta Oscillation
Assembly discretization occurs during SPW-Rs in our model. However, gamma oscillations (of varying frequency) are often nested with theta and implicated in memory tasks [30,31]. These gamma oscillations are also found during SPW-Rs [3, 32, 33]. Thus, we investigated if high-frequency assemblies would persist in theta states (Supplementary Figure 14).
Providing a stronger background current to SHOT-CA3E results in RO-CA1E firing as phase precessing assemblies nested on theta. This occurs after learning with the local Fourier rule (as in Figure 5). Thus, theta phase segregation and Hebbian learning rules can generate theta-nested assemblies. This is similar to the operations of a working memory buffer (suggested in [34]). Here, the number of assemblies/items in this buffer is The period of the INP-MS oscillation (θ M S ) −1 acts as the buffer's temporal resolution while the IGTS period (θ IN T − θ M S ) −1 determines the capacity of the buffer (Materials and Methods). Our model parameters yield a maximum of 16 items in the buffer.
Here, we consider a more controlled and phenomenological model of stochastic SPW-R activation, rather than recurrent excitation (considered in Figure 2). We disabled INP-MS post-learning and stochastically turned on extra excitation to initiate SPW-Rs ( Figure 5K). This phenomenological model of SPW-Rs initiation caused random replays of the learned memory (seen in Figure 5K). We found that depending on when the SPW-R was initiated, replays could display an error we term phase distortion. In phase distortion, the replay is fragmented into disordered pieces. Phase distortion occurs when SPW-Rs are initiated in the wrong phase of θ IN T ( Figure  5K). To reduce phase distortions, SHOT-CA3I can bias the probability of SPW initiation to the correct phase range of θ IN T . We considered the probability p(t) of initiating a SPW-R as where ψ i is the phase that minimizes phase distortion. The parameter p 0 serves as the background SPW activation rate and p A determines how strong the SPW generation probability oscillates. The oscillatory component is decoded from SHOT-CA3I (see Materials and Methods). A relative refractory period was also incorporated (see Materials and Methods). This strategy to minimize phase distortions has an observable effect on inter-SPW-R-interval distributions. The inter-SPW-R-intervals now display a strong θ IN T modulation ( Figure 5L). This is qualitatively similar to distributions recently reported in [36] which are multi-modal with peaks at 110 ms, 220 ms, 330 ms, and 440 ms. The inter-SPW-R-intervals from [36] may imply some theta modulation in SPW-R initiation with a potential θ IN T ≈ 9.1 Hz. Thus, our results support the idea that inhibition can also influence SPW-R initiation and inter-SPW-R-interval distributions to minimize replay fragmentation.

Discussion
Our memory systems have to handle multiple simultaneous constraints to record, store, and catalog important events. One constraint is that stimuli are often presented only once. Thus, we investigated whether singletrial learning was possible utilizing background spiking activity as a backbone for learning. The background spiking activity took the form of the Internally Generated Theta Sequence (IGTS) [6,7], generated with a dual oscillator model. We interpreted one oscillator as intra-hippocampal that controlled spike timing during SPWs. The second oscillator was externally applied and interpreted as the Medial Septum (MS). We found that this was sufficient to create stable firing-fields in pyramidal, but not interneurons. Removal of the INP-MS oscillator can trigger compression of the firing-fields as ordered bursts. A dedicated population of interneurons can also induce reverse replays. We derived a local learning rule that could successfully bind to-be-learned sequences that were generated by external inputs onto the IGTS. Forward or reverse replays could be induced by removing the INP-MS. Surprisingly, through a combination of the local Fourier rule and theta oscillations during learning, highfrequency assemblies emerged during replays. These assemblies also persisted in the presence of MS inhibition, where they nested on the theta oscillation. These assemblies triggered ripple-like oscillations in RO-CA1. The RO-CA1I were the dominant component of the population activity and displayed phase locking to the ripples in the population activity. The ripples served to segregate learned assemblies and prevent assembly overlap. Finally, we demonstrated that SHOT-CA3I can prevent memory fragmentation due to phase distortions during replay by modulating the probability of SPW-R initiation. We predicted that this would have an observable effect on inter-SPW-R interval distributions that has been observed experimentally [36]. Our results support the role of interneurons in generating, compressing, reversing, and refining spike sequences for stable learning.

A Link Between Theta-Phase Compression and Sharp-Wave-Ripple Compression
Our results suggest that the mechanism for SPW-R compression and compression within a theta cycle may be one and the same. The timing of spikes during both events can be coordinated by a single population of interneurons, operating in two modes, which are delineated by the presence of MS inhibition. When INP-MS was present, the sequence compression within a theta cycle is due to the SPW regulating current, θ IN T .
The hypothesis that these compression mechanisms are related is supported by two lines of experimental evidence. First, the SPW-R compression ratio is only slightly larger (30%) than theta-phase compression [3,37]. Second, in extended replay, SPW-Rs often occur in clusters containing multiple SPW-Rs [36,38]. The individual SPW-Rs within a cluster have a 100-150 ms gap, or a theta period, between them [3, 36,38]. Strikingly, in [36], the authors find theta-harmonics in the inter-SPW-R-interval as peaks in the distribution at 110, 220, 330, and 440 ms. It is important to note that SPW-R clusters also form due to phase locking in 7-10 Hz sleep-spindles [39]. However, SPW-R clusters are reliably evoked during quite awake, when an animal is exploring longer mazes [36,38].

Experimental Predictions from Our Model
We suggest a novel experiment to assess whether dual oscillators support IGTSs. First, employing the left/right alternation task as in [6] should yield stable firing-fields during wheel running [6, 7] ( Figure 6A). We suggest using a closed loop feedback system to [40][41][42] optogenetically drive the MS to oscillate at different frequencies during wheel running ( Figure 6A). For a sufficiently large driving frequency, θ D > θ IN T , interference theory predicts sequence reversal in the IGTS and phase recession ( Figure 6A). For further model predictions, see supplementary discussion (Supplementary Note S1) and Table 2.

Limitations of the Model
While our model displays many of the intrinsic rhythms and behaviours of the hippocampus, like all models, it is a simplification. Here, we consider discrepancies between our model and experimental data. We list the following four limitations: 1) interneurons display limited or no phase precession, 2) the presence of attractor network dynamics in CA3, 3) the gamma/high-frequency oscillation that we predict due to CA3 is faster than observed and 4) the aperiodicity of observed IGTSs. We discuss the limitations of our model more in the supplementary discussion (Supplementary Note S1).
While pyramidal cells precess in phase under a variety of conditions [6, 10], interneuron phase precession is less clear. Some interneurons in CA1 display phase precession due to strong, monosynaptic connectivity with a phase precessing pyramidal neuron [43]. Other measurements of interneuron phase precession demonstrate a greater variability in phase precession slopes [19], with some interneurons even displaying bouts of phase recession (see Figure 3,4 in [19]). Interestingly, the authors ([19]) note that interneurons fire wider bursts which might make phase precession more difficult to assess. However, as the leaky-integrate-and-fire model has limitations in producing the rich dynamics of real neurons, interneuron phase precession in our network model may be eliminated by considering other neuron models.
Our model uses no recurrent excitatory connectivity to create sequences. The recurrent activity is used only to activate interneuron networks. Both phase precessing and SPW spike sequences are internally created in our model by inhibitory weights. While interneuron spike sequences have been documented [44], both experimental work [45] and theoretical modelling [11,12,46] have demonstrated that excitatory connectivity is an important factor for sequence generation and compression. It would be interesting to combine our model with attractor-like dynamics in the future.

Other Models of Replay and Theta Sequences
Dual oscillator models have been suggested as a mechanism for phase precession long before this study [10]. Other studies have considered somato-dendritic interference as a source of phase precession [47]. Recent studies have demonstrated how multiple intra-hippocampal oscillators with uniform phase clustering at discrete values could yield a slower oscillation in the global population activity [11,12]. This mirrors our result with the continuous uniform distribution endowing INP-MS control over population activities. Here, we consider the novel hypothesis that phase precession arises from the need to dilate pre-existing SPW sequences for singletrial-learning of new sequences.
Multiple-oscillator models were also considered as a mechanism of grid cell firing [9]. In these models, the theta frequency is linked to animal velocity in velocity controlled oscillator's (VCO's). By using three VCO's with phase-shifts separated by 2π 3 , hexagonal grid patterns emerge in space. However, recent experiments have cast doubt on this mechanism as grid cells still occur in bats without theta oscillations [48]. Further, in vivo patch recordings demonstrate that an underlying ramp of depolarization encodes significantly more spatial information than the theta oscillation envelope [49,50]. Here, we considered a different hypothesis: the interference arises between an oscillator controlling the spike-times during SPWs, and the medial septum. This creates a compressible IGTS for single-trial learning.
Finally, other models have either linked theta sequences to SPW-Rs [46] or considered the problem of fast-learning in the hippocampus. We summarize these in the supplementary discussion (Supplementary Note S1).

Conclusion
A long-standing hypothesis is that interneurons control spike timing. Here, we demonstrate how populations of interneurons can 1) generate theta oscillations in the hippocampus through recurrent and septal inhibition 2) create phase precessing theta sequences, 3) trigger compressed replay of these sequences in forward or reversed modes 4) facilitate single-trial learning of spike trains as high-frequency-assemblies and 6) refine said assemblies by minimizing their overlap in time, 7) nest these assemblies on a theta oscillation after learning and finally 8) modulate the inter-SPW-R-intervals to prevent the replay fragmentation. This intricate control over spiking in our model was facilitated with little more than a pair of oscillators, interneuron/pyramidal reciprocal connectivity, and Hebbian plasticity.

Bibliography
[1] Buzsáki, G. Two-stage model of memory trace formation: a role for "noisy" brain states. Neuroscience 31, 551-570 (1989). [21] Gerstner, W., Kempter, R., van Hemmen, J. L. & Wagner, H. A neuronal learning rule for sub-millisecond temporal coding. Nature 383, 76 (1996  [31] Lisman, J. The theta/gamma discrete phase code occuring during the hippocampal phase precession may be a more general brain coding scheme. Hippocampus 15, 913-922 (2005).             We suggest running the same task, only optogenetically stimulating interneurons in the medial septum with an oscillatory input of driving frequency of θ D . The driving frequency θ D is varied from wheel running trial to trial. When θ D > θ IN T , the firing field orientation should reverse. Further, we predict that the resulting firing fields will display similar transitions as in our sweep protocol in Supplementary Figure 8.

Data Availability Statement
The data that support the findings in this study are available from the corresponding author upon request.

Code Availability Statement
The relevant code to generate the results of this paper can be found on modelDB [51], under accession number 243447.

Leaky Integrate and Fire Network
The Septal Hippocampal Oscillatory Theta (SHOT-CA3) network consists of coupled leaky integrate-and-fire neurons: where CA3E, CA3I denote the excitatory and inhibitory populations of SHOT-CA3, respectively while REV denotes the reversion interneurons (see Materials and Methods Section 4.9). The neurons receive a constant background current I α for α = CA3E, CA3I, REV current set at or above the threshold (see Table 1 for parameters). The parameter R = 1 · 10 9 Ω serves as the resistance. When the voltage variables reach a threshold value, v threshold , they are immediately reset to v reset followed by an absolute refractory period, τ ref during which the neuronal dynamics are quenched at the reset value. The membrane time constant, τ m determines the degree to which the synaptic currents are filtered by the voltage dynamics. The spikes themselves are filtered by the double exponential synapse: where τ r is the synaptic rise time, and τ d is the decay time t jk is the kth spike fired by the jth neuron. The weight matrices ω CA3I,CA3I , ω CA3E,CA3I are described in further detail in the section FORCE Training Internal Sequences. All weight matrices that we consider are unitless with the units of current (pA) carried by the synaptically filtered spike trains r(t) (see Equation (14)). The SHOT-CA3E also receive a current from the reversion interneuron population (in Figures 2 and 4), I CA3E,REV syn (t), which is described in greater detail in the section Reversion Interneuron Population. In all cases, when we refer to the filtered spike train, we mean r α (t) as given by Equation (13). However, the filtering time constants differ for different figures and different populations of neurons (see Table of Parameters). Note that the only connections in Equations (10)-(12) are from the SHOT-CA3I (r CA3I (t)) or the reversion interneurons (r REV (t)).
For Figure 4, the RO-CA1 layer was included in the simulation: where the weights ω CA1E,CA3E ij were trained using the local Fourier Rule. A series of external input currents I CA1,EC syn (t) were transiently applied to elicit bursting. Further details are described in the supplementary notes. Finally, we added another population of interneurons to RO-CA1 in Figure 5: where ζ i (t) is an independent white noise term with mean 0 and standard deviation σ = 0.2 pA. This noise term was added to the interneurons to prevent pathological synchrony in the RO-CA1Is. The weight matrices ω CA1E,CA1I , ω CA1E,CA1I are untrained and discussed further in the section Further Methods for Figure 5. The inputs from the entorhinal cortex to the RO-CA1Es, I CA1E,EC i (t) are also described in the section Further Methods for Figure 5. Finally, the weight matrix ω CA1E,CA3E is trained by the local Fourier rule, which is discussed further in the Supplementary Materials Section S1. Finally, we note that the synaptic currents from SHOT-CA3I, I CA3E,CA3I syn (t), I CA3I,CA3I syn (t), display characteristics of an interference pattern between θ IN T and θ M S oscillations. These currents can be non-dimensionalized as follows: where b 1 i carries the amplitude (and dimensions) of the current, b 0 i carries the baseline amount, and z i (t) is the fluctuating component that displays the interference between oscillations. We consider all analytical derivations on the dimensionless current, z i (t). The parameters b 1 i and b 0 i are auxiliary-parameters in the subsequent derivations and need not be explicitly computed.

FORCE Training Inhibitory Connections in a Balanced E/I Network
The weight matrices ω αβ = S αβ + L αβ , α, β = CA3E, CA3I are decomposed into static terms (S αβ ) that are chosen to form an initial balanced spiking network, and learned terms (L αβ ) that are FORCE trained [8]. The static term is given by: With N CA3E excitatory neurons and N CA3I interneurons, each neuron receives precisely C CA3E = pN CA3E excitatory connections and C CA3I = pN CA3I inhibitory connections from the rest of the network where p is the degree of sparsity in the network. The variable ι αβ ij = 1 if a connection is present between neuron j in population β (presynaptic) and neuron i in population α (postsynaptic), and 0 otherwise. The variable g controls the coupling strength while γ = C CA3I /C CA3E is used to balance the weight matrix in cases where C CA3I = C CA3E . Here however, we do not train a balanced E/I network, but rather a balanced I network (see Materials and Methods: FORCE Training a Balanced I Network) [54]. We include this section for completion.
In normal FORCE training, the learned term is given by the following: where q determines the amount of learned recurrence the network receives and η α i and φ β j are referred to as the neural encoders and decoders, respectively [8]. The encoders help determine the tuning preferences of the neurons. The decoders are trained using Recursive Least Squares (RLS), an online L 2 minimization scheme described in further detail in the section Recursive Least Squares. The network approximant of the intended dynamics is given by the followingx Here we want to pre-specify which neuron is excitatory and which is inhibitory (although that would not be necessary for hippocampus [3]). To that end, it is sufficient (but not necessary) to constrain the learned component of the weight matrix: While implementing the boundary condition (Equation (21)) on L αβ ij is sufficient, it is time consuming during learning as it requires O(N 2 ) computations to bound every weight in the weight matrix. We opt to impose a condition on the encoders and decoders separately. Consider the following operations: where the operation (s) ± sets s to 0 if its negative (positive) and retains its value otherwise. First, we decompose the decoders and encoders as follows: Note that the terms on line (23) are positive while the terms on line (24) are negative. Thus, we can consider the following weights which implements an alternate boundary condition to (21). Furthermore, implementing (25) only requires O(N ) computations as only encoders and decoders are bounded by their signs. This derivation, among other results, is originally owed to [52].

Recursive Least Squares
The decoders are determined dynamically to minimize the squared error between the approximant and intended dynamics, e(t) =x(t) − x(t). The Recursive Least Squares (RLS) technique updates the decoders accordingly: and r(t) = (r CA3E (t), r CA3I (t)) T . RLS and FORCE training is described in greater detail in [8,53]. The network is initialized with φ(0) = 0, P (0) = I N µ −1 , where I N is an N -dimensional identity matrix, and µ controls the learning rate of RLS.

FORCE Training a Balanced I Network
In addition to the procedure described above, a balanced network can be generated with recurrent inhibition alone [54]. In particular, we can consider the same network equations as before with the constraints: where the interneuron network receives a super-threshold I CA3I > I threshold = −40 pA background current that is being balanced by the recurrent inhibition. The excitatory neurons in the population also receive inhibitory connections and thus have a super-threshold background current after FORCE training is concluded. However, during FORCE training, the current I CA3E is clamped to below the threshold (I CA3E < I threshold ) to prevent firing in the excitatory neurons thereby eliminating both EE and IE weights. This occurs through the dependence of Equation (26) on r CA3E (t), and φ(0) = 0.

Internally Generated Theta Sequence
To construct a SHOT-CA3 network with a recurrently generated theta sequence of firing, we FORCE train the network to learn the following supervisor: Here, θ IN T is the frequency of the oscillation that is embedded in the recurrent weights of the network and χ i is a phase shift for the m components of the supervisor. The m components of the supervisor yield an N × m dimensional matrix of encoders and decoders with m = 100 components. For each neuron, exactly one encoder element is non-zero and is set to 1. The phase shift is randomly selected and uniformly distributed over [0, 2π]. When trained in this fashion, the inhibitory currents and INP-MS combine to generate an interference pattern in the excitatory neurons (see Supplementary Figure 1-2): We interpret z i (t) to be the fluctuating, dimensionless, component of the current an excitatory SHOT-CA3 neuron receives. Note that we have written ψ i as opposed to χ i as the phases of the excitatory neurons are only indirectly related to the phases in the supervisor (Equation (28)). This is due to the operations applied in Equation (25).

Interference Based Control of the Population Activity by INP-MS
With the currents arriving at each neuron given by Equation (29), we consider the conditions under which INP-MS influences the population activity of SHOT-CA3. In particular, we will set out to provide sufficient conditions under which the population activity ρ(t) is a periodic function in time with period θ −1 M S . Proceeding generally, we assume that a network transforms its currents via some unspecified differentiable non-linearity F (z) where F (z) ≥ 0 into spike rates. This implies that the population activity is given by the following: The function ρ(t) is the population activity for a network of neurons with identical tuning curves F (z) receiving heterogeneous currents, z i (t). The heterogeneity is in the phases of the currents, z i (t) through the variableψ i . We assume that this parameter comes from a static distribution with density function ρ ψ (ψ). Now, consider ρ(t + θ −1 M S ): If the phase distribution is uniform (ρ ψ (ψ) = 1 2π ), then: where line (35) is justified by the fact that G(ω, t) = F [cos(2πθ IN T t + ω) + cos(2πθ M S t)] is periodic in ω with period 2π and phase shifts (2πθ IN T θ −1 M S in line (34)) do not alter integrals of periodic functions. This implies that a uniform distribution is a sufficient condition for θ M S control of population activities. Finally, we consider ρ(t) when INP-MS is removed: where F z denotes the derivative of the function F (z), F (z). If we again assume that ψ is uniformly distributed, then we arrive at the following:

Interference Based SPW-R and Theta Phase Compression of Spike Sequences
Multiple-oscillator models are well studied in the literature [9, 10,55,56]. However, here we demonstrate that reversible sharp-wave ripple compression is also well explained by interference theory, and that the phase precession based compression and sharp-wave based compression share the same mechanism in interference theory. Considering the current arriving at the SHOT-CA3E, z i (t), again, = 2 cos 2πt Equation (38)

Reverse Compressed Replay with Dual Oscillator Based Mechanisms
If INP-MS is removed, the interference pattern collapses and the neurons receive the input where F OR is short for forward compression. Removal of INP-MS triggers a compression provided that the neurons also receive an additional super-threshold constant current. Thus, the default compression mode is a forward compression, where bursts occur in the same order as the firing-fields when INP-MS is present. This compression can be reversed in two ways: 1) reverse the order of the firing fields z i (t) and maintain z F OR i (t) or 2) reverse z F OR i (t) and maintain the order of firing fields. In mirror reversion, INP-MS is present and satisfies the frequency relationship θ M S > θ IN T . The firing fields are reversed in order relative to the compression phase. However, this is not a valid mechanism for reverse compression as θ M S > θ IN T , both excitatory and interneuron populations recess in phase. To our knowledge, phase recession has only been observed once [19]. In harmonic reversion, the firing-fields reverse near harmonics of θ IN T (see Supplementary Figure 8). For θ M S just slightly larger than θ IN T /2, the firing fields reverse and display phase precession instead of phase recession.
Finally, directly reversing z F OR i (t) also triggers a reverse compressed replay. An additional population of interneurons, the reversion interneurons, is sufficient for this type of reversion. Reverse compression corresponds to SHOT-CA3E receiving the current: Thus, reversion interneurons should be inhibited by SHOT-CA3I (through sin(2πθ IN T t)) while subsequently inhibiting the excitatory neurons. The amount of inhibition an excitatory neuron in SHOT-CA3 receives from the reversion population is controlled by the amplitude function AM P (ψ i ) = 2 sin(ψ i ).

Reversion Interneuron Population
The reversion population is a population of N REV interneurons with LIF dynamics that receives inhibition from SHOT-CA3I. For simplicity, we take N REV = N CA3I and ω REV,CA3I = −ω REV,CA3I I N CA3I , where I N CA3I is the N CA3I ×N CA3I identity matrix and ω REV,CA3I is a scalar quantity with (ω REV,CA3I = 1) that determines the amount of inhibition the reversion population receives. The only parameter difference between the reversion population and the other neurons is the background current, I REV , which varies when INP-MS is present or absent (see Figure 3, Figure 4 Captions).
From the preceding section, we require the reversion population to transmit the signal 2 sin(ψ i ) sin(2πθ IN T t) in a biologically plausible manner. First, the reversion interneurons must be locked into the internally generated theta oscillation through the term sin(2πθ IN T t). Thus, we train a decoder φ REV that can decode sin(2πθ IN T t) from the activity of the reversion population: where φ REV is trained using the decoded internal oscillators of the SHOT-CA3 network as a supervisor. The training is not performed online as in RLS but is learned immediately through L 2 optimization: where δ = 50 acts as a regularization parameter and T = 25 s was used to compute φ REV . Finally note that the current:  Here, we will derive the relationship between, θ M S , the interneuron recurrent oscillation frequency, θ IN T , and the high-frequency assemblies during compressed replay (without INP-MS). We refer to the frequency of these assemblies as Γ. As in experimental data, these assemblies nest on both theta and SPW states [30-32, [57][58][59] The quantity n x denotes the number of cycles of the oscillation x (either x = M S or x = IN T ) in a specific time period while τ M S = θ −1 M S , τ IN T = θ −1 IN T , and τ = (θ IN T − θ M S ) −1 is the period of the firing field. First, as the number of theta cycles (θ M S ), determines the number of assemblies ( Figure 5F), then the total number of theta cycles is then the frequency of assemblies during SPW-R compression, which we define as Γ = n Γ τ IN T is given by which in our case yields = 0.5 Hz, θ M S = 8 Hz, yields an assembly frequency of Γ = 136 Hz. This is in good agreement with the results of Figure 4 and Figure 5. In Figure 4, we see 13 assemblies over an 88 ms interval, yielding Γ ≈ 147.7 Hz (forward replay), 13 assemblies over 98 ms, yielding Γ ≈ 132.65 (reverse replay), while in Figure 5 we see 8 assemblies over a 70 ms interval yielding Γ ≈ 114 Hz.

Statistical Methods
Due to the deterministic nature of the model, the only statistical tests performed in this manuscript were in Supplementary Figure 9, to determine if network responses to current or amplitude ramps were non-zero. An F-test was used in all cases, as implemented in the MATLAB fitlm model. All residuals were assumed to be normal, although this was not formally tested. For small current ramps, the data in Supplementary  Figure 9 showed some skewness inconsistent with normality. We opted not to exclude this data as this effect vanishes for larger values of I CA3I . Additionally, some of the data tested in Supplementary Figure 9 displayed a minor amount of nonlinearity, inconsistent with a linear regression assumption. However, we opted for a simple first order (linear fit), rather than more elaborate methods. No statistical methods were used to predetermine sample sizes, due to the underlying system being deterministic. See Life Science Reporting Summary for additional details.