Abstract
Genetic risk for autism spectrum disorder (ASD) is associated with hundreds of genes spanning a wide range of biological functions1,2,3,4,5,6. The alterations in the human brain resulting from mutations in these genes remain unclear. Furthermore, their phenotypic manifestation varies across individuals7,8. Here we used organoid models of the human cerebral cortex to identify cell-type-specific developmental abnormalities that result from haploinsufficiency in three ASD risk genes—SUV420H1 (also known as KMT5B), ARID1B and CHD8—in multiple cell lines from different donors, using single-cell RNA-sequencing (scRNA-seq) analysis of more than 745,000 cells and proteomic analysis of individual organoids, to identify phenotypic convergence. Each of the three mutations confers asynchronous development of two main cortical neuronal lineages—γ-aminobutyric-acid-releasing (GABAergic) neurons and deep-layer excitatory projection neurons—but acts through largely distinct molecular pathways. Although these phenotypes are consistent across cell lines, their expressivity is influenced by the individual genomic context, in a manner that is dependent on both the risk gene and the developmental defect. Calcium imaging in intact organoids shows that these early-stage developmental changes are followed by abnormal circuit activity. This research uncovers cell-type-specific neurodevelopmental abnormalities that are shared across ASD risk genes and are finely modulated by human genomic context, finding convergence in the neurobiological basis of how different risk genes contribute to ASD pathology.
This is a preview of subscription content, access via your institution
Access options
Access Nature and 54 other Nature Portfolio journals
Get Nature+, our best-value online-access subscription
$29.99 / 30 days
cancel any time
Subscribe to this journal
Receive 51 print issues and online access
$199.00 per year
only $3.90 per issue
Buy this article
- Purchase on Springer Link
- Instant access to full article PDF
Prices may be subject to local taxes which are calculated during checkout
Similar content being viewed by others
Data availability
Read-level data from scRNA-seq and scATAC-seq, along with proteomics data, supporting the findings of this study have been deposited in a controlled access repository at https://www.synapse.org with accession number project ID syn26346373 Count-level data and metadata have been deposited at the Single Cell Portal (https://singlecell.broadinstitute.org/single_cell/study/SCP1129). The electrophysiology materials and data are available from the corresponding authors on request. Public data used in this paper include the GRCh38 human reference genome and the EnsDb.Hsapiens.v86 annotation package.
Code availability
The code used for data analysis is available at GitHub (https://github.com/AmandaKedaigle/mutated-brain-organoids).
References
Lord, C. et al. Autism spectrum disorder. Nat. Rev. Dis. Primers 6, 5 (2020).
Rosenberg, R. E. et al. Characteristics and concordance of autism spectrum disorders among 277 twin pairs. Arch. Pediatr. Adolesc. Med. 163, 907–914 (2009).
Sanders, S. J. et al. De novo mutations revealed by whole-exome sequencing are strongly associated with autism. Nature 485, 237–241 (2012).
Ruzzo, E. K. et al. Inherited and de novo genetic risk for autism impacts shared networks. Cell 178, 850–866 (2019).
Grove, J. et al. Identification of common genetic risk variants for autism spectrum disorder. Nat. Genet. 51, 431–444 (2019).
Satterstrom, F. K. et al. Large-scale exome sequencing study implicates both developmental and functional changes in the neurobiology of autism. Cell 180, 568–584 (2020).
Cooper, D., Krawczak, M., Polychronakos, C., Tyler-Smith, C. & Kehrer-Sawatzki, H. Where genotype is not predictive of phenotype: towards an understanding of the molecular basis of reduced penetrance in human inherited disease. Hum. Genet. 132, 1077–1130 (2013).
Zlotogora, J. Penetrance and expressivity in the molecular age. Genet. Med. 5, 347–352 (2003).
Velasco, S. et al. Individual brain organoids reproducibly form cell diversity of the human cerebral cortex. Nature 570, 523–527 (2019).
de Rubeis, S. et al. Synaptic, transcriptional and chromatin genes disrupted in autism. Nature 515, 209–215 (2014).
Stessman, H. A. F. et al. Targeted sequencing identifies 91 neurodevelopmental-disorder risk genes with autism and developmental-disability biases. Nat. Genet. 49, 515–526 (2017).
Yuen, R. K. C. et al. Whole genome sequencing resource identifies 18 new candidate genes for autism spectrum disorder. Nat. Neurosci. 20, 602–611 (2017).
Sanders, S. J. et al. Insights into autism spectrum disorder genomic architecture and biology from 71 risk loci. Neuron 87, 1215–1233 (2015).
Bernier, R. et al. Disruptive CHD8 mutations define a subtype of autism early in development. Cell 158, 263–276 (2014).
Faundes, V. et al. Histone lysine methylases and demethylases in the landscape of human developmental disorders. Am. J. Hum. Genet. 102, 175–187 (2018).
Vals, M. et al. Coffin-Siris syndrome with obesity, macrocephaly, hepatomegaly and hyperinsulinism caused by a mutation in the ARID1B gene. Eur. J. Hum. Genet. 22, 1327–1329 (2014).
Lodato, S. & Arlotta, P. Generating neuronal diversity in the mammalian cerebral cortex. Annu. Rev. Cell Dev. Biol. 31, 699–720 (2015).
Greig, L. C., Woodworth, M. B., Galazo, M. J., Padmanabhan, H. & Macklis, J. D. Molecular logic of neocortical projection neuron specification, development and diversity. Nat. Rev. Neurosci. 14, 755–769 (2013).
Langfelder, P. & Horvath, S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinform. 9, 559 (2008).
Wickramasekara, R. N. & Stessman, H. A. F. Histone 4 lysine 20 methylation: a case for neurodevelopmental disease. Biology 8, 11 (2019).
Garaschuk, O., Linn, J., Eilers, J. & Konnerth, A. Large-scale oscillatory calcium waves in the immature cortex. Nat. Neurosci. 3, 452–459 (2000).
Adelsberger, H., Garaschuk, O. & Konnerth, A. Cortical calcium waves in resting newborn mice. Nat. Neurosci. 8, 988–990 (2005).
Wang, Z.-J. et al. Autism risk gene KMT5B deficiency in prefrontal cortex induces synaptic dysfunction and social deficits via alterations of DNA repair and gene transcription. Neuropsychopharmacology 46, 1617–1626 (2021).
Villa, C. E. et al. CHD8 haploinsufficiency alters the developmental trajectories of human excitatory and inhibitory neurons linking autism phenotypes with transient cellular defects. Preprint at bioRxiv https://doi.org/10.1101/2020.11.26.399469 (2020).
Wang, P. et al. CRISPR/Cas9-mediated heterozygous knockout of the autism gene CHD8 and characterization of its transcriptional networks in cerebral organoids derived from iPS cells. Mol. Autism 8, 11 (2017).
Tuncbag, N. et al. Network-based interpretation of diverse high-throughput datasets through the Omics Integrator software package. PLoS Comput. Biol. 12, e1004879 (2016).
Szklarczyk, D. et al. STRING v11: protein-protein association networks with increased coverage, supporting functional discovery in genome-wide experimental datasets. Nucleic Acids Res. 47, D607–D613 (2019).
Rubenstein, J. L. R. & Merzenich, M. M. Model of autism: increased ratio of excitation/inhibition in key neural systems. Genes Brain Behav. 2, 255–267 (2003).
Gogolla, N. et al. Common circuit defect of excitatory-inhibitory balance in mouse models of autism. J. Neurodev. Disord. 1, 172–181 (2009).
Dani, V. S. et al. Reduced cortical activity due to a shift in the balance between excitation and inhibition in a mouse model of Rett syndrome. Proc. Natl Acad. Sci. USA 102, 12560–12565 (2005).
Mariani, J. et al. FOXG1-dependent dysregulation of GABA/glutamate neuron differentiation in autism spectrum disorders. Cell 162, 375–390 (2015).
Marchetto, M. C. et al. Altered proliferation and networks in neural cells derived from idiopathic autistic individuals. Mol. Psychiatry 22, 820–835 (2017).
Adhya, D. et al. Atypical neurogenesis in induced pluripotent stem cells from autistic individuals. Biol. Psychiatry 89, 486–496 (2020).
Wade, A. A., Lim, K., Catta-Preta, R. & Nord, A. S. Common CHD8 genomic targets contrast with model-specific transcriptional impacts of CHD8 haploinsufficiency. Front. Mol. Neurosci. 11, 481 (2019).
Moffat, J. J., Smith, A. L., Jung, E. M., Ka, M. & Kim, W. Y. Neurobiology of ARID1B haploinsufficiency related to neurodevelopmental and psychiatric disorders. Mol. Psychiatry https://doi.org/10.1038/s41380-021-01060-x (2021).
Velmeshev, D. et al. Single-cell genomics identifies cell type-specific molecular changes in autism. Science 364, 685–689 (2019).
Willsey, A. J. et al. Coexpression networks implicate human midfetal deep cortical projection neurons in the pathogenesis of autism. Cell 155, 997–1007 (2013).
Bourgeron, T. From the genetic architecture to synaptic plasticity in autism spectrum disorder. Nat. Rev. Neurosci. 16, 551–563 (2015).
Chen, A. E. et al. Optimal timing of inner cell mass isolation increases the efficiency of human embryonic stem cell derivation and allows generation of sibling cell lines. Cell Stem Cell 4, 103–106 (2009).
Church, G. M. The personal genome project. Mol. Syst. Biol. 1, 2005.0030 (2005).
Velasco, S., Paulsen, B. & Arlotta, P. Highly reproducible human brain organoids recapitulate cerebral cortex cellular diversity. Protoc. Exchange https://doi.org/10.21203/rs.2.9542/v1 (2019).
Lovell-Badge, R. et al. ISSCR guidelines for stem cell research and clinical translation: the 2021 update. Stem Cell Rep. 16, 1398–1408 (2021).
Doench, J. G. et al. Rational design of highly active sgRNAs for CRISPR-Cas9–mediated gene inactivation. Nat. Biotechnol. 32, 1262–1267 (2014).
Hsu, P. D. et al. DNA targeting specificity of RNA-guided Cas9 nucleases. Nat. Biotechnol. 31, 827–832 (2013).
Mangeot, P. E. et al. Genome editing in primary cells and in vivo using viral-derived nanoblades loaded with Cas9-sgRNA ribonucleoproteins. Nat. Commun. 10, 45 (2019).
Park, Y.-G. et al. Protection of tissue physicochemical properties using polyfunctional crosslinkers. Nat. Biotechnol. 37, 73–83 (2019).
Schneider, C. A., Rasband, W. S. & Eliceiri, K. W. NIH Image to ImageJ: 25 years of image analysis. Nat. Methods 9, 671–675 (2012).
Ohgane, K. Quantification of gel bands by an Image J macro, band/peak quantification tool. protocols.io https://doi.org/10.17504/protocols.io.7vghn3w (2019).
Pachitariu, M. et al. Suite2p: beyond 10,000 neurons with standard two-photon microscopy. Preprint at bioRxiv https://doi.org/10.1101/061507 (2017).
Müller, J. et al. High-resolution CMOS MEA platform to study neurons at subcellular, cellular, and network levels. Lab Chip 15, 2767–2780 (2015).
Ostasiewicz, P., Zielinska, D. F., Mann, M. & Wisniewski, J. R. Proteome, phosphoproteome, and N-glycoproteome are quantitatively preserved in formalin-fixed paraffin embedded tissue and analyzable by high-resolution mass spectrometry. J. Proteome Res. 9, 3688–3700 (2010).
Wiśniewski, J. R. Quantitative evaluation of filter aided sample preparation (FASP) and multienzyme digestion FASP protocols. Anal. Chem. 88, 5438–5443 (2016).
Bairoch, A. & Apweiler, R. The SWISS-PROT protein sequence data bank and its supplement TrEMBL in 1999. Nucleic Acids Res. 27, 49–54 (1999).
Consortium, T. U. UniProt: a worldwide hub of protein knowledge. Nucleic Acids Res. 47, D506–D515 (2018).
Käll, L., Storey, J. D., MacCoss, M. J. & Noble, W. S. Posterior error probabilities and false discovery rates: two sides of the same coin. J. Proteome Res. 7, 40–44 (2008).
Zhang, X. et al. Proteome-wide identification of ubiquitin interactions using UbIA-MS. Nat. Protoc. 13, 530–550 (2018).
Ritchie, M. E. et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43, e47 (2015).
Subramanian, A. et al. Gene set enrichment analysis: a knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl Acad. Sci. USA 102, 15545–15550 (2005).
Tuncbag, N. et al. Simultaneous reconstruction of multiple signaling pathways via the prize-collecting steiner forest problem. J. Comput. Biol. 20, 124–136 (2013).
Akhmedov, M. et al. PCSF: an R-package for network-based interpretation of high-throughput data. PLoS Comput. Biol. 13, e1005694 (2017).
Shannon, P. et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 13, 2498–2504 (2003).
Yoon, S. et al. GScluster: network-weighted gene-set clustering analysis. BMC Genom. 20, 352 (2019).
Quadrato, G., Sherwood, J. L. & Arlotta, P. Long term culture and electrophysiological characterization of human brain organoids. Protoc. Exchange https://doi.org/10.1038/protex.2017.049 (2017).
Zheng, G. X. Y. et al. Massively parallel digital transcriptional profiling of single cells. Nat. Commun. 8, 14049 (2017).
Stuart, T. et al. Comprehensive integration of single-cell data. Cell 177, 1888–1902 (2019).
McCarthy, D. J. et al. Cardelino: computational integration of somatic clonal substructure and single-cell transcriptomes. Nat. Methods 17, 414–421 (2020).
Huang, Y., McCarthy, D. J. & Stegle, O. Vireo: Bayesian demultiplexing of pooled single-cell RNA-seq data without genotype reference. Genome Biol. 20, 273 (2019).
Korsunsky, I. et al. Fast, sensitive and accurate integration of single-cell data with Harmony. Nat. Methods 16, 1289–1296 (2019).
Pardy, C. mpmi: mixed-pair mutual information estimators (2020).
Bates, D., Mächler, M., Bolker, B. & Walker, S. Fitting linear mixed-effects models using lme4. J. Stat. Softw. 67, 1–48 (2015).
Fonseka, C. Y. et al. Mixed-effects association of single cells identifies an expanded effector CD4+ T cell subset in rheumatoid arthritis. Sci. Transl. Med. 10, eaaq0305 (2018).
Cao, J. et al. The single-cell transcriptional landscape of mammalian organogenesis. Nature 566, 496–502 (2019).
Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550 (2014).
Lun, A. T. L. & Marioni, J. C. Overcoming confounding plate effects in differential expression analyses of single-cell RNA-seq data. Biostatistics 18, 451–464 (2017).
Yu, G., Wang, L.-G., Han, Y. & He, Q.-Y. clusterProfiler: an R package for comparing biological themes among gene clusters. Omics 16, 284–287 (2012).
Nuclei Isolation from Mouse Brain Tissue for Single Cell ATAC Sequencing Rev B (10x Genomics, 2019).
Corces, M. R. et al. An improved ATAC-seq protocol reduces background and enables interrogation of frozen tissues. Nat. Methods 14, 959–962 (2017).
Stuart, T., Srivastava, A., Madad, S., Lareau, C. A. & Satija, R. Single-cell chromatin state analysis with Signac. Nat. Methods 18, 1333–1341 (2021).
Rainer, J. EnsDb.Hsapiens.v86: ensembl based annotation package (2017).
Heinz, S. et al. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities. Mol. Cell 38, 576–589 (2010).
Acknowledgements
We thank J. R. Brown (from the P.A. laboratory) for input and assistance in editing the manuscript and all of the members of the Arlotta laboratory for discussions; M. Daly, E. Robinson and B. Neale for discussions on autism genetics; T. Nguyen (from the G.Q. laboratory) for support with organoid generation; X. Jin (from the P.A. laboratory) for helping with designing and sequencing edited lines; A. Shetty (from the P.A. laboratory) for help with scRNA-seq cell type classification; N. Haywood for help with scRNA-seq experiments; D. Di Bella (from the P.A. laboratory) for help with final edits of the manuscript; S. Andreadis and S. Getz (from the P.A. laboratory) for help with editing images; F. Zhang and J. Pan for supporting the creation of the HUES66 CHD8-mutant line; members of the Talkowski laboratory for the GM08330 line; members of the Cohen laboratory for the Mito210 line; members of the Ricci laboratory for providing nanoblades for the generation of SUV420H1 edited cell lines; L. M. Daheron at the Harvard Stem Cell Facility for expanding edited lines; and B. Budnik at the Harvard Center for Mass Spectrometry for assisting with proteomics experiments. This work was supported by grants from the Stanley Center for Psychiatric Research, the Broad Institute of MIT and Harvard, the National Institutes of Health (R01-MH112940 to P.A. and J.Z.L.; P50MH094271, U01MH115727 and 1RF1MH123977 to P.A.), the Klarman Cell Observatory to J.Z.L. and A.R., and the Howard Hughes Medical Institute to A.R. A.R. was a Howard Hughes Medical Institute and a Koch Institute extramural member while conducting this study. The HUES66 CHD8-mutant line was created with support from the Simons Foundation (346073 to F. Zhang) and the National Institutes of Health (MH099448 to J. Pan).
Author information
Authors and Affiliations
Contributions
P.A., B.P., S.V., A.J.K., M.P., and G.Q. conceived the experiments. A.J.D. designed SUV420H1 and ARID1B gRNAs and generated the Mito210 and Mito294 SUV420H1-edited lines with B.P.; L.B. generated the ARID1B-edited lines. N.E.S. and X.S. designed CHD8 gRNAs and generated the HUES66 CHD8 line. A.J.D. and B.P. generated the Mito210 and Mito294 CHD8-edited lines. S.V., B.P., M.P., R.S., C.A., A.T and S.N.S. generated, cultured and characterized all of the organoids used in this study and P.A. supervised their work. X.A. performed scRNA-seq experiments with help from B.P., S.V., M.P., R.S. and G.Q. under the supervision of P.A. and J.Z.L.; A.J.K., K.K., S.K.S. and J.Z.L. performed scRNA-seq analysis and J.Z.L. and A.R. supervised their work. S.V., B.P., M.P., R.S., A.U., G.Q. and A.J.K. worked on cell type assignments and data analysis. K.T., M.P. and A.J.K. performed proteomic analysis, supervised by K.L.; S.M.Y., P.S. and A.P. performed the calcium imaging experiments and analysis, supervised by E.S.B. and P.A.; S.M.Y. and R.S. performed the MEA recordings and analysis supervised by P.A.. A.A. performed whole-organoid imaging under the supervision of K.C. P.A., B.P., S.V., A.J.K. and M.P. wrote the manuscript with contributions from all of the authors. All of the authors read and approved the final manuscript.
Corresponding authors
Ethics declarations
Competing interests
P.A. is a SAB member at Herophilus, Foresite Labs, and Rumi Scientific, a consultant for the New York Stem cell Foundation, and is a co-founder of esalius Therapeutics. A.R. is a founder and equity holder of Celsius Therapeutics, an equity holder in Immunitas Therapeutics and, until 31 August 2020, was a SAB member of Syros Pharmaceuticals, Neogene Therapeutics, Asimov and Thermo Fisher Scientific. From 1 August 2020, A.R. is an employee of Genentech. N.E.S. is an advisor to Vertex and Qiagen.
Peer review information
Nature thanks Prisca Liberali, Annie Vogel-Ciernia and Irina Voineagu for their contribution to the peer review of this work.
Additional information
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Extended data figures and tables
Extended Data Fig. 1 Cortical organoids cultured for one, three and six months generate the cellular diversity of the human cerebral cortex with high organoid-to-organoid reproducibility.
a, scRNA-seq and immunohistochemistry analysis of organoids cultured for one month (32 d.i.v.), three months (98 d.i.v.), and six months (190 d.i.v.). Left, t-SNE plots (n = 3 organoids per timepoint, co-clustered). Cells are coloured by cell-type. Right, immunohistochemistry for specific markers. Neural progenitor marker SOX2 (magenta) and postmitotic neuronal marker TBR1 (green) are shown at one month. CPN marker SATB2 (magenta) and CFuPN marker CTIP2 (green) are shown at three months. The astroglia markers S100B (magenta) and GFAP (green) are shown at six months. Below, schematic images of brain organoids in each timepoint. Scale bars are 100 μm. b, Immunohistochemistry for neuronal (MAP2), dorsal forebrain neural progenitor (EMX1, SOX2), CFuPN (CTIP2), and CPN (SATB2) markers in GM08330 organoids at one, three, and six months. Scale bars: whole organoids (leftmost column), 200 μm; others, 50 μm. c, Immunohistochemistry for cell-type specific markers in Mito210 organoids, as in b. d, Top, t-SNE plots of the scRNA-seq data from individual replicates from three organoids at one month, three organoids at three months, and three organoids at six months from the GM08330 cell line shown in b. Bottom, bar charts showing the cell-type composition of each individual organoid. On top of the bar charts, mutual information (MI) scores between cell-type proportions and organoid identities are displayed. A MI score of 0 would indicate identical cell type proportions between organoids, while a score of 1 would indicate completely divergent profiles. In previous work, MI scores for endogenous brain datasets were reported to range from 0.008 to 0.0649. e, scRNA-seq data of organoids from the Mito210 cell line at one month (35 d.i.v.), three months (92 d.i.v.), and six months (178 d.i.v.), as in d. Organoids for the one and three month timepoints are the same as the control organoids in Extended Data Fig. 4f and Extended Data Fig. 5b. f, Expression of selected marker genes used in cell-type identification. Violin plots show distribution of normalized expression in cells from GM08330 organoids at one, three and six months (n = 3 individual organoids per timepoint). g, Expression of marker genes in Mito210 organoids, as in f. Number of organoids used for each analysis can be found in the Methods under “Statistics and reproducibility”. aRG, apical radial glia; DL, deep layer; UL, upper layer; PN, projection neurons; oRG, outer radial glia; IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; GABA INP, GABAergic interneuron progenitors; GABA IN, GABAergic interneurons.
Extended Data Fig. 2 Expression of selected ASD risk genes in cortical organoids cultured for one, three, and six months.
a, t-SNE plots of 58,568 cells from nine organoids from the GM08330 cell line, shown in Extended Data Fig. 1d, after Harmony batch correction. Cells are coloured according to cell-type (left) and timepoint (right). b, Gene set expression scores for a set of 102 genes associated with ASD risk6 across cell-types, in cells from a. Scores above 0 indicate enriched expression over similar sets of randomly chosen genes. c, t-SNE plots showing normalized expression of selected ASD risk genes in cells from a. d, Average expression of 102 genes associated with ASD risk across cell-types and timepoints in the GM08330 cell line. e, t-SNE plots of nine organoids from the Mito210 cell line, shown in Extended Data Fig. 1e, after Harmony batch correction. Cells are coloured according to cell-type (left) or timepoint (right). f, Gene set scores for the set of ASD risk genes as in b, in cells from e. Scores above 0 indicate higher expression than similar modules of randomly chosen genes. g, t-SNE plots showing normalized expression of selected ASD risk genes in cells from e. h, Expression of 102 genes associated with ASD risk across cell-types and timepoints in Mito210 cell line. RG, radial glia (aRG, oRG, and oRG/Astroglia), IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; EN, Excitatory neurons (CPN, CFuPN and PN); GABA IN, GABAergic interneurons.
Extended Data Fig. 3 Generation and characterization of SUV420H1, ARID1B, and CHD8 mutant organoids.
a, Protein domain structure of SUV420H1. Arrow indicates the region (N-domain) mutated in the Mito210, PGP1 and Mito294 parental lines (bottom). b, Protein domain structure of ARID1B. Arrow indicates the region before the ARID domain mutated in the Mito210 and Mito294 parental lines (bottom). c, Protein domain structure of CHD8. Arrows indicates the helicase C-terminal (HELC) domain mutated in the HUES66, H1, GM08330, Mito294 and Mito210 lines (bottom). d–f, Western blot analysis showing presence of SUV420H1 (d), ARID1B (e) and CHD8 (f) protein expression in control lines, and its reduction in the mutant lines. Molecular weight in kDa is shown on the left of the gel. H4K20me3, a hallmark of SUV420H1 activity, and total levels of histone H4 were also detected in control and in SUV420H1+/− lines (d). ARID1B was not detectable in Mito294 even with a longer exposure of the blotted membrane (e, right). Asterisks indicate the bands used for quantification. Bottom, protein levels in control and mutant lines were quantified and normalized for housekeeping genes β-Actin or GAPDH. For gel source data, see Supplementary Fig. 1. g, Immunohistochemistry for neuronal (MAP2), dorsal forebrain neural progenitor (EMX1, SOX2) and CFuPN (CTIP2) markers in organoids at 35 d.i.v. derived from the Mito210 SUV420H1+/−, Mito210 ARID1B+/− and HUES66 CHD8+/− and isogenic control cell lines. Scale bar, 300 μm. h–j, Size quantification of control and SUV420H1+/− (h), ARID1B1+/− (i) and CHD8+/− (j) organoids across lines and at different timepoints. The ratio of organoid size compared to the average of control organoids in each batch is plotted. Differentiation batch (b.) is indicated by the shade of the points. Lower and upper hinges correspond to the first and third quartiles (the 25th and 75th percentiles) and middle hinge is the median (50th). Both whiskers extends from the hinge to the largest or smallest value no further than 1.5 * IQR from the hinge (where IQR is the inter-quartile range, or distance between the first and third quartiles). P-values from a two-sided t-test, after Bonferroni adjustment within each mutation. Number of organoids and differentiations used for the measurement are summarized in Supplementary Table 2 and in the Methods under “Statistics and reproducibility”.
Extended Data Fig. 4 Cell-type composition of SUV420H1+/− and isogenic control organoids.
a, Immunohistochemistry of Mito210 SUV420H1+/− and control organoids cultured for one month (35 d.i.v.). Optical section from the middle of whole-organoid dataset. Scale bars are 500 μm. SOX2, a marker of neuronal progenitors, is shown in red, and nuclei (Syto16) are shown in blue. b, Immunohistochemistry for the postmitotic excitatory neuronal marker TBR1 and GABAergic marker DLX2 in Mito294 control and SUV420H1+/− organoids at one month (35 d.i.v.). Scale bars: 200 μm. c–e, scRNA-seq data from one month (Mito294 35d.i.v. (c), PGP1 35 d.i.v. (d) and Mito210 28 d.i.v., batch I (e)) control and SUV420H1+/− organoids. Bar charts show the percentage of cells for all the cell populations in each control and mutant organoid. Adjusted p-values for a difference in cell type proportions between control and mutant, based on logistic mixed models (see Methods) are shown. f, scRNA-seq data from Mito210 35 d.i.v. (batch II) control and SUV420H1+/− organoids. Left top shows combined t-SNE plots of control and mutant organoids (n = 3 single organoids per genotype). Cells are coloured by cell-type, and the total number of cells per plot is indicated. Left bottom, t-SNE plots for control and mutant individual organoids. Immature deep-layer projection neuron populations are highlighted in colour. Right, bar charts show the percentage of cells for all the cell populations in each control and mutant organoid, as in c–e. g, Enriched GO terms for genes upregulated and downregulated in SUV420H1+/− vs. control across lines. Genes were calculated using cells from the partition of interest. The top 5 most significant terms per dataset are shown. Size of dot indicates the proportion of genes belonging to each term found in the list of dysregulated genes (“GeneRatio”). Colour indicates enrichment adjusted p-value. Numbers in parentheses along the y axis indicate the number of DEGs in that dataset. As control, we calculated GO term enrichment for 100 random gene sets of the same size sampled from genes expressed in each dataset, and found no significant enrichment of these terms (see Methods). Number of organoids used for each analysis can be found in the Methods under “Statistics and reproducibility”. aRG, apical radial glia; DL, deep layer; UL, upper layer; PN, projection neurons; oRG, outer radial glia; IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; GABA N, GABAergic neurons.
Extended Data Fig. 5 Cell-type composition, full pseudotime trajectories, and gene modules in SUV420H1+/− and isogenic control organoids.
a–c, scRNA-seq data from three month Mito294 89 d.i.v. (a), Mito210 92 d.i.v. batch I (b), and 90 d.i.v. batch II (c) control and SUV420H1+/− organoids. Left top shows combined t-SNE plots of control and mutant organoids (n = 3 single organoids per genotype). Cells are coloured by cell type, and the number of cells per plot is indicated. Left bottom, t-SNE plots for control and mutant individual organoids. Cell-types of interest are highlighted in colour. Right, bar charts show the percentage of cells for all the cell populations in each control and mutant organoid. Adjusted p-values for a difference in cell type proportions between control and mutant, based on logistic mixed models (see Methods) are shown. d, Pseudotime trajectory from the full dataset of Mito210 SUV420H1+/− 35 d.i.v. (batch II) and control organoids, calculated with Monocle3. The partition highlighted by a box was subsetted and the trajectory is shown in Fig. 1d. e, Module scores (top) and their distribution across mutant and control cells (bottom) for all modules resulting from WGCNA analysis of the partition of interest from Mito210 SUV420H1+/− and control organoids at 35 d.i.v. (batch II). Cells were downsampled to have an equal number of cells per organoid. Names were assigned to each module based on the known functions of the genes included in each one. Horizontal bars show median scores, and dots show average score per organoid. Adjusted p-values show differences between control and mutant based on linear mixed models (see Methods). aRG, apical radial glia; DL, deep layer; UL, upper layer; PN, projection neurons; CP/CH, Choroid Plexus/Cortical Hem; oRG, outer radial glia; IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; GABA INP, GABAergic interneuron progenitors; GABA IN, GABAergic interneurons; GABA N, GABAergic neurons.
Extended Data Fig. 6 scATAC-seq analysis in SUV420H1+/− and isogenic control organoids.
a, UMAPs of scATAC-seq data in Mito210 SUV420H1+/− and control organoids at one month (31 d.i.v., upper left) and three months (93 d.i.v., upper right), and coembedded UMAPs with scRNA-seq in Mito210 SUV420H1+/− and control organoids at one month (28 d.i.v., middle bottom left) and three months (90 d.i.v., middle bottom right). Number of nuclei per plot is indicated. b, Enriched GO terms for the nearest genes to regions with increased and decreased accessibility in SUV420H1+/− compared to control organoids. c, Genome tracks showing read coverage for representative regions with increased accessibility between SUV420H1+/− and control organoids. Scales for the y axes (normalized counts) are displayed on the top right. d, Genome tracks showing read coverage for representative regions with increased accessibility between Mito210 SUV420H1+/− and control organoids, split by cell-type. Scales for the y axes (normalized counts) are displayed on the top right. e, Top 5 de novo motifs enriched in the regions with altered accessibility in Mito210 SUV420H1+/− compared to control organoids at one month (31 d.i.v.) and three months (93 d.i.v.), as calculated with HOMER (see Methods). Regions that showed increased accessibility in mutant compared to control organoids are in green and purple, while those with decreased accessibility are in red and blue. Transcription factors with known binding sites resembling the discovered motifs are shown.
Extended Data Fig. 7 Neuronal spontaneous activity in SUV420H1+/− and isogenic control organoids.
a, Left, Representative image of a PGP1 SUV420H1 organoid infected with SomaGCaMP6f2. Right, ΔF/F signal at the peak of a network burst. Scale bar: 100 μm. b, Heat map of calcium signal from individual cells (rows), showing the effect of 2 μM TTX. c, Top, representative trace of spontaneous calcium signal (corresponding to cell #3 in Fig. 1g). Bottom, high magnification traces of calcium transients, displaying the difference in amplitude between the isolated event and the network burst (top), and normalized traces (bottom) showing their kinetics and the multiple peaks of the burst signal. d, MEA recordings. Top, Spatial configuration of recording electrodes. Middle, example raw traces for the numbered electrodes shown at the top, and the effect upon 2 μM TTX application. Yellow columns indicate the network bursts. Right, individual (grey) and average (colour) spike waveforms for each electrode. High magnification of the trace #4 showing the individual spikes (asterisk) during a burst event. Bottom, average spike waveform (right) in a unit of electrodes (left), extracted at the time points determined by the spikes in electrode #4. e, f, Synchronous network activity in human brain organoids. Heat map of cross-correlation index (e) and cross-correlogram against a reference signal (cell #135) for a representative recording. As a control, the signal of 10 cells were circularly shifted by a random phase and the cross-correlation was then calculated. In f, the average value was plotted, and the synchronous activations as well as the periodic bursting can be seen (“All cells” in red). g, Effect of NBQX on neuronal activity. Representative traces for individual cells were normalized (3 traces for SUV420H1+/− are superimposed) and post-NBQX residual/isolated calcium transients are indicated by asterisks. h, Effect of NBQX on calcium signal. Heat map of ΔF/F signal for 15 representative cells in control (top) and SUV420H1+/− (bottom) organoids. i, j, Inter-spike interval (ISI) analysis for the network bursting. Relative frequency (top) and cumulative frequency distribution (bottom) of ISI for control and SUV420H1+/− organoids. In j, two-sided Kolmogorov-Smirnov test (n = 5 organoids per genotype). Number of organoids used for each analysis can be found in the Methods under “Statistics and reproducibility”.
Extended Data Fig. 8 Cell-type composition, full pseudotime trajectories, and gene modules of ARID1B+/− and isogenic control organoids.
a, Immunohistochemistry for the postmitotic excitatory neuronal marker TBR1 (magenta) and GABAergic marker DLX2 (green) in Mito210 control and ARID1B+/− organoids at one month (35 d.i.v.). Scale bars: 200 μm. b, c, scRNA-seq data from Mito210 one month (35 d.i.v. batch I in b, batch II in c) control and ARID1B+/− organoids. Bar charts show the percentage of cells for all the cell populations in each control and mutant organoid. Adjusted p-values for a difference in cell-type proportions between control and mutant, based on logistic mixed models (see Methods) are shown. d, Immunohistochemistry for TBR1 (magenta) and DLX2 (green) in Mito210 control and ARID1B+/− organoids at three months (90 d.i.v.). Scale bars: 100 μm. e, scRNA-seq data from Mito210 three months (90 d.i.v.) control and ARID1B+/− organoids. Left top shows combined t-SNE plots of control and mutant organoids (n = 3 single organoids per genotype). Cells are coloured by cell-type, and the number of cells per plot is indicated. Left bottom, t-SNE plots for control and mutant individual organoids. GABAergic interneurons are highlighted in colour. Left, bar charts show the percentage of cells for all the cell populations in each control and mutant organoid, as shown in b, c. Two out of three mutant organoids show an increase in GABAergic interneurons, vs. zero out of three controls (adjusted p = 0.19, logistic mixed models). f, Immunohistochemistry for the postmitotic excitatory neuronal marker TBR1 (magenta) and GABAergic marker DLX2 (green) in Mito210 control and ARID1B+/− organoids at three months (90 d.i.v.). Three out of four mutant organoids contain DLX2-positive cells, while no DLX2 expression is detected in the four controls. Scale bars: 500 μm. g, scRNA-seq data from Mito294 one month (35 d.i.v.) ARID1B+/− and control organoids. Left top shows combined t-SNE plots of control and mutant organoids (n = 3 single organoids per genotype). Cells are coloured by cell type, and the number of cells per plot is indicated. Left bottom, t-SNE plots for control and mutant individual organoids. GABAergic neurons, newborn deep-layer projection neurons and immature deep-layer projection neuron populations are highlighted in colour. Right, bar charts show the percentage of cells for all the cell populations in each control and mutant organoid, as in b, c, e. h, Pseudotime trajectories from the full dataset of Mito210 ARID1B+/− 35 d.i.v. batch II and control organoids, calculated with Monocle3. The partition highlighted by a box was subsetted and the trajectory is shown in Fig. 2c. i, Module scores (top) and their distribution across mutant and control cells (bottom) for all modules resulting from WGCNA analysis of the partition of interest from Mito210 ARID1B1+/− and control organoids at 35 d.i.v. Cells were downsampled to have an equal number of cells per organoid. Names were assigned to each module based on the known functions of the genes included in each one. Horizontal bars show median scores, and dots show average score per organoid. Adjusted p-values show differences between control and mutant based on linear mixed models (see Methods). Number of organoids used for each analysis can be found in the Methods under “Statistics and reproducibility”. aRG, apical radial glia; DL, deep layer; UL, upper layer; PN, projection neurons; CP/CH, Choroid Plexus/Cortical Hem; oRG, outer radial glia; IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; GABA NP, GABAergic neuron progenitors; GABA N, GABAergic neurons; GABA INP; GABAergic interneuron progenitors; GABA IN, GABAergic interneurons.
Extended Data Fig. 9 Cell-type composition, immunohistochemistry, and full pseudotime trajectories and gene modules of CHD8+/− and isogenic control HUES66 organoids.
a, b, scRNA-seq data from HUES66 3.5-month (109 d.i.v. (a), batch I and 107 d.i.v. (b). batch II) CHD8+/− and control organoids. Bar charts show the percentage of cells for all the cell populations in each control and mutant organoid. Adjusted p-values for a difference in cell-type proportions between control and mutant, based on logistic mixed models (see Methods) are shown. c, Immunohistochemistry for the postmitotic excitatory neuronal marker TBR1 (magenta) and GABAergic marker DLX2 (green) in HUES66 control and CHD8+/− organoids at 3.5 months (107 d.i.v.). Scale bars: 100 μm. d, Immunohistochemistry for neuronal (MAP2), dorsal forebrain neural progenitor (EMX1, SOX2) and CFuPN (CTIP2) markers in HUES66 CHD8+/− and control organoids at 3.5 months (107 d.i.v., top), and six months (190 d.i.v., bottom). Scale bars: whole organoids, 500 μm; others, 100 μm. e, scRNA-seq data from HUES66 CHD8+/− and control organoids at six months (190 d.i.v.). Top left shows combined t-SNE plots of control and mutant organoids (n = 3 single organoids per genotype). Cells are coloured by cell-type, and the number of cells per plot is indicated. Top right, t-SNE plots for control and mutant individual organoids. GABAergic interneurons are highlighted in colour. Bottom, bar charts show the percentage of cells for all the cell populations in each control and mutant organoid, as in a, b. f, Immunohistochemistry for the post mitotic neuronal marker TBR1 (magenta) and GABAergic marker DLX2 (green) in HUES66 control and CHD8+/− organoids at six months (190 d.i.v.). Scale bars: 100 μm. g, Pseudotime trajectories from the full dataset of HUES66 batch I CHD8+/− and control organoids at 109 d.i.v., calculated with Monocle3. The partition highlighted by a box was subsetted and the trajectory is shown in Fig. 3c. h, Module scores (top) and their distribution across mutant and control cells (bottom) for all modules resulting from WGCNA analysis of the partition of interest from HUES66 CHD8+/− and control organoids at 109 d.i.v. Cells were downsampled to have an equal number of cells per organoid. Names were assigned to each module based on the known functions of the genes included in each one. Horizontal bars show median scores, and dots show average score per organoid. Adjusted p-values show differences between control and mutant based on linear mixed models (see Methods). Number of organoids used for each analysis can be found in the Methods under “Statistics and reproducibility”. aRG, apical radial glia; DL, deep layer; UL, upper layer; PN, projection neurons; oRG, outer radial glia; IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; GABA INP, GABAergic interneuron progenitors; GABA IN, GABAergic interneurons.
Extended Data Fig. 10 Bulk RNA-seq and scRNA-seq of CHD8+/− and isogenic control organoids from multiple cell lines.
a, Bulk RNA-seq data from HUES66, GM83330 and H1 35 d.i.v. organoids. Enriched GO terms for genes differentially expressed in CHD8+/− vs. control organoids. The top 5 most significant terms per dataset are shown. Size of dot indicates the proportion of genes belonging to each term found in the list of dysregulated genes (“GeneRatio”). Colour indicates enrichment adjusted p-value. Numbers in parentheses along the y axis indicate the number of DEGs in that dataset. b–d, scRNA-seq data from control and CHD8+/− organoids at 3.5 months (GM83330 108 d.i.v., batch I (b), GM83330 108 d.i.v., batch II (c) and H1 105 d.i.v. (d)). Left top shows combined t-SNE plots of control and mutant organoids (n = 3 single organoids per genotype). Cells are coloured by cell type, and the number of cells per plot is indicated. Left bottom, t-SNE plots for control and mutant individual organoids. GABAergic interneurons are highlighted in colour. Right, bar charts show the percentage of cells for all the cell populations in each control and mutant organoid. Adjusted p-values for a difference in cell-type proportions between control and mutant, based on logistic mixed models (see Methods) are shown. aRG, apical radial glia; DL, deep layer; UL, upper layer; PN, projection neurons; CP/CH, Choroid Plexus/Cortical Hem; oRG, outer radial glia; IPC, intermediate progenitor cells; CPN, callosal projection neurons; CFuPN, corticofugal projection neurons; GABA INP, GABAergic interneuron progenitors; GABA IN, GABAergic interneurons; GABA N, GABAergic neurons.
Extended Data Fig. 11 Convergent differential expressed genes for the three mutations.
a, Log fold change of all genes which showed significant change (adjusted p < 0.05) in all three of the 1 month datasets: Mito210 SUV420H1+/− 35 d.i.v., Mito210 ARID1B+/− 35 d.i.v., and HUES66 CHD8+/− 35 d.i.v. DEGs were calculated using all cells as a pseudobulk for Mito210 SUV420H1+/− and Mito210 ARID1B+/−. b, Differential expression of all 102 genes associated with ASD risk6 in the three datasets Mito210 SUV420H1+/− 35 d.i.v., Mito210 ARID1B+/− 35 d.i.v. and in HUES66 CHD8+/− 35 d.i.v. compared to relative controls. Expression of risk genes was calculated using all cells (pseudobulk) for Mito210 SUV420H1+/− and Mito210 ARID1B+/−. Boxes are coloured according to -log10(adjusted p value) according to whether they are upregulated (purple), or downregulated (turquoise) in mutant vs. control. Genes are ordered according to hierarchical clustering (using Euclidean distance) of those values. c, d, Enriched GO terms for genes upregulated (c) and downregulated (d) in mutant vs. control. Genes were calculated using the cells as in a, b. The top 5 most significant terms per dataset are shown. Size of dot indicates the proportion of genes belonging to each term found in the list of dysregulated genes (“GeneRatio”). Colour indicates enrichment adjusted p-value. Numbers in parentheses along the x axis indicate the number of DEGs in that dataset.
Extended Data Fig. 12 Convergent differentially expressed proteins for the three mutations.
a–c, Volcano plot showing fold change versus adjusted p-value of measured proteins in MS experiments on Mito210 SUV420H1+/− (a), Mito210 ARID1B+/− (b), and HUES66 CHD8+/− (c) vs. control organoids at 35 d.i.v. (n = 4 single organoids per genotype for SUV420H1, 4 controls and 5 mutants for ARID1B, and n = 3 single organoids per genotype for CHD8). To detect statistically significant differential protein abundance between conditions a moderated t-test was performed (see Methods, FDR threshold of 0.1). Significant DEPs are shown in red (FDR < 0.1). d–f, Selected enriched GO terms for DEPs in Mito210 SUV420H1+/− (d), Mito210 ARID1B+/− (e), and HUES66 CHD8+/− (f) vs. control organoids cultured for 35 d.i.v. GO terms and KEGG pathways were calculated using the GSEA software (see Methods) and FDR q-values < 0.05 were considered statistically significant. g, Protein-protein interaction network using the top 50 DEPs from the three sets of mutant versus control organoids, created using the prize-collecting Steiner forest algorithm (see Methods). Protein nodes are coloured by the mutant in which they were differentially expressed. Gray nodes indicate “Steiner nodes”, proteins that did not result from any screen but were included by the algorithm to connect DEPs. Lines between nodes indicate physical protein-protein interactions from the STRING database, where line thickness correlates with interaction confidence. Subclusters of the network and significantly enriched terms for those subclusters are highlighted with gray rectangles and black text. h, Protein set distances between pairs of differentially expressed protein sets. For each pair of mutations, a PPI-weighted protein set distance was calculated between all significant DEPs (FDR < 0.1, pink diamond). To determine if this distance was smaller than would be expected by chance, size-matched sets were randomly chosen from the proteins detected in each experiment, and distance between these random sets was calculated 1000 times per pair. P-values were assigned by counting the fraction of times this random distance was less than the actual distance value between differential sets. Lower and upper hinges correspond to the first and third quartiles (the 25th and 75th percentiles) and middle hinge is the median (50th). Both whiskers extends from the hinge to the largest or smallest value no further than 1.5 * IQR from the hinge (where IQR is the inter-quartile range, or distance between the first and third quartiles). i, Protein set distances between the top 50 DEPs per mutation. For each pair of mutations, a PPI-weighted protein set distance was calculated as in h. Number of organoids used for the analyses are summarized in the Methods under “Statistics and reproducibility”. DEPs: differentially expressed proteins. MS: mass spectrometry.
Extended Data Fig. 13 Mutations in ASD risk genes in human brain organoids converge on asynchronous development of shared neuronal classes.
Conceptual schematics highlighting main results.
Supplementary information
Supplementary Information
Supplementary Notes and Supplementary References, and legends for Supplementary Fig. 1, Supplementary Tables 1–10 and Supplementary Video 1.
Supplementary Fig. 1
Uncropped gel images of western blots shown in this study.
Supplementary Video 1
Network bursting in a human brain organoid.
Rights and permissions
About this article
Cite this article
Paulsen, B., Velasco, S., Kedaigle, A.J. et al. Autism genes converge on asynchronous development of shared neuron classes. Nature 602, 268–273 (2022). https://doi.org/10.1038/s41586-021-04358-6
Received:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1038/s41586-021-04358-6
This article is cited by
-
Insights for disease modeling from single-cell transcriptomics of iPSC-derived Ngn2-induced neurons and astrocytes across differentiation time and co-culture
BMC Biology (2024)
-
Genetics of human brain development
Nature Reviews Genetics (2024)
-
iPSC-derived models of PACS1 syndrome reveal transcriptional and functional deficits in neuron activity
Nature Communications (2024)
-
Mechanisms guiding the slow pace of maturation in human neurons uncovered
Nature (2024)
-
Humanized brain organoids-on-chip integrated with sensors for screening neuronal activity and neurotoxicity
Microchimica Acta (2024)
Comments
By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.
