Orexin neurons inhibit sleep to promote arousal

Humans and animals lacking orexin neurons exhibit daytime sleepiness, sleep attacks, and state instability. While the circuit basis by which orexin neurons contribute to consolidated wakefulness remains unclear, existing models posit that orexin neurons provide their wake-stabilizing influence by exerting excitatory tone on other brain arousal nodes. Here we show using in vivo optogenetics, in vitro optogenetic-based circuit mapping, and single-cell transcriptomics that orexin neurons also contribute to arousal maintenance through indirect inhibition of sleep-promoting neurons of the ventrolateral preoptic nucleus. Activation of this subcortical circuit rapidly drives wakefulness from sleep by differentially modulating the activity of ventrolateral preoptic neurons. We further identify and characterize a feedforward circuit through which orexin (and co-released glutamate) acts to indirectly target and inhibit sleep-promoting ventrolateral preoptic neurons to produce arousal. This revealed circuitry provides an alternate framework for understanding how orexin neurons contribute to the maintenance of consolidated wakefulness and stabilize behavioral state.

O rexin (Ox; also called hypocretin) neurons in the lateral hypothalamus play an indispensable role in stabilizing and maintaining wakefulness. For example, loss of Oxproducing neurons in people with narcolepsy, results in debilitating excessive daytime sleepiness, sleep attacks, and wakeinstability 1 . The narcoleptic phenotype is nearly fully recapitulated in Ox-deficient animal models [2][3][4] . While existing circuit models hold that Ox neurons stabilize wake by exciting other brain arousal nodes, -for example, Ox neurons increase wake by exciting noradrenergic locus coeruleus neurons 5 , we hypothesized here that Ox neurons may instead (or also) produce their wakepromoting and wake-stabilizing effects through inhibition of neurons of the ventrolateral preoptic nucleus (VLPO), a region that is essential for initiating and maintaining sleep. The cellular VLPO contains neurons that are active during sleep, lesions of the VLPO result in profound and persisting insomnia, and acute stimulation of the VLPO, in particular resident GABA/galanincontaining VLPO neurons (VLPO GABA/Gal ), produces sleep [6][7][8][9][10] . Local administration of Ox into the VLPO also rouses animals from sleep, although how Ox acts on its postsynaptic cellular targets, and the identity of these targets, within the VLPO remains unresolved 11 .
Given that Ox is an excitatory peptide, and that Ox neurons co-release glutamate, we predicted that Ox/glutamate input would not produce arousal by directly inhibiting sleep-promoting VLPO GABA/Gal neurons, as some other arousal inputs to the VLPO appear to do 12,13 , but rather would inhibit VLPO GABA/Gal neurons indirectly. Ox might, for example, inhibit VLPO GABA/Gal neurons by activating a GABAergic afferent input to the VLPO, or do so via an additional, inhibitory synaptic relay. We hypothesized that Ox regulation of arousal depends upon a polysynaptic pathway from the lateral hypothalamus to an intra-VLPO circuit that comprises Ox-responsive GABAergic (VLPO GABA ) interneurons whose collaterals inhibit sleeppromoting VLPO GABA/Gal neurons to produce arousal.
To test our hypothesis, we first sought to determine if acute activation of presynaptic Ox inputs (Ox terminals) within the VLPO could trigger behavioral and electroencephalographic arousal from sleep. We then employed in vitro pharmacology, optogenetically-based circuit mapping, combined with single-cell transcriptomics, to determine the neurochemical and molecular phenotypes of the postsynaptic cellular targets of Ox input within the VLPO. In our final experimental step, we sought to identify and characterize the intra-VLPO circuitry and synaptic mechanisms through which Ox inputs to the VLPO might promote and support arousal. In this work, our findings identify a polysynaptic circuit, including an intra-VLPO inhibitory feedforward circuit, through which Ox signaling can rapidly produce, and maintain, arousal.

