Divergent transcriptional regulation of astrocyte reactivity across disorders

Abstract

Astrocytes respond to injury and disease in the central nervous system with reactive changes that influence the outcome of the disorder^1-4. These changes include differentially expressed genes (DEGs) whose contextual diversity and regulation are poorly understood. Here we combined biological and informatic analyses, including RNA sequencing, protein detection, assay for transposase-accessible chromatin with high-throughput sequencing (ATAC-seq) and conditional gene deletion, to predict transcriptional regulators that differentially control more than 12,000 DEGs that are potentially associated with astrocyte reactivity across diverse central nervous system disorders in mice and humans. DEGs associated with astrocyte reactivity exhibited pronounced heterogeneity across disorders. Transcriptional regulators also exhibited disorder-specific differences, but a core group of 61 transcriptional regulators was identified as common across multiple disorders in both species. We show experimentally that DEG diversity is determined by combinatorial, context-specific interactions between transcriptional regulators. Notably, the same reactivity transcriptional regulators can regulate markedly different DEG cohorts in different disorders; changes in the access of transcriptional regulators to DNA-binding motifs differ markedly across disorders; and DEG changes can crucially require multiple reactivity transcriptional regulators. We show that, by modulating reactivity, transcriptional regulators can substantially alter disorder outcome, implicating them as therapeutic targets. We provide searchable resources of disorder-related reactive astrocyte DEGs and their predicted transcriptional regulators. Our findings show that transcriptional changes associated with astrocyte reactivity are highly heterogeneous and are customized from vast numbers of potential DEGs through context-specific combinatorial transcriptional-regulator interactions.

Extended Data Fig. 1.. RNA sequencing, Transcriptional Regulator Enrichment Analysis (TREA) and previously published TRs.

a. Flow diagram of RNAseq procedure. To minimize technical differences, we used mice of similar age and genetic background for all experimental disorder models and examined the same anatomical region (thoracic spinal cord). Spinal cord tissue from different experiments were frozen until all experiments were completed. All tissue was then processed at the same time to limit the potential for technical variations. The RiboTag procedure was used to harvest ribosome-associated RNA transcripts specifically from reactive astrocytes. RiboTag hemagglutinin (HA) was transgenically-targeted specifically to astrocytes by using well-characterized mGFAP-Cre. RNA-sequencing and analysis were conducted under identical conditions. b. Specificity of RiboTag-HA targeting to astrocytes and not microglia or other cells. Two sets of orthogonal (3D) scans from uninjured spinal cord or after spinal cord injury (SCI) with multichannel immunofluorescence for HA targeted to astrocytes plus Sox9 and Iba1 as markers of astrocytes or microglia respectively. The same areas are shown with different fluorescent wavelength filters and different orthogonal slices that demonstrate 3D staining associated with astrocytes or microglia. HA is robustly present in Sox9-positive astrocytes but is not detectable in Iba1-postive microglia in either uninjured cord or after SCI. Absence of HA-targeting to neurons or oligodendrocytes has been demonstrated previously. These immunohistochemical comparisons were repeated independently three times with similar results. c. Venn diagrams show that the relative proportions of shared astrocyte reactivity DEGs and TRs identified in spinal cord astrocytes after EAE, SCI or LPS do not detectably differ when using thresholds of FDR<0.1 or FDR<0.05. Table shows PCA analysis of DEGs identified in spinal cord astrocytes after EAE, SCI or LPS using thresholds of FDR<0.1 or FDR<0.05. The relative locations of the three disorders in PC space when compared to non-reactive astrocytes do not detectably differ when using thresholds of FDR<0.1 or FDR<0.05 as reflected in the percent of total vector length and the angles between vectors. d. Flow diagram of Transcriptional Regulator Enrichment Analysis (TREA) procedure for TR identification by upstream analysis of DEGs in reactive astrocytes. To identify TRs of astrocyte gene expression, we applied a conservative, multi-step algorithm that draws on both computationally- and biologically-derived regulator-target gene interaction data from multiple resource databases: i) ChEA: transcription factor regulation inferred from integrating genome-wide ChIP-X experiments, ii) JASPAR and iii) TRANSFAC transcription factor DNA-binding preferences as position weight matrices, and iv) Ingenuity Pathway Analysis Upstream Regulator Analytic (IPA^®, Qiagen, Valencia, CA). Using these resource databases, TREA identifies TRs implicated in regulating DEGs by interrogating multiple forms of TR-target gene regulatory interactions that include findings from experimental studies involving techniques such as chromatin immunoprecipitation and genetic loss-of-function studies, as well as well-validated predictive computational DNA binding ‘motif’ analytics. In this manner, TR-target gene interactions considered by TREA include traditional direct TR-DNA binding mechanisms, as well as indirect forms of gene expression regulation wherein a TR may act through different types of intermediaries to effect expression of downstream target genes, including chromatin modifiers and other forms of epigenetic regulators. TRs were included if they met either of two criteria: (1) convergence across resource databases; (2) differential gene expression of the TR plus convergence with at least one resource database. Together, these databases allow for interrogation of gene expression datasets for enrichment of downstream targets for approximately 1350 TRs. Resource database output files containing statistically enriched TRs and their downstream astrocyte target gene IDs were processed for TREA using in-house python scripts available at GitHub repository (https://github.com/burdalab/TREA). Final TREA libraries containing significantly enriched TRs and associated astrocyte target genes were then generated for each condition’s DEG dataset. TREA libraries were used for all comparisons of astrocyte TRs and gene expression profiles across disorders and experimental conditions and to generate a resource database of reactive astrocyte TRs and the DEGs that they regulate across a broad spectrum of CNS disorders and conditions. This database can be accessed via an open-source website http://tr.astrocytereactivity.com and has multiple search parameters according to TR, DEG or condition. e. Published astrocyte reactivity TRs plus literature references. * 4 of 62 published TRs not predicted in EAE, LPS or SCI.

Extended Data Fig. 2.. Single nucleus ATAC sequencing.

a. Experimental models: All ATACseq experiments examining LPS or SCI treatments or healthy controls used wild-type or transgenic mice of 657Bl6 background strain. Transgenic mice expressing mGFAP-Cre were used for astrocyte-specific deletion (cKO) of Smarca4 or Stat3. b. The same region of thoracic spinal cord at T9-T10 was harvested for all evaluations. Spinal cord tissue from different experiments were frozen until all experiments were completed. All tissue was then processed at the same time to limit technical variations. c. Two spinal cords from the same experimental group were pooled to prepare suspensions of nuclei. These suspensions were enriched for astrocyte nuclei with a Sox9 antibody and magnetic beads precipitation. Two such suspensions were prepared from a total of n=4 mice per experimental condition. d. Two biological replicates, each consisting of nuclei from two mice, were used for single nucleus ATACseq. All biological replicates from all experimental conditions were sequenced at the same time to avoid batch effects. e. Box plots compare the distribution of the per-nucleus averages of unique DNA fragments per nucleus or TSS enrichment per nucleus for the seven experimental conditions listed, with each dot representing the average value of all nuclei collected from each experimental group (whiskers show range, box encompass 25–75% quartiles and the centre line indicates the median; n=8 WT healthy mice and n=4 mice for all other experimental conditions). f. UMAP clustering based on differential ATAC peaks across 145,973 high quality nuclei isolated across all conditions showed separation of nuclei into multiple distinct clusters. g. UMAP distribution of examples of specific DAGs used to identify the dominant cell types in different clusters. h. Heatmap of DAGs used to identify specific cell types. Each line represents the per-nucleus z-score for a given gene averaged across all nuclei in a cluster. i. UMAP and stacked bar graph show relative contribution of biological replicates for different experimental conditions. Because experiments for LPS and SCI were conducted at different times, groups of healthy control mice were collected for each experiment. All biological replicates showed similar contributions to their respective clusters, confirming that separation of treatment groups with essentially no overlap of LPS with either healthy or SCI was due to biological variation and not due to technical artifacts. * each biological replicate consisted of nuclei from two mice.

Extended Data Fig. 3.. Identification of TR proteins in reactive astrocytes by immunohistochemistry.

For immunohistochemical (IHC) detection of TR protein we focused on previously unpublished, TREA-predicted astrocyte reactivity TRs for which sensitive antibodies were commercially available whose specificity was supported by western blots. a. Newly identified TRs co-localized by IHC in reactive astrocytes in thoracic spinal cord (T9-T-10) after LPS treatment. TRs are show in alphabetical sequence. b. Summary of all (newly identified and previously published) astrocyte reactivity TRs after LPS treatment identified here by at least two experimental procedures, either prediction from DEGs by TREA, or prediction based on significant change in motif access determined by ATACseq or by IHC or by all three. c. Newly identified TRs co-localized by IHC in reactive astrocytes in thoracic spinal cord (T9-T-10) after SCI (continued in Extended data figures 4,5). TRs are show in alphabetical sequence. Each immunohistochemical evaluation was repeated at least three times with similar results.

Extended Data Fig. 4. . Immunohistochemistry and motif analysis of TR proteins in reactive astrocytes after SCI.

a. Newly identified TRs co-localized by IHC in reactive astrocytes in thoracic spinal cord (T9-T-10) after SCI (see also Extended data figures 3,5). Each immunohistochemical evaluation was repeated at least three times with similar results. b. Summary and heatmap of all (newly identified and previously published) astrocyte reactivity TRs after SCI that have detectable DNA-binding motifs and were identified here as exhibiting significant change in motif access determined by ATACseq.

Extended Data Fig. 5.. Immunohistochemistry of TR proteins in reactive astrocytes after SCI.

a. Newly identified TRs co-localized by IHC in reactive astrocytes in thoracic spinal cord (T9-T-10) after SCI (see also Extended data figures 3,4). b. Summary of all (newly identified and previously published) astrocyte reactivity TRs after LPS treatment identified here by at least two experimental procedures, either prediction from DEGs by TREA, or prediction based on significant change in motif access determined by ATACseq or by IHC or by all three. c. Nine examples of previously published TRs co-localized by IHC in reactive astrocytes in thoracic spinal cord (T9-T-10) after SCI show that these exhibit similar staining patterns to TRs newly identified here. Each immunohistochemical evaluation was repeated at least three times with similar results.

Extended Data Fig. 6.. Astrocyte-specific TR deletion (cKO).

a. Experimental models: Transgenic mice expressing mGFAP-Cre were used for astrocyte-specific deletion (cKO) of Smarca4 or Stat3. Smarca4 was selected as a newly identified TREA-predicted reactivity TR that lacks DNA binding motifs and acts as a chromatin regulator via protein-protein interactions. Stat3 is a well-established astrocyte reactivity TR that acts via DNA-binding motifs^, and that was predicted by both TREA and ATACseq motif analysis as a reactivity TR in both LPS and SCI. b. Heatmaps show that Smarca4-astro-cKO had minimal effects on overall gene expression or on the expression of highly enriched astrocyte genes under basal conditions in untreated mice, shown as mean FPKM (fragments per kilobase of transcript sequence per million mapped fragments). c. Immunohistochemistry images and graphs (mean ± sem) of various staining parameters show that Smarca4-astro-cKO had no visibly detectable or quantifiably significant effects on the appearance or number of astrocytes, neurons or microglia in untreated mice. Unpaired t tests, ns non-significant WT versus Smarca4-cKO, n=4 mice per group. d. Effects of Smarca4-astro-cKO on RNAseq reads of Gfap and seven predicted Smarca4-regulated DEGs that are not expressed in WT untreated and are upregulated by LPS in WT but not Smarca4-astro-cKO mice. e. In situ hybridization shows predicted Smarca4-regulated gene, Cnr1, expresed in WT, but not Smarca4-astro-cKO mice after LPS. f. Multichannel immunofluorescence demonstration of Stat3 and Srebf1, or of Smarca4 and Zfp36, in the same reactive astrocytes after SCI. g. Heatmaps compare changes from healthy in differential gene expression (DEG) or differential chromatin accessibility (DAG) across the same genes after LPS or SCI in WT, Smarca4-cKO or Stat3-cKO mice, and graphs show an 86 to 96% congruence between changes in gene expression and changes in chromatin accessibility in both Smarca4- and Stat3-regulated DEGs after LPS or SCI respectively. h. Heatmaps show that Stat3-cKO significantly alters the changes from healthy normally observed after SCI in both differential gene expression (DEG) and differential chromatin accessibility (DAG) across the same 31 genes that lack Stat3-binding motifs. i. Heatmaps show that Stat3-cKO significantly alters the changes from healthy normally observed after SCI in differential gene expression (DEG) of 71 chromatin regulators. j. Multichannel immunofluorescence of Irf9 and Cxcl10 in the same reactive astrocyte after SCI. k. Graphs (mean ± sem) show effects of Smarca4-astro-cKO and Stat3-astro-cKO on various microglial histopathological responses to LPS, P values are cKO+LPS versus WT+LPS, one-way ANOVA with Bonferroni’s test, ns nonsignificant; for all graphs WT with no LPS n = 18 mice, all other conditions n = 6 mice. PCA shows composite microglia histopathology score derived from histopathological quantifications in the four graphs. l. Immunohistochemistry and graphs of mean ± sem (n ≥ 6) cell counts of the neuronal marker, NeuN, shows that high dose systemic LPS sufficient to cause pronounced microglial activation and behavioral effects did not lead to detectable neuronal loss in either spinal cord or brain. One-way ANOVA with Bonferroni’s test, ns nonsignificant; for both graphs WT with no LPS n = 18 mice, all other conditions n = 4 mice. m. Distance biplot of PCA for the effects of Smarca4-astro-cKO and Stat3-astro-cKO on composite locomotor scores after SCI. The locations of values for each individual locomotor parameter, from either the longitudinal observer scored Open Field (OF) evaluations or ladder walk testing at 28 days after SCI, indicates their contributions to defining the PC space. Graph (mean ± sem) shows composite SCI locomotor score derived from PCA of all OF and ladderwalk locomotor parameters recorded over 28 days of recovery; P values are cKO versus WT, one-way ANOVA with Bonferroni’s test, n = 11 mice per group. n. Immunofluorescence images and graph (mean ± sem) of staining intensity show effects of Smarca4-astro-cKO and Stat3-astro-cKO on CD68 and Gfap 28 days after SCI. o. Bar graph (mean ± sem) shows composite histology scores derived from PCA of Gfap, CD13, MBP and CD68 quantifications at 28 days after SCI; P values are cKO versus WT, one-way ANOVA with Bonferroni’s test, WT n = 8, Smarca4cKO n = 7, Stat3cKO n = 9 mice. Line graph shows correlation analysis of composite locomotor and histological scores. p. DEG-associated functional signaling pathway analysis shows that deletion of either Smarca4 or Stat3 alters many functions normally associated with WT astrocyte reactivity after LPS or SCI. Each immunohistochemical or in situ hybridization evaluation shown (c,e,f,j,l,n) was repeated at least three times with similar results.

Extended Data Fig. 7.. DEGs and TRs in reactive astrocytes compared with astrocytes in healthy tissue.

a,b. Venn diagrams and graph compare proportions of astrocyte reactivity DEGs that derive from genes that either are, or are not, expressed in healthy astrocytes in EAE, LPS or SCI. c-f. Heatmap, Venn digrams and graphs compare proportion of TREA-predicted reactivity TRs that derive from TRs that either are, or are not, predicted in healthy astrocytes in EAE, LPS or SCI. In all three disorders much higher proportions of TRs than DEGs were recruited from those not detectably active in healthy tissue. Notably, of the TRs shared by all three disorders 87% were derived from TRs already active in healthy states, whereas 70%, 76% and 54% of disorder unique TRs respectively in EAE, SCI or LPS, were derived from TRs not detectably active in healthy states.

Extended Data Fig. 8.. Comparison of DEGs and TRs across disorders.

8. a. Heatmap and graph show that in every disorder examined, TRs overlapped more with other disorders than did DEGs. b. Bar graphs show numbers of astroctye DEGs or TRs across all 15 disorders and conditions examined here. c. 3D distance plots of PCA for DEGs and TRs from disorders with WT astrocytes. Graphs on left show PCA of all 8 disorders, and on right show PCA of 5 disorders with three most divergent disorders removed. d. TRs predicted in three or more of six neurodegenerative disorders with genetic mutations or polymorphisms, compared with their predictions in disorders with WT astrocytes. No TRs are unique to multiple neurodegenerative disorders. e. DEGs upregulated in three or more of six neurodegenerative disorders with genetic mutations or polymorphisms, compared with their expression levels in disorders with WT astrocytes. No DEGs are unique to multiple neurodegenerative disorders.

Extended Data Fig. 9.. Schematic of a working model of astrocyte reactivity transcriptional regulation.

Divergent non-cell autonomous, disorder-selective reactivity triggers lead to context-dependent and combinatorial TR interactions in which a core set of TRs (TR_A & TR_B) is active across many if not most forms of astrocyte reactivity. These core TRs can nevertheless regulate different cohorts of DEGs in different disorders and contexts via complex and interdependent combinatorial interactions that also involve disorder- or context-selective TRs (TR_C - TR_n). These TR interactions can be influenced by astrocyte cell autonomous factors such as mutations, polymorphisms, or by differing basal starting conditions that can vary with regional or local astrocyte heterogeneity or with other factors such as exposure to previous insults. These TR interactions give rise to the exquisitely heterogenous DEG profiles associated with different astrocyte reactivity states in different disorders and different contexts.

Fig. 1.. Astrocyte reactivity DEGs and TRs vary across disorders.

a. TR prediction schematic. b. Heatmap, bar graph and violin plot comparing astrocyte reactivity DEGs that either are, or are not (H⊘), expressed in healthy astrocytes. c. Venn diagrams and PCA of DEGs or predicted TRs in SCI, EAE or LPS compared with healthy (Hlth). Red lines show vectors of TR contributions to disorders. d. Immunohistochemistry comparing astrocytes in Healthy, LPS and SCI. e. Astrocyte ATACseq schematic. f. UMAP clustering of astrocyte nuclei based on differential ATAC peaks. g. Volcano plots and hierarchical clustering of ATAC peaks in LPS or SCI relative to healthy. h. Schematic and Venn diagrams comparing DEGs and DAGs in LPS or SCI. i. Correlation of congruence between DEGs and DAGs. j. Heatmaps and ATAC peak opening or closing in two examples of up- or down-regulated DEGs. k. Heatmaps comparing chromatin accessibility of the same genes in LPS and SCI.

Fig. 2.. Astrocyte reactivity TRs identified by multiple approaches: DEG analysis, chromatin accessibility changes at DNA-binding motifs and immunohistochemistry (IHC).

a. Schematic and pie charts of TREA-identified TRs with significant changes in motif access (orange). b. UMPA clusters of representative TRs with motif access differences in healthy, LPS and SCI. c. Venn diagram of TREA-TRs with significant changes in motif access. d. Heatmaps of TREA-predicted TRs with differentially accessible motif (DAM) z-scores. All such TRs are shown for LPS or both, with 25 selected TRs for SCI. For all SCI TRs see Extended Data figure 4b. e. TREA-predicted TRs without known DNA-binding motifs confirmed by IHC of TR-protein. f,g. IHC of representative new TRs in reactive astrocytes; for individual channels of all IHC+ TRs see Extended Data figures 3–5.

Fig. 3.. TR regulation of astrocyte reactivity DEGs, chromatin accessibility and disorder outcome.

a. Experimental models. b. Smarca4 immunohistochemistry. c. Smarca4-astro-cKO or Stat3-astro-cKO effects on Smarca4- or Stat3-regulated genes. d. TREA-predicted Smarca4- or Stat3-regulated DEGs in both LPS and SCI, and astro-cKO effects. e,f. DEGs regulated in the same (e) or in opposite (f) directions by Smarca4 or Stat3. g. Stat3 and Smarca4 in same reactive astrocyte. h. TREA-predicted Smarca4- and Stat3-regulated DEGs and astro-cKO effects on those DEGs. i. TRs significantly regulated as DEGs; Smarca4 and Pitx1 in same reactive astrocyte. j. UMAP clustering of astrocyte nuclei from WT, Smarca4-astro-cKO or Stat3-astro-cKO mice based on chromatin accessibility after LPS or SCI. k. Stat3-cKO effects on DEG and DAG. l. Smarca4-astro-cKO or Stat3-astro-cKO effects on Serpine1 DEG, DAG; and ATAC peaks in Serpine1 TSS. m. Stat3-astro-cKO effects on DEG and DAG of 10 genes without Stat3 DNA-binding motifs. n. Stat3-astro-cKO effects on 10 chromatin regulators. o. Cxcl10 DEG and DAG; ATAC peaks in Cxcl10 TSS and Irf9 DNA-binding motif; Irf9 and Cxcl10 in same reactive astrocyte. p. Timp1 DEG and DAG; ATAC peaks in Timp1 TSS and DNA-binding motifs of 7 predicted TRs; Runx1 and Srebf1 in same reactive astrocyte. q. Sickness behavior after LPS. r. Microglia immunohistochemistry after LPS. s. PCA of composite histopathology after LPS. t. Open field (OF) locomotor recovery after SCI. u,v. Images and graph of mean±sem staining intensity for Gfap, CD13 and myelin basic protein (MBP). w. ETS family TRs DEG and differential motif access (DAM) after SCI. x. PU.1 (Spi1) in reactive astrocytes. y. Fcgr2b DEG and DAG; ATAC peaks in Fcgr2b TSS and PU.1 DNA-binding motif. Quantitative values in q,s,t are mean ± sem analyzed by one-way (q,s) or two-way (t) ANOVA with Bonferroni’s test, P values are cKO versus WT LPS or SCI, NS non-significant. n = mice.

Fig. 4.. Astrocyte reactivity DEGs and TRs diverge across CNS disorders in mice and humans.

a. Experimental procedures. b. Heatmaps of pairwise overlap of DEGs or TRs. c. Unsupervised clustering of disorders based on DEG and TR PCA vector lengths. d. DEG-associated functional pathways in disorders with WT astrocytes or genetic mutations and polymorphisms. e. Heatmap comparing DEGs upregulated in disorders with WT astrocytes. f. Venn diagram of TRs predicted in FAD^sn and AD^M. g. Astrocyte reactivity TRs predicted in all eight (bold) or seven of eight mouse or human disorders with WT astrocytes, plus confirmation of motif access changes or protein immunohistochemistry in LPS or SCI.