Single-nucleus sequencing reveals enriched expression of genetic risk factors in extratelencephalic neurons sensitive to degeneration in ALS

Abstract

Amyotrophic lateral sclerosis (ALS) is a neurodegenerative disorder characterized by a progressive loss of motor function linked to degenerating extratelencephalic neurons/Betz cells (ETNs). The reasons why these neurons are selectively affected remain unclear. Here, to understand the unique molecular properties that may sensitize ETNs to ALS, we performed RNA sequencing of 79,169 single nuclei from cortices of patients and controls. In both patients and unaffected individuals, we found significantly higher expression of ALS risk genes in THY1⁺ ETNs, regardless of diagnosis. In patients, this was accompanied by the induction of genes involved in protein homeostasis and stress responses that were significantly induced in a wide collection of ETNs. Examination of oligodendroglial and microglial nuclei revealed patient-specific downregulation of myelinating genes in oligodendrocytes and upregulation of an endolysosomal reactive state in microglia. Our findings suggest that selective vulnerability of extratelencephalic neurons is partly connected to their intrinsic molecular properties sensitizing them to genetics and mechanisms of degeneration.

Fig. 1. Cellular susceptibility to ALS–FTD in the human postmortem cortex.

a, Diagram of the workflow for isolation of nuclei from cortices of patients with ALS and age-matched controls followed by snRNAseq and assessment of expression of gene modules associated with neurodegenerative diseases. Exc, excitatory neurons; Inh, inhibitory neurons; Oligo, oligodendrocytes; OPCs, oligodendrocyte progenitor cells; Micro, microglia; Astro, astrocytes; Endo, endothelial cells. b–d, tSNE projections and violin plots of z-scores for expression of genes associated with ALS–FTD (b), AD (c) and MS (d) in the different cell types identified (the bars denote median for each side of the violin and the symbols indicate the average score per individual). e–g, tSNE projections and violin plots of z-scores for expression of genes associated with the ALS–FTD (e), AD (f) and MS (g) in the different subtypes of excitatory neurons (the bars denote median for each side of the violin and the symbols indicate the average score per individual) (in a–g, n = 3 control individuals and n = 5 patients with sALS).

Fig. 2. ALS–FTD susceptible neurons are layer V ETNs.

a, tSNE projection of presumptive layer V neurons. b, A dot plot representing expression of layer V markers c, A dot plot for markers of LVb ETNs of human motor cortex. d, A dot plot representing expression of von Economo markers. e, A dot plot representing expression LR-SCPN markers. f, Representative violin plot for markers of layer V ETNs of human motor cortex (geometric box plots for median and interquantile ranges). g, Visual depiction of layers identified by Maynard et al. 2021 (n = 7, g–j, publicly available) (L, layer; WM, white matter). h, A spot plot depicting expression of layer Vb motor cortex marker, STMN2, identified as enriched in THY1-Exc1, with corresponding box plot quantification. i,j, Box plots (j) and corresponding spot plots (i) for the expression of the top five ALS–FTD-associated genes expressed in Exc1 (box plots for mean and interquantile ranges). Source data

Fig. 3. ALS excitatory neurons present increased expression of stress-related pathways.

a, A schematic of DGE analysis. b, A dot plot representing scores for genes upregulated in each subgroup of Exc neurons (DGE0–6) and globally upregulated in all Exc (DGEall). c, Comparison of genes globally upregulated in ALS (DGEall) with genes upregulated in classes of L5-ETNs (genes expressed by >10% of cells, >2 FC and adjusted P value <0.05). d, Violin plots of z-scores for genes globally upregulated in all excitatory cells (DGEall) in all excitatory neurons (n = 3 controls, n = 5 patients with sALS; geometric box plots represent median and interquantile ranges and symbols indicate average score per individual). e, Violin plots of z-scores for genes upregulated in classes of L5 ETNs (DGE1, DGE5 and DGE6) in the three groups (geometric box plots represent median and interquantile ranges, symbols indicate average score per individual). f, GO analysis for genes upregulated in L5 ETNs classes (DGE1, DGE5 and DGE6); the highlighted terms are shared between the three. PD, Parkinson’s Disease; AD, Alzheimer’s Disease; HD, Huntington’s disease; CC, cellular components; KEGG, Kyoto Encyclopedia of Genes and Genomes. g–h, Western blot quantification of ubiquitin accumulation and 20S proteasome (prot.) subunit from motor cortices of separate cohort of patients with ALS (n = 6) and age-matched controls (n = 7) (t-test).