Results
In vivo stimulation of orexin terminals in the VLPO rapidly triggers arousals from both NREM and REM sleep. Orexin neurons innervate the preoptic region including the VLPO [14][15][16] (Fig. 1a). We first sought to determine whether activation of Ox terminals in the VLPO could produce or otherwise trigger arousal from sleep. To do so, we injected AAV-DIO-ChR2-mCherry or AAV-DIO-ChR2-eYFP into the Ox field of Ox-IRES-Cre mice and then optogenetically stimulated ChR2(+) Ox terminals in the VLPO in vivo, during both NREM and REM sleep. (Fig. 1). Selective expression of mCherry or YFP in Ox neurons was confirmed by double immunofluorescence labeling ( Fig. 1d and Supplementary Table 1). On average, photostimulation of the Ox terminals in the VLPO consistently produced short-latency arousals during both NREM ( Fig. 1e and g) and REM sleep ( Fig. 1f and h), an effect not seen in sham trials (n = 5; paired bootstrapping, p > 0.05). When stimulated during NREM sleep, we found that arousal probability in the ChR2-mCherry group (n = 5) was elevated above that in the wildtype group (WT; n = 3) for all stimulation levels (significant main effect of genotype F (1, 24) = 22.92, p = 0.003; non-significant interaction between genotype and stimulation frequency F (4, 24) = 1.66, p = 0.192, two-way ANOVA; Fig. 1i, top). However, the significant increases in arousal probability from NREM sleep were observed only during 10 and 20 Hz stimulations (t = 2.98 and 3.73, p = 0.019 and 0.004, respectively; post-hoc Holm-Šidák t-tests), not at the 1 and 5 Hz stimulation levels (t = 0.20 and 2.20, p = 0.842 and 0.073, respectively; post-hoc Holm-Šidák t-tests; n = 5). Comparison of the bootstrap confidence intervals and effect sizes at different stimulation frequencies revealed that the effect of 5 Hz stimulation was consistent with the effects of 10 and 20 Hz stimulations but inconsistent with the non-effect at the 1 Hz level (Fig. 1i, bottom). When stimulation was applied during REM sleep (interaction between genotype and stimulation frequency F( 2, 12 ) = 7.02, p = 0.010, two-way ANOVA one-repeated-measure; Fig. 1j, top), arousal probability increased significantly in association with both the 1 and 10 Hz stimulations (t = 3.82 and 5.16, p = 0.002 and <0.001, respectively; sham versus 1 and 10 Hz conditions within the ChR2-mCherry group, post-hoc Holm-Šidák ttests; n = 5; Fig. 1j, bottom). Interestingly, the arousal response during REM sleep appears to be specific to the Ox input to the VLPO given that similar optogenetic stimulation of GABAergic terminals, also arising from the lateral hypothalamic region, in the VLPO rouses mice from NREM but not REM sleep 13 .
These results demonstrate that activation of Ox input to the VLPO produces rapid arousal from both NREM and REM sleep, and in a frequency-dependent manner. The arousals induced by photostimulation were typically brief. We also found that stimulation during NREM sleep, as compared with stimulations in REM sleep, doubled the length of the evoked arousals ( Supplementary Fig. 1). Hence while activation of Ox terminals in the VLPO can evoke immediate arousal from NREM or REM sleep, only stimulation during NREM sleep was associated with longer-term effects, including an immediate post-activation suppression of sleep propensity. This finding might suggest that Ox input has downstream effects on separate NREM and REM sleep regulatory components, but these could both be located within the VLPO itself. This is consistent with previous pharmacological work showing that injections of Ox directly into the VLPO produces arousal by reducing both NREM and REM sleep and also confirming that Ox has local effects in the VLPO 11 .
Effects of orexin on the cellular VLPO neurons. Hence, we next studied the effects of orexin-A (Ox-A) in VLPO neurons in brain slices of WT mice (Fig. 2a-c). We found a dual and opposite effect of Ox-A on VLPO neurons. Ox-A (0.3-1 µM) excited 35% of the neurons in the VLPO (firing frequency: +173.29 ± 41.29%; n = 21; resting membrane potential: +2.80 ± 0.44 mV; n = 26; Fig. 2d), but unexpectedly inhibited 60% of the VLPO neurons (firing frequency: −84.37 ± 4.66%; n = 37; resting membrane potential: −4.53 ± 0.47 mV; n = 44; Fig. 2e). Ox-A had no effects on the remaining 5% of neurons tested. We also found that while the Ox-A mediated inhibition was blocked by tetrodotoxin (TTX) or bicuculline, neither compound had an effect on Ox-A mediated excitation ( Supplementary Fig. 2). These results demonstrate that while Ox directly excites some VLPO neurons, Ox also inhibits other VLPO neurons, but this inhibition must occur indirectly and, we predict, involves synaptic GABA signaling. We also found that the VLPO neurons that were excited by Ox-A were also excited by the cholinergic agonist carbachol and that the VLPO neurons that were inhibited by Ox-A were, instead, inhibited by carbachol ( Supplementary Fig. 3). These findings suggest that the VLPO contains at least two distinct neuronal populations that respond in opposite manners to wake-promoting signals, suggesting that these two subpopulations differentially regulate sleep and arousal.
VLPO GABAergic neurons show a dual response to orexin. The VLPO is largely composed of GABAergic neurons, including a subpopulation that express the peptide galanin (VLPO GABA/Gal ). VLPO GABA/Gal neurons are active during sleep and promote NREM sleep 6,9,17 . Currently, little is known about the role in sleep-wake regulation of VLPO GABAergic neurons that lack galanin (VLPO GABA ). To study these neurons, we first focused on the effects of orexin (Ox-A) on the VLPO GABA population. We recorded from fluorescently labeled VLPO GABA neurons using Vgat-IRES-Cre mice into which we had placed microinjections of AAV-DIO-TdTomato or AAV- Fig. 1 Optogenetic stimulation of orexin terminals in the VLPO rapidly triggers arousals from both NREM and REM sleep. a Orexin innervation of the POA and VLPO (scale bars: 1 and 0.5 mm). b AAV-DIO-ChR2-mCherry was injected bilaterally into the Ox field of Ox-IRES-Cre mice and optical fibers bilaterally placed above the VLPO. c Orexin neurons transduced with ChR2-mCherry (red; left) and optical fiber (*) above the VLPO (right; scale bars: 500 µm). Optical fiber placements (red * ChR2-mCherry and grey * WT; bottom). d ChR2-eYFP(+) neurons double-labeled for Ox-A (in red; scale bars: 200; 100, and 20 µm). Laser-induced arousals from NREM (e) and REM (f) sleep (10 s-long periods before and after 10 s light pulses, top). EMG and EEG recordings and continuous EEG wavelet transforms (bottom). Trial-averaged responses to laser stimulation in NREM (g) and REM (h) sleep. EEG wavelet power for trials containing arousals (top; n = 5). Black contour lines: 5% significance level against sham trials (paired bootstrap confidence interval of the mean wavelet power difference). EMG mean responses (grey region: 95% bootstrap confidence interval; bottom). Mean differences in arousal probability, from NREM (i) and REM sleep (j) in ChR2-mCherry (n = 5) and WT (n = 3; top). On the bottom, the mean differences in arousal probability (stimulation vs. sham) within the ChR2-mCherry group (NREM: effect of genotype F (1, 24) = 22.92, p = 0.003; interaction between genotype and stimulation frequency F (4, 24) = 1.66, p = 0.192, two-way ANOVA one-repeated-measure. REM: interaction between genotype and stimulation frequency F (2, 12) = 7.02, p = 0.010, two-way ANOVA one-repeated-measure). Mean difference in arousal probability from NREM (i, bottom; t = 2.98 and 3.73, p = 0.019 and 0.004 for 10 and 20 Hz stimulations) and REM sleep (j, bottom; t = 3.82 and 5.16, p = 0.002 and 0.0007 for 1 and 10 Hz stimulations; post-hoc Holm-Šidák t-tests, two-way ANOVA one-repeated-measure; n = 5; Source Data file). Error bars: SEM; ‡: significant for factor genotype (twoway ANOVA); #, significant interaction between the factors genotype and stimulation frequency (two-way ANOVA). Mean differences (stimulations vs sham within the ChR2-mcherry group) plotted as bootstrap sampling distributions (bottom). Dots: mean difference in arousal probability; vertical error bars: 95% confidence intervals; *p < 0.05 post-hoc Holm-Šidák t-test. All testing was two-tailed, and n refers to the number of independent animals. 3V 3 rd ventricle, Opt optic tract, f fornix, Ac anterior commissure, CPu corpus striatum, BF basal forebrain and LH lateral hypothalamus.  Fig. 4a,b). Both AAVs showed selective expression of fluorescent protein in Vgatexpressing neurons (Supplementary Fig. 5 and Supplementary  Table 1). We found that Ox-A excited 42.8% of VLPO GABA neurons (firing frequency: +225.61 ± 95.01%; n = 6 and resting membrane potential: +2.95 ± 1.06 mV; n = 6) and inhibited 50% of VLPO GABA neurons (firing frequency: −72.92 ± 14.26%; n = 7 and resting membrane potential: −3.40 ± 1.13 mV; n = 7). Ox-A was without effect on the remaining 7% of the VLPO GABA neurons. These results indicate that, within the VLPO, and similar to our initial findings, there are at least two populations of VLPO GABA neurons, one that is excited by Ox and one that is instead inhibited.
Previous studies showed that the VLPO GABA/Gal sleeppromoting neurons are inhibited by noradrenaline (NA) 12,18 . Accordingly, we found that the VLPO GABA neurons that were inhibited by Ox-A were also inhibited by NA (n = 12) and that they expressed Gad1, Gad2, and galanin, indicating these are the VLPO GABA/Gal neurons. In contrast, we found that VLPO GABA neurons that were excited by Ox-A were also excited by NA (n = 24) and that they expressed Gad1, and Gad2, but not galanin, i.e., were VLPO GABA neurons ( Fig. 3 and Supplementary  Fig. 4). These results confirm and extend previous studies demonstrating that the VLPO contains at least two subgroups of GABAergic neurons. The first of these subgroups, we predict, are VLPO GABA/Gal neurons, which express galanin, are inhibited by wake-promoting signals such as NA and Ox, and comprise the subpopulation of VLPO neurons that promote sleep 9,12,17 . The second subgroup, VLPO GABA neurons, lack galanin, are excited by NA and Ox and, we also predict, could be responsible for Oxmediated feedforward inhibition of VLPO GABA/Gal neurons that would drive arousal.
Orexin inhibits VLPO GABA/Gal neurons by increasing GABAergic afferent input. Our next experiment focused on sleep-promoting VLPO GABA/Gal neurons. We recorded from TdTomato-labelled VLPO GABA/Gal neurons in Gal-IRES-Cre mice that had received injections of AAV-DIO-TdTomato into the VLPO. TdTomato-labeled VLPO GABA/Gal neurons expressed GABAergic markers, and as predicted, they were inhibited by NA. The majority of VLPO GABA/Gal neurons exhibited low-threshold spikes (LTS) (88.5%; n = 26), and importantly, they were inhibited by Ox-A (76.5%; action potential frequency: −84.40 ± 6.05%; n = 12; and membrane potential: −4.60 ± 0.93 mV; n = 13; Supplementary Fig. 6). Given again that Ox is an excitatory peptide, the finding that Ox-A inhibits VLPO neurons suggests that the observed inhibition by Ox-A must be an indirect effect and, moreover, might imply at least one additional synaptic relay and possibly the involvement of synaptic GABA. We, therefore, tested the effects of Ox-A on GABAergic synaptic input to VLPO GABA/Gal neurons. Ox-A increased the frequency of spontaneous inhibitory postsynaptic currents (sIPSCs) in 75% of the VLPO GABA/Gal neurons (+58.76 ± 13.29%; n = 9), indicating that Ox likely inhibits VLPO GABA/Gal neurons by increasing the GABAergic afferent input ( Supplementary Fig. 7).
Orexin mediates feedforward inhibition of VLPO GABA/Gal neurons via VLPO GABA neurons. Given that most VLPO GABA neurons are directly excited by Ox, we next sought to determine whether VLPO GABA neurons could mediate feedforward inhibition of VLPO GABA/Gal neurons. We first tested whether VLPO GABA neurons inhibit VLPO GABA/Gal neurons (VLPO GABA → VLPO GABA/Gal input) and then whether this input was positively regulated by Ox (Fig. 4). Using an intersectional approach, we placed microinjections of Flp-dependent AAV-fDIO-ChR2-eYFP and Cre-dependent AAV-DIO-TdTomato into the VLPO of Vgat-Flp::Gal-IRES-Cre mice. This approach produced expression of ChR2 in VLPO neurons that express the vesicular GABA transporter (Vgat; VLPO Vgat ) and TdTomato in VLPO neurons that express galanin. RNA scope in situ hybridization confirmed selective expression of YFP and TdTomato in VLPO Vgat and VLPO GABA/Gal neurons, respectively (Supplementary Fig. 8 and Supplementary  Table 1). Histological assessment of the AAV injections further confirmed that ChR2 expression was restricted to the VLPO region (Fig. 4b). We recorded from the fluorescently labeled VLPO GABA/Gal neurons at the reversal potential of the ChR2-mediated current 19,20 while photostimulating VLPO Vgat neurons. Photostimulation of the VLPO Vgat neurons produced short-latency opto-evoked IPSCs (oIPSCs) in VLPO GABA/Gal neurons (Fig. 4).
As both VLPO GABA and VLPO GABA/Gal neurons express Vgat we then tested whether collateral projections from VLPO GABA/Gal neurons (VLPO GABA/Gal → VLPO GABA/Gal ) could contribute to the oIPSCs evoked in the VLPO GABA/Gal neurons when photostimulating the VLPO GABA neurons. To explore this possibility, we placed injections of AAV-DIO-ChR2-mCherry into the VLPO of Gal-IRES-Cre mice and recorded from mCherry-labeled VLPO GABA/Gal neurons. While photostimulation of VLPO Vgat input triggered oIPSCs in 92.7% (n = 55) of VLPO GABA/Gal neurons, photostimulation of VLPO GABA/Gal input failed to produce oIPSCs in any of the VLPO GABA/Gal neurons tested (n = 6; Fig. 4c,d). These results demonstrate an absence of reciprocal connectivity between VLPO-GABA/Gal neurons, and more importantly show that activation of VLPO Vgat neurons by Ox is the mechanism by which the VLPO GABA → VLPO GABA/Gal circuit is activated.
We next applied peak-scaled non-stationary fluctuation analysis 21-23 to assess whether Ox-A enhances the VLPO GABA → VLPO GABA/Gal input by increasing the unitary GABA A current or the number of activated GABA A receptors or both. For each cell, we obtained parabolic variance vs current amplitude curves. We estimated the ionic channel current (i) and the number of activated channels (N) open at the peak of the oIPSCs, in control and during Ox-A application. We found that Ox-A increased the oIPSC amplitude by enhancing the number of activated GABA A channels (+89.33 ± 33.23 %; n = 8, p = 0.0091, paired t-test) without affecting the GABA A unitary current (Fig. 4j). These  Table 3). We then compared the molecular identity of GABAergic neuronal groups in VLPO that either express or lack expression of galanin, by performing differential gene expression. The VLPO GABA/Gal and VLPO GABA groups comprised 3 clusters (515 neurons) and 12 clusters (2,911 neurons), respectively. We found that the VLPO-GABA/Gal and VLPO GABA groups have very distinct transcriptional profiles ( Fig. 5a-b and Supplementary Table 4), but also share commonalities. For example, both cell populations expressed  markers for GABA synthesis and transmission (Slc32a1, Gad1, Gad2) and neither expressed the Vglut2 gene (Slc17a6). In agreement with our electrophysiology data, orexin/hypocretin receptors (Hcrtrs) were not expressed in any of the three VLPO GABA/Gal clusters, whereas cluster #1 of the VLPO GABA group (Fig. 5c, Tables 6,7). Hcrtr1 RNA (Ox1R) was not detected in any of the VLPO clusters, as previously reported by in situ hybridization 25 (Fig. 5b).
We confirmed these bioinformatic results using RNA scope in situ hybridization. Our in situ analysis confirmed an absence of expression of Hcrtr2 mRNA in VLPO neurons that expressed both Slc32a1 (Vgat) and Gal mRNAs (i.e., VLPO GABA/Gal ). In addition, Hcrtr2 mRNA was expressed in 28.5% of VLPO neurons that express Slc32a1 but not Gal mRNAs (i.e., VLPO GABA ) ( Fig. 5e; Supplementary Table 8). These results confirm our in vitro recording findings that only VLPO GABA neurons, or a subset thereof, are capable of responding directly to Ox and hence represent a likely 'gateway' for Ox signaling into the cellular VLPO network.
Inputs from orexin neurons excite VLPO GABA neurons. We next verified functional synaptic connectivity between Ox and VLPO neurons using in vitro CRACM recordings. We expressed ChR2 in Ox neurons by injecting AAV-DIO-ChR2-eYFP or AAV-DIO-ChR2-mCherry into the lateral hypothalamus of Ox-IRES-Cre mice. We then recorded from VLPO neurons in brain slices while photostimulating the Ox axon terminals (Fig. 6a). We confirmed selective expression of YFP or mCherry in Ox neurons by immunolabeling ( Fig. 1d and Supplementary Table 1). Histological assessment of the AAV injections confirmed that expression of ChR2 was restricted to the Ox field (Fig. 6b). Photostimulation of the Ox → VLPO input produced optoevoked excitatory postsynaptic currents (oEPSCs) in 24% of recorded VLPO neurons (n = 63). DNQX (20-200 µM; n = 4) abolished the oEPSCs, indicating release of glutamate and AMPA-mediated signaling in VLPO neurons (Fig. 6c-e). Neurons that responded to photostimulations with glutamatemediated oEPSCs displayed no changes in holding currents nor in resting membrane potentials during or after 60 s-long photostimulation trains indicating no detectable release of Ox (Fig. 6f).
To identify the postsynaptic target neurons in VLPO of this Ox input, we tested the VLPO neurons that responded to photostimulation for their response to NA and/or for expression of Gal by scRT-sqPCR (n = 7). We found that 6 of the recorded neurons were excited by NA (Fig. 6g) whereas 1 was inhibited. All VLPO neurons that responded to photostimulation of the Ox → VLPO input lacked Gal mRNA (Fig. 6h). Altogether, these results demonstrate that Ox nerve terminals exclusively contact and activate VLPO GABA neurons (Fig. 7).

Discussion
Here we demonstrate the existence of a functional polysynaptic circuit between the Ox neurons and sleep-promoting VLPO GABA/ Gal neurons. We further demonstrate that selective activation of the Ox → VLPO GABA → VLPO GABA/Gal circuit is potently wakepromoting in vivo. And our in vitro results demonstrate that Ox directly activates Hcrtr2-expressing VLPO GABA neurons, but not VLPO GABA/Gal neurons. These results identify a clearly delineated polysynaptic circuit by which Ox neurons can effectively "turn off" sleep-promoting VLPO GABA/Gal neurons to promote arousal in vivo. Our findings specifically support a circuit model in which the Ox input (and co-released glutamate) exerts its excitatory influence on VLPO GABA interneurons, including in the case of Hcrtr2-expressing VLPO GABA interneurons, which in turn inhibit sleep-promoting VLPO GABA/Gal neurons to promote and stabilize wakefulness. In this verified circuit model, loss of Ox, as occurs in narcolepsy, would result in a reduction in excitatory tone on VLPO GABA neurons, thereby reducing feedforward inhibition of VLPO GABA/Gal neurons and biasing their sustained activity and, hence, the sleep state.
A reciprocal inhibitory relationship between the VLPO and components of the brain's arousal system, including hypothalamic Ox neurons 26,27 , is thought to support both the consolidation of behavioral state and the ability to rapidly and completely transition between behavioral states, e.g., wake to sleep 28 . This circuit arrangement has been conceptualized as being analogous to an electronic flip-flop switch, which possesses the desirable properties of being both self-stabilizing and highly stable in their 'on' or 'off' states. Hence, under normal neurobiological conditions, this flip-flop circuit arrangement ensures strong state boundary control, i.e., unstable intermediate states are avoided. However, when this circuitry is damaged, as occurs in various parasomnias and sleep-wake disorders such as narcolepsy, the 'weakened' switch operates closer to the boundary zone, leading to a greater percentage of time spent in unstable transition or intermediate states. In the case of narcolepsy, these frequent and unwanted transitions into the boundary zone result in the manifestation of debilitating physical symptoms including fragmentation of sleep or wake, cataplexy, sleep paralysis, hallucinations, and sleep-onset rapid-eye-movement sleep events (SOREMs). A similar phenomenon occurs during aging, when cell loss in the VLPO leads to weakening of the "sleep side" of the switch, in turn resulting in sleep fragmentation and daytime napping, both of which are frequent complaints in the elderly 29,30 .
During wakefulness, sleep-promoting VLPO GABA/Gal neurons are inhibited by afferent inputs from many sources 16 , including the Ox neurons. Intermingled however with VLPO GABA/Gal neurons are additional GABAergic neurons (VLPO GABA ) that lack galanin and appear to be active during wakefulness 12,18,[31][32][33] . In this regard, they are similar to the GABA neurons in the more caudal lateral hypothalamus that also promote wakefulness by inhibiting the VLPO GABA/Gal neurons 13 . Here, we show that Ox directly excites a subset of Hcrtr2-expressing VLPO GABA neurons, and in turn through their local collaterals, these neurons inhibit the sleep-promoting VLPO GABA/Gal cell population. Findings from the present study, therefore, support a circuit model in which local VLPO GABA neurons function as an interface between afferent inputs to the VLPO and sleep-promoting VLPO GABA/Gal neurons. In other words, Ox, and possibly other wake-promoting signals, via feedforward inhibition, inhibit VLPO GABA/Gal neurons to produce arousal, and sleep-promoting signals, in turn via disinhibition, activate VLPO GABA/Gal neurons to produce sleep. This intra-VLPO circuit mechanism for promoting state transitions could be shared by most, if not all, wake-and sleep-promoting afferences that modulate VLPO activity 32,34-37 . Our model further predicts that VLPO GABA neurons are wake-active and, possibly, wake-promoting. In fact, recordings in the VLPO have uncovered neurons that are active in wakefulness 8,33 and optogenetic activation of Gad2-expressing VLPO neurons [which likely activates both VLPO GABA/Gal and VLPO GABA neurons] produces arousal, suggesting that co-activation of VLPO GABA and VLPO GABA/Gal neurons can override the ability of VLPO-GABA/Gal neurons to produce sleep 38 .
A better understanding of the synaptic, cellular, and circuit bases by which Ox neurons stabilize and maintain wakefulness has important implications for treating patients with neurological disorders associated with arousal dysfunction. For instance, patients suffering from narcolepsy exhibit debilitating sleepiness, sleep attacks, and wake-instability, which take a tremendous toll on quality of life. At present, these symptoms are largely addressed clinically through the provision of stimulants, which themselves possess unwanted side effects. We show here that Ox exerts its wake-promoting effects, at least in part, through indirect inhibition of sleep-promoting VLPO GABA/Gal neurons, and our electrophysiological and single-cell transcriptomics findings uncover and characterize two molecularly-distinct VLPO cell populations. There exists considerable redundancy in the brain's arousal circuitry, with only a few 'nodes' having been established as truly necessary for arousal maintenance 39-42 , i.e., loss or disruption produces a chronic reduction in arousal level. Similarly, Ox neurons likely produce arousal (and contribute to state stabilization) through multiple projections 14 , wherein each postsynaptic 'target' of Ox neurons likely contributes to arousal control and stabilization, yet any given individual 'target' in isolation is unlikely necessary for arousal maintenance. Taken together, our findings suggest that VLPO GABA neurons, or a subset thereof, may represent a 'common point of entry' for a wide range of inputs into the intra-VLPO cellular network. These findings also suggest the possibility of new 'targets' for the development of more selective pharmacologic strategies for treating the crippling inability to maintain consolidated To activate the input from Ox neurons to VLPO (in vivo and in vitro recordings), we expressed ChR2-mCherry or ChR2-eYFP in the Ox neurons using two Cre-dependent AAVs: the AAV-DIO-ChR2-mCherry (rAAV8-Ef1a-DIO-hChR2(H134R)-mCherry; 1.5 × 10 13 virus molecules/ml; UNC Gene Therapy Center) and the AAV-DIO-ChR2-eYFP (AAV10-Ef1a-DIO-hChR2(H134R)-eYFP; 6.6×10 12 virus molecules/ml; provided by Dr. P.M.F.). As previously described 45  To study the VLPO internal circuit by in vitro CRACM recordings, we used two Cre-dependent AAVs: one AAV that codes for ChR2 (AAV-DIO-ChR2-mCherry) and one that codes for TdTomato (AAV-DIO-TdTomato). We also used a Flpdependent AAV that codes for ChR2-eYFP, the AAV-fDIO-ChR2-eYFP (AAVDJ-Ef1a-fDIO-ChR2-eYFP; 2 × 10 12 virus molecules/ml; UNC Gene Therapy Center). Specifically, to activate collateral inputs between VLPO galanin neurons, we unilaterally injected AAV-DIO-ChR2-mCherry Generation of Ox-IRES-Cre knock-in mice. These mice were generated, validated, and kindly provided by Drs. D.K., T.M., B.B.L., and T.E.S. The Ox-IRES-Cre mice express the Cre-recombinase under the control of the endogenous orexin/hypocretin (Ox/Hcrt) gene locus. Briefly, the generation of the gene targeting construct was performed through the help of a modified BAC clone 46 that contained the Ox/ Hcrt genomic sequence, as previously reported 47 . This BAC-derived construct spans the region 5 kb upstream of the IRES site and 2 kb downstream of the Cre site inserted into 3′-end genomic region. Then, it was electroporated into mouse embryonic stem (ES) cells (W4/129S6, Taconic, Germantown, NY). Correctly targeted ES cells were identified and then injected into blastocysts to generate chimeras. The male chimeras mice were bred to mice bearing a Flp-recombinase transgene, in order to remove the selection marker neomycin. The Crerecombinase activity was verified in the hypothalamus by crossing male chimeras to Ai14 Cre reporter line (obtained by The Jackson Laboratory: Jax B6;129S6-Gt(ROSA)26Sor tm14(CAG-tdTomato)Hze /J, Cat. #007908) 47 .
EEG/EMG and optical fiber implants. We used Ox-IRES-Cre (n = 5 mice) and WT (n = 3 mice as controls). Following bilateral injections of AAV-DIO-ChR2-mCherry into the lateral hypothalamus mice underwent a second surgery to implant a headstage for EEG and EMG recordings and optical fibers for bilateral optogenetic stimulation of the Ox input to VLPO. The headstages were custommade in-house by assembling a 6-pin connector (Heilind Electronics, catalog #MMX853-43-006-10-001000), 4 EEG screws (Pinnacle, catalog# 8403), and 2 flexible EMG wire electrodes (Plastics One, catalog #E363). During surgery, two burr holes were drilled in the skull of the anesthetized mice immediately above the VLPO (AP = + 0.4 mm; ML = ± 0.7 mm) for placement of the optical fibers. Each optical fiber was stereotaxically guided into position (DV = −5.0 mm). We drilled 4 additional burr holes (0.7 mm diameter) for EEG electrode placement. The EMG electrodes were then guided down the back of the neck underneath the trapezius muscle. EEG and EMG electrodes and optical fibers were glued into place using a mixture of dental cement and cyanoacrylate glue to provide insulation and structural stability. The optical fibers were manufactured in-house by assembling a fiber optic cable (200 μm outer diameter, ThorLabs, Cat. #FT200EMT) inserted into a ferrule (ceramic, 200 µm internal diameter, ThorLabs, Cat. CFLC230-10) and epoxied into place (Precision Fiber Products, Chula Vista, CA, Cat. #353ND). The optical fiber was cut to size (5.5 mm) and then flat cleaved and polished at the mating end of the ferrule.
Sleep-wake monitoring, signal processing, optogenetic experiments and immunohistochemistry. Sleep-wake recordings were conducted 2 weeks after surgeries. For recordings the mice were housed individually in transparent barrels in an insulated sound-proofed recording chamber maintained at an ambient temperature of 22 ± 1°C and on a 12 h light/dark cycle (lights-on at 6 AM, Zeitgeber time: ZT0) with food and water available ad libitum. Mice were habituated for at least 3 days before commencing polysomnographic recordings. Cortical EEG (ipsilateral frontoparietal leads) and EMG signals were amplified x5000 (Model 3600; AM Systems, Sequim, WA) and digitized with a resolution of 500 Hz using a Micro 1401-3 (Cambridge Electronic Design; Cambridge, UK).
Offline, EEG signals were digitally filtered using a 2 nd order Butterworth, bandpass, zero-phase filter with cut-off frequencies 0.5 Hz and 225 Hz filtfilt function, MATLAB (version R2020B; Mathworks, Natick, MA). EMG signals were likewise filtered using cut-off frequencies 100 Hz and 225 Hz. Also, for the trial-averaged analysis of EMG responses to light stimulation, EMG time series were integrated, smoothed (MATLAB; smooth function; 500-sample window), and normalized according to the pre-stimulus mean. Continuous wavelet transform of EEG signals were performed using MATLAB based tools developed by Grinsted and colleges (Available for download at: https://github.com/grinsted/wavelet-coherence).
For optogenetic experiments, experiments were carried out between ZT3-ZT9 (a period of low activity in nocturnal mice). We modified an open-source Online Sleep Detection script (Spike 2, Cambridge Electronic Design, source code available from the Cambridge Electronic Design website) to detect when mice had been in NREM sleep or REM sleep for at least 20 s, to trigger a 12-15 mW, 473 nm blue light (R471003GX, LaserGlow, Toronto, Canada) stimulation to the VLPO. The stimulation consisted of light pulses (5 ms) delivered at frequencies 1, 5, 10, 20 Hz over a 10 s period, delivered in random order throughout the recording session. A refractory period of at least 3 min was allowed between each stimulation, no mouse received more than 200 stimulations over any recording session and each recording session was separated by at least two days.
Sleep scoring, analysis and statistics. We used a custom-written MATLAB script for scoring arousals within stimulation periods on a trial-by-trial basis. The scorer was presented with a 30 s epoch of EEG and EMG recording centered on the stimulation period. The scorer was blinded to the animal genotype and stimulation condition (i.e., sham versus 1, 5, 10, or 20 Hz). The scorer determined, through visual inspection, if arousal occurred within the 10s-long stimulation period. Scored arousals required an EEG change lasting at least three seconds. We excluded trials where the 10 s pre-stimulus period contained a mixture of states (e.g., microarousals, muscle twitching activity, NREM-to-REM transitionary states). Sham stimulation trials were interleaved amongst the real laser stimulations.
To examine the time course of evoked EEG and EMG responses to light stimulation, we computed the continuous wavelet transform of the EEG and the integral of the EMG for each stimulation trial, inclusive of the 10 s before and the 20 s after the stimulation period. Mean EEG wavelet power and integrated EMG signals from each animal were subsequently averaged (n = 5). For the mean EMG signal we calculated the bootstrap confidence interval (MATLAB; bootci function). For average EEG power plots, confidence interval contours were computed as follows. For each animal we calculated the difference in wavelet power (real minus sham) across the time-frequency plane ((−10 to +30 s) x (0.5-30 Hz)). The data was smoothed in the time-domain (500-sample-wide moving average) and were subsequently down-sampled by a factor of 10. We used bootstrap resampling, with replacement, of the paired difference plots (all 3,125 permutations) to calculate the 95% bootstrap confidence interval of the mean power difference at each time-frequency coordinate. 5% significance contours were drawn around regions in the spectral plane where effect size confidence intervals excluded zero. Very small clusters of significance were removed for clarity and statistical conservatism.
All parametric statistical testing were carried out using Sigmastat (version 4.0, Systat Software Inc. San Jose, CA). Paired bootstrap resampling was performed using MATLAB; in each case, 5000 bootstrap samples were taken, and confidence intervals were bias-corrected and accelerated. The details of the statistical tests used for each experiment may be found in the results. Effect sizes were considered statistically significant for p < 0.05. All experimental data were subject to histological validation. Data were excluded if the conditions of the histological validation were not met (i.e., cases in which there was not adequate bilateral transduction of the viral vector or optical fiber placement was not correctly positioned). All behavioral recordings were scored by an investigator that was blinded to the recording conditions. We represented data as mean ± SEM and n refers to the number of animals per group.
Brain slice preparations and in vitro electrophysiological recordings. We deeply anesthetized mice with isoflurane via inhalation (5% in O2) and transcardially perfused them with ice-cold modified ACSF (N-Methyl-D-glucamin-(NMDG)-based ACSF solution). We quickly removed the mouse brains and sectioned them in coronal slices (250 µm-thickness) in ice-cold NMDG-based ACSF using a vibrating microtome (VT1200S, Leica, Bannockburn, IL). We incubated the recording slices first for 5 min at 37°C in NMDG-based ACSF, then in normal ACSF (Na-based solution) for an additional 10 min at 37°C, and then we let them return to room temperature (RT) for at least 1 h.
We recorded VLPO neurons in brain slices submerged and perfused (1-1.5 ml/ min) with Na-based ACSF. We recorded under infrared differential interference contrast (IR-DIC) visualization and we recorded fluorescently labeled neurons using a combination of fluorescence and IR-DIC microscopy. We used a fixedstage upright microscope (BX51WI, Olympus America Inc.) equipped with a Nomarski water immersion lens (Olympus 40X/0.8NAW) and an IR-sensitive CCD camera (ORCA-ER, Hamamatsu, Bridgewater, NJ). Real-time images were acquired using MATLAB (MathWorks) script software. We recorded VLPO neurons in whole-cell and cell-attached configurations using a Multiclamp 700B amplifier (Molecular Devices, Foster City, CA), a Digidata 1322 A interface, and Clampex 9.0 software (Molecular Devices). Neurons showing changes in input resistance of more than 10% over time, were excluded from the analysis. We recorded action potential firing in Na-based ACSF in cell-attached (voltage-clamp mode; Vh = 0 mV) or in whole-cell (current-clamp mode) configurations using the K-gluconate-based solution. For cells that were not spontaneously active we raised the ACSF KCl concentration from 2.5 to 6.3 mM.
We recorded spontaneous inhibitory postsynaptic currents (sIPSCs) and optoevoked IPSCs (oIPSCs) in Na-based ACSF using a Cs-methane-sulfonate-based pipette solution at Vh = 0 mV. We recorded VLPO GABA/Gal neurons while photostimulating the local VLPO GABA/Gal →VLPO GABA/Gal and VLPO GABA → VLPO GABA/Gal inputs at the reversal potential of the ChR2mediated current 19,20 determined for each cell (Vh = −5 to −15 mV) . We recorded opto-evoked excitatory postsynaptic currents (oEPSCs) in Na-based ACSF ( V h = −70 mV) using the K-gluconate-based pipette solution. In all the recordings we added 0.5% biocytin in the pipette solutions to mark the recorded neurons.
To ensure unbiased detection of the synaptic events, IPSCs and EPSCs were detected and analyzed automatically using MiniAnalysis. We measured synaptic event frequency and amplitude in control ACSF (last 5 min before drug applications), drug application (last 5 min of 10-15 min of drug applications), and in washout (last 5 min of 30-40 min washouts).
We used the nonparametric Kolmogorov-Smirnov test (K-S test; MiniAnalysis) to evaluate the responses of drugs on sIPSC frequency or amplitude (statistical significance for p < 0.05, K-S test). We statistically compared mean sIPSC interevent interval cumulative distributions using two-way repeated-measures (RM) ANOVA followed by Bonferroni's multiple comparisons post-hoc test. We considered VLPO neurons to be responsive to photostimulation if the oEPSC or oIPSC probability during the first 50 ms that follows the light pulses was greater than baseline EPSC/IPSC probability + five times the SEM (baseline EPSC probability = 22.43 ± 1.84%, n = 61; baseline IPSC probability = 9.08 ± 1.83, n = 14) 23 . We calculated the latency of the oEPSCs and oIPSCs as the time difference between the start of the light pulse and the 5% rise point of the first synaptic event 48 .
We calculated firing frequency and membrane potential changes by comparing values in control ACSF (last 2 min before drug applications), drug application (last 2 min of the 4-6 min of drug applications), and in washout (last 2 min of the 15-20 min washout).
We represented data as mean ± SEM and n refers to the number of cells per group, unless otherwise specified. We compared group means using one-way or two-way RM ANOVA followed by Bonferroni's multiple comparisons post-hoc test (adjusted-p-value, adj-p) or paired t-tests. Values showing p < 0.05 were considered significant.
Peak-scaled non-stationary fluctuation analysis. We used peak-scaled nonstationary fluctuation analysis to determine changes in GABA A single-channel current (i) or in the number of GABA A channels activated (N) in response to Ox-A 21-23,49 . As previously described 23 , we low-pass filtered and then aligned the peaks of photo-evoked IPSCs. We scaled the averaged mean-current waveform to the peak amplitude of individual oIPSCs, squared the difference, and then sampled this variance time series in 30 bins of equal current decrement from peak to baseline. The binned variance, was plotted against the mean-current amplitude (Fig. 4j). We estimated the N and the i values by least-squares fitting of the peakscaled variance and mean-current curve to the equation: δ 2 = iI-I 2 /N + b where δ 2 is the variance, I is the mean-current, and b is baseline variance. For each cell, we selected a minimum of 20 oIPSCs that had no overlapping spontaneous IPSCs. All the analysis was done using software written in Python 3.
To visualize the Ox innervation of VLPO we processed for Ox immunolabeling sections from WT mice (n = 4 mice). We intracardially perfused the mice with formalin ( overnight) and then in Alexa Fluor-555-conjugated donkey anti-rabbit secondary antibodies (1:500, Cat. #: A31572; Lot #: 2286312; Invitrogen; 2 hr). Sections were washed first in TBS and then in PBS, mounted, and coverslipped in deionized water. We imaged and photographed the sections using a confocal microscope (Zeiss LSM 5 Pascal). The Vgat probe was provided by Dr. Shigefumi Yokota (University School of Medicine, Izumo, Japan) 23 .
We combined RNA scope in situ hybridization for Vgat mRNA (in Vgat-Flp mice) and for Gal mRNA (in Gal-IRES-Cre mice) with immunolabeling for YFP and TdTomato (Supplementary Fig. 8 and Supplementary Table 1). Mice were injected into the VLPO with AAV-fDIO-ChR2-eYFP (n = 2 Vgat-Flp mice) or with AAV-DIO-TdTomato (n = 2 Gal-IRES-Cre mice). Four weeks after AAV injections, we intracardially perfused the mice with formalin (10% buffered solution) under deep anesthesia with isoflurane (5% in O 2 ). Brains were post-fixed in formalin (10% buffered solution, overnight), cryoprotected (20% sucrose) and cut into 30 μm, and processed for RNA scope in situ hybridization for Vgat mRNA following the protocol described above. For the sections from Vgat-Flp mice injected with AAV-fDIO-ChR2-eYFP the Vgat mRNA was labeled with Cy3 (in red) and then the sections were immunolabeled for GFP (in green) with chicken anti-GFP primary antibodies We imaged and photographed the covered sections (Vectashield mounting medium; Vector Laboratories) with a confocal microscope (Leica Stellaris 5) at final magnification of 20X and 63X. We used a z-stack at 3-5 µm intervals to image throughout the slice after tile acquisition. We viewed stacks of images using image J software to identify the region containing the VLPO and the fluorescent signals for RNA scope and immunohistochemistry at cellular level.
used Primer 3web 4.1.0 software (http://primer3.wi.mit.edu/) to verify the Gal, Gapdh, Gad1, and Gad2 primers. The sequencing of the second amplification samples that contain the Gal, Gapdh, Gad2 and Gad1 amplicons were conducted by Quintara Bioscience (Boston, MA). We verified the sequence results to the Gal, Gapdh, Gad2 and Gad1 sequences (NCBI GenBank) using Genetyx software. Molecular identification was based on the results from the second reaction. We purchased all the custom primers from Sigma-Aldrich, St. Louis, MO.
For scRT-sqPCR, the cDNAs (n = 38 neurons) were produced as described in the section above. We performed TaqMan™ based probe assay for semi-quantitative gene expressions using the 7500 Fast Real-Time PCR (Applied Biosystems; Foster City, CA). We added cDNA template (2.5-3 µl) enzyme master mix (12.5 µl; TaqMan Universal PCR Master Mix, Applied Biosystems) and TaqMan Gene Expression assays for the amplification of Gal and Gapdh cDNAs (1.25 µl; ThermoFisher Scientific). Volumes for the reactions were adjusted to 25 µl with H 2 O. The Gal cDNA probe was tagged to FAM/MGB and the Gapdh cDNA to VIC/MGB. TaqMan TM Gene Expression assays for the Gal cDNA (TaqMan Assay Reagents, Gal, 20X, Cat. #4331182) and for Gapdh cDNA (Pre-Developed Gapdh TaqMan Assay Reagents, 20X, Cat. #4352339E). The thermal cycling conditions were as follows: enzyme activation (95°C; 10 min), followed by DNA denaturation (95°C; 15 s; 45-50 cycles) and annealing/extension (60°C; 1 min). Relative gene expression between samples was calculated using the 2 -Ct method 53,54 . Data were normalized to the Gapdh expression levels. As control group, we used TdTomato labeled VLPO GABA/Gal neurons. We used values of 50 for non-detects Ct. Each sample was run in duplicates (2 times, in 2 separate reactions). We considered scRT-PCR-and scRT-sqPCR-negative, samples in which Gapdh expression was undetectable. We represented data as mean ± SEM and n refers to the number of cells per group and we compared group means using Mann-Whitney unpaired t-test.
Single cell RNA sequencing data analysis. Single cell RNA-sequencing (scRNAseq) data of the POA 24 were retrieved from the Gene Expression Omnibus (GEO) repository (GSE113576). Data were available as unique Digital Gene Expression (DGE) matrix. Briefly, the DGE matrix was imported into R software (v.3.6.3) and all the downstream processing was performed using functionalities provided by the R library Seurat (Seurat v.3.1.5) [55][56][57] . Expression data were filtered according to the following criteria: (i) cells expressing >200 genes and with a mitochondrial gene expression rate <10% were retained; (ii) genes detected in >2 cells were retained. Seurat3 was used to cluster all POA cells and to perform Differential Gene Expression between clusters. Genes expressed only or expressed predominantly in a specific cluster were identified as marker genes. (Supplementary Fig. 9, Supplementary Table 2). Afterwards, clusters corresponding to glial cell types were discarded and neuronal clusters were retained for downstream analyses. According to the spatial information gathered from Multiplexed Error-Robust Fluorescence In Situ Hybridization (Multiplexed Error-Robust Fluorescence In Situ Hybridization) 24 and in situ hybridization from the Allen Brain Atlas, we isolated only neuronal clusters located in the VLPO (#i5, i8, i10, i20, i26, i29, i35, i37, i39, e3, e7, e11, e20 and e24). Briefly, the following Seurat v.3.1.5 functions were used: (1) SCTransform () to normalize and scale the data, to identify the top 3000 variable genes, and to correct for number of unique molecular identifiers (UMIs), percentage of mitochondrial gene expression, and difference between S and G2M cell cycle scores. Cell cycle scores were inferred using CellCycleScoring () function 58 ; (2) runPCA was used to calculate the first 50 Principal Components (PCs) of the 3000 most variable features; (3) FindNeighbors () to construct a Shared Nearest Neighbor (SNN) Graph. This method first determines the k-nearest neighbors (knn) of each cell and then uses this knn graph to construct the SNN graph by calculating the neighborhood overlap (Jaccard index) between every cell and its k.param nearest neighbors; (4) FindClusters to define the cell clusters and their granularity using the original Louvain algorithm at resolution of 0.6 for all cells and 1.2 for neurons; (5) runTSNE for t-Distributed Stochastic Neighbor Embedding (t-SNE) dimensional reduction method, to visualize cell clusters in 2 dimensions. (6) FindAllMarkers to perform a differential gene expression between each cluster versus all the other cells of the dataset using the non-parametric Wilcoxon Rank Sum statistics. A gene was defined as differentially expressed if absolute logFC was >0.25 and adjusted p-value (adj-p) (Bonferroni corrected) was <0.05. (7) finally, we used FindMarkers to compare two different cell clusters/groups. The abovementioned parameters were also applied to this function. Cells and neuronal clustering were performed using the same bioinformatic pipeline.
Statistics and Reproducibility. Figure 1a: an example of Ox innervation of the POA and VLPO by immunolabeling for Ox-A. The immunolabeling was conducted in n = 4 WT mice. Figure 1c an example of immunolabeling for ChR2-mCherry of n = 5 Ox-IRES-Cre mice used for whole-animal optogenetic studies. Figure 1d: an example of double labelling for ChR2-mCherry and Ox-A immunoreactivity to test selectivity of the AAV-DIO-ChR2-eYFP (and AAV-DIO-ChR2-mCherry) for Ox neurons. Immunolabeling conducted in n = 4 Ox-IRES-Cre mice injected with AAV-DIO-ChR2-eYFP (and AAV-DIO-ChR2-mCherry) into the Ox field. Figure 2a and b: map of n = 29 recorded neurons distributed in n = 19 slices from n = 14 WT mice. Figure 4a: an example of VLPO GABA neurons (in green) and VLPO GABA/Gal neurons (double-labeled in red and green) following the injection of a mix of AAV-fDIO-ChR2-eYFP and AAV-DIO-TdTomato into the

Data availability
The data generated in this study are presented within this paper or its supplementary materials as well as source data are provided with this paper. These include all individual data points and average values showed in both figures and supplementary information. The raw data for in vivo optogenetic stimulation, in vitro electrophysiology and imaging experiments are available from the corresponding authors, upon reasonable request. For scRNA-seq analysis we provided expression count matrix, barcodes and gene IDs for POA and VLPO and figures at the following link: https://doi.org/10.5281/zenodo. 6570978 as source data. For scRNA-seq raw data, the access code on the GEO repository is GSE113576. Source data are provided with this paper.

Code availability
All custom codes used in this manuscript are available at the following repository: https://doi.org/10.5281/zenodo.6570978 (version 3) or by the corresponding authors upon request. Other custom codes used in this study were previously published and are available in the relevant citations.