Abstract
Brain nuclei are traditionally defined by their anatomy, activity, and expression of specific markers. The hypothalamus contains discrete neuronal populations that coordinate fundamental behavioral functions, including sleep and wakefulness, in all vertebrates. Particularly, the diverse roles of hypocretin/orexin (Hcrt)-releasing neurons suggest functional heterogeneity among Hcrt neurons. Using single-cell RNA sequencing (scRNA-seq) and high-resolution imaging of the adult male and female zebrafish hypothalamic periventricular zone, we identified 21 glutamatergic and 28 GABAergic cell types. Integration of zebrafish and mouse scRNA-seq revealed evolutionary conserved and divergent hypothalamic cell types. The expression of specific genes, including npvf, which encodes a sleep-regulating neuropeptide, was enriched in subsets of glutamatergic Hcrt neurons in both larval and adult zebrafish. The genetic profile, activity, and neurite processing of the neuronal subpopulation that coexpresses both Hcrt and Npvf (Hcrt+Npvf+) differ from other Hcrt neurons. These interspecies findings provide a unified annotation of hypothalamic cell types and suggest that the heterogeneity of Hcrt neurons enables multifunctionality, such as consolidation of both wake and sleep by the Hcrt- and Npvf-releasing neuronal subpopulation.
Significance Statement
The study reveals the intricate heterogeneity within the hypothalamic periventricular zone of zebrafish, identifying 21 glutamatergic and 28 GABAergic cell types through single-cell RNA sequencing (scRNA-seq) and high-resolution imaging. Comparative analysis with mouse scRNA-seq data revealed conserved and divergent cell types, transcriptional regulatory mechanisms, and neuropeptide localization. Notably, we identified a unique neuronal subpopulation coexpressing both hypocretin/orexin and neuropeptide VF neuropeptides in zebrafish. The distinct genetic profiles, activity patterns, and neurite processing of this subpopulation suggest a role in regulating both sleep and wakefulness.
Introduction
The hypothalamus, an evolutionary conserved mediobasal region of the brain, regulates fundamental aspects of physiology and behavior, including arousal, sleep, feeding, energy balance, stress, reward, reproduction, and motivated behavior (Machluf et al., 2011; Saper and Lowell, 2014). The hypothalamus modulates these diverse functions via specific neuronal populations and corelease of neurotransmitters and neuropeptides, which are conserved from fish to mammals (Xie and Dorsky, 2017). Historically, hypothalamic neuronal networks have been characterized based on biochemical, electrophysiological, and histological assays (Guillemin, 1978; Flament-Durand, 1980; Stern, 2001) and traditionally categorized based on their neuropeptide or neurotransmitter identity. However, recent technological advancements in single-cell transcriptomics enabled the identification of intricate diversity within hypothalamic neuronal populations, highlighting their complexity (Aldridge and Teichmann, 2020; Fong et al., 2023). Single-cell RNA sequencing (scRNA-seq) has been used to shape cell-type–specific transcriptional signatures of the hypothalamus of zebrafish, mice, and humans (Romanov et al., 2017; Shafer et al., 2022; Wake et al., 2022; Berkhout et al., 2024; Guo et al., 2024). Despite these advancements, a combination of comprehensive genetic, anatomical, and functional assays is required to understand the complex roles of each hypothalamic neuronal subtype.
In fish, hypothalamic neuronal cell types, located in the medial–ventral diencephalon, contain a lower number of cells but are functionally analogous to their mammalian counterparts. The transparent zebrafish is a powerful genetic model to study hypothalamic development and functions (Lovett-Barron et al., 2020). For example, studies in zebrafish revealed conserved functional roles of Orthopedia (Otp), a pivotal homeodomain transcription factor, in regulating the development of hypothalamic nuclei across vertebrates (Amir-Zilberstein et al., 2012; Fernandes et al., 2013). Thus, studying the genetic profile, anatomy, and activity within discrete hypothalamic regions in simple vertebrate can dissect the development and roles of neuronal subpopulations.
Underlying the functional diversity of the hypothalamus is a small population of neurons releasing hypocretin/orexin (Hcrt) neuropeptides, which regulate an array of physiological and behavioral states, including the transition between sleep and wake (Sakurai, 2007; Appelbaum et al., 2009; Elbaz et al., 2012; Adamantidis and de Lecea, 2023). Loss of Hcrt neurons causes narcolepsy, which is characterized by fragmented sleep and wake, and cataplexy, i.e., sudden muscle weakness triggered by emotional stimuli (Liblau et al., 2023). In rodents and zebrafish, Hcrt neurons are located in the lateral hypothalamus area (LHA) and hypothalamic periventricular zone (PVZ), respectively. The zebrafish PVZ encompasses several regions, including the ventral zone (Hv) and dorsal zone (Hd) of the periventricular hypothalamus, zona limitans (ZL), lateral hypothalamic nucleus (LH), anterior tuberal nucleus (ATN), and the paraventricular organ (PVO). Based on their Hcrt peptide expression, Hcrt neurons are considered a homogeneous population. However, recent studies in mammals and zebrafish suggested that Hcrt neurons can be divided into functionally and genetically distinct subpopulations, which express heterogeneous molecular actors that mediate the diverse functions of Hcrt neurons (Sagi et al., 2021).
Here, we aimed to identify neuronal subpopulations within the PVZ. Using scRNA-seq of the zebrafish PVZ and cross-species integration with mouse transcriptomic hypothalamic datasets, as well as whole-mount hybridization chain reaction fluorescent in situ hybridization (HCR-FISH), we profiled the cellular repertoire and found divergent and conserved evolutionary landscape of cell types in this hypothalamic region. In zebrafish, we identified a subpopulation of Hcrt neurons with distinct activity, which express Hcrt and additional sleep-related neuropeptide and project to specific brain regions.
Materials and Methods
Zebrafish husbandry and transgenic lines
Adult WT and tg(hcrt:EGFP) zebrafish (Appelbaum et al., 2009) were raised and maintained in automated zebrafish housing systems (Aquazone; temperature 28 ± 0.5°C, conductivity 500 μS), pH 7.0, subjected to a light cycle of 14 h of light followed by 10 h of darkness, and were fed twice daily. Before decapitation and brain extraction, adult zebrafish were anesthetized using 0.016% tricaine and immediately put in ice-cold water. Embryos were obtained through natural spawning and raised in E3 medium containing methylene blue (0.3 ppm) in a light-controlled incubator at 28 ± 0.5°C, adhering to the 14 h light/10 h dark cycle. All procedures involving fish were evaluated and granted approval by the Bar-Ilan University Bioethics Committee.
Adult zebrafish hypothalamus collection and scRNA-seq
A total of 28 male and female tg(hcrt:EGFP) zebrafish were used for scRNA-seq analysis. Six samples overall were run with 10x Genomics—10X_39_1 and 10X_39_2 included 9 brains, 10X_47_1 and 10X_47_2 included 9 brains, and 10X_49_1 and 10X_49_2 included 10 brains. All zebrafish used were ∼1 year of age. To initiate the dissociation process, zebrafish brains were extracted and embedded in molds containing 1.5% low-melting-temperature agarose. Once the agarose solidified, the brains were cooled and placed on a Leica VT1200S Automated Vibrating Microtome, where they were sliced into coronal sections measuring 500 μm in thickness. The PVZ brain area was dissected from vibratome sections, guided by the EGFP signal of Hcrt neurons, and immediately microdissected in cold aCSF. To facilitate tissue digestion, we subjected the microdissected tissue pieces to a solution containing papain and DNase. Specifically, a papain digest solution was prepared by reconstituting Vial 2 of the Worthington Papain system in 5 ml of aCSF. Additionally, 5% DNase was prepared by reconstituting Vial 3 of the Worthington Papain system in 500 μl of aCSF. The tissue pieces were incubated in 800–1,000 μl of the papain digest solution along with 5% DNase for 20–25 min at a temperature of 34°C. During this time, mechanical trituration using a wide-diameter fire–polished glass pipette allowed for the effective separation of most of the tissue. Any remaining undigested pieces were eliminated by filtering the suspension through a 30 μm cell strainer that had been equilibrated with aCSF. The filtered suspension was collected in a microcentrifuge tube coated with bovine serum albumin (BSA). To pellet the cells, the tube was centrifuged at 200 × g for 5 min at 4°C, and the resulting cell pellet was resuspended in 200 μl of aCSF containing 2.5% DNase I, which was prepared by reconstituting Vial 3 of the Worthington Papain system in 500 μl of aCSF.
To remove myelin and debris, the cell suspension was layered on top of 1 ml of 5% OptiPrep (Sigma-Aldrich) in aCSF within a BSA-coated microcentrifuge tube. The tube was then centrifuged at 150 × g and 4°C for 6 min using a slow ramping speed. This centrifugation step helped separate the cell pellet from myelin and debris. The resulting cell pellet was resuspended in a minimal volume of aCSF and examined in a Burker counting chamber to ensure that the cells had intact morphologies, high viability, and successful removal of debris. Throughout the entire dissociation and cell suspension preparation process, from perfusion to the final single-cell suspension, the zebrafish brain tissue or cells were maintained in ice-cold carboxygenated aCSF (95% O2, 5% CO2), except during the papain digestion step where the temperature was maintained at 34°C. After the dissociation and preparation of single-cell suspensions, the cells were diluted to a concentration of 500–1,000 cells/μl. The 10X Chromium-v3 GEM was used for scRNA-seq following the manufacturer's instructions. The aim was to capture 5,000–6,000 cells per sample, each sample containing 9–10 tg(hcrt:EGFP) zebrafish PVZs. The resulting sequencing libraries were multiplexed and sequenced on Illumina NextSeq or NovaSeq NGS platforms, with a target depth of >35–40 K reads per cell.
aCSF solution components were as follows: NaCl (82 mM; MW, 58.4; g/L, 4.79; g/250 ml 8× stock, 9.58), KCl (2.5 mM; MW, 74.6; g/L, 0.19; g/250 ml 8× stock, 0.37), NaH2PO4*H2O (1.25 mM; MW, 138.0; g/L, 0.17; g/250 ml 8× stock, 0.35), NaHCO3 (26 mM; MW, 84.0; g/L, 2.18; g/250 ml 8× stock, 4.37), sucrose (20 mM; MW, 342.3; g/L: 6.85; g/250 ml 8× stock, 13.69), glucose (20 mM; MW, 180.2; g/L, 3.60; g/250 ml 8× stock, 7.21), HEPES (5 mM; MW, 238.3; g/L, 1.19; g/250 ml 8× stock, 2.38), sodium ascorbate (5 mM; MW, 198.0; g/L, 0.99; g/250 ml); thiourea (2 mM; MW, 76.1; g/L, 0.15; g/250 ml as a 1× solution), sodium pyruvate (3 mM; MW, 110.0; g/L, 0.33; g/250 ml), MgSO4.7H2O (10 mM; MW, 246.5; prepared as a 10 ml stock; 2.5 ml used), and CaCl2.2H2O (0.5 mM; MW, 147.0; prepared as a 0.5 ml stock; 0.125 ml used), with an estimated osmolarity of 307.5 mOsm.
Computational preprocessing and filtering of sequencing data
The scRNA-seq data obtained from tg(hcrt:EGFP) zebrafish were aligned to the reference genome and transcriptome (NCBI RefSeq assembly GCF_000002035.5, appended with GFP and markers). Using the 10x Genomics Cell Ranger software (version 5.0.1), the mRNA molecules were quantified, and all sequencing runs were consolidated into a single database containing metadata for each cell. This comprehensive dataset comprised 8,700 valid cells, after filtering. To address paralog genes, their expression levels were aggregated into a single gene, resulting in a final count of 21,427 genes. To ensure data normalization and reduce noise, several steps were implemented. First, mean-centering was applied to the dataset. Next, normalization was performed by adjusting the counts to a common molecule count. Standardization was then carried out by dividing the counts by their respective standard deviations. Finally, a log transformation was applied to the data. Careful selection of features and dimensionality reduction are crucial for effective clustering. Highly variable genes were chosen based on the coefficient of variance (CV) as a function of the mean. Genes detected in fewer than 20 cells or in >60% of all cells were deemed invalid and excluded from the analysis of highly variable genes. Additionally, immediate early genes were excluded. The transcriptome of cells with a minimum of 300 genes and 500 mRNA molecules were analyzed. In addition, non-neuronal cells expressing the markers myelin protein zero, S100 calcium-binding protein beta (s100b), fatty acid binding protein 7 a (fabp7a), CD74 (cd74a), and oligodendrocyte lineage transcription factor 2 (olig2) were excluded from the analysis. The selection process involved calculating the log2 of the mean and log2 of the CV for each gene. A linear regression was performed, and genes displaying the greatest deviation from the fitted curve, indicative of higher-than-expected variance, were selected. The number of genes selected for analysis was determined using an “elbow plot,” a common heuristic method for identifying an appropriate cutoff point. For dimensionality reduction, principal component analysis (PCA) was employed using highly variable genes. The number of principal components used for downstream analysis was determined using an “elbow plot,” assisting in the selection of the most informative components.
Clustering and classification of cell types
Clustering and classification analyses were performed to identify and categorize different cell types within the dataset. A multilevel clustering approach was employed. The first step involved splitting cells into major classes and subsequently subdividing these major classes into subclasses. Cell types were annotated based on known markers specific to each cell type. The selection of features for each level was based on the entire cell set and projected using PCA. To account for batch effects and technical variations, the Harmony algorithm (Korsunsky et al., 2019) was applied to the PCA, generating a shared embedding that grouped cells based on cell type rather than dataset-specific conditions. This resulted in a new PCA representation. A 2D embedding was then performed using t-distributed stochastic neighbor embedding (t-SNE) on the PCA space. Subsequently, cells were clustered based on their 2D distances using the density-based spatial clustering of applications with noise (DBSCAN) nonparametric algorithm. This allowed for the identification of cell types with similar characteristics.
The analysis was carried out on three levels. In Level 1, samples were pooled by the tissue, and various procedures, including preprocessing, clustering, classification, gene enrichment, and marker gene detection, were performed. Level 2 involved the splitting of cells based on major class, and manual annotation was conducted to indicate major cell classes, such as glutamatergic excitatory neurons and GABAergic inhibitory neurons. Additionally, cell types that represented clear doublets between majors or exhibited poor quality were identified and removed. For Level 3 analysis, cells within the major classes were further split into subclasses using the same analysis steps employed in Levels 1 and 2. To create a consolidated dataset, thorough annotation and naming were performed for each cluster based on enriched genes and known markers specific to cell types. The Level 3 analysis served as the foundation for downstream analysis and subsequent investigations.
Zebrafish and mouse data integration and comparison
To compare cell types between zebrafish and mouse, we integrated our zebrafish dataset with published mouse hypothalamus scRNA-seq dataset (Romanov et al., 2017; Zeisel et al., 2018). We employed a data integration approach using the Harmony algorithm (Korsunsky et al., 2019). First, we combined the datasets of both species, normalized them, and excluded genes that were not present in both datasets. Highly variable genes were selected individually for each species and then combined to create a merged gene list for feature selection. Next, we normalized the merged datasets and performed dimensionality reduction using PCA to mitigate batch effects and technical variations. The Harmony algorithm (Korsunsky et al., 2019) was applied to the PCA space, resulting in a shared embedding where cells were grouped based on cell type rather than dataset-specific conditions. This new PCA accounted for experimental factors but not biological ones.
To facilitate the comparative analysis of zebrafish and mouse cell types, we utilized two methods: DBSCAN clustering (Fig. 2) and K-nearest neighbors (KNN, Fig. 2). In the KNN approach, zebrafish cells were classified based on their k-nearest mouse cell neighbors in the PCA space. We calculated the fraction of each specific cell type according to the shared mouse cell type they were classified into. Various k values (10, 25, 50, 100) were tested, but the classification fractions remained stable with a k value of 25. In the clustering approach, the new PCA space was further embedded using t-SNE and then clustered using the DBSCAN algorithm. This clustering method allowed for the identification of similar cell types from both species, as they were positioned close to their counterparts in the 2D t-SNE embedding. Finally, we calculated the percentage of each cell type within the integrated cell types and sorted them based on their percentage values. Higher percentages for both mouse cell type and its corresponding zebrafish cell type in an integrated cluster indicated a high similarity or analogy between the two. This cross-species comparison and integration approach enabled us to gain insights into the relationships and similarities between cell types in zebrafish and mice, providing information for understanding conserved and divergent cellular characteristics across species.
Microinjections and transient expression assays
In transient mosaic expression assays, the pT2-hcrt:EGFP vector (Faraco et al., 2006) was diluted to 45 ng/μl and microinjected into one-cell-stage eggs using a micromanipulator and PV830 Pneumatic Pico Pump (World Precision Instruments). The embryos were raised in Petri dishes under a 14/10 h light/dark cycle at 28°C. EGFP expression in Hcrt neurons was monitored throughout development using an epifluorescent stereomicroscope (model M165FC, Leica). At 6 dpf, larvae expressing EGFP in single Hcrt neurons in at least one hemisphere were sorted out and prepared for immunohistochemistry assay.
Immunohistochemistry assay
Larvae were fixed with 4% paraformaldehyde (PFA) in phosphate-buffered saline (PBS) with 0.3% Triton (PBT) for 18 h at 4°C. Immunohistochemistry was performed based on a previously described protocol (Corradi et al., 2022). Briefly, larvae were washed with PBT and incubated 5 min at room temperature (RT) in 150 mM Tris–HCl, pH 9, followed by 15 min at 70°C, and then washed in PBT and incubated 40 min on ice in Trypsin EDTA (Sigma-Aldrich T4299, diluted 1:50 in PBT). Next, larvae were blocked 1 h at RT in blocking solution (5% goat serum, 1% BSA, and 1% DMSO in PBT). The larvae were incubated for 96 h at 4°C in blocking solution containing the primary antibodies (1:500 dilution): total ERK (Cell Signaling Technology; catalog #4696; RRID, AB_390780), phosphorylated ERK (Cell Signaling Technology; catalog #4370; RRID, AB_2315112), and NPVF/RFRP (Abcam, catalog #AB228496). Following three PBT washes, larvae were incubated for 72 h at 4°C in the following secondary antibodies (1:300 dilution): AF647-α-rabbit (Thermo Fisher Scientific; catalog #A-21244; RRID, AB_2535812), AF405-α-mouse (Thermo Fisher Scientific; catalog #A-31553; RRID: AB_221604), and AF564-α-guinea pig (Thermo Fisher Scientific; catalog #A-11075; RRID, AB_2534119).
Whole-mount HCR-FISH assay
Adult zebrafish brains or 6–7 dpf larvae were fixed overnight in 4% PFA in 1× PBS at 4°C. After washing (3 times in 1× PBT, 5 min each), larvae were permeabilized for 10 min in 100% (v/v) methanol at −20°C and then washed in SSCx2 with 0.1% Triton (SSCTx2). Hybridization with split probes was performed overnight in hybridization buffer (Molecular Instruments) at a probe concentration of 4 nM at 37°C. Probes were designed according to the split initiator approach of third-generation HCR-FISH v.3.0 (Choi et al., 2018). Briefly, even and odd 22–25-nucleotide DNA antisense oligonucleotide pairs carrying split B1, B3, or B5 initiation sequences were tiled across the length of the mRNA transcript, synthesized by Integrated DNA Technologies, and used without further purification (for HCR probe sets, see Extended Data Table 3-1). The next day, larvae were washed in wash buffer (Molecular Instruments) at 37°C for 30 min, followed by two washes in SSCx5 with 0.1% Triton (SSCTx5) at RT for 20 min each and then incubated in amplification buffer (Molecular Instruments). During this time, dye-conjugated hairpins (Molecular Instruments) were heated at 95°C for 90 s and cooled to RT in the dark for 30 min. Hairpin amplification was performed by incubation in amplification buffer with B1, B3, and B5 probes at concentrations of 240 nM overnight in the dark. Finally, samples were washed three times with SSCTx5 for 20 min each, then mounted in 1.5% low-melting-point agarose, covered in SSCx5, and immediately imaged using a confocal microscope.
RNAscope on adult mouse brain sections
We prepared aCSF-perfused brains by cryomold-embedding them in cryoprotective OCT (Tissue-Tek), followed by flash-freezing in isopentane on dry ice and storing them at −80°C. Coronal cryosections (10 μm) were collected on Superfrost slides (Thermo Fisher Scientific) or 2% APTES silanized glass slides. After quick postfixation in 4% PFA, dehydration, and storage in 70% ethanol, sections were pretreated with protease 4 and subjected to multiplexed fluorescence in situ hybridization using the RNAScope Fluorescent Multiplex (3-plex) Reagent Kit (ACDBio) according to the manufacturer's instructions. Mouse probes targeting specific genes were combined for 3-plexing in alternating channels; probes specific for HCRT and NPVF genes were designed to ensure high specificity and sensitivity, using RNAscope NPR form. Imaging was performed on a Nikon Eclipse Ti2 epifluorescence microscope, and image processing was conducted using the NIS Elements software (Nikon).
Imaging and image analysis
An epifluorescent stereomicroscope (M165FC, Leica) and a Leica Application Suite imaging software (version 4.12; Leica) were used image brain sections. A confocal microscope (LSM710, Zeiss) equipped with a 20× objective (W Plan-Apochromat 20×/1.0 DIC VIS-IR, Zeiss) was used to image HCR-FISH and immunohistochemistry experiments. In the HCR-FISH experiments, larvae and brain sections were imaged, and colocalization of cells labeled with HCR probe sets for hcrt, npvf, crhb, qrfp, anxa13l, hmx3a, and lhx9 was determined using one optical plane and analyzed using the ImageJ software (National Institutes of Health). In the immunohistochemistry experiments, every two channels were simultaneously imaged in unidirectional resonant scanning mode (4× line average). The relative fluorescent intensities of pERK and tERK immunostainings were quantified using ImageJ in GFP-positive cells of tg(hcrt:EGFP) larvae. In the axonal projection experiments, the two channels were imaged simultaneously in unidirectional resonant scanning mode (8× line average).
For RNAscope imaging on mouse brains, the Nikon Eclipse Ti2 epifluorescence microscope was employed. This setup allowed for high-resolution visualization of the fluorescent signals from the multiplexed in situ hybridization. Image acquisition and processing were carried out using NIS Elements software (Nikon), ensuring precise detection and analysis of the specific gene expressions targeted by the probes. This comprehensive imaging approach facilitated a detailed examination of the spatial distribution and colocalization of HCRT and NPVF within the hypothalamic sections.
Experimental design and statistical analysis
Statistical significance for neurite projection experiments was determined using Fisher's exact test to test for differences in neurite projections proportions of Hcrt neuron subpopulations. P values were adjusted for multiple comparisons with the Benjamini–Hochberg (FDR) procedure. Differences were considered statistically significant if the adjusted p value is <0.05.
Statistical significance for neuronal activity experiments was determined using ordinary two-way ANOVA with main effects only, following Tukey's multiple-comparison test, with individual variances computed for each comparison, using GraphPad Prism (GraphPad, version 9.2.0). Differences were considered statistically significant if the adjusted p value is <0.05. A description of the number of animals used for each experiment and the number of neurons considered for each analysis can be found in the figure legends.
Data and code accessibility
The datasets produced in this study are available in the following databases: https://www.ncbi.nlm.nih.gov/bioproject/PRJNA1019490. To facilitate ease of reference, we have compiled the comprehensive single-cell expression dataset, including annotations and metadata, in a user-friendly tabular format to allow for convenient exploration and analysis of the data. The datasets for all the samples utilized in this study are accessible on the Figshare website: https://figshare.com/s/5e48b85dcf1d4320e7b4. We provide an online resource that enables exploration of genes and cell types. The website is currently accessible through the link provided as follows: https://storage.googleapis.com/www_zeisellab/zebrafish_hypothalamus/singlegene_violin_png_230910/VIP_alltypes_violinSingle.png. We have uploaded the customized codes that we made using MATLAB to a GitHub repository. These resources are readily accessible to the research community at https://github.com/muhammadtibi/Zebrafish_SC_hypothalamus. The code files provided in the GitHub link above are written in MATLAB programming language. The code utilizes various libraries and packages for data analysis and visualization. It is recommended to use MATLAB Version 9.11 (R2021b) to run the codes. Please refer to the individual code files for specific instructions on running the code and any additional dependencies. Any additional information required to reanalyze the data reported in this paper is available from the lead contact upon request.
Results
Single-cell transcriptomic profiling of the PVZ in adult zebrafish revealed distinct cell types
Neuronal populations, which release neurotransmitters and neuropeptides, are organized in complex intermingled compartments in the PVZ (Levitas-Djerbi et al., 2015; Lovett-Barron et al., 2020; Kenney et al., 2021). To characterize the gene expression profile of cell types in the PVZ of adult zebrafish, we performed scRNA-seq. Since the structural anatomy of hypothalamic nuclei is indefinite, and in order to yield maximum tissue-specific RNA reads, we microdissected the adjacent region surrounding fluorescent Hcrt neurons in adult tg(hcrt:EGFP) zebrafish brain sections (Fig. 1A). In adult zebrafish, Hcrt neurons are located in the dorsomedial part of the PVZ, a layer of neurons surrounding the diencephalic ventricle (Kaslin, 2004; Appelbaum et al., 2009, 2010), and encompass several regions, including the Hv, Hd, ATN, ZL, and PVO (Wullimann et al., 1996; Rink and Wullimann, 2001; Kenney et al., 2021; Fig. 1B,C). Following enzymatic dissociation of the brain samples, 25,006 cells were suspended from the PVZ of a total of 28 male and female tg(hcrt:EGFP) zebrafish. Among these cells, 8,700 neurons were retained after quality assessments (Extended Data Table 1-1). A total of 3,538 neurons were classified as excitatory by expressing vesicular glutamate transporter 2 (vglut2, also named slc17a6), and 5,162 neurons were classified as inhibitory by expressing vesicular GABA transporter (vgat, also named slc32a1, Fig. 1D). To resolve neuronal cell populations in unsupervised clustering of single cells, iterative DBSCAN clustering approach was performed. GABAergic neurons were subdivided into 28 distinct cell types, while glutamatergic neurons were subdivided into 21 distinct cell types (Fig. 1E,F). The neuronal types were manually annotated based on neuropeptidergic markers, transcription factors, or enriched marker genes (Fig. 1E,F).
scRNA-seq of the PVZ in the adult zebrafish brain. A, A schematic diagram showing workflow for microdissection of the PVZ, single-cell isolation, and scRNA-seq. Red dots represent Hcrt neurons. The yellow area represents the PVZ. B, Lateral view of the adult zebrafish brain (Kenney et al., 2021) and the PVZ, containing the ZL, PVO, and the Hv of the periventricular hypothalamus. C, The dissected brain region. Coronal view of two adjacent brain sections (100 μm each) of tg(hcrt:EGFP) adult zebrafish. Schematic and enlarged views are shown at the bottom panel. D, Neurons were classified as GABAergic or glutamatergic based on the expression of slc32a1 or slc17a6a. Unsupervised DBSCAN clustering of neuronal cell types represented in a t-SNE plot (n = 8,700 cells). Cell types are color-coded. E, F, Unsupervised clustering of GABAergic (E, n = 5,162 cells) and glutamatergic (F, n = 3,538 cells) neuronal cell types in the PVZ is represented in a t-SNE plot. Cell-type clusters are color-coded. G, t-SNE plots of neuropeptides and transcription factors with specific cluster expression (hcrt, gal, neurod1, lhx8a) or wide expression (vgf, bdnf, dlx6a, lhx9). Red, high expression; gray, no expression. H, I, Dot plots representing the normalized gene expression distribution [gene unique molecular identifiers (UMIs)/total gene UMIs] of neuronal neurotransmitter and neuropeptide markers (H) and transcription factors (I) in GABAergic or glutamatergic neuronal cell types. Quality assessment of scRNA sequencing data can be found in Extended Data Table 1-1.
Table 1-1
Quality assessment of scRNA sequencing data. The median UMI per cell, the median genes per sample, the UMI molecule count, and the total number of cells, before and after processing the data, in each one of the 6 samples and in total. Download Table 1-1, XLSX file.
Among all transcriptomes, some genes exhibited unique expression in specific clusters (e.g., hcrt, gal, neurod1, and lhx8a), while other genes were expressed across multiple cell types (e.g., vgf, bdnf, dlx6a, and lhx9, Fig. 1G). For example, the GABA19kiss2-gal cell type is defined by the unique expression of kisspeptin2 (kiss2) and galanin (gal), but it also has a high expression of proenkephalin a (penka), which shares expression across multiple cell types (Fig. 1H). In contrast, specific cell types were defined by shared expression of transcription factors but were distinguished by their specific expression of neuropeptides (Fig. 1H,I). For example, GABA15lhx9-qrfp, GABA16lhx9-bdnf, GABA17lhx9-id2a, and Glut14lhx9-hcrt-npvf cell types were defined by the shared expression of the transcription factor LIM homeobox 9 (lhx9, Fig. 1I) but were distinguished by their specific expression of neuropeptides, such as pyroglutamate RFamide peptide (qrfp) for GABA15lhx9-qrfp, neuropeptide VF precursor (npvf), and hcrt for Glut14lhx9-hcrt-npvf and brain-derived neurotrophic factor (bdnf) for GABA16lhx9-bdnf (Fig. 1H,I).
The neuronal cell types corresponded to specific gene expression. For instance, the expression of glutamate decarboxylase 2 (gad2) is mostly linked to slc32a1+ cell types, while the expression of adenylate cyclase-activating polypeptide 1b (adcyap1b) is primarily associated with slc17a6+ cell types (Fig. 1H,I). Moreover, we identified a variety of neuropeptidergic genes that were highly expressed in nearly all of the glutamatergic cell types, including bdnf (Fig. 1G,H) and tachykinin 1 (tac1), as is the case in mice (Mickelsen et al., 2019) or across all of the GABAergic cell types, like vgf nerve growth factor inducible (vgf; Fig. 1G,H) and proenkephalin b (penkb). In contrast, we recognized several GABAergic and glutamatergic populations that exhibit unique expression of neuropeptidergic markers. For example, corticotropin-releasing hormone binding protein (crhbp, GABA2crhbp)-, vasoactive intestinal peptide (vip, GABA12vip)-, qrfp (GABA 15lhx9-qrfp)-, and gal and kiss2 (GABA19kiss2-gal)-expressing neurons were transcriptionally distinct GABAergic populations (Fig. 1H,I). Moreover, hcrt- and npvf-expressing neurons represent a transcriptionally distinct glutamatergic population with low abundance of glutamate marker (Glut14lhx9-hcrt-npvf). In addition, we identified two molecularly distinct populations of arginine vasopressin (avp)-expressing neurons (Glut15avp-dusp2 and Glut19avp-nr5a2), two populations of corticotropin-releasing hormone b (crhb)-expressing neurons (Glut7otpa-nts and Glut14lhx9-hcrt-npvf), and thyrotropin-releasing hormone (trh)-expressing neurons (Glut7otpa-nts and Glut8otpa-trh), suggesting functional subpopulations among the PVZ neuronal networks (Fig. 1H,I). These results provide a comprehensive single-cell transcriptomic view of the zebrafish PVZ and show that subpopulations in the PVZ are regulated by distinct or multiple transcription factors and express unique or several neuropeptides.
Comparative neuronal characterization of the hypothalamus in zebrafish and mice identified divergent and conserved cell types
Cell types are considered “evolutionary units” with the capability of undergoing independent evolutionary changes. Throughout animal evolution, cell types that emerged from common ancestral lineages demonstrate shared core regulatory transcription factors (Arendt et al., 2019). Profiling and comparing cell-type identities between vertebrate species could benefit anatomical and functional studies of the hypothalamus. This comparative analysis is made possible by the advantages of scRNA-seq technologies (Butler et al., 2018; Ding et al., 2019; Shafer, 2019). To better understand the evolutionary similarities and differences between fish and mammals, we compared and integrated the zebrafish scRNA-seq data of the PVZ with the scRNA-seq dataset of the mouse hypothalamus (Romanov et al., 2017; Zeisel et al., 2018; Fig. 2A). We merged the datasets of GABAergic neurons and glutamatergic neurons from both species and integrated the data until the cell atlases of the two species converged, using the Harmony algorithm (Korsunsky et al., 2019; Fig. 2B). Subsequently, we employed two approaches for comparative analysis: (1) clustering of the two-species integrated–cell types in a group of cell types using DBSCAN algorithm (Fig. 2C, Extended Data Table 2-1) and (2) a KNN classifier that assigned zebrafish to mouse cells based on their similarity in integrated PCA space (Fig. 2D, Extended Data Table 2-2).
Cross-species integration and coclustering of the mouse and zebrafish hypothalamus. A, Adult zebrafish and adult mouse coronal brain sections. Thin fragmented black lines indicate areas that were dissected and sequenced: zebrafish PVZ and mouse hypothalamus. B, A schematic diagram showing workflow for cross-species interrogation. C, Upper panel, DBSCAN-based clustering and representative integrated cell types shared between zebrafish and mouse. Bottom, the numbers in the blue box represent the proportion of the original species-specific cluster in the new combined cluster (Extended Data Table 2-1). D, Upper panel, Cross-species comparison using KNN classifier. Bottom, selected shared cell types between zebrafish and mouse. The numbers in the blue box represent the similarity between the two cell types (Extended Data Table 2-2). The names of the original mouse cell types are mentioned as they appear in the original data sets (Romanov et al., 2017; Zeisel et al., 2018). In the t-SNE plots, each cell is colored based on the selected gene expression and its origin (gray, no expression; green, expression in zebrafish; pink, expression in mouse).
Table 2-1
DBSCAN analysis. Comparative analysis of mouse and zebrafish hypothalamus using density-based spatial clustering of applications with noise (DBSCAN) clustering to mixed-species pseudo–cell types. Download Table 2-1, XLSX file.
Table 2-2
KNN analysis. Comparative analysis of mouse and zebrafish hypothalamus using KNN (k-nearest neighbor) classifier of zebrafish cells to mouse cell types (K=25). Download Table 2-2, XLSX file.
The zebrafish and mice integrated data revealed 49 cell types (Fig. 2C, Extended Data Table 2-1), including six cell types that shared >35% of both zebrafish and mouse original cell types and coexpress common transcription factors/neuropeptide markers (Fig. 2C). Intriguingly, both KNN and integrated cluster analysis found an evolutionary conserved or diverged correlation between the expression of transcription factors and neuropeptides in those cell populations. For example, integrated cluster 28 includes 90% of zebrafish Glut8otpa-trh cells and 35–87% of mouse Trh cell types. This shared cell type is cohesive with the KNN analysis, which revealed 88% proximity among zebrafish Glut8otpa-trh cells and mouse Trh-expressing cells (Fig. 2D). Indeed, Otp is required for the development of Trh neurons (Acampora et al., 1999; Fernandes et al., 2013), suggesting a conserved Otp-dependent mechanism of differentiation of Trh neurons in mouse and zebrafish. Moreover, Otp is required for somatostatin (Sst) neuron differentiation in mouse and zebrafish (Acampora et al., 1999; Fernandes et al., 2013), and these findings were confirmed in integrated cluster 36 that combine 41–78% of mouse and zebrafish cell types that express both genes. The control of both Trh and Sst expression by Otp in mouse and zebrafish emphasizes the idea that a conserved array of regulatory factors is required to drive differentiation and development of a specific hypothalamic neuronal population in vertebrates. Future gain and loss-of-function experiments to these candidate transcription factors could resolve the mechanism that regulates the differentiation and development of hypothalamic neurons across evolution.
Mice and zebrafish neuronal clustering revealed neuronal subpopulation that coexpresses both hcrt and npvf (Hcrt+Npvf+) in zebrafish
The hypothalamus includes several neuropeptide-producing neurons that modulate the sleep–wake cycle. Traditionally, these hypothalamic neuronal populations are defined as sleep- or wake-promoting neurons based on the expression and function of specific neuropeptides. We identified the two-species integrated cluster 29, which is composed of 96% of zebrafish Glut14lhx9-hcrt-npvf cells, 86% of mouse Hcrt cells, and 67% of mouse Npvf cells (Fig. 2C). The KNN analysis also confirmed 80% proximity among the mouse Hcrt cell type and the zebrafish Glut14lhx9-hcrt-npvf cell type (Fig. 2D). In zebrafish, Hcrt neurons consolidate sleep–wake transitions (Elbaz et al., 2017), while neurons releasing the RFamides neuropeptide Npvf (also named Rfrp or GnIH) has been implicated in sleep regulation in Caenorhabditis elegans, Drosophila melanogaster, and zebrafish (Chartrel et al., 2011; Lenz et al., 2015; D. A. Lee et al., 2017, 2020). In mammals, Hcrt neurons stabilize wakefulness (Adamantidis et al., 2007; Oesch and Adamantidis, 2021), while the role of Npvf in mammalian sleep remains unexplored. The genetic heterogeneity within integrated cluster 29 suggests a similar regulatory mechanism of differentiation of Hcrt and Npvf neurons and mutual promoter regulation of both hcrt and npvf genes. In accordance, the expressions of the regulatory factor X 4 (Rfx4) and Lhx9 transcription factors were conserved in both species (Fig. 3A,B), although the neuron-specific localization of the neuropeptides was diverged between zebrafish and mice (Fig. 3C,D).
Characterization of Glut14 lhx9-npvf-hcrt cell type in zebrafish and mouse. A, t-SNE of zebrafish and mouse shared cell types; each cell is color-coded by a unique shared cluster. B, Zoom in on shared Cluster Number 29. Each cell is color-coded by gene expression and its origin (gray, no expression; green, expression in zebrafish; pink, expression in mouse). C, Log2 expression of Hcrt and Npvf in zebrafish (green) and mouse (pink). In zebrafish, hcrt and npvf are coexpressed in the majority of the cells, while in mouse, there is no coexpression. D, Venn diagram of Npvf and Hcrt mRNA expression in zebrafish and mouse. In zebrafish, hcrt is coexpressed in a subpopulation of npvf-expressing cells. In mouse, Hcrt and Npvf are two distinct populations. The size of the circle represents the relative number of cells in the dataset. E, Images of hcrt (cyan) and npvf (yellow) mRNA expression using FISH in tg(hcrt:EGFP) zebrafish and WT mouse brains (for zebrafish probe sets, see Extended Data Table 3-1). Mouse brain was stained with DAPI (gray). Scale bar in zebrafish brain, 15 μm. Scale bar in mouse brain, 250 mm. F, Bar graph representing the top nine genes mostly enriched in the Glut14 lhx9-hcrt-npvf cluster in zebrafish, excluding hcrt. G, A bar graph representing the top nine genes highly correlated with hcrt gene in zebrafish. H, I, t-SNE plots of selected genes with high correlation with hcrt and high enrichment score in the cluster containing hcrt, colored by expression level. Red, high expression; gray, no expression. Black box marks the location of the high magnification view (I) on Glut14lhx9-npvf-hcrt cell type.
Table 3-1
HCR probes. Sequences of all probes used to target genes in hybridization chain reaction fluorescent in situ hybridization (HCR-FISH), see Methods. Download Table 3-1, XLSX file.
The sequencing findings suggest that the same individual neuron coexpresses both sleep and wake-regulatory neuropeptides in zebrafish, although Npvf and Hcrt neurons are considered independent hypothalamic neuronal types. Previous anatomical assays did not identify colocalization of both neuropeptides in mammals and zebrafish (Legagneux et al., 2009; Yelin-Bekerman et al., 2015; Madelaine et al., 2017). However, transcriptome profiling of Hcrt neurons in zebrafish using RNA-seq showed enriched expression of npvf (Yelin-Bekerman et al., 2015), and scRNA-seq of zebrafish and Mexican tetra showed cell population of hcrt+ subcluster, which also expresses npvf (Shafer et al., 2022). In contrast, in mouse, RNA-seq and clustering analyses identified Npvf+ and Hcrt+ as two separate cell types (Romanov et al., 2017; Fig. 3C,D). Nevertheless, more recently, Npvf expression was observed in a small number of Hcrt neurons (Mickelsen et al., 2019). The discrepancy between anatomical and sequencing data raised the possibility of low sensitivity and specificity of traditional ISH assays. In order to stain mRNA of npvf and hcrt at high resolution, we used HCR-FISH and RNAscope in adult zebrafish and mice brains, respectively (Fig. 3E). In zebrafish, both Hcrt neurons and Npvf neurons are localized in the PVZ. In mice, Hcrt neurons are primarily found in the LHA, while Npvf neurons are predominantly in the dorsomedial hypothalamus, with some cells also present in the LHA (Fig. 3E). Consistent with the scRNA-seq results (Fig. 3A–D), the expression of npvf and hcrt partially colocalize in zebrafish, but not in mouse (Fig. 3E). These results identified Hcrt+Npvf+ subpopulation of neurons in zebrafish, which is associated with both sleep and wake regulation, and underline the evolutionary divergence with distinct hypothalamic populations in mice.
Molecular heterogeneity of Hcrt neurons in larva and adult zebrafish
Studying the subpopulation of all Hcrt neurons is challenging in mammals, which contain ∼1,000 and 70,000 cells in mice and humans, respectively (Thannickal et al., 2000, 2018; McGregor et al., 2017). The zebrafish Hcrt neuronal network is relatively simple, comprised of ∼20 and 60 neurons in larvae and adults, respectively (Kaslin, 2004; Faraco et al., 2006; Prober et al., 2006; Yokogawa et al., 2007). The zebrafish Hcrt system includes a single hcrt gene, which encodes to Hcrt1 and Hcrt2 peptides (Faraco et al., 2006), and only one Hcrt receptor, which is structurally similar to the mammalian HcrtR2 (Yokogawa et al., 2007). Similar to mammals, zebrafish Hcrt neurons project to the telencephalon, diencephalon, mesencephalon, rhombencephalon, habenula, and pineal gland toward the noradrenergic, dopaminergic, serotonergic, cholinergic, histaminergic, and melatonergic nuclei (Kaslin, 2004; Faraco et al., 2006; Prober et al., 2006; Yokogawa et al., 2007; Appelbaum et al., 2009, 2010; Panula, 2010; Elbaz et al., 2017).
The neuropeptide Hcrt plays a crucial role in vertebrates, controlling various functions, such as sleep, metabolism, feeding, anxiety, reward, and addiction (Sakurai, 2007; Eban-Rothschild et al., 2018). However, it is likely that Hcrt may not be the sole regulator in these multifunctional neurons, and additional neurotransmitters and neuropeptides may also be involved (Sagi et al., 2021). The identification of the sleep-promoting npvf gene in a specific subpopulation of Hcrt neurons prompted us to investigate the genetic profile of these neurons in single-cell resolution. The scRNA-seq results revealed diverse gene expression in the glutamatergic Hcrt neurons. We found a significant enrichment of npvf, crhb, lhx9, adrenomedullin (adm), annexin A13 like (anxa13l), DENN/MADD domain containing 1B (dennd1b), inhibin subunit beta Aa (inhbaa), somatostatin 6 (sst6), and ETS variant transcription factor 1 (etv1) in the Glut14lhx9-hcrt-npvf cell type (Fig. 3F). Moreover, we found a high correlation of lhx9, npvf, dennd1b, vestigial-like family member 2a (vgll2a), rfx4, rapunzel 5 (rpz5), H6 family homeobox 3a (hmx3a), myelin transcription factor 1b (myt1b), and anxa13l with the hcrt gene (Fig. 3G). These genes can have wider expression across the PVZ (npvf, crhb, anxa13l, dennd1b; Fig. 3H) or specific expression in the Glut14lhx9-hcrt-npvf cell type (hcrt, hmx3a; Fig. 3H). Intriguingly, hcrt and crhb coexpressed with npvf in approximately half of the cells in the Glut14lhx9-hcrt-npvf cell type, but their expression does not overlap (Fig. 3I).
The finding of the novel Hcrt+Npvf+ neuronal subpopulation (Fig. 3) and the scRNA-seq data in adult zebrafish suggest that Hcrt neurons comprise several genetic heterogeneous subpopulations of neurons in both larvae and adult zebrafish. To validate the scRNA-seq data (Fig. 4A) and determine the spatial expression pattern of the genes in 6 dpf larvae and adults, HCR-FISH was performed using probes for hcrt and six key transcripts: crhb, npvf, qrfp, anxa13l, hmx3a, and lhx9 (Fig. 4B–E). Adult coronal brain sections and 6 dpf whole WT larvae were stained using multicolored mRNA antisense probe sets. Confocal imaging confirmed the expression of lhx9 in Hcrt neurons in both larvae and adults (100 and 91%, respectively; Fig. 4B–E). Similarly, colocalization of hmx3a and hcrt mRNA was found in the majority of Hcrt neurons in larvae (99%) and adults (93%). Moreover, anxa13l mRNA, a member of the annexin calcium-binding protein family, was expressed in 52% of Hcrt neurons in adults but only in 24% of Hcrt neurons in larvae. We also observed a minimal overlap of the sleep-related qrfp mRNA with Hcrt neurons in both adults (4%) and larvae (1%), consistent with the finding in rodents (Romanov et al., 2017; Seifinejad et al., 2019; T. M. Takahashi et al., 2020), suggesting that Npvf is not the sole sleep-related RF-amide peptide forming a subpopulation within Hcrt neurons and both subpopulations might play a role in the regulation of sleep and wakefulness. Lastly, as predicted in the scRNA-seq data, limited expression of crhb was found in Hcrt neurons of adults (5%), and colocalization was not observed in larvae. As we described in adults (Fig. 3), hcrt and npvf coexpression was detected in 57% of Hcrt neurons of adults and 18% in larvae (Fig. 4B–E). The expansion of Hcrt+Npvf+ neuronal subpopulation along development suggests age-dependent changes in the mechanism by which Hcrt and Npvf regulate sleep and wake.
Spatial expression of markers in hcrt-expressing and adjacent neurons. A, Coexpression based on scRNA-seq. B, C, HCR-FISH assays show coexpression of hcrt mRNA with selected markers in adult brains (B) and 6 dpf larvae (C). Data are represented as mean ± SEM. D, One plane images of coronal sections of adult brains. E, Ventral view of the hypothalamus of 6 dpf larvae. Yellow arrowheads indicate coexpression. Scale bar, 15 μm.
Unique activity and projection patterns of Hcrt+Npvf+ neurons in zebrafish
In zebrafish larvae, the whole population of Hcrt neurons project to the anterior hypothalamus (AH), lateral forebrain bundle (LFB), posterior tuberculum (PT), and intermediate hypothalamus (IH). These neurons also send long ipsilateral descending projections toward posterior regions, including the caudal hypothalamus (CH), raphe nucleus (RN), locus ceruleus (LC), and spinal cord (SC). Moreover, they bilaterally project to anterior regions with projections ascending to the subpallium (SP), dorsal pallium (DP), olfactory bulb, habenula, pineal gland, and optic tectum (Fig. 5A; Kaslin, 2004; Faraco et al., 2006; Prober et al., 2006; Appelbaum et al., 2009; Elbaz et al., 2017; Yasmin et al., 2023). Similarly, the Npvf neuronal population projects bilaterally to the SP, DP, LFB, AH, IH, CH, LC, and RN (Legagneux et al., 2009; Madelaine et al., 2017; Fig. 5B).
Hcrt+Npvf+ neuronal projections and activity. A–C, Ventral view (head point to the top) of 6 dpf tg(hcrt:EGFP) larva immunostained with Npvf antibody. hcrt- and Npvf-expressing neurons in the hypothalamus project to the SP, DP, PT, LFB, AH, IH, CH, RN, and LC. D, E, Ventral view (head point to the top) of neurite projections of representative Hcrt+Npvf+ neuron (D) and Hcrt+Npvf− neuron (E) in pT2-hcrt:EGFP-injected 6 dpf larva immunostained with Npvf antibody. Dorsal views of single-plane optical section of the neuron cell body (D′, E′) and axon projection to the CH (D′, E′). Yellow arrowheads indicate brain regions with significant difference in neurite projection between Hcrt+Npvf+ neurons and Hcrt+Npvf− neurons. F, Quantification of the percentage of neurons projecting to specific brain region out of the total number of either Hcrt+Npvf+ or Hcrt+Npvf− neuronal subpopulations. *p < 0.05; Fisher's exact test and p values were adjusted for multiple comparisons with the Benjamini–Hochberg (FDR) procedure. n = 7 Hcrt+Npvf+ neurons, n = 25 Hcrt+Npvf− neurons. G, 6 dpf tg(hcrt:EGFP) fish immunostained with antibodies against Npvf (red), tERK (blue), and pERK (white). Dashed box marks the location of the high magnification view. Yellow arrowheads indicate Hcrt+Npvf+ neurons. H, The activity (pERK/tERK) of Hcrt+Npvf−, Npvf+Hcrt−, and Hcrt+Npvf+ neurons in 6 dpf tg(hcrt:EGFP) larvae during day (ZT4), night (ZT18), and following nighttime SD (ZT18-SD) represented as violin plots. Hcrt+Npvf− neurons, nday = 226; nnight = 248; nnight SD = 237. Npvf+Hcrt− neurons, nday = 156; nnight = 179; nnight; SD = 151. Hcrt+Npvf+ neurons, nday = 80; nnight = 73; nnight; SD = 74. ****p < 0.0001; two-way ANOVA test followed by Tukey's multiple-comparison test.
The genetically unique Hcrt+Npvf+ subpopulation is restricted to the anterior–medial–dorsal region of the PVZ related to other Hcrt neurons (Fig. 5A–C,G). To determine whether Hcrt+Npvf+ neurite projections are distinct, we microinjected the pT2-hcrt:EGFP construct into one-cell-stage WT embryos. Taking advantage of mosaic expression, 6 dpf embryos expressing EGFP in single Hcrt neurons were selected, and Npvf immunostaining was performed. Imaging of single Hcrt neurons showed that Hcrt+Npvf+ neurons projected to the SP, DP, PT, LFB, AH, IH, CH, RN, LC, and SC (Fig. 5D,D′). In contrast, Hcrt+Npvf− neurons projected to the SP, DP, PT, LFB, AH, IH, CH, LC, and SC (Fig. 5E,E′). Notably, we found increased neurite projections of Hcrt+Npvf+ neurons compared with Hcrt+Npvf− neurons in the SP, DP, and RN (Fig. 5D–F; p < 0.05). These data provide spatial characteristics of Hcrt+ subpopulations, suggesting potential functional distinctions among these subpopulations in different brain regions.
To further understand the unique role of Hcrt+Npvf+ subpopulation, we turned to investigate the activity of Hcrt subpopulations. The activity of Hcrt neurons during day and night is complex and dynamic. In mammals, Hcrt neurons are most active during active wakefulness, virtually silent during slow-wave sleep and rapid eye movement (REM) sleep, and their activity increases during transitions from sleep to wakefulness (Estabrooke et al., 2001; M. G. Lee et al., 2005; Mileykovskiy et al., 2005; K. Takahashi et al., 2008; Duffet et al., 2022; Zhou et al., 2022). In contrast, a subpopulation of Hcrt neurons in mice becomes more active during REM sleep (Feng et al., 2020). In zebrafish, the activity of Hcrt neurons is associated with periods of increased locomotor activity (Naumann et al., 2010). However, the activity of subpopulations of Hcrt neurons was not studied in zebrafish. To determine the activity of Hcrt+Npvf− and Hcrt+Npvf+ neuronal subpopulations, 6 dpf tg(hcrt:EGFP) larvae were sampled during the day (zeitgeber time, ZT, 4), night (ZT 18), and following 4 h of sleep deprivation (SD, ZT 18). The localization of Hcrt+Npvf−, Hcrt−Npvf+, and Hcrt+Npvf+ neurons and the relative neuronal activity was quantified using Npvf antibody and the phosphorylated extracellular signal-regulated kinase (pERK) marker (Fig. 5G). ERK is rapidly phosphorylated following neuronal activation (Rosen et al., 1994), and the intensity of pERK immunofluorescence staining is a reliable proxy of the activity of neurons in zebrafish and other vertebrates (Randlett et al., 2015; Corradi et al., 2022). Notably, the activity of Hcrt+npvf− neurons is consistently lower than the activity of both Hcrt−Npvf+ and Hcrt+Npvf+ neurons during both sleep and wakefulness (Fig. 5H; p < 0.0001). These results demonstrate that the activity of Hcrt+Npvf+ subpopulation is distinct from other Hcrt neurons and is similar to the activity of Npvf neurons. Future studies aiming to understand the dynamic release of Hcrt and Npvf neuropeptides during day and night could correlate these changes in neuronal activity with sleep–wake regulation. Altogether, these genetic, anatomical, and functional experiments identified distinct Hcrt neuronal subpopulations in zebrafish.
Discussion
In this study, we used scRNA-seq and high-resolution brain imaging to classify and characterize neuronal types and subpopulations in zebrafish PVZ. The results will provide a resource to study the evolutionary landscape of neuronal types across vertebrates and facilitate the comparison with other scRNA-seq datasets (Zeisel et al., 2018; Hain et al., 2022; Wei et al., 2022; Woych et al., 2022; Tibi et al., 2023; Anneser et al., 2024). We identified 49 GABAergic and glutamatergic neuron types in the zebrafish PVZ, as well as conserved and divergent evolution of cell types between zebrafish and mice. Specifically, we identified the Hcrt+Npvf+ neuronal subpopulation in zebrafish, which demonstrated distinct genetic profile, activity levels, and axonal projection targets. This work provides a basis for investigating the mechanisms by which the hypothalamus regulates diverse functions in vertebrates.
The identification and classification of distinct neuronal cell types within the zebrafish PVZ provide valuable insights into the complexity of this brain region. The diverse gene expression profiles and shared features among GABAergic and glutamatergic populations offer a foundation for understanding the functional characteristics of hypothalamic subpopulations. Specifically, we found that multiple cell types (GABA15lhx9-qrfp, GABA16lhx9-bdnf, GABA17lhx9-id2a, and Glut14lhx9-hcrt-npvf) were defined by the expression of Lhx9, a transcription factor belongs to the LIM homeobox family, but each cell type was distinguished by specific neuropeptide expression (Fig. 1H,I). Indeed, Lhx9 is conserved from invertebrates to mammals and is necessary and sufficient for the specification of Qrfp and Hcrt neurons in zebrafish and mouse (Dalal et al., 2013; Liu et al., 2015). Each cell type expressed additional transcription factors, suggesting that specific combinations of transcription factors can interact to induce the specification and differentiation of specific cell types. For example, the expression of the transcription factors lhx9 and etv1 can induce the specification of Glut14lhx9-hcrt-npvf cell type, while a combination of lhx9 and orthopedia homeobox b (otpb) can induce the specification of GABA15lhx9-qrfp cell type. The identification of distinct neuropeptidergic markers and transcription factor combinations provides a roadmap for potential applications, such as generating hypothalamic cell types from induced pluripotent stem cells. This comprehensive single-cell transcriptomic view of the zebrafish PVZ opens avenues for future research at both the single-cell and whole circuit levels.
The comparative analysis of zebrafish and mouse hypothalamic cell types provides valuable insights into the evolutionary conservation and divergence of regulatory mechanisms. Shared transcription factors and neuropeptides across species suggest a common genetic toolkit for the development and differentiation of specific cell populations. Intriguingly, our results revealed six conserved cell types (Fig. 2C), which express common transcription factors and neuropeptidergic markers in both zebrafish and mice. A prominent transcription factor found in our analysis is Otp, which frequently been described as a key player in the development of the hypothalamic neuroendocrine system of vertebrates, and in Otp knock-out mice, oxytocin (Oxt), Avp, corticotropin-releasing hormone (Crh), Trh, and Sst-expressing cells are lost (Acampora et al., 1999; Wang and Lufkin, 2000). In contrast to mice, zebrafish mutants for otpa display normal development of Oxt, Avp, and Crh (Amir-Zilberstein et al., 2012). Indeed, our analysis identifies shared cell types involving the genes, oxt, avp, or crh. The identification of Otp-dependent mechanisms in both zebrafish and mouse emphasizes the evolutionary significance of this transcription factor in driving the differentiation of specific hypothalamic neuronal populations. Furthermore, Integrated Cluster 2 comprises genes in the Crh system, which includes crh and crhbp (Fig. 2C). While Crhbp and Otp can modulate Crh activity (Amir-Zilberstein et al., 2012; Kalin, 2018), the transcription factors that regulate the development of crh- and crhbp-expressing neurons remain unclear. Our results propose candidate transcription factors, including regulator of G-protein signaling 4 (rgs4), LIM homeobox 6 (lhx6), and CCAAT enhancer binding protein delta (cebpd), and future gain and loss of function experiments to these candidate transcription factors could resolve the mechanism that regulates the differentiation and development of Crh neurons across evolution. Overall, this comparative analysis provides a valuable resource for advancing our understanding of the evolution and multiple functions of the hypothalamus.
The integration of zebrafish and mouse scRNA-seq data has revealed a noteworthy convergence in transcription factors and divergence in neuropeptides within Integrated Cluster 29, associated with sleep–wake regulation (Fig. 2C). This cluster includes a significant proportion of zebrafish Glut14lhx9-hcrt-npvf cell type, along with a substantial representation of mouse Hcrt and Npvf cell types. The observed proximity in the KNN analysis further supports the evolutionary conservation of these cell types (Fig. 2D). In zebrafish, Hcrt neurons play a role in regulating sleep–wake transitions, while Npvf neurons are implicated in promoting sleep. Oppose to previous findings (Yelin-Bekerman et al., 2015; Madelaine et al., 2017), our anatomical results confirmed partial colocalization of npvf and hcrt in zebrafish and validated that individual neurons coexpress both sleep- and wake-regulatory neuropeptides (Figs. 3E, 4D,E). In mice, RNA-seq clustering analyses and RNAscope assay identified Npvf+ and Hcrt+ as two separate cell types (Fig. 3E), suggesting evolutionary divergence of Hcrt and Npvf neurons. In contrast, the conserved expression of regulatory factors (Fig. 3B), such as Rfx4 and Lhx9 transcription factors, suggests a shared regulatory mechanism for the differentiation of Hcrt and Npvf neurons across species. The finding of coexpression of the transcription factors lhx9, rfx4, and hmx3a in Hcrt neurons was previously described in zebrafish, rodents, and human cells (Dalal et al., 2013; Liu et al., 2015; Yelin-Bekerman et al., 2015; Luo et al., 2021). These transcription factors may drive the expression, differentiation, and maturation of an ensemble of Hcrt neuron-specific genes, which regulate the diverse functions of Hcrt and Npvf neurons across evolution.
The comprehensive examination of glutamatergic Hcrt neurons through scRNA-seq has uncovered a rich genetic landscape, shedding light on the intricate molecular architecture underlying the diverse functions of Hcrt in vertebrates. The expression patterns of hcrt, npvf, and crhb within the Glut14lhx9-hcrt-npvf cell type hint that npvf-expressing neurons can split into two subpopulations, each expressing either hcrt or crhb (Fig. 3I). Utilizing imaging of single Hcrt neuronal axons, we observed that the genetically unique Hcrt+Npvf+ subpopulation demonstrates increased projections to key brain regions, including the SP, DP, and RN, indicating possible functional differences among Hcrt neurons in various brain regions (Fig. 5D–F). Moreover, the projection to the RN was uniquely attributed to the Hcrt+Npvf+ subpopulation. The serotoninergic RN are among the densest regions containing Hcrt receptor (Hcrtr) mRNA expression in mammals (Marcus et al., 2001). The wake-active nature of the RN suggested that it control wakefulness (Weber and Dan, 2016; Saper and Fuller, 2017; Scammell et al., 2017); however, recent studies suggested that the RN play a role in the initiation and maintenance of sleep in zebrafish and mice (Iwasaki et al., 2018; Zhang et al., 2018; Oikonomou et al., 2019; Venner et al., 2020). The dual functionality of neurons in the RN may be explained by the unique projection of Hcrt+Npvf+ neurons, possibly releasing both neuropeptides. Furthermore, the activity of Hcrt+Npvf+ neurons, characterized by increased pERK expression during both sleep and wakefulness, stands out as distinct from other Hcrt neurons and aligns more closely with Npvf neurons. These findings emphasize the functional diversity within Hcrt neurons and highlight the need for further exploration into the dynamic release of Hcrt and Npvf neuropeptides during the sleep–wake cycle. Understanding the conserved role of these systems in regulating sleep and wakefulness in different animal models will provide insights into the evolutionary aspects of sleep regulation. Moreover, it has the potential to reveal the mechanism by which the Hcrt system stabilizes and consolidates both sleep and wakefulness.
This study combined scRNA-Seq and HCR-FISH techniques in zebrafish and characterized the cellular and molecular landscape of the PVZ in larvae and adult zebrafish. Transcriptomic comparative analysis of specific cell types in zebrafish and mammals provides multispecies cell types and evolutionary divergent and conserved neuronal populations in both vertebrates. Among these findings are a subpopulation of Hcrt neurons expressing sleep- and wake-associated neuropeptides, which may explain how Hcrt neurons consolidate both sleep and wakefulness. Furthermore, the identification of multiple Hcrt neuron-specific factors may lead to the discovery of autoantigens that trigger narcolepsy. Future functional studies in zebrafish and mammals aiming to manipulate neuronal activity and the expression of specific cell-type markers can identify the molecular mechanism and neuronal circuits that regulate diverse hypothalamic functions across evolution. These studies could provide insights into the function of the human hypothalamus in both health and disease.
Footnotes
We thank the members of the Appelbaum and Zeisel Laboratories for their technical assistance. We thank Jennifer Benichou-Israel-Cohen for the statistical analysis and Yael Laure for English editing. The research in the Appelbaum Laboratory was supported by the Israel Science Foundation (ISF; 961/19), the United States–Israel Binational Science Foundation (BSF; 2021177), and the National Institutes of Health (NIH; R01 MH116470-01). Research in the Zeisel Lab was supported by the European Research Council (TYPEWIRE-852786), Human Frontiers Science Program (CDA-0039/2019-C), ISF (559/20), and the Swedish Brain Foundation (Hjärnfonden; PS2020-0026).
The authors declare no competing financial interests.
- Correspondence should be addressed to Lior Appelbaum at lior.appelbaum{at}biu.ac.il.
