Single-Nucleus RNA Sequencing Data Collection and Processing

A targeted search was conducted on PubMed for research articles until 2022 that contained the terms 'Alzheimer's disease' and 'Single-nucleus RNA sequencing'. Each article retrieved was subjected to an in-depth review. Our selection criteria were as follows: 1) Studies must use human brain tissue samples without any pre-selection for specific cell populations; 2) The sample size should be greater than 10; 3) The study must include both AD and control groups; 4) The research should be supported by peer-reviewed published literature. Finally, three snRNA-seq datasets of AD were included in this study, including the data GSE138852 (Grubman et al. 2019), GSE157827 (Lau et al. 2020) and GSE174367 (Morabito et al. 2021) (Table 1).

Table 1 Dataset characteristics

The dataset GSE138852 included postmortem entorhinal cortex (EC) samples from six AD patients and six sex- and age-matched controls. The GSE157827 dataset comprised 21 prefrontal cortex (PFC) tissue samples, including 12 from AD patients and 9 from control subjects. The GSE174367 included PFC data from 20 postmortem brains, with 12 cases diagnosed with AD. A total of 53 individuals were included in the study, comprising 20 females and 33 males. Among the AD patients, the age ranged from 60 to 95 years, while the control group had an age range of 67 to 94 years.

The processed data of each dataset was separately loaded into R using Seurat (version 4.0.5). The quality of cells was assessed based on five metrics: (1) the number of total UMI count per cell (library size) was below 30,000; (2) the number of detected genes was above 500 and below 6000; (3) the percentage of mitochondrial genes was below 50%. (4) Cells with > 200 genes detected were selected for further analyses. And (5) cells with > 20% of transcripts derived from mitochondria were considered apoptotic, and thus excluded. Then, doublets were removed using DoubletFinder with default settings for each sample.

After discarding low-quality cells, the data was processed using the Harmony package in R for data integration with default parameters. The integration involved aligning cells from different batches or experimental conditions based on the "orig.ident" variable. The convergence plot was utilized to visualize the convergence during the Harmony algorithm's execution. The k-means algorithm was initialized with 20 initial cluster assignments and allowed for a maximum of 10,000 iterations to converge. The Harmony algorithm itself was set to iterate a maximum of 200 times for convergence. Following the completion of the Harmony run, further dimensionality reduction and clustering were performed on the integrated data to identify cellular subpopulations. Cell clusters in each sample were identified by examining the top 30 principal components (PCs) across highly variable genes (HVGs). The markers for each cell cluster were identified using the Seurat 'FindAllMarkers' function, and singleR was applied for cell cluster annotations. To categorize the major cell lineages, a set of marker genes based on previous studies was utilized. For astrocytes, markers included GFAP, FGFR3, GJA1, AQP4, ALDH1L1; for excitatory neurons, SYN3, RBFOX3, CAMK2A; for inhibitory neurons, ERBB4, NXPH1, GAD1, GAD2; for microglia, HLA-DRA, CX3CR1, C1QB, CSF1R; for oligodendrocytes, MOBP, MBP, PLP1; for oligodendrocyte progenitor cells (OPCs), PCDH15, MEGF11; and for endothelial cells, FLT1, CLDN5 (Chen et al. 2023; Darmanis et al. 2015).

After determining the major cell lineages, subclustering of astrocytes was conducted using Seurat. Different resolution parameters (from 0.1 to 2.0) were tested. For each parameter, the marker genes (Log fold change > 0.25 and adjusted p-value < 0.05) of each subcluster were checked. To obtain as many subclusters as possible, the highest resolution parameter was chosen when the marker genes of any subcluster showed less than 30% overlap with those of any other subclusters. Additionally, some manual adjustments were made to avoid over-clustering. Signature scores were then computed using the Seurat 'AddModuleScore' function, utilizing the gene signature of interest (DAAs marker and A1-specific transcripts) to identify the neurotoxic astrocytes subcluster.

Gene Ontology (GO) Enrichment and Kyoto Encyclopedia of Genes and Genomes (KEGG) Pathway Analysis of the Subcluster-Specific Genes

GO enrichment analyses were performed in R using the clusterProfiler package. Metascape was used to perform the KEGG pathway analysis. Functional and pathway enrichment analyses were conducted for subcluster-specific genes of neurotoxic astrocytes. In this analysis, an adjusted p-value < 0.05 was considered significant for the screening of significant GO terms.

Bulk RNA-Sequencing Data Collection and Co-Differentially Expressed Genes Identification

Bulk RNA-sequencing data of AD patients with complete follow-up information were collected from Gene Expression Omnibus (GEO) databases. According to the following criteria, datasets were considered eligible for our analysis: (1) datasets with AD samples; (2) datasets supported by peer-reviewed PubMed-indexed publications; (3) datasets with brain tissue from the EC or PFC. Two datasets related to AD, GSE33000 and GSE48350, were selected for co-differentially expressed genes analysis. Each dataset was isolated to minimize batch effects and inter-study variability. A comprehensive preprocessing routine was then applied to each dataset. The impute package was used to supplement missing data (Troyanskaya et al. 2001). Next, the normalize Between Arrays function in the limma package was used to normalize gene expression. DEGs in each dataset were identified using the limma package in the R computing environment, determined by an absolute log2 fold change (| Log 2 FC|) greater than 0 and an adjusted p-value less than 0.05 (Ritchie et al. 2015). Volcano plots were generated using ggplot2 in R. Finally, the Venn Diagram, a web-based tool, was employed to identify key genes from upregulated DEGs and subcluster-specific genes of neurotoxic astrocytes. Interactions among key genes and neurotoxic astrocyte-specific genes overlapping with upregulated DEGs were investigated using STRING (selecting high confidence 0.7) for protein–protein interactions (PPI) network analysis, and the results were visualized in Cytoscape (version 3.10.1).

Clinical Correlation Analysis and Gene Set Enrichment Analysis

We enrolled the dataset GSE106241 and compared the expression level of key genes in different disease stages. We investigated their associations with alpha-secretase activity, gamma-secretase activity, β-secretase activity, braak stage, and Aβ42 levels in samples from GSE106241 using the Spearman correlation analysis. The Agora’s Gene Comparison Tool on Synapse.org was used for key gene validation. Next, Gene Set Enrichment Analysis (GSEA) was conducted on the GSE48350, GSE33000 and GSE106241 datasets to identify Biological Process (BP) and GO terms associated with the key genes. The analysis was performed using the clusterProfiler package in R. As reference gene sets, the c5.bp.v7.0.symbols.gmt datasets from the MsigDB Version 6.2 database were utilized. Gene sets that showed a significantly enriched result with an adjusted p-value < 0.05 after 1,000 permutations were considered as significant findings.

Animals

The 5 × FAD model mice and wild-type (wt) mice, six in each group, were utilized for the study. The care, feeding, and handling of all animals were performed in adherence to the Guide for the Care and Use of Laboratory Animals from Chongqing Medical University. The experimental design was conducted randomly and in a double-blind fashion with adult mice. All experiments were designed to minimize animal suffering and to limit the number of animals used.

Primary Perinatal Astrocyte Culture and Isolation

Primary astrocyte cultures were prepared from the cerebral neocortex of P0-P3 neonatal mice, as previously described with slight modifications (Huang et al. 2016). Briefly, brain tissue was isolated from the skull, meninges were gently removed, and the tissue was minced and passed through a 70 μm filter to generate a single-cell suspension. The mixed glial cells were then cultivated in T75 flasks with 10% heat-inactivated fetal bovine serum (CTCC-002–071,Meisen, Zhejiang, China) and 1% penicillin/streptomycin (BL505A, Biosharp, Hefei, China) in DMEM (21,068,028, Gibco, Massachusetts, USA) for 12–14 days.

Preparation of Aggregated Aβ42

The toxic oligomers of Aβ42 (107,761–42-2, Sigma, USA). were prepared as described (Stine et al. 2011). Briefly, peptide was dissolved in hexafluoroisopropanol, dried under vacuum, and stored at − 20 °C. Immediately prior to use, the peptide residue was reconstituted in DMEM media to make a stock solution at 0.1 mM and incubated at 4 °C for 24 h to form diffusible oligomers. Aβ42 oligomers at a final concentration of 2 μM were assayed for astrocytes.

Transfection

WW domain-containing transcription regulator 1 (WWTR1) overexpression plasmid was constructed by TSINGKE (Wuhan, China). Primary astrocytes were transfected with WWTR1 overexpression plasmids. Cells were cultured at 60%-70% confluence in 10 cm plates and were transfected using Lipofectamine 3000 (L3000001, Invitrogen, Boston, USA). The empty vector plasmids were transfected in the same way. At 48 h post-transfection, the efficacy of transfection was identified by RT-PCR.

Western Blotting

Western blotting (WB) was performed according to the description mentioned earlier. Brain tissues from 9-month-old wt and 5 × FAD mice were collected. Tissue and cell proteins were extracted using RIPA lysis buffer (P0013B, Beyotime, Shanghai, China) and phenylmethanesulfonyl fluoride (PMSF, ST506, Beyotime, Shanghai, China). After centrifugation at 16,000 rpm (4 °C), the supernatant was collected and stored at − 80 °C. The protein concentration was determined using the BCA protein assay kit (P0010, Beyotime, Shanghai, China). The supernatant samples were separated on a 10% SDS-PAGE gel and transferred onto a 0.22 μm polyvinylidene fluoride (PVDF) membrane (GVHP29325, Millipore, Billerica, MA, USA). The membrane was blocked with 5% skim milk in tris-buffered saline with Tween-20 (TBST) at room temperature for 2 h. Then, it was incubated with the primary antibody (diluted in TBST buffer at 1:3000) specific for WWTR1 (66,500–1-Ig, Proteintech, Wuhan, China, RRID: AB_2881864) overnight at 4 °C. The next day, the membrane was washed and then incubated with the PBST secondary antibody (1:10,000, SA00001-1, Proteintech, Wuhan, China, RRID: AB_2722565) at room temperature for 1 h. The bands were visualized using an enhanced chemiluminescence reagent (WBKLS0100, Thermo, Marina, CA, USA) and an image analysis system (Bio-Rad, USA). The membrane was then stripped using a stripping buffer (P0025, Beyotime, Shanghai, China), re-blocked and incubated with the HRP-conjugated β-actin antibody (1:10,000, HRP-60008, Proteintech, Wuhan, China, RRID: AB_2819183). Relative quantification analysis was performed using ImageJ software and reference proteins.

Immunostaining

The procedure involves fixing brain slices with cold 4% paraformaldehyde, followed by washing and permeabilization using a 0.5% Triton X-100 PBS solution, and blocking with 5% BSA. Subsequently, the samples are incubated with specific antibodies, including mouse anti-WWTR1 (1:500, 66,500–1-Ig, Proteintech, Wuhan, China, RRID: AB_2881864), rabbit anti-GFAP (1:500, 16,825–1-AP, Proteintech, Wuhan, China, RRID: AB_2109646), rabbit anti-NeuN (1:500, 26,975–1-AP, Proteintech, Wuhan, China, RRID: AB_2880708), and rabbit anti-Iba1 (1:500, 10,904–1-AP, Proteintech, Wuhan, China, RRID: AB_2224377). The secondary antibodies used in this experiment are DyLight 488 Goat Anti-Mouse IgG (1:500, A23210, Abbkine Inc., Wuhan, China) and Dylight 594 Goat Anti-Rabbit IgG (1:500, A23420, Abbkine Inc., Wuhan, China). DAPI (1:1000, C1002, Beyotime, Shanghai, China) is used for nuclear visualization. The primary antibodies used in this study were validated and widely used in previous publications, including WWTR1 (Hu et al. 2021), GFAP (Boivin et al. 2021), IBA1(Xu et al. 2021a), NEUN(Li et al. 2021), HRP-conjugated β-actin (Shiraishi et al. 2023). To clarify the cellular localization of WWTR1, three random images were taken at a 60X magnification and then analyzed using ImageJ software. To analyze the differences in WWTR1 expression within the brain sections, the stained samples were scanned by using an automatic slide scanner (VS200, Olympus, Tokyo, Japan). The digital images were analyzed using the Visiopharm software module (Visiopharm, Hørsholm, Denmark). The whole cortical and hippocampal regions of the brain sections were selected separately for this purpose. The software automatically identified all nuclei expressing GFAP, and these were subsequently auto-labeled. Nuclei co-expressing GFAP and WWTR1 were marked in red, while those expressing GFAP but not WWTR1 were marked in white. The software's AI capabilities were utilized to calculate the proportion of WWTR1-positive cells in the WT and 5 × FAD groups.

Quantitative Real-­Time PCR (qRT-PCR)

Total RNA from the mouse cortex was extracted using an RNA-Easy™ Isolation Reagent Vazyme Cat (RC112-01, Vazyme, Nanjing, China). HiScript II Q RT SuperMix for qPCR (+ gDNA wiper) (R223-01,Vazyme, Nanjing, China) was used for reverse transcription. Quantitative PCR was performed using Universal SYBR Green Fast qPCR Mix was acquired from ABclonal (RK21203, Wuhan, China). The fold change of gene expression was calculated using the 2−ΔΔCt method. The primers were synthesized by tsingke Biotech and presented in Table S1.

Statistical Analysis

Sample size calculations utilized the formula E = total number of animals per group—number of groups. Each experiment ensured a value of E > 10(Charan and Kantharia 2013). Following previous research guidelines and setting a significance criterion at α = 0.05 and power at 0.80, the minimum required sample size with this effect size for a Wilcoxon-Mann Whitney test between two groups was determined to be N = 3 (Reid et al. 2023). In order to bolster data robustness, 5–6 animals were included in each experimental group. According to existing literature (Reid et al. 2023), the least sample size required is N = 4. Nevertheless, factoring in the diversity in effect size and data, it is postulated that a sample size of N = 6–7 per group would be adequate. The absence of a priori sample size calculation emerged as a potential limitation in this method. Parametric data was presented as mean ± SD and non-parametric data as medians ± interquartile range. Normality was evaluated using the Shapiro–Wilk test and variance homogeneity was assessed by the Levene test. For comparing two groups, when the data was normally distributed and exhibited homogeneity of variances, we utilized a two-tailed unpaired t-test. In cases where the data was normally distributed, but variances were not equal, Welch's t-test was employed. All statistical analyses were conducted using GraphPad Prism software, with p-values of < 0.05 deemed statistically significant. To preserve the objectivity of the research, the authors were blinded to the experimental protocol and did not have access to statistical calculations during the execution of experiments. The detailed information on statistical tests was in the Table S2-5.