Microglia Impede Oligodendrocyte Generation in Aged Brain

Introduction

Extensive brain atrophy is a characteristic manifestation of aged brain.^1,2 White matter degeneration is accompanied by encephalatrophy, leading to irreversible neurological and cognitive impairments.^3,4 Remyelination is a natural protective and regenerative process that will be initiated in response to the degeneration of white matter.⁵ It is a complex process involving oligodendrocyte precursor cell (OPC) activation, migration, differentiation, and maturation to be oligodendrocytes (OLGs).⁵ Accumulating evidence indicates that remyelination is disturbed in the aged brain.⁶ However, a clear understanding of the effect of brain aging on the process of remyelination is still lacking.

Brain aging heightened the neuroinflammatory profile of the cerebral microenvironment, which activated microglia by attributing it to aging-related changes in remyelination.⁷ Previous studies have reported that activated microglia secrete a heterogeneity array of signaling molecules, including nitric oxide, reactive oxygen species, Il-6, Il-1β, and tumor necrosis factor (TNF), contributing to myelin damage and hindering proliferation or differentiation of OPCs.^8–10 However, other studies indicated that activated microglia could promote OLG survival or OPC differentiation by releasing several regenerative factors including Igf-1,¹¹ Igf-2,¹² galectin-3,¹³ activin-A,¹⁴ and Il-1β,¹⁵ and clearing myelin debris.¹⁶ These diverse results may be driven by the heterogeneous subpopulation of microglia, OLGs, and OPCs.¹⁷ Specific subgroups of OLGs and OPCs in different states may respond inconsistently to the stimulus of activated microglia.⁵ In addition, our previous study revealed that there were six subgroups of microglia with divergent functions in aged brain.¹⁸ A unique type of highly activated microglia was observed in aged mice only, with functional implications in immuno-inflammatory response.¹⁸ It remains poorly understood how this specific age-related subgroup of microglia effect on remyelination in aged brain. The exact mechanism of aged microglia in regulation of distinct subgroups of OLGs and OPCs needs to be determined.

In the present study, we aimed to explore the subclusters of OLGs and OPCs by analyzing the single-cell RNA sequence (scRNAseq) data of both young and aged brains. Then, cell–cell interaction analysis was performed to examine the effect of aged microglia on the OLG survival and OPC differentiation. Finally, we expected to screen several factors that could transform senescent microglia into tissue repair type to promote oligodendrocyte generation.

Materials and Methods

Animals

Young male C57BL/6 mice (8–10 weeks old) were purchased from SLAC Laboratory Animal Company Limited (Shanghai, China). Aged male C57BL/6 mice (18 months old) were purchased from Beijing Vital River Laboratory Animal Technology Co. Ltd (Beijing, China). Mice were maintained in a humidity- and temperature-controlled animal facility, under a controlled photoperiod (12-h light and 12-h dark). The animal experiment protocol was approved by the Institutional Ethics Committee of the Second Affiliated Hospital, Zhejiang University School of Medicine. All procedures followed the Guide for the Care and Use of Laboratory Animals (National Institutes of Health). All efforts were made to minimize animal suffering and the number of animals used.

Diffusion Tensor Imaging (DTI) Preprocess and Tractography

