Dataset acquisition and study participants
The processed MDD single-nucleus RNA sequencing data from GBM samples (GSE144136) were obtained from the Gene Expression Omnibus (GEO, http://www.ncbi.nlm.nih.gov/geo/) database, and the GWAS summary statistics for MDD patients (170,756 cases and 329,443 controls) were obtained from the IEU OpenGWAS (https://gwas.mrcieu.ac.uk/) database with a GWAS ID of ieu-b-102. The published MDD cohort48 consisted of 17 male patients and controls. Cases met the diagnostic criteria for MDD and died by suicide, whereas individuals in the control group died by accident (n = 6) or natural causes (n = 11). Postmortem brain samples were dissected from Brodmann area 9 (dlPFC), and snRNA sequencing was performed (Supplementary Data 1).
Single-cell data processing and cell clustering
We downloaded the processed snRNA sequencing data from the GEO dataset, which had been aligned and quantified using the Cell Ranger Single-Cell Software Suite (10X Genomics Cellranger v2.1.0). The gene expression data were also mapped to the human genome (GRCh38-1.2.0). For the quality check procedure, we first removed the cluster cells originally labeled “Mix” from the dataset. Next, any cells were removed if they expressed fewer than 201 genes or if more than 20% of the UMIs were linked to mitochondrial genes74,75. To eliminate batch effects, the R package Harmony (version 1.2.0) was used to remove batch correction before clustering analysis. Other steps, including Seurat object creation, data normalization, data scaling, identification of variable genes, dimensional reduction, clustering, and UMAP projection, were performed following the standard pipeline of the Seurat R package (version 4.4.0, https://satijalab.org/seurat).
Cluster annotation
All the nuclei were annotated as distinct major cell types on the basis of the average gene expression of known canonical marker genes, including excitatory neurons (markers: SATB2, SLC17A7, and CAMK2A), inhibitory neurons (GAD1 and GAD2), astrocytes (GFAP, AQP4, GLUL, SOX9, GJA1, NDRG2, and ALDH1A1), oligodendrocytes (MBP, PLP1, MOBP, and MAG), oligodendrocyte precursor cells (PCDH15, PDGFRA, OLIG1, OLIG2, and PTGDS), endothelial cells (CLDN5 and VTN), and microglia/macrophage (CX3CR1, MRC1, SPI1, and TMEM119) cells. In addition, SNAP25, STMN2, and RBFOX3 were used as common neuron markers across all the neuronal populations.
We further identified subclusters and annotated them as specific cell subtypes within excitatory and inhibitory neurons on the basis of the average expression of the corresponding gene sets. The marker genes for each subcluster within the major cell types were identified using the FindAllMarkers function in the Seurat package. For excitatory neurons, subclusters were annotated according to layer-specific gene markers as follows: layers II–IV (CUX2, RASGRF2, and PVRL3), layer IV/V (RORB), layer V (SULF2, PCP4, and HTR2C), layers V/VI (TOX, ETV1, RXFP1, and FOXP2), and layer VI (NR4A2, SYNPR, TLE4, and NTNG2). For inhibitory neurons, subclusters were annotated on the basis of well-known markers, including CCK, CALB, VIP, SST, PVALB, SLC32A1, and DLX.
Differential expression analysis
Analysis of differentially expressed genes (DEGs) for key ligands and receptors between the suicide and control groups was performed using two-sided Wilcoxon tests within each responding cluster, with Benjamini–Hochberg correction applied to adjust P values for multiple testing. A BH-adjusted P value less than 0.05 was considered statistically significant. Additionally, the jjVolcano function from the scRNAtoolVis package (version 0.0.7) was used to visualize DEGs in a Manhattan plot between the suicide and control groups across all clusters.
Group preference of each cell type
To quantify the cell type enrichment across groups, we compared the observed to expected cell numbers in each type according to the following formula as we described previously76,77:
$${Ro}/e=\frac{{observed}}{{expected}}$$
(1)
where the expected cell numbers for each combination of cell types and groups are obtained via the chi-square test. Ro/e is the ratio of observed to expected cell numbers in each cell type. Ro/e > 1 suggests a preference for this cell type in this group, and Ro/e < 1 suggests that cells of the given cell type are observed with less frequency than random expectations in the specific group.
Cell‒cell interaction analysis
CellChat
To compare differences in cell‒cell communication networks across different sources, we utilized CellChat (version 1.6.1, https://github.com/sqjin/CellChat)78. The normalized expression data, preprocessed using Seurat, were input into CellChat for downstream analysis. We applied standard preprocessing functions with default parameters, including identifyOverExpressedGenes, identifyOverExpressedInteractions, and projectData, to identify significantly overexpressed genes and interactions and project the data onto inferred interaction networks. We subsequently ran the core functions computeCommunProb, computeCommunProbPathway, and aggregateNet using default settings to compute the communication probabilities and aggregate the interaction networks. The MDD and control groups were merged for further analysis using computeNetSimilarityPairwise, allowing pairwise comparisons of network differences. To visualize the communication patterns, the functions compareInteractions, rankNet, and netVisual_bubble were used.
Initially, we performed CellChat analysis of all the nuclei and determined that excitatory neurons were most regulated by astrocytes. We then subset the relevant nuclei from excitatory neurons and astrocytes for a more focused CellChat analysis followed by the abovementioned steps.
NicheNet
We used NicheNet (version 2.1.7, https://github.com/saeyslab/nichenetr)79 to infer interactions between astrocytes and excitatory neurons, designating astrocytes as “Sender” cells and excitatory neurons as “Receiver” cells. Only genes expressed in more than 10% of the cells within clusters were considered for further ligand‒receptor interaction analysis. We extracted ligands with upregulated activity that overlapped with the CellChat results from the control vs. suicide comparison (i.e., ligands with downregulated activity in the suicide group compared with the control group) and the top 1000 targets from the differentially expressed genes of the “Sender” and “Receiver” cells for paired ligand‒receptor activity analysis. The ligand_activity_target_heatmap function was applied to visualize ligand regulatory activity, with activity scores ranging from 0 to 1. In addition, a heatmap was generated to illustrate the expression of differentially expressed ligands and receptors, which was calculated by averaging gene expression in the relevant cell populations and scaling across the indicated clusters.
Pathway enrichment analysis
Gene ontology (GO) analysis
To perform GO analysis, the potential target genes of excitatory neurons regulated by astrocytes through the PTN pathway, as predicted by NicheNet, were analyzed using the online tool Metascape80 (http://metascape.org/gp/index.html) to identify the enriched pathways. Only gene sets in the “GO Biological Processes”, “GO Molecular Functions” and “GO Cellular Components” categories were considered, and a P value less than 0.05 was considered statistically significant.
Gene set enrichment analysis (GSEA)
To explore the potential downstream signaling pathways affected by PTN signaling, we conducted GSEA using GSEABase (version 1.64.0) and gene sets from the Molecular Signatures Database (https://www.gsea-msigdb.org/gsea/msigdb/human/collections.jsp)81. The DEGs in excitatory neurons between the suicide and control groups were subjected to GSEA. Genes were ranked by log2-fold change from the differential expression analysis, GSEA was conducted using hallmark (H) and curated (C2) gene sets with 1000 permutations, and a P adjusted value less than 0.05 was considered statistically significant.
Pathway-based polygenic regression (scPagwas) to identify MDD-relevant clusters
To identify the genetic associations of cell clusters with MDD, we used scPagwas (https://dengchunyu.github.io/)82, which uses an optimized polygenic regression model to integrate snRNA-seq and GWAS data. For GWAS data preparation, we downloaded the GWAS summary statistics for MDD from the IEU OpenGWAS database and conducted data pruning following the scPagwas guidelines, including GWAS data reprocessing, SNP extraction, plink, filtering of LD information and integration of the result files for chromosomes 1–22 with default parameters. For snRNA-seq preparation, we subset the normalized and scaled expression data from the suicide group. The prepared GWAS and snRNA-seq data were then input into scPagwas for downstream analysis, with 1000 top genes and 200 iterations specified.
Animals
Male C57/BL6 mice (8–10 weeks) were purchased from the Institute of Experimental Animals of Guangdong Medicine Experimental Animal Center. Male Aldh1l1-Cre/ERT2 (Cat# 031008) mice (8–10 weeks) were purchased from the Jackson Laboratory. All the mice were housed on a standard 12-h light/dark cycle (light from 7 a.m. to 7 p.m.) at a constant temperature of 25 ± 1 °C and 50% humidity, with food and water available at all times. Ten mice per group were used in behavioral tests, electrophysiological recording, and dendritic spine analysis. Six mice per group were used in RT-PCR, western blotting (WB), and immunohistochemistry (IHC). All experimental procedures were approved by the Use Committee of Sun Yat-sen University Cancer Center and the Animal Care Committee (No. L102012024120D, No. L102012024228P) and were conducted in accordance with the guidelines of the National Institutes of Health (NIH).
Chronic restraint stress (CRS) model
The CRS model was established on the basis of procedures reported previously83 with slight modifications. The mice were placed in well-ventilated 50-mL plastic tubes without food or water from 9:00 AM to 15:00 PM for 14 consecutive days. The plastic tubes prevented the mice from moving freely or turning around but did not squeeze the animals or cause any pain or physical injury. After being restrained, the mice were sent back to their home cages immediately. The mice that were not subjected to restraint were kept in their usual home cages without any restrictions for 14 days.
Chronic social defeat stress (CSDS) model and social interaction test (SIT)
CSDS was performed as followed84. Adult male CD-1 mice (4–6 months old) were used as aggressors. Before each defeat, aggressors were screened for aggressive behavior for 3 consecutive days. Two days before the start of defeat, the CD-1 mice were housed on one side of a perforated Plexiglass partition. During 10 consecutive days of CSDS, experimental mice (7–8 weeks old) were subjected to direct physical interaction with a CD-1 mouse for 10 min per day, and for the rest of the day were placed on the other side of the Plexiglass divider, allowing for sensory but not direct physical contact. The experimental mice were exposed to a new CD-1 aggressor every day for 10 days.
Before the SIT, the mice were placed in the behavioral suite and allowed to adapt for 1 h. After each test, the behavioral suite was cleaned with 75% alcohol to prevent odor interference. SIT was performed 24 h after CSDS modeling in a new CD-1 mouse. The test lasted for 5 min. The experimental mice were allowed to explore freely in the field (44 × 44 cm) in the first 2.5 min, where a rectangular container (10 × 6 cm) was placed. The CD-1 mice were placed in the container 2.5 min later and observed in the “interaction area” (14 × 26 cm). The statistical formula for the social interaction ratio was as follows: time in the interaction area when there was a CD-1 target mouse/time in the interaction area when there was no CD-1 target mouse × 100%.
Virus microinjection
Under continuous isoflurane inhalation anesthesia, the mice were placed in a stereotaxic frame (RWD Life Science Co., Ltd.), and the body temperature was maintained at 36 °C with a heating pad. The mice were then treated with a bilateral stereotaxic injection of the virus into the dmPFC (anteroposterior (AP), +1.7 mm; mediolateral (ML), ±0.4 mm; and dorsoventral (DV), −2.3 mm relative to bregma). The virus was infused bilaterally through a microinjector with a 33 G needle. The virus mixture (150 nl) was infused for 10 min. After infusion, the needle was retained at the injection site for an additional 10 minutes before being withdrawn slowly. The mice were allowed to recover for 3 weeks to allow stable transgene expression.
The shRNA target sequences used in the study are listed in Supplementary Table 1. CMV-PTN-shRNA was constructed to knock down the expression of PTN, and CMV-PTPRZ1-shRNA was constructed to knock down the expression of PTPRZ1. AAV-DIO-PTN-shRNA was used to specifically knock down the expression of astrocytic PTN in the dmPFC in Aldh1l1-Cre/ERT2 mice, with the assistance of intraperitoneally administered tamoxifen (Sigma‒Aldrich, dissolved in corn oil, 75 mg/kg/day for 5 consecutive days) to ensure the activation of Cre recombinase. CaMKIIα-Cre along with AAV-DIO-PTPRZ1-shRNA was applied to specifically knock down the expression of PTPRZ1 in dmPFC excitatory neurons. We constructed a DIO-PTN-Overexpression (DIO-PTN OE) to overexpress the astrocytic PTN in the dmPFC in Aldh1l1-Cre/ERT2 mice. All Aldh1l1-Cre/ERT2 mice received a 10-day rest for recovery after the consecutive administration of tamoxifen.
For sparse labeling of CaMKIIα+ excitatory neurons in the dmPFC, the mice were bilaterally injected with a total of 400 nL of viral cocktail (1:1) of rAAV-CaMKIIα-FLP-WPRE-pA and rAAV-nEf1α-FDIO-EYFP-WPRE-pA.
Exogenous chemicals administration
For infusions into the dmPFC, mice received 1 µl of 5 µg/ml PTN (MedChemExpress). For microinjections of AKT activator SC79 or AKT inhibitor MK2206, 1 µl of 1 mg/ml SC79 (MedChemExpress; dissolved in dimethyl sulfoxide in saline) or 1 µl of 1 mg/ml MK2206 (MedChemExpress; dissolved in dimethyl sulfoxide in saline) was administered.
Behavioral tests
Open field test (OFT)
The mice were placed in an open chamber (50 × 50 × 40 cm), which was made of gray polyvinyl chloride. The mice were allowed to freely explore the apparatus. The mice were gently placed in the center, and movement was recorded for 5 min. The total distance traveled and time spent in the center (25 × 25 cm) were recorded by a video tracking system and analyzed with software (Shanghai Jiliang Software Technology, Co., Ltd.)
Elevated plus maze (EPM)
The apparatus consisted of two opposing open arms (35 × 5 cm) and two opposing enclosed arms (35 × 5 × 15 cm), which were connected by a central platform (5 × 5 cm) and positioned 50 cm above the ground. The activity of the experimental subject within the apparatus was tracked for 5 min via an overhead digital video camera. The time spent in the open arms was recorded and analyzed.
Forced swimming test (FST)
The FST was performed in a glass cylinder (height 30 cm, diameter 12 cm), which was filled to 25 cm with water (22–24 °C). The mice were placed in the cylinder for 6 min. The first 2 min were for the mice to acclimate, and the duration of immobility was recorded during the last 4 min. The immobility time was recorded and analyzed with software (Shanghai Jiliang Software Technology, Co., Ltd.).
Tail-suspension test (TST)
The mice were suspended by the tail 50 cm above the floor. The activity was automatically monitored during the last 4 minutes of the 6-min test with a threshold defining immobility behavior. The latency to the first immobilization was also recorded.
Sucrose preference test (SPT)
The mice were acclimated to two bottles (50-ml tubes with fitted ball-point sipper tubes) for 24 h. One bottle was filled with 1% sucrose, and the other was filled with drinking water. Every 12 h, the bottle positions were swapped to prevent position bias. After another 24 h, the amounts of sucrose and water consumed were recorded, and sucrose preference (%) was calculated as sucrose consumption/(sucrose + water consumption) × 100%.
RT-PCR
TRIzol was used to extract total RNA, and reverse transcription was performed following the protocol of the polymerase chain reaction (PCR) production kit (Accurate Biology, AG 11706). The primer sequences of the investigated mRNAs for the PCR assay are shown in Supplementary Table 2. The reaction cycle conditions were as follows: initial denaturation at 95 °C for 3 min, followed by 40 cycles of 10 s at 95 °C, 20 s at 58 °C, and 10 s at 72 °C. The ratio of mRNA expression in the dmPFC tissues was analyzed via the 2−ΔΔCT method.
Western blotting (WB) analysis
Mouse brain tissues were lysed in ice-cold lysis buffer (RIPA buffer with proteinase inhibitor) for 30 min and centrifuged at 10,000 × g for 20 min at 4 °C. Protein samples were separated by SDS‒PAGE and transferred onto polyvinylidene difluoride (PVDF) membranes (Millipore). The membranes were blocked with 5% milk in TBST for 1 h at room temperature and then incubated with primary antibodies against S100β (1:1000, Proteintech, Cat# 66616), SOX9 (1:1000, Proteintech, Cat# 55152), GFAP (1:1000, Proteintech, Cat# 60190), PTN (1:1000, Proteintech, Cat# 27117), PTPRZ1 (1:1000, Proteintech, Cat# 55125), PSD95 (1:1000, Affinity, Cat# AF5283-50), AKT (1:1000, CST, Cat# 4685), p-AKT (Ser473) (1:1000, CST, Cat# 4060) and GAPDH (1:1000, Fdbio Science, Cat# FD0063) overnight at 4 °C. Information on the antibodies used in the study can be found in Supplementary Table 3. On the next day, the membranes were washed with TBST three times and incubated with horseradish peroxidase (HRP)-conjugated secondary antibodies (Fdbio Science) for 1 h at room temperature. Protein abundance was quantified by analyzing the western blot bands using ImageJ software. Quantified band intensities were normalized to GAPDH levels. Unprocessed scans of the blots are provided in the Source Data file.
Immunohistochemistry (IHC)
The mice were anesthetized with phenobarbital and perfused with 4% PFA in 0.1 M PBS. The brain tissues were postfixed overnight in 4% PFA at 4 °C and transferred to 30% sucrose in 0.1 M PBS. The brain tissues were cut into 20 μm-thick sections on a freezing microtome. Then, the sections were washed with PBS three times and incubated first in a blocking buffer containing 3% bovine serum albumin in 0.2% Triton X-100/PBS for 2 h at room temperature and then with primary antibodies against PTN (1:200, Proteintech, Cat# 27117), PTPRZ1 (1:100, Proteintech, Cat# 55125), S100β (1:200, Proteintech, Cat# 66616), SOX9 (1:200, Invitrogen, Cat# 14-9765-80), and GFAP (1:200, Proteintech, Cat# 60190), in blocking buffer overnight at 4 °C. After that, the sections were incubated with secondary antibody (Cy3 and Alexa fluor 488) at 37 °C for 60 min. The stained sections were examined using with a Nikon confocal microscope equipped, and images were captured with a Nikon DS-Qi2 camera.
Isolation of astrocytes and neurons from the mouse brain
Astrocytes and neurons were isolated from the dmPFC of mice using magnetic-activated cell sorting (MACS)23. In brief, after the mice were anesthetized with pentobarbital sodium and perfused with ice-cold sterile PBS, the dmPFC tissue was dissociated at 37 °C for 15 min using the Adult Brain Dissociation Kit (Miltenyi Biotec, Cat# 130-107-677) with a gentleMACS Dissociator (Miltenyi Biotec, Cat# 130-093-235).
To isolate astrocytes, the cells were incubated with FcR Blocking Reagent after the cell pellet was passed through 70-μm nylon mesh. Next, the cells were incubated with Anti-ACSA-2 MicroBeads (Miltenyi Biotec, Cat# 130-097-679) and then subjected to MACS through the LS column.
To isolate neurons, after the brain tissue was dissociated, the cell pellet was treated with cold debris removal solution and red blood cell removal solution. Following the instructions of the Adult Neuron Isolation Kit (Miltenyi Biotec, Cat# 130-126-603), the cells were incubated with the Adult Non-Neuronal Cell Biotin-Antibody Cocktail and Anti-Biotin MicroBeads and then subjected to MACS through the LS column. The specificity of the Anti-ACSA-2 MicroBeads kit and Adult Neuron Isolation Kit was validated and is shown in Supplementary Fig. 9E, F.
Electrophysiological recording
The mice were anesthetized with isoflurane, and the whole brain was quickly dissociated into prechilled and oxygenated dissection fluid containing (in mM) 213 sucrose, 10 glucose, 26 NaHCO3, 3 KCl, 1 NaH2PO4·2H2O, 10 MgCl2, and 0.5 CaCl2. Acute brain slices (300 µm) containing the dmPFC were acquired in chilled dissection fluid using a microtome (VT1000S, Leica). The sections were transferred to the incubation chamber and immersed in artificial cerebrospinal fluid (in mM; 125 NaCl, 26 NaHCO3, 5 KCl, 1.2 NaH2PO4, 2.6 CaCl2, 1.3 MgCl2, and 10 glucose) at 30 °C for 1 h. After incubation, the sections were placed in the slice chamber for electrophysiological recording with continuous perfusion of artificial cerebrospinal fluid (saturated with 95% O2/5% CO2).
The dmPFC excitatory neurons were identified by CaMKIIα-EGFP and visualized for recording using an infrared-differential interference contrast microscope (Olympus, Japan). We used patch pipettes (4–8 MΩ, WPI, USA) for the whole-cell patch clamp recordings. The signal was amplified by a MultiClamp 700B amplifier (Molecular Devices, USA). The miniature excitability postsynaptic current (mEPSC) was sampled in the presence of tetrodotoxin (1 μM) and picrotoxin (100 μM) at –70 mV. The patch pipettes were filled with an intracellular solution containing (in mM) 122 potassium-gluconate, 5 NaCl, 2 MgCl2, 0.3 CaCl2, 10 HEPES, 5 EGTA, 4 Mg-ATP, and 0.3 Na3-GTP. The action potentials were recorded using the injected current (from 0 pA to 120 pA in 20 pA increments, 500 ms duration) without any synaptic transmission blockers. We calculated the number of action potentials induced by each injected current (shown as firing rate–injected current (f–I) curves), frequency (f–I), and amplitude of the mEPSCs using pCLAMP10.7.
Dendritic spine analysis
The mice were bilaterally injected with a total of 400 nL of a viral cocktail (1:1) of rAAV-CaMKIIα-FLP-WPRE-PA or rAAV-nEf1α-FDIO-EYFP-WPRE-PA to label CaMKIIα+ excitatory neurons in the dmPFC. We used laser confocal photography to photograph dendritic spines, acquired multilayer Z axes, and reconstructed the dendritic spine image. The spines were classified into one of four morphological subtypes: filopodial, thin, stub, or mushroom-shaped. ImageJ was used (http://rsbweb.nih.gov/ij) to calculate the density of thin, stub, and mushroom-shaped dendritic spines. Approximately ten randomly selected neurons were analyzed per condition across two coverslips. The density of spines was scored in dendritic segments 10 μm in length. Finally, we counted and used the number of dendritic spines per 10 μm to describe the density of the dendritic spines.
Statistical analysis
The statistical methods used for bioinformatics analysis in this paper are indicated in the legend of each figure. For experimental validation, SPSS 21.0 (IBM) was used to analyze the data, which are presented as the means ± SEMs. The normality of the data was examined via the Shapiro‒Wilk test. The normally distributed data were tested by a two-sided unpaired t test, one- or two-way analysis of variance (ANOVA). The nonnormally distributed data were analyzed by the Wilcoxon test. Statistical significance was set at *P < 0.05, **P < 0.01, ***P < 0.001, N.S. not significant.
Reporting summary
Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.