Fig. 4. Proteostatic stress in hPS cell-derived neurons resembles changes in excitatory neurons from brain of patients with ALS.

a, Diagram of neuronal differentiation from PSCs and treatment with proteasome inhibitors for bulk RNA sequencing. b, Immunofluorescence of TDP-43 localization after treatment. c, Venn diagram depicting shared upregulated genes between treated hPS cell-derived neurons and excitatory neurons from patients with ALS. d. GO analysis for shared genes in c, highlighted terms involved in protein folding and neurodegenerative diseases. CC, cellular components. Source data

Fig. 5. In ALS, oligodendroglial cells decrease their myelinating machinery in favor of a neuro-engaged state.

a, tSNE projection of OPCs and oligodendrocytes markers. b, tSNE projection of oligodendroglia (ALS n = 8,372 nuclei and control n = 11,168 nuclei). c, tSNE projection of subclusters within oligodendroglia (Wilcoxon–Mann–Whitney). d, Distribution of subclusters by diagnosis) (mean ± s.e.m. e, GO analysis for genes characteristic of control-enriched oliglia0 highlighted terms involved in myelination. CC, cellular components. f, GO analysis for genes characteristic of ALS-enriched oliglia1 highlighted terms involved in neuro-engaged functions. g, Violin plots of representative genes for neuro-supportive functions (left) and myelination (right) (geometric box plots for median and interquantile ranges; symbols indicate log₂(average expression) per individual (fraction of cell expressing). h, Volcano plot of DEGs in oligodendroglia. Highlighted genes identified in GO terms related to myelination (orange) and neuro-engaged functions (green). i, Violin plots representing z-score for selected GO terms and related tSNE projection (boxplot representing median and interquantile ranges; symbols indicate average score per individual). j,k, Western blot (j) and quantification (k) of CNPase and MBP from motor cortices of patients with ALS and age-matched controls (t-test). l. Diagram illustrates shift of oligodendrocytes states (t-test) (for a–l, n = 3 control and n = 5 patients with sALS). nUMI - normalized Unique Molecular Identifier.

Fig. 6. DAM signature in ALS.

a, tSNE projection of microglia (ALS n = 759 nuclei and control n = 693 nuclei). b,c, Volcano plot of genes upregulated in microglia from ALS. Genes identified in GO terms for endocytosis and exocytosis (b) and genes associated with neurodegenerative diseases (c). d, Violin plots of representative genes upregulated in patients with ALS associated with reactive microglia (geometric box plots represent median and interquantile ranges; symbols indicate log₂(average expression) per individual) (fraction of cell expressing). e, A dot plot representing expression of genes associated with ALS–FTD pathogenesis. f, Violin plots of representative ALS–FTD genes upregulated in ALS (geometric box plots represent median and interquantile ranges; symbols indicate log₂(average expression) per individual) (fraction of cell expressing). g, GO analysis for genes upregulated in ALS microglia, highlighted terms involved in myeloid cells biology and/or pathogenesis of ALS (WP, WikiPathways). h, Violin plots representing z-scores for selected, statistically significant GO terms from f (geometric box plots represent median and interquantile ranges; symbols indicate average score per individual). i, Comparison of genes upregulated in microglia from ALS with genes upregulated in microglia in other neurodegenerative diseases (for a–i, n = 3 control individuals and n = 5 patients with sALS). nUMI, normalized Unique Molecular Identifier.

Fig. 7. Apoptotic neurons upregulate lysosomal genes in microglia.

a, Diagram of microglia and neuronal differentiation from PSCs and induction of apoptosis neurons and feeding to iMGLs (piNs, patterned induced neurons). b, Brightfield images of untreated day 40 iMGLs and day 40 iMGLs that were fed apoptotic neurons for 24 h. c. RT–qPCR quantification of homeostatic and DAM genes after feeding (AN, apoptotic neurons). d, RT–qPCR quantification of selected ALS–FTD-associated and lysosomal genes 24 h after feeding iMGLs with apoptotic neurons (AN, apoptotic neurons). (t-test, *P < 0.05, **P < 0.01, ***P < 0.001, n = 3 biological replicates). piNS, induced neurons from hPS cells.

Extended Data Fig. 1. Technical parameters of snRNAseq and cell-type distribution across individuals.

a. Schematic diagram of workflow for isolation of nuclei from cortices of ALS patients and age-matched controls followed by single-cell RNA sequencing by DropSeq, library generation and Quality Controls for analysis with Seurat 3.0.2 b. Frozen tissue from one of the individuals. c. Staining for TDP-43 in one of the patient sample, note neuron with skein-like inclusions and faint nuclear staining (scale bar 25μm, n = 3 patients analysed). d. Quality controls post-filtering (FC – Frontal Cortex): number of total nuclei detected (barcodes), average number of genes per nucleus (nFeatures), and average number of UMIs (Unique Molecular Identifiers) per nucleus (nCounts). e. t-SNE projections of the whole cohort with expression of broad cell type markers. f. Dotplot representing percentage of cells expressing additional cell type specific markers. g. t-SNE distribution of whole cohort with annotated cell types split by diagnosis (ALS patients n = 5, age-matched Controls n = 3, n = 79,169 total nuclei). h. Quality controls post-filtering (FC – Frontal Cortex): number of total nuclei detected (barcodes), average number of genes per nucleus (nFeatures), and average number of UMIs (Unique Molecular Identifiers) per nucleus (nCounts). i. Fraction of each cell types identified in whole cohort split by diagnosis (mean ± SEM, ALS patients n = 5, age-matched Controls n = 3).

Extended Data Fig. 2. Expression of ALS-FTD associated genes in different cellular subtypes and excitatory neurons subtypes.

a. Dotplot representing expression of gene associated with the ALS-FTD spectrum in each cell type identified in the whole cortex split by diagnosis. b. t-SNE projection of excitatory neurons clusters (ALS n = 15,227 nuclei, Control n = 17,583 nuclei). c. t-SNE projection of subclusters identified in excitatory neurons represents different, biologically relevant neuronal layers (FindNeighbor(res=0.2)). d. Dotplot representing percentage of cells expressing broad markers for different cortical layers. e. Distribution of excitatory neurons in subclusters by diagnosis (mean ± SEM). f. t-SNE projection of excitatory neurons by clusters and by diagnosis. g. Gene Set Enrichment Analysis for the ALS-FTD associated genes in Exc1 excitatory neuron subtype. h. Heatmap representing expression of gene associated with the ALS-FTD spectrum in each excitatory neurons identified split by diagnosis.

Extended Data Fig. 3. L5-ETNs/CSMNs-like neurons express higher levels of ALS-FTD related genes.

a,b, Spotplot and corresponding boxplot from Maynard et al. for the expression of layer Vb Motor Cortex marker, SNCG and THY1, identified as enriched in Exc1 (boxplots for mean and interquantile ranges + -SD). c,d. Dotplot and representative Violin plots for markers of L5 ExtraTelencephalic neurons of human Motor Cortex in Schirmer et al. e,f. Dotplot and representative Violin plots for markers of L5 ExtraTelencephalic neurons of human Motor Cortex in Velmeshev et al. (geometric boxplots for median and interquantile ranges) g-i. Dotplot representing expression of Layer V markers (d), von Economo markers (e), LR-SCPN markers (f) in Schirmer et al. and Velmeshev et al. l. Violin plots and corresponding Gene Set Enrichment Analysis of z-scores for expression of ALS-FTD-associated genes in THY1-neurons identified by Schimer et al. (bars denote median). m. Violin plots and corresponding Gene Set Enrichment Analysis of z-scores for expression of ALS-FTD-associated genes in THY1-neurons identified by Velmeshev et al. (bars denote median). n. Boxplot from Maynard et al. for the expression of additional top ALS/FTD associated genes identified in Exc1 (geometric boxplots for mean and interquantile ranges + -SD) (n=multiple individuals as reported in publicly available datasets)

Extended Data Fig. 4. Classes of L5-ETNs express higher levels of stress pathways.

a. DGE-based PCA plots showing separation of individuals by diagnosis in Exc1. b. Violin plots of z-scores for genes upregulated in each class of L5-ETNs (DGE1, DGE5, DGE6) in the three groups split by diagnosis (bars denote median – symbols: average score per individual). c, Comparison of DGEs signature (>10% cells, >2-FC, adj.p < 0.05) and redoDGE signatures identified in each subgroup with random seeding of equal cell numbers per diagnosis (>10% cells, >2-FC, adj.p < 0.05, max.cells.per.ident=lowest.number.per.group). d. Violin plots of z-scores for redoDGE signatures (median and interquantile ranges – symbols: average score per individual). e. Correlation plot of Log Fold Changes in ALLminus1 comparisons. f. Venn diagram comparing genes shared between ALLminus1 analyses. g. Venn diagram comparing genes in common between ALLminus1 analysis and ExcAll. h. Violin plots of z-scores for genes globally downregulated in Excall (geometric boxplots represent median and interquantile ranges – symbols: average score per individual). i. Comparison of genes globally upregulated in ALS (DGEall) with genes upregulated in CUX2-exc0 (>10% of cells, >2-FC, adj.p < 0.05) j. Gene Ontology analysis of terms for DGE0 highlighted terms involved in synaptic biology (CC=Cellular Components). k. Gene Ontology analysis of terms for DGEall highlighted terms involved in synaptic biology (CC=Cellular Components). l,m. Selected GO terms for DGE0, DGE1 and DGEall. n. Venn Diagram for shared upregulated genes between other excitatory neurons and global signature. o,p. Venn Diagrams for shared upregulated genes between treated MS excitatory neurons from Schirmer et al. and AD excitatory neurons from Mathys et al. with changes found in this study.

Extended Data Fig. 5. Global protein-protein interaction network for genes upregulated in ALS excitatory neurons.

Color-coding derived from MCL clustering to identified closely related groups of proteins.

Extended Data Fig. 6. Proteostatic stress in hPSC-derived neurons resembles changes in excitatory neurons from brain of ALS patients.

a. Diagram of neuronal differentiation from Pluripotent Stem Cells and treatment with proteasome inhibitors for bulk RNA-sequencing. b. Quantification of proteasome inhibition (mean ± SD, n = 3 biological replicates, n = 9 technical total, (t-test, ****p < 0.001). c. Principle Component Analysis plot showing strong effect of treatments compared to untreated controls (n = 3 biological replicates). d. Protein-protein interaction network of shared genes. e. Venn Diagram depicting shared upregulated genes in treated hPSC-derived neurons (MG132), excitatory neurons from ALS patients (DGEall) and genes misregulated in human neurons after TDP-43 siRNA from Klim et al. (proportional test) f-g. Western Blot quantification of soluble and insoluble TDP-43 from motor cortex of ALS patients and age-matched controls (t-test). Source data

Extended Data Fig. 7. Oligodendrocytes polarize between myelinating and neuro-engaged states.

a,b. t-SNE projection and Violin plot of markers of Oligodendrocyte Progenitor Cells (OPCs) and mature oligodendrocytes (geometric boxplots represent median and interquantile ranges – symbols: log2(AverageExpression) per individual) (fraction of cell expressing). c,d. t-SNE projection of markers of actively myelinating oligodendrocytes and violin plots representing z-score for selected GO terms by cluster (geometric boxplots represent median and interquantile ranges – symbols: average score per individual). e,f. t-SNE projection of markers of neuronally-engaged oligodendrocytes and violin plots representing z-score for selected GO terms by cluster (geometric boxplots represent median and interquantile ranges – symbols: average score per individual). g. DGE-based PCA plot showing separation of individuals by diagnosis. h. Violin plot for z-score for genes upregulated and downregulated in ALS patients (geometric boxplots represent median and interquantile ranges – symbols: average score per individual). i. Correlation plot of Log Fold Changes in ALLminus1 comparisons. j. Venn diagram comparing genes shared between ALLminus1 analyses. k. Venn diagram comparing genes in common between ALLminus1 analysis and OligoAll. l. Violin plot of markers of actively myelinating oligodendrocytes (geometric boxplots represent median and interquantile ranges – symbols: log2(AverageExpression) per individual) (fraction of cell expressing). m. Violin plot of markers of neuronally-engaged oligodendrocytes (geometric boxplots represent median and interquantile ranges – symbols: log2(AverageExpression) per individual) (fraction of cell expressing). n,o. GO analysis for genes downregulated and upregulated in ALS oligodendrocytes, highlighted terms involved in myelination (CC=Cellular Component). p. Dotplot representing genes characteristic of maturation and development of OPCs in myelinating oligodendrocytes in each subcluster split by diagnosis.

Extended Data Fig. 8. Comparison of ALS-driven changes with study with similar signature disrupted in disease (MS).

a,b. t-SNE projection and violin plot representing z-score for highly myelinating, OPALIN⁺ oligodendrocytes in Jäkel et al (geometric boxplots represent median and interquantile ranges – symbols: average per individual). c. Comparison of genes downregulated in oligodendroglia from ALS patients with genes characteristic of highly myelinating, OPALIN⁺ subtypes identified by this study (oliglia0) and by Jäkel et al (Jäkel6), highlighted genes are shared with GO terms shown in figures. d,e. t-SNE projection and violin plot representing z-score for genes of mature, not-actively myelinating oligodendrocytes in Jäkel et al (geometric boxplots represent median and interquantile ranges – symbols: average per individual). f,g. Comparison of genes upregulated in oligodendroglia from ALS patients with genes characteristic of mature, lowly myelinating groups in this study (oliglia1 and 4) and by Jäkel et al (Jäkel1), highlighted genes are shared with GO terms shown in figures. h. Dotplot representing z-scores for the genetic signatures identified in the actively myelinating cells, the mature lowly myelinating cells and DEGs identified in this study. (for whole panel n = 3 Control individuals, n = 5 sALS patients)

Extended Data Fig. 9. Shared features between ALS-driven changes and reactive subcluster of microglia.

a. Correlation plot of Log Fold Changes in ALLminus1 comparisons. b. Venn diagram comparing genes shared between ALLminus1 analyses. c. Venn diagram comparing genes in common between ALLminus1 analysis and MicroAll. d. DGE-based PCA plot showing separation of individuals by diagnosis. e. Violin plot for z-score for genes upregulated and downregulated in ALS patients (geometric boxplots represent median and interquantile ranges – symbols: average score per individual). f. t-SNE projection of subclusters identified within microglia (Micro0 = Homeo = homeostatic, Micro1 = DAMs = Disease-associated microglia, Micro2 = Cycling cells)). g. Distribution of microglia within clusters by diagnosis. h. Dotplot representing genes identified as characteristic of Homeostatic microglia and DAMs by subcluster. i. Dotplot representing genes identified as characteristic of Homeostatic microglia and DAMs by diagnosis. j. Volcano plot of statistically significant differentially expressed genes between Control and ALS microglia (top ten upregulated and top ten downregulated genes highlighted). k. Violin plots of representative DEGs downregulated in ALS patients of genes associated with homeostatic microglia (geometric boxplots represent median and interquantile ranges – symbols: log2(AverageExpression) per individual) (fraction of cell expressing). l. Gene Ontology analysis of terms associated with genes characteristic of DAMs microglia, highlighted terms playing important role in microglial biology and/or pathogenesis of the disease. m. t-SNE projections representing z-score for selected, statistically significant GO terms.

Extended Data Fig. 10. Apoptotic neurons upregulate lysosomal genes in microglia.

a. Schematic of workflow and results from the Connectivity Map project for the genes upregulated in ALS microglia. Heatmap shows what cellular signature is most closely related to the query. b. Diagram of microglia and neuronal differentiation from Pluripotent Stem Cells, induction of apoptosis neurons and feeding to iMGLs. c. RT-qPCR quantification of cell cycle-associated genes after feeding (t-test, *p < 0.05, **p < 0.01, ***p < 0.001; n = 3 biological replicates).