DTI imaging data of 100 volunteers were obtained from “The Human Connectome Project (HCP)-A Mapping the Human Connectome During Typical Aging” database (https://nda.nih.gov/), which included 53 volunteers aged 20–60 years old and 47 volunteers aged >60 years. The DTI preprocess and tractography were performed by DSI studio (http://dsi-studio.labsolver.org). A multi-shell diffusion scheme was used, and the b-values were 1500 and 3000 s/mm². The number of diffusion sampling directions was 46. The in-plane resolution was 1.5 mm. The slice thickness was 1.5 mm. The diffusion tensor was calculated including fractional anisotropy (FA), mean diffusivity (MD), axial diffusivity (AD), and radial diffusivity (RD). A deterministic fiber tracking algorithm was used.¹⁹ A seeding region was placed at the whole brain. The anisotropy threshold was 0.2. The angular threshold was 55 degrees. Tracks with lengths shorter than 20 or longer than 200 mm were discarded. A total of 5000 seeds were placed.

Immunostaining of Brain Sections

Mice were deeply anaesthetized and transcardially perfused with 25 mL ice-cold phosphate buffered saline (PBS). The whole brain was immersed in 4% paraformaldehyde for 24 h, then dehydrated in serial 15% and 30% sucrose solutions at 4°C. The brain samples were then cut into 20-μm-thick coronal slices and collected in cryoprotectant (40% PBS, 30% glycerol, 30% ethylene glycol). The slices were kept at −20°C until they were immunostained. Sections were washed with PBS twice, and then 5% BSA and 0.2% Triton-X in PBS were used to block the sections for 2 h at room temperature. The sections were incubated at 4°C overnight with anti-MBP (Santa Cruz, sc66064, 1:200), anti-IBA1 (Abcam, ab5076, 1:250), anti-phospho-STAT1 (Tyr701) (absin, abs130925, 1:200), anti-PDGFRA (Abcam, ab203491, 1:200), anti-CC1 (Calbiochem, OP80, 1:200). Donkey secondary antibody conjugated with Alexa Fluor 488 or 594 (Invitrogen, 1:500) Fluoro were utilized to incubate the sections at room temperature in the dark for 1 h. For OPCs and OLGs, TUNEL Assay Kit (Abcam, ab66110) was used for TUNEL staining. Briefly, before incubated with secondary antibody, the sections were labeled with Br-dUTP for 1 hour at 37°C, and incubated with anti–BrdU-Red antibody for 30 minutes at room temperature. The sections were then washed three times with PBS and mounted on glass slides and coverslipped with Fluoromount-G containing DAPI (Yeasen Biotech).

Image Acquisition and Quantification

Sections were observed and analyzed using a Leica TCS SP8 confocal microscope (Leica Microsystems). For the statistical analysis, the area of layer I–III motor cortex was acquired from each sample for MBP intensity as described by Wang et al.²⁰ Quantification of MBP intensity was performed using Fiji 2.1.0/1.53c. The areas of cortex, corpus callosum, external capsule and striatum were acquired from each sample for quantification of OPCs, OLGs and microglia. Numerical counting was performed manually by an investigator blinded to experimental groups. Brightness and contrast of images were adjusted using Fiji 2.1.0/1.53c.

Transmission Electron Microscopy

Transmission electron microscopy (TEM) was performed to measure the myelin thickness in the external capsule (EC) area. Briefly, mice were deeply anaesthetized and transcardially perfused with ice-cold saline, followed by 4% PFA and 2.5% glutaraldehyde. EC area tissue was then dissected into 1 mm³ blocks. Following overnight fixation in 2.5% glutaraldehyde, tissue was washed 3 times in PBS and then post-fixed in 1% OsO4 for 1.5 h. The tissue was then washed 3 times in PBS and stained in 2% uranyl acetate for 0.5 h. A graded series of 50–100% ethanol, 100% propylene oxide was used to dehydrate the tissue. After that, the tissue was infiltrated in 1:3 mixture of propylene oxide: epoxy resin overnight at room temperature. Tissue was embedded in molds after several changes of 100% resin for 24 h. Ultrathin 60 nm sections were cut on a Leica UC7 ultramicrotome. Sections were imaged using a transmission electron microscope at 100 kV. Two sections from each animal were analyzed. Five to 12 images were randomly selected within the EC area from each section. The magnification was 4200x. The results were analyzed by a researcher blinded to experimental groups with ImageJ. At least 40 random axons per animal were analyzed. G-ratios are defined as the ratio of the inner axonal diameter to the total outer diameter.

Basic Processing and Clustering Analysis of Single-Cell Transcriptome Data

ScRNAseq dataset in the gene expression omnibus (GEO) database containing the expression matrix of the whole brain in young and aged mice (GSE129788) was used for downstream analysis. Seurat package (v.3.2.2) was used for basic processing and visualization of the scRNAseq data in R (v.4.1.0).²¹ Cells with genes less than 200 or more than 2500 or containing >10% of mitochondrial genes were excluded from analysis. The data were then normalized by the total expression before they were log-transformed. FindVariableFeatures function with default parameters was performed to detect highly variable genes. Linear scaling was subsequently conducted. Unique molecular identifier (UMI) counts and fraction of mitochondrial reads were removed using the ScaleData function. Principal component analysis (PCA) was performed on the scaled data for dimensional reduction, and top 20 principal components were included. Clusters were identified with the FindClusters() function by use of the shared nearest neighbor modularity optimization with a clustering resolution set. Resolution sets between 0.1 and 1.2 were tested for clustering of whole brain cells, OPCs and OLGs. Heatmap of differently expressed genes with each resolution was performed to visualize the clusters. Finally, resolution set 0.1 for whole brain cells and 0.3 for OLGs generally yielded good clustering. Resolution 0.4 divided OPCs into 6 groups. Three out of the 6 groups showed no significant functional differences with GO analysis so they were merged into the OPC C1 cluster. Microglia were clustered as mentioned in Jin et al.¹⁸ Results were visualized with non-linear dimensional reduction methods including uniform manifold approximation and projection (UMAP) and t-distributed stochastic neighbor embedding (t-SNE). Average scores on different gene sets were evaluated with AddModuleScore function in Seurat. CellCycleScoring function was used to set the cell cycle phase of each cell based on its expression of G2/M and S phase markers.

Differentially Expressed Gene Expression

After initial quality-control preprocessing and determination of cellular identities, differentially expressed genes (DEGs) between distinct subgroups were found using the FindAllMarkers (or FindMarkers) function that ran Wilcoxon rank sum tests with Bonferroni adjust P value. In order to detect significance in very small changes in transcription, we used the MAST package (v.1.6.1) in R (v.4.1.0) to perform DEGs analysis between aged and young group.²² MAST generated False Discovery Rate (FDR) adjusted P values and logFC (based on natural log of the FCs), using a hurdle model with normalized nUMI as a covariate. Due to the Bayes shrinkage approach leveraged by MAST, significance of subtle transcriptional alterations could be detected. We defined significant DEGs as those adjusted P-value <0.05 and an absolute value of log fold change exceeding 0.1.

Gene Ontology Enrichment Analysis

Gene Ontology (GO) enrichment analysis was performed with the online tool Metascape (http://metascape.org).²³ All genes in the mouse genome are used as the enrichment background. Briefly, a set of differentially expressed genes were submitted to Metascape, and a list of statistically GO enriched terms would be returned. For the predicted activation state of the GO terms, the activation Z-score was calculated for each term by the R package GOplot (v.1.0.2).²⁴ A term was predicted to be significantly changed with an absolute value of Z score ≥2 and a P-value <0.05. The results were visualized with bubble plot or bar plot using GOplot (v.1.0.2) or prism software (v.8.4.3, GraphPad Software Inc.).

Pseudo-Time Analysis

The Monocle 2 R package (v.2.14.0) was applied to determine the developmental pseudo-time of oligodendrocytes and OPCs.^25–27 Briefly, all the significantly upregulated DEGs were defined as ordering genes. Discriminative dimensionality reduction with trees (DDRTree) algorithm was performed with reduceDimension function in Monocle to reduce the dimensionality of data. The cellular trajectory was generated using plot_cell_trajectory. Gene expression was plotted across pseudo-time using the function plot_genes_in_pseudotime.

Cell–Cell Interaction Analysis

IntercellDB package developed by us before was used for cell–cell interaction analysis of scRNA-seq data. Briefly, the significant DEGs in each cluster calculated as described above were mapped to the mouse gene pairs reference database, which contained a list of data sources, including GO, STRING, UniProt, NCBI-gene, COMPARTMENTS, etc. Derived gene pairs were visualized with the sankeyplot or dotplot using R package IntercellDB or networkD3 (v.0.4).

Network Analysis of Protein–Protein Interactions

Gene symbols obtained from IntercellDB or transcription factors inferred from TRRUST v2 (https://www.grnpedia.org/trrust/) were used to investigate the direct and indirect interactions with target proteins.²⁸ CytoScape (v.3.8.0) and String database (https://string-db.org/) were used to generate the interaction network and analyze the network topologies.

Processing and Analysis of Other Datasets

Dataset from Marques et al was downloaded from GSE75330.¹⁷ Following normalization and natural log transformed, the expression of marker genes for oligodendrocytes and OPCs were extracted. Pearson’s correlation coefficients were calculated between gene expression lists of distinct subgroups from Marques et al and the brain data we used. Datasets and white matter-associated microglia (WAM) signature genes from Safaiyan et al were available through their supplementary material.²⁹ Pearson’s correlation coefficients were calculated on normalized gene expression lists similar to previous process.

SCENIC Analysis

SCENIC analysis was applied to evaluate the activities of different transcription factors (regulonAUC) and conducted as described previously.³⁰ We used the pySCENIC package (v.0.11.2), a lightning-fast python implementation of the SCENIC pipeline. GRNBoost was applied to run the co-expression modules. Gene-motif rankings (10 kb around the transcription start site (TSS)) were used to determine the search space around the TSS for RcisTarget. Heatmap of normalized regulonAUC was plotted in R (v.4.1.0).

Statistical Analysis

Wilcoxon rank sum test was used as the statistical method for single-cell RNA-seq data analysis. For statistical comparison of the means between the two groups, the student’s t test or the Wilcoxon Rank Sum test (both two-tailed) was carried out. Other data were presented as mean ± standard deviation. All statistical analyses were performed with R (v.4.1.0) or GraphPad Prism (v.8.4.3). Statistical significance was defined as a P-value <0.05.

Results

Brain Aging Leads to Degeneration of Myelin and Decrease of OPCs

DTI was utilized for the assessment of white matter integrity in aged and young brain (Figure 1A–C). A total of 53 young (48.24±7.32 years) and 47 aged (72.91±8.98 years) adults underwent DTI examination. Representative DTI tractography showed obvious atrophy fiber in aged brain compared with young brain (Figure 1A). The diffusion tensor of the whole brain was calculated in both young and aged subjects. Consistent with previous studies,³¹ significant decreased intravoxel FA accompanied by increases in MD, AD, and RD were observed (P-value <0.0001, Figure 1B). Furthermore, the linear regression revealed a positive correlation between RD and age (P-value <0.0001, Figure 1C). These results indicated age-related loss of myelin sheath in the human brain. To further explore myelin dynamics during aging, myelin basic protein (MBP) expression was examined in layers I–III of the motor cortex. MBP intensity was significantly decreased in aged versus young mice (Figure 1D). Transmission electron microscopic (TEM) analysis was performed to detect the change in myelin integrity. Many axons with thin myelin sheaths were observed in aged mice (Figure 1E). Increases in the g ratio were observed in aged versus young mice (Figure 1E), indicating a decrease in myelin thickness. Scatterplots of g ratio as a function of axon diameter also revealed higher g ratios in aged mice (Figure 1E). Together, these data suggested that degeneration of myelin occurred in both aged human and mice. We further examined PDGFRA-positive OPCs and CC1-positive mature OLGs in 8–10 week-old and 18-month-old brains (Figure 1F). We found the number of OPCs decreased in aged brain, while OLGs remained constant (Figure 1F). To investigate the alterations of OLGs and OPCs during normal aging, we analyzed the scRNAseq data of whole brain cells from young and aged mice (Supplementary Figure 1A). On the basis of quality control, 37,069 cells from whole brain of eight young and eight aged mice were used for unsupervised clustering analysis. A total of 12,360 OLGs and 2,221 OPCs were extracted for downstream analysis (Figure 1G and Supplementary Figure 1B). Differential expression analysis showed that OPCs highly expressed Pdgfra, Cspg5, and Ptprz1 as previously reported,^17,32 and OLGs specifically expressed classical myelin related genes including Mal, Mog, and Cldn11 (Supplementary Figure 1C).^17,33 Consistent with the results from immunofluorescence, the proportion of OPCs in whole brain cells was decreased, while the proportion of OLGs remains constant comparing aged brain to young brain (Supplementary Figure 1D). The ratio of OPCs and OLGs decreased significantly in aged brain (Supplementary Figure 1E). These data indicated that brain aging might damage the potential ability of remyelination by reducing the number of OPCs. To further elucidate age-related molecular changes of OLGs and OPCs, we identified DEGs (adjusted P-value <0.05) in OPCs and OLGs respectively. The functional implications of these DEGs were further explored by GO enrichment analysis. Several crucial functions in cell division, microtubule cytoskeleton organization, and ATP synthesis were down-regulated in aged OPCs compared with those in young mice (Figure 1H and Supplementary Figure 1F, Supplementary Table 1). Other functions including central nervous system myelination and cholesterol metabolic process were downregulated in aged oligodendrocytes (Figure 1H and Supplementary Figure 1G, Supplementary Table 1). Downregulation of these biological processes will lead to inadequate remyelination.^34,35 We further applied the IntercellDB toolkit to determine the effect of brain resident cells on OPCs and OLGs. Microglia and endothelial cells were the main resident cells which impacted on the oligodendrocyte lineage cells (Figure 1I). These results provided evidence that age-associated alterations in cerebral microenvironment may influence the status of OPCs and OLGs thus hinder remyelination process.

Figure 1 Age-related changes of myelin sheath and oligodendrocyte lineage cells in CNS. (A) Representative DTI-based white matter mapping images from young and aged adults. With the presented color-coding, red represents fibers crossing from left to right, green crossing in the posterior anterior direction, and blue crossing in the inferior–superior direction. (B) Comparison of DTI derived metrics between young and aged adults (****P-value < 0.0001 by student’s t test). (C) Linear regression plot showing the change in radial diffusivity with age (P-value < 0.0001). (D) Representative confocal microscopic images and quantification of mean MBP intensity (**P-value < 0.01 by student’s t test, n = 5 samples per condition). Scale bar: 50μm. (E) Representative TEM images in the external capsule for young and aged mice. Red asterisks and red arrows showing myelin thickness changes of axons with similar diameters, respectively. Scale bar: 2μm. Violin plots showing the g ratio of myelinated axons in aged mice versus young mice (****P-value < 0.0001 by student’s t test. The numbers of axons from four mice are indicated). Scatterplots depicting the g ratio as a function of axon diameters. (F) Representative confocal microscopic images in young and old cortex staining for PDGFRA and CC1. Scale bar: 100μm. Bar plots depicting the quantification of PDGFRA+ and CC1+ cells (** P-value < 0.01 by student’s t test. n = 5 samples per condition). (G) T-distributed stochastic neighbor embedding (t-SNE) projection of OPCs and OLGs from young and aged mice. (H) Gene ontology (GO) enrichment network plot was achieved using Metascape. Two gene lists of DEGs between aged and young OPCs or OLGs were submitted. The significantly overexpressed GO terms were grouped into several clusters based on their membership similarity. The contribution percentage of different cells was shown in each GO node. (I) Sankey diagram showing the interactions between resident cells in brain and OPCs or OLGs on GO term immune response (GO:0006955). The proportional flow represented the number of gene pairs. Abbreviations: CNS, central nervous system; FA, fractional anisotropy; MD, mean diffusivity; AD, axial diffusivity; RD, radial diffusivity; OPCs, oligodendrocyte precursor cells; OLGs, oligodendrocytes; MG, microglia; EC, endothelial cells; MAC+ImmN, macrophages and immature neurons; ASC, astrocytes; EPC, ependymocytes; mNEUR, mature neurons; MBP, Myelin basic protein; WB, whole brain; DTI, diffusion tensor imaging; TEM: Transmission electron microscopy.

Distinct Age-Related Transcriptomic Changes of OLGs Involving in Cell Survival and Remyelination

Myelin sheath was mainly constituted by mature OLGs, which were reported to be impaired in aged brain. Unsupervised sub-clustering analysis was applied to explore the heterogeneity of OLGs in both young and aged brains. The t-distributed stochastic neighbor embedding (t-SNE) map revealed 6 distinct subclusters of OLGs (Figure 2A and Supplementary Figure 2A) based on the enriched genes for each cluster (Supplementary Figure 2B). All six subclusters expressed OLG markers including Mbp, Plp1 and Cd9 (Supplementary Figure 2C).³³

Figure 2 Subgroup division and age-related transcriptomic alterations of OLGs. (A) T-distributed stochastic neighbor embedding (t-SNE) projection of six OLG subclusters from young and aged mice. (B) Dot plot showing the proportion of each subcluster in young and aged OLGs (**P-value < 0.01, *P-value < 0.05 by Wilcoxon Rank Sum test, n = 8 samples per condition). (C) Cellular trajectory of all OLGs subclusters, generated by Monocle DDRTree dimensionality reduction algorithm. (D) Pseudo-time trajectory plots showing the expression of representative well-known immature (Opalin) and mature (Klk6) oligodendrocyte marker genes. (E) Heatmap of Pearson correlation between each OLG subclusters with the corresponding populations mentioned by Marques et al (F) Heatmap depicting the transcriptomic alterations of well-known myelin protein or myelination regulator genes in aged OLGs compared to young ones. (G) Heatmap showing part of significant (P-value < 0.05) aging-related GO enrichment terms in oligodendrocyte subclusters. Numbers in legend represented Z-score, positive values indicated upregulation and negative ones indicated downregulation. (H) Heatmaps depicting the average gene expression on aging, apoptosis, cholesterol synthesis and remyelination module in subclusters from young and aged OLGs, respectively. (I) Representative confocal microscopic images in young and old cortex staining for CC1 and TUNEL. Scale bar: 50μm. Abbreviations: OLGs, oligodendrocytes; MOL, mature oligodendrocyte; MFOL, myelin forming oligodendrocyte; GO, gene oncology.

Four out of six subsets showed significant alterations in cell numbers (Figure 2B). Specifically, the cell number of the subset C2, C4, and C5 was decreased, while the cell number of the subset C3 dramatically increased in aged brain (Figure 2B). The subset C0 occupied about 50% of the total cell number with comparable cell numbers between young and aged mice (Figure 2B). Those DEGs between C0 cluster and the other subclusters of OLGs were included in GO enrichment analysis, which showed downregulation of various biological processes (Supplementary Figure 2D). Hence, we considered the C0 subcluster as OLGs in a resting state. With function features from GO enrichment analysis using DEGs between other subsets compared to the C0 subset, we defined the C1 subset as ‘Axon-related OLGs’ with marked protein Clasp2 involving in morphogenesis associated microtubule-binding, the C2 subcluster as ‘Mature oligodendrocytes’ with mature marker Klk6, the C3 subcluster as “Stress-related oligodendrocytes” with stress associated protein Hspa1a, the C4 subset as “Transmembrane transport-related oligodendrocytes” with Na⁺/K⁺ ATPase subunit Atp1a2, and the C5 subset as “Premature oligodendrocytes” with Tmem141 (Supplementary Figure 2C and E, Supplementary Table 2).¹⁷ We further performed monocle pseudo-time analysis revealing a gradual transition from immature OLGs from the C5 cluster (marked by Opalin) to mature OLGs from the C2 cluster (marked by Klk6) (Figure 2C and D). A strong correlation was observed between our OLG subclusters and myelin-forming cells mentioned by Marques et al (Figure 2E).¹⁷ Consistent with our findings of monocle pseudo-time analysis, the C5 subset was myelin forming OLGs (MFOL) which was newly formed from OPCs, and the other subsets were mature OLGs (Figure 2E).

To investigate the biological alterations of OLG subsets by brain aging, we calculated DEGs for each subset of OLGs between aged and young mice. Genes encoding myelin proteins or myelination regulators were down-regulated in all subsets of OLGs (Figure 2F). GO enrichment analysis was further performed to explore the biological function changes of OLG subsets influenced by brain aging, which identified three major functional categories including ribosome biogenesis, proteometabolism, and myelin biosynthetic (Figure 2G). Specifically, aged OLGs from the C2 subset showed down-regulation of biosynthetic process in glycosylceramide, galactolipid, and galactosylceramide, which were all essential steps in myelination (Figure 2G).³⁶ Furthermore, we calculated scores on several gene sets of aging, apoptosis, cholesterol synthesis, and remyelination to explore the age-related changes of each OLG cluster (Figure 2H and Supplementary Table 3). Aged OLGs from all six subsets got higher aging score compared to their corresponding young status (Figure 2H). Aged OLGs from the C3 subset showed an obvious tendency to apoptosis (Figure 2H). Immunofluorescence staining also revealed a few apoptotic OLGs in aged brain (Figure 2I). Cholesterol synthesis is one of the important steps in myelination and remyelination.³⁵ Genes encoding cholesterol synthesis-related proteins including Hmgcs1, Hmgcr, and Mvk, were highly expressed in OLGs from the C2 and C5 subsets and decreased with age (Figure 2H).³⁵ In addition, OLGs from the C2 cluster also highly expressed those genes involving in remyelination (Figure 2H).³⁵ Brain aging reduced the ability of the OLGs from the C2 cluster to evoke remyelination (Figure 2H).

Together, these data indicated that all OLGs were impacted by brain aging. The ability of premature OLGs to generate mature OLGs was impaired in aged mice. And a fraction of aged OLGs (the C3 subset) were at the edge of apoptosis.

Aged Microglia Were Activated by OLGs to Promote Remyelination

To further investigate the influence of aged OLGs on brain microenvironment, we performed cell–cell interaction analysis among brain cells in aged mice (Supplementary Figure 3A). The results showed that aged OLGs have impacted on microglia (Supplementary Figure 3A). Our previous study has identified six distinct clusters of microglia in both young and aged mice.¹⁸ Those microglia from the cluster 3 was almost only existed in aged mice, which showed a phenotype of inflammation-active state (Supplementary Figure 3B). Cell-cell interaction analysis between subsets of aged OLGs and subclusters of aged microglia revealed a strong interaction between mature OLGs (the C2 subset) and activated microglia (the C3 subset) (Supplementary Figure 3C). These data indicated that aging OLGs might activate microglia in aged brain.

To elucidate the impact of activated microglia on all the OLGs in aged brain, we performed cell–cell interaction analysis. The results indicated that activated microglia released a bundle of cytokines to influence the C2 subset of OLGs (Figure 3A and B, Supplementary Table 4). These cytokines included Spp1, Lgals3, Igf1 and Lpl, which were found to have a positive correlation with the myelin regenerative capacity.^37–39 In addition, Safaiyan et al identified white matter-associated microglia (WAM) which were engaged in clearing degenerated myelin to promote remyelination.²⁹ Similar transcriptomic profiles between our activated microglia (the C3 subset) and WAM were observed (Figure 3C), which up-regulated genes involving in phagosome and lipid metabolism, and down-regulated microglia homeostatic genes (Figure 3D). Several gene pairs including Spp1~Itgb5, Spp1~Itgav, Spp1~Gsn, and Lpl~Sorl1 were involved in communications between mature OLGs (the C2 subset) and activated microglia (the C3 subset) (Figure 3E and Supplementary Figure 3D). All these data indicated that the activated microglia have positive effect on mature OLGs in aged brain.

Figure 3 Cell–cell interactions between microglia and OLGs in aged CNS. (A) Dot plot depicting the interactions between different subclusters in microglia and OLGs regarding GO term immune response (GO:0006955), using DEGs for each subcluster in aged brain. Node color and size represented the strength and the count of interaction pairs, respectively. (B) Volcano plot showing the cytokines acting on OLG C2 subset released by microglia C3 subset. DEGs compared with other microglia subclusters (|log2 fold change| > 0.5 and Bonferroni adjusted P-value < 0.05) were colored, red for upregulation and blue for downregulation. (C) Heatmap of the average expression of white matter associated microglia (WAM) signature genes (Geneset 1–4) for each microglia subcluster, along with the populations mentioned by Safaiyan et al. Values were normalized per gene. (D) Expression profiles of phagosome and lipid metabolism related genes and homeostatic genes in microglia were shown with Uniform manifold approximation and projection (UMAP) visualization approach. (E) Dot plot depicting gene pairs between microglia and OLGs, along with all their co-occurring cell–cell interactions. DEGs for each subcluster in aged were used. Dot size represented the strength of gene pairs using the sum of log2 fold change of related genes and dot color meant the aggregated confidence for gene pairs. Abbreviations: CNS, central nervous system; MG, microglia; OLGs, oligodendrocytes; DEGs, differentially expressed genes; GO, gene oncology.

Four Subsets of OPCs Showed Dramatic Alterations in Aged Mice

To explore the reserve capacity of remyelination, we performed a deep analysis of OPCs in both young and aged mice (Figure 4A and Supplementary Figure 4A). A total of 2,221 OPCs were divided into 4 distinct cell types (Supplementary Figure 4A) based on the enriched genes for each cluster (Supplementary Figure 4B). Specific markers involving in various biological progresses were used to define the function of OPC subsets (Supplementary Figure 4C). Mitosis related genes Mki67 and Top2a were highly expressed in the C0 subset (Supplementary Figure 4C). No unique gene markers were expressed in the C1 subset, which only highly expressed those well-known markers for OPCs (Cspg4 and Pdgfra) (Supplementary Figure 4C). OPCs from the C2 subset was marked by Bmp4 and Tns3 involving in cell migration (Supplementary Figure 4C). Several myelin related genes, Mag, Cldn11, Mal, and Plp1, were highly expressed in the OPCs from the C3 subset (Supplementary Figure 4C). These genes normally were expressed in newly formed oligodendrocytes or myelinating oligodendrocytes. Heatmap of Pearson correlation revealed OPC C3 had a close relationship to OLGs (Supplementary Figure 4D). We thus considered OPC C3 as OPC/OLG for an intermediate transition state between OPCs and OLGs. Furthermore, all biological functions from the GO enrichment of DEGs between the C1 subset and the other subsets were down-regulated, which indicated that the OPCs from the C1 subset were in a resting state (Supplementary Figure 4E).

Figure 4 Subgroup division and age-related transcriptomic alteration of OPCs (A) T-distributed stochastic neighbor embedding (t-SNE) projection of four OPC subclusters from young and aged mice. (B) Dot plot showing the proportion of each subcluster in young and aged OPCs (**P-value < 0.01, *P-value < 0.05 by Wilcoxon Rank Sum test, n = 8 samples per condition). (C) Cellular trajectory of all OPC subclusters, generated by Monocle DDRTree dimensionality reduction algorithm. (D) Heatmap showing part of significant (P-value < 0.05) GO enrichment terms using DEGs between the C0, C2, C3 and C1 subset of OPCs. Numbers in legend represented Z-score, positive values indicated upregulation. (E) Heatmap of Pearson correlation between OPCs subclusters with the related populations mentioned by Marques et al (F) Heatmaps depicting the average gene expression on mitosis, cholesterol synthesis, morphogenesis and apoptosis regulation module in four subclusters from young and aged OPCs, respectively. (G) Plot showing part of GO terms predicted to be inhibited significantly (Z-score ≤ −2, P-value < 0.05) along with aging in OPC C1. These terms were divided into two functional clusters (ribosome and mitochondrial respirasome) according to the biologic functions they were involved in. (H) Representative confocal microscopic images in young and old cortex staining for PDGFRA and TUNEL. Scale bar: 50μm. Abbreviations: OPCs, oligodendrocyte precursor cells; COP, committed OPC; NFOL, newly formed oligodendrocyte; DEGs, differentially expressed genes; GO, gene oncology.

Three out of four subsets of OPCs showed significant alterations in cell numbers (Figure 4B). Specifically, the cell number of the subset C0 and C2 was decreased, while the cell number of the subset C1 was increased in aged brain (Figure 4B). Pseudo-time analysis indicated an underlying differentiation process from the C0 to C1 subset, and finally formed the C2 subset (Figure 4C and Supplementary Figure 4F). The major GO enrichment term of the C0 subset was vigorous mitosis, indicating that the proliferation of OPCs was intact in aged brain (Figure 4D). The C2 subset was those OPCs initiating differentiation which highly expressed those genes involving in oligodendrocyte differentiation and ensheathment of neurons (Figure 4D). Furthermore, Pearson correlation analysis revealed that the C2 subset was close to committed OPCs and newly formed OLGs (Figure 4E). Together, these data indicated that brain aging prevented the OPCs differentiation to newly OLGs.

To further examine the effect of brain aging on each subcluster of OPCs, we calculated scores on several gene sets of mitosis, morphogenesis, cholesterol synthesis, and apoptosis (Figure 4F and Supplementary Table 3). The C0 subset of aged OPCs showed a significant decrease in mitosis and cholesterol synthesis (Figure 4F). Cell-cycle phase analysis also revealed a significant reduction of OPCs in G2/M phase in the elderly (Supplementary Figure 4G). These results indicated that brain aging reduced OPC proliferation capacity. The C2 subset of aged OPCs reduced the function of cholesterol synthesis (Figure 4F). The majority of OPCs (the C1 subset) have no significant changes in mitosis, morphogenesis, cholesterol synthesis, and apoptosis. We further performed GO enrichment analysis of DEGs (adjusted P-value <0.05) of OPCs from the C1 subset between aged and young mice (Figure 4G). Functions regarding mitochondrial respirasome, peptide biosynthetic, and ribosome were evidently suppressed in aged OPCs from the C1 subset (Figure 4G). These functions were reported to be important in OPCs differentiation.⁴⁰ All the four subsets of aged OPCs did not show obvious alterations of apoptosis-related genes (Figure 4F). Immunofluorescence staining further suggested that there were no more apoptotic OPCs in the aged brain compared with young ones (Figure 4H). Together, these findings suggest that insufficient remyelination in the elderly may result from impeded differentiation of OPCs due to dysfunctional energy metabolism.

Aged Microglia Might Impede the Differentiation of OPCs

We have revealed that aged OLGs evoke the activation of microglia in aged brain (Supplementary Figure 3A). To evaluate the effect of activated microglia on OPCs, we performed cell-cell interaction between subsets of microglia and subclusters of OPCs in aged brain (Supplementary Figure 5A). We observed a strong interaction between the activated microglia (the C3 subset) and the main subcluster of OPCs (the C1 subset) (Supplementary Figure 5A). To further screen potential candidate factors involving in this interaction, we calculated all the significant changed receptors of the OPCs from the C1 subset compared with the other subclusters of OPCs (Figure 5A). A total of 19 significant up-regulated receptors were identified for downstream analysis (Figure 5A). And only 12 out of these receptor genes were up-regulated in aged OPCs compared to those in young mice (Figure 5B). Furthermore, we observed that only Epha5, Lrp1, Dner, Ramp1, Gpr37l1, and Pdgfra were highly expressed in the OPCs from the C1 subset, but not in the OLGs (Figure 5C and Supplementary Figure 5B). Subsequently, GO enrichment analysis indicated that Lrp1, Pdgfra, Epha5, Gpr37l1, and Dner were broadly involved in OPCs differentiation related functional terms including glial cell differentiation, regulation of cell projection organization, and actin cytoskeleton organization (Figure 5D). Previous studies have reported that Lrp1, Pdgfra, and Gpr37l1 have negative effects on OPC differentiation.^41–43 Finally, we explored potential ligands of Lrp1, Pdgfra, Epha5, Gpr37l1, and Dner encoding receptors, which found that cytokines such as App, Apoe, Cxcl16, Ccl4, and Ccl9 might be candidate factors releasing from the activated microglia to impact on OPCs (Figure 5E and Supplementary Figure 5C and D). The interactions among App, Apoe and Lrp1 attracted our attention. Remarkable mitochondrial dysfunction in aged C1 subset from OPCs is one of the prominent features of Alzheimer’s Disease patients, in which App and Apoe are the most frequently concerned risk genes.⁴⁴ Moreover, App and its product Amyloid beta (Aβ) carried into the cells by Apoe have been reported to cause toxicity and disrupt the electron transport chain.^45–47 Together, these data provide evidence that the abnormal microenvironment caused by age-activated microglia might impede the differentiation potential of OPCs in homeostasis.

Figure 5 Cell–cell interactions between aged microglia and OPCs. (A) Dot plot depicting differentially expressed receptor genes in OPC C1 compared with other subclusters. Receptors (|log2 fold change| > 0.1 and Bonferroni adjusted P-value < 0.05) were colored, red for upregulation and green for downregulation. (B) Heatmap showing the aging-related expression changes of those receptors marked by red in panel (A). (C) Violin plot depicting expression of receptors which were upregulated along with aging in panel (B) and differentially expressed in OPC C1 subset compared with other OPC subsets and OLGs. (D) Heatmap of part of gene ontology (GO) terms related to myelination in which the receptors in panel (C) were involved. (E) Network plot showing the interactions between extracellular ligands from microglia C3 subset and receptors in panel (C). Red for receptors and blue for ligands. Abbreviations: OPCs, oligodendrocyte precursor cells; OLGs, oligodendrocytes.

Regulation of Stat1 Might Reverse Aged Microglia to Promote Remyelination

To clarify the overall effect of activated microglia on OPCs, we performed Gene Set Enrichment Analysis (GSEA) to fetch the enrichment scores of GO term related to cell differentiation (Figure 6A). The results showed that oligodendrocyte differentiation was evidently suppressed by the C3 subset of microglia, indicating a more negative other than pro-remyelination role of this subset. To explore the regulation of microglia to promote remyelination, according to the findings of previous studies, we constructed two categories of genes. Spp1, Lgals3 and Lpl highly expressed in the activated microglia in aged mice were included in the myelin promotion group (Figure 6B).^37,43 Apoe, App, Pdgfa, Pdgfb and Sdc4 were included in the myelin impediment group (Figure 6B) due to their potential hindrance in OPC differentiation.^48,49 Beneficial cytokines Vegfb and Igf1 did not show significantly high expression in activated microglia (Figure 6B). Subsequently, we attempted to identify several regulatory factors that can upregulate those factors in the myelin promotion group and down-regulate those factors in the myelin impediment group. The TRRUST database was used to construct the network between regulatory factors and cytokines (Figure 6C). Jun, Stat1, Ep300, Runx2 exhibited extensive regulation weight in the regulatory network (Figure 6C). Among those regulatory factors, Jun formed dimeric transcription factor AP-1 was involved in Spp1 overexpression and could induce inhibition of Apoe expression in macrophages.^50,51 Runx2 was previously confirmed as critical transcriptional up-regulators of Spp1 and Lgal3 but Stat1 functions as a cytoplasmic attenuator of Runx2.^52–54 Beyond that, macrophage-specific up-regulation of Apoe gene expression by Stat1 was achieved previously.⁵⁵ Stat1 could also downregulate other positive cytokines such as Lpl and Igf1.^56,57 Furthermore, regulon activity scores of those regulatory factors were calculated with pySCENIC package to clarify candidate regulatory targets (Figure 6D). Both Jun and Stat1 showed active transcription activity in age-activated microglia compared to microglia in homeostasis. All these data indicated that Stat1 might be candidate regulatory factor of microglia to promote remyelination. Immunofluorescence further suggested a higher percent of phosphorylated STAT1-positive microglia in aged versus young brains (Figure 6E). The potential mechanism of Stat1 regulating microglia to promote OPC differentiation was shown in Figure 6F.

Figure 6 Microglia might impede remyelination via Stat1 pathway in aged CNS. (A) Heatmap of enrichment scores derived from GSEA showing the C3 subset of microglia suppressed OPC differentiation. Numbers in legend represented scaled enrichment scores. (B) Violin plot depicting expression of ligands identified to promote remyelination (first line) or impede OPC differentiation (second line) in aged microglia cluster. The third line showed that aged microglia did not highly express Vegfb, Igf1, Tnf or Il1b. (C) Network plot showing interactions between cytokines in panel (B) and their possible transcription factors inferred from TRRUST, yellow for cytokines and blue for transcription factors. (D) Heatmap derived from pySCENIC package showing the activity scores of different transcription factors in microglia C3 compared with the C0 subset. Values were normalized per gene. (E) Representative confocal microscopic images in young and old cortex staining for IBA1 and pSTAT1. DAPI labels cell nuclei. Scale bar: 50μm or 25 μm. Bar plots depicting the percent of pSTAT1+ microglia (**** P-value < 0.0001 by student’s t test. n = 5 samples per condition). (F) Schematic diagram depicting that Stat1 might regulate the production of different cytokines (upregulate Apoe and downregulate Spp1, Lgals3, Igf1 and Lpl) in microglia. Those cytokines acted on OPCs and OLGs alternatively, evoking OLG remyelination and impeding the differentiation process of OPCs. Abbreviations: CNS, central nervous system; GSEA, Gene Set Enrichment Analysis; ES, enrichment score; pSTAT1, phosphorylated STAT1; OPCs, oligodendrocyte precursor cells; OLGs, oligodendrocytes; MG, microglia.

Discussion

This study revealed that brain aging leads to degeneration of myelin with the features of OLG apoptosis and OPC reduction. Aged OLGs evoked the activation of microglia, which released cytokines to promote OLG survival and hinder OPC differentiation. Regulation of Stat1 might transform the phenotype of the aged microglia to promote remyelination. Our findings indicated that microglia might be a new therapeutic target for age-related remyelinating diseases.

OLGs have been considered as a functionally homogeneous population in the CNS. Here, we identified six OLG subpopulations in different functional states, although they shared the expression of many myelin-related genes, such as Mbp and Plp1. All subgroups of aged OLGs have a downward trend in myelin-related protein. Pseudo-time analysis revealed that the C2 and C5 located at the two ends of the cellular trajectory, representing OLGs in premature and terminal differentiation status, respectively. This result was verified through comparison with data from Marques et al.¹⁷ The C5 got a lower score on aging module and overrepresented cholesterol synthetic and cytoskeleton organization, indicating the formation of multiple processes to wrap the axons.³⁴ Klk6 highly expressed by the C2 subset was a typical marker of mature OLGs.¹⁷ Kallikrein 6 (Klk6) is a secreted serine protease preferentially expressed by OLGs in CNS white matter and plays a functional role in the maintenance of OLG number and the expression of myelin proteins.^58,59 It is currently thought that myelin can not only be repaired by the generation of new OLGs from OPCs but the remaining mature OLGs are also able to participate.^35,60 Overrepresented cholesterol synthetic capacity and a high degree of overlap with up-regulated genes in remyelination further implicated the role of C2 in myelination repair. Consistent with previous reports that Klk6+ oligodendrocytes were sensitive to injury,⁶¹ we observed a decline in C2 proportion in aged brain, suggesting an age-related decreased repair ability. The remaining C0, C1, C3, and C4 subset of OLGs were all similar to MOL5/6, the main constituent of OLGs mentioned by Marque et al. Although they did not show large discrepancy in cellular trajectory, we still detected some functional differences among these four groups as mentioned above.

The CNS contains an abundant population of OPCs that generate OLGs and repair myelin after demyelination. Immunofluorescence revealed that OPCs decreased in aged brain, which support previous studies that there was a decline in the number of OPCs during natural aging.^62,63 However, some researchers pointed out that OPC numbers were unaltered in aged brain.^20,62 Mice of various ages and the diverse areas used for cell counting may account for the difference in results. It was reported that the velocity of OPC aging depends on brain regions.⁶⁴ What we focused on here was OPCs in the whole brain. The areas of cortex, corpus callosum, external capsule, and striatum were all acquired for quantification, which were different from the motor cortex or somatosensory cortex concerned in previous studies.^20,62 OPCs were divided into three groups, namely OPCs in mitosis, in basic state and in differentiated stage. We found OPCs in aged underwent an obstructed differentiation, leading to reduced OPC C2 and OLG C5 and thus inadequate compensatory for the myelin loss. Furthermore, OPC proliferation was reported to be enhanced following demyelination.⁶⁵ Although obvious demyelination was observed in aged, the proportion of OPCs exhibiting a remarkable capacity for self-renewal was found to be decreased along with aging. This may result from the fact that the increase in OPCs in resting state could reject proliferation in a self-repulsion way to maintain cell density and control their uniform distribution.⁶⁶ OPC differentiation into mature myelinating OLGs depends on high metabolic and mitochondrial demand.^40,67 Glucose and drugs that improve mitochondrial function could fuel OPC survival and promotes myelination.⁶⁷ Nevertheless, oxidative phosphorylation and respiratory electron transport-related pathway were evidently suppressed in OPC C1, which was considered as a hallmark of aging and linked to restricted differentiation.⁶⁸ This indicated that OPCs in the basic pool were more susceptible to aging. Furthermore, such age-related changes are more likely caused by extrinsic factors instead of intrinsic ones since aged OPCs, which in their normal environment rarely underwent division, into the neonatal CNS, and found that they proliferated and differentiated to the same extent as the endogenous neonatal OPCs.⁶⁹ This provided evidence that an enhanced remyelination could be achieved through an improved microenvironment.

Gene expression studies have shown aged microglia develop a more inflammatory phenotype, particular in white matter.^4,18 Aged microglia we focused here exhibited close relationship between WAM identified in Safaiyan et al,²⁹ which represent a protective metabolic response required to clear the increasing amounts of myelin debris that accumulate during aging. Moreover, Lgals3 and Axl suggesting an enhanced phagocytosis, Lpl which was associated with hydrolysis of lipoproteins, Ch25h suggesting that microglia were actively suppressing their endogenous cholesterol biosynthetic enzymes were all distinctly upregulated in aged microglia (Supplementary Table 4), providing evidence that they were dealing with lipid-rich and difficult-to-digest myelin debris. Further, Igf1 was reported to make trophic support for the remyelinating OLGs and Lgals3 was considered to drive OLGs differentiation to control myelin integrity and function.^13,38 These results elucidated the side of age-related microglia trying to promote myelin repair. Cell-cell analysis revealed that Spp1 secreted by microglia played an important role in this process. Osteopontin (Spp1) is a kind of secreted glycoprotein with chemotactic, cytokine-like and anti-apoptotic properties. It contains an Arg-Gly-Asp (RGD) cell adhesion motif that mediates interactions with several integrins including Itgb5 and Itgav. Osteopontin upregulated after peripheral nerve injury has been showed to promote proliferation and inhibit apoptosis of Schwann cells.⁷⁰ Increased myelin basic protein synthesis and enhanced myelin formation were also observed in vitro after the Spp1 treatment.⁷¹ The interaction between aged microglia and OLG C2 (such as Spp1~Itgb5, Spp1~Itgav) may provide a therapeutic target for promoting remyelination in aged. As the main lipid transporter in the brain, Apoe is synthesized and secreted primarily by microglia and some other glia cells.⁷² It is well described that Apoe-mediated clearance of toxic Aβ peptides can occur at different locations in the brain.^72,73 Such endocytosis was efficiently facilitated by LDL receptor-related protein 1 (LRP1) receptor which was abundantly expressed on OPCs. Consistent with our findings, disruption of OPCs was seemed as an early sign of pathology in the mouse model of Alzheimer’s disease characterized by Aβ accumulation.⁷⁴ Moreover, it was shown that hypertrophic OPCs surrounded Aβ plaques and was involved in clearance of Aβ.^74,75 With the accumulation of Aβ during aging, intracellular Aβ peptide could further potentiate mitochondrial dysfunction through localizing to mitochondrial membranes, blocking the transport of nuclear-encoded mitochondrial proteins, disrupting the electron transport chain, and ultimately inducing mitochondrial damage,⁴⁶ which was highly correlated with aging phenotypes of OPC C1. At last, we proposed that inhibiting Stat1 pathway in microglia could benefit remyelination process by up-regulating regeneration factors, down-regulating factors that hinder differentiation. Signal transducer and activator of transcription-1 (Stat1) has an important role in inflammation and was critically involved in microglia activation, representing a suitable target for the prevention and treatment of neurodegeneration.⁷⁶ It was reported that Stat1 negatively regulated the formation of spatial memory and induced the memory-impairing effect of Aβ.⁷⁷ As the main upstream signal ligand in Stat1 pathway, immune cytokine interferon-γ (IFN-γ) was proven to play a vital role in inhibition of developmental myelination or remyelination in demyelinated lesions in CNS and Stat1 deletion diminished oligodendrocyte loss.⁷⁸ Nevertheless, these studies are based on the inhibition of Stat1 pathway in all somatic cells, the influence of Stat1 pathway regulation on remyelination specifically in microglia needs further study. Although the factors such as Pdgfa and Pdgfb secreted by microglia and receptor Pdgfra, Gpr37l1 on OPCs can hinder the differentiation of OPCs,^43,48 functional changes revealed by GO analysis in aged OPCs showed no hint of their roles. However, our analysis did not exclude their possible influences and more efforts are needed for further investigation. IL-1β and TNF-α released by activated microglia, which were broadly explored previously, were shown to be able to enhance remyelination or induce cell death and demyelination. Nevertheless, both of these two cytokines were not highly expressed by aged microglia (Figure 6B). This may be because microglia consist of a spectrum rather than a polarization of activation states. Aged microglia we identified here represented a progressive activation state in aging which was different from those activated by lipopolysaccharides (LPS) or some other cytokines.

Several limitations existing in the present study should be noted. First, our findings were mainly based on the analysis of single-cell transcriptomic data, while changes in the transcriptome do not always dictate molecular alterations at the protein or functional level. Nevertheless, our study was to unmask the transcriptome alteration of the OLG lineage cells and provided source to screen for potential targets to modulate microglial functions to promote remyelination in aged subjects. Immunofluorescence STAT1 staining provided weak evidence that Stat1 might be a target to modulate microglial functions for regulating OLG regeneration. Further studies are needed to verify the function of STAT1 in aged microglia. Second, the present research only focused on the interactions between microglia and oligodendrocyte lineage cells. Endothelial cells were also shown to interact a lot with the OLG lineage cells, which need to be further studied.

Conclusion

In conclusion, these data exhibited that brain aging caused the degeneration of myelin. Activated microglia by aged OLGs may be a candidate target for promoting remyelination. Our findings will help to advance a variety of efforts towards understanding and modulating the aging process and exploring molecular and cellular therapeutic targets for aging-related demyelination diseases.