Skip to main content

Thank you for visiting nature.com. You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

  • Article
  • Published:

Integrating genetics with single-cell multiomic measurements across disease states identifies mechanisms of beta cell dysfunction in type 2 diabetes

Nature Genetics volume 55, pages 984–994 (2023)Cite this article

Subjects

Abstract

Dysfunctional pancreatic islet beta cells are a hallmark of type 2 diabetes (T2D), but a comprehensive understanding of the underlying mechanisms, including gene dysregulation, is lacking. Here we integrate information from measurements of chromatin accessibility, gene expression and function in single beta cells with genetic association data to nominate disease-causal gene regulatory changes in T2D. Using machine learning on chromatin accessibility data from 34 nondiabetic, pre-T2D and T2D donors, we identify two transcriptionally and functionally distinct beta cell subtypes that undergo an abundance shift during T2D progression. Subtype-defining accessible chromatin is enriched for T2D risk variants, suggesting a causal contribution of subtype identity to T2D. Both beta cell subtypes exhibit activation of a stress-response transcriptional program and functional impairment in T2D, which is probably induced by the T2D-associated metabolic environment. Our findings demonstrate the power of multimodal single-cell measurements combined with machine learning for characterizing mechanisms of complex diseases.

This is a preview of subscription content, access via your institution

Access options

Buy this article

  • Purchase on Springer Link
  • Instant access to full article PDF
Buy now

Prices may be subject to local taxes which are calculated during checkout

Fig. 1: Beta cells exhibit changes in chromatin activity in T2D.
Fig. 2: Machine learning identifies two beta cell subtypes with differential abundance in T2D.
Fig. 3: The two beta cell subtypes are distinguished by chromatin accessibility, gene expression and function.
Fig. 4: GRNs defining the two beta cell subtypes.
Fig. 5: Beta cell functional and gene regulatory changes in T2D.
Fig. 6: T2D risk variants affect beta cell subtype chromatin accessibility.

Similar content being viewed by others

Data availability

snATAC-seq data and processed data are available through the Gene Expression Omnibus (GEO) under accession number GSE169453, single-nucleus multiome data are available under accession number GSE200044, and genotyping data are available under accession number GSE170763. UCSC Genome Browser Sessions of aggregated snATAC-seq data are available at https://genome.ucsc.edu/s/gaowei/hg19_cell_type and https://genome.ucsc.edu/s/gaowei/hg19_beta_cell. Previously published16,17 Patch-seq data are available as raw sequencing reads in GEO under accession numbers GSE124742 and GSE164875. Additional Patch-seq data are accessible at the HPAP database (https://hpap.pmacs.upenn.edu). Blacklisted regions from ENCODE are found at https://www.encodeproject.org/annotations/ENCSR636HFF. Source data are provided with this paper.

Code availability

Custom codes for the main analysis used in this study have been deposited on GitHub at https://github.com/gaoweiwang/Islet_snATACseq.

References

  1. Noguchi, G. M. & Huising, M. O. Integrating the inputs that shape pancreatic islet hormone release. Nat. Metab. 1, 1189–1201 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  2. Wojtusciszyn, A., Armanet, M., Morel, P., Berney, T. & Bosco, D. Insulin secretion from human beta cells is heterogeneous and dependent on cell-to-cell contacts. Diabetologia 51, 1843–1852 (2008).

    CAS  PubMed  Google Scholar 

  3. Dominguez-Gutierrez, G., Xin, Y. & Gromada, J. Heterogeneity of human pancreatic β-cells. Mol. Metab. 27, S7–S14 (2019).

    CAS  PubMed Central  Google Scholar 

  4. Benninger, R. K. P. & Kravets, V. The physiological role of β-cell heterogeneity in pancreatic islet function. Nat. Rev. Endocrinol. 18, 9–22 (2022).

    CAS  PubMed  Google Scholar 

  5. Segerstolpe, Å. et al. Single-cell transcriptome profiling of human pancreatic islets in health and type 2 diabetes. Cell Metab. 24, 593–607 (2016).

    CAS  PubMed  PubMed Central  Google Scholar 

  6. Chiou, J. et al. Single-cell chromatin accessibility identifies pancreatic islet cell type– and state-specific regulatory programs of diabetes risk. Nat. Genet. 53, 455–466 (2021).

    CAS  PubMed  PubMed Central  Google Scholar 

  7. Cohrs, C. M. et al. Dysfunction of persisting β cells is a key feature of early type 2 diabetes pathogenesis. Cell Rep. 31, 107469 (2020).

    CAS  PubMed  Google Scholar 

  8. Chen, C., Cohrs, C. M., Stertmann, J., Bozsak, R. & Speier, S. Human beta cell mass and function in diabetes: recent advances in knowledge and technologies to understand disease pathogenesis. Mol. Metab. 6, 943–957 (2017).

    CAS  PubMed  PubMed Central  Google Scholar 

  9. Fadista, J. et al. Global genomic and transcriptomic analysis of human pancreatic islets reveals novel genes influencing glucose metabolism. Proc. Natl Acad. Sci. USA 111, 13924–13929 (2014).

    CAS  PubMed  PubMed Central  Google Scholar 

  10. Wigger, L. et al. Multi-omics profiling of living human pancreatic islet donors reveals heterogeneous beta cell trajectories towards type 2 diabetes. Nat. Metab. 3, 1017–1031 (2021).

    CAS  PubMed  Google Scholar 

  11. Xin, Y. et al. RNA sequencing of single human islet cells reveals type 2 diabetes genes. Cell Metab. 24, 608–615 (2016).

    CAS  PubMed  Google Scholar 

  12. Lawlor, N. et al. Single-cell transcriptomes identify human islet cell signatures and reveal cell-type-specific expression changes in type 2 diabetes. Genome Res. 27, 208–222 (2017).

    CAS  PubMed  PubMed Central  Google Scholar 

  13. Fang, Z. et al. Single-cell heterogeneity analysis and CRISPR screen identify key β-cell-specific disease genes. Cell Rep. 26, 3132–3144.e7 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  14. Wang, Y. J. & Kaestner, K. H. Single-cell RNA-seq of the pancreatic islets—a promise not yet fulfilled? Cell Metab. 29, 539–544 (2019).

    CAS  PubMed  Google Scholar 

  15. Dorrell, C. et al. Human islets contain four distinct subtypes of β cells. Nat. Commun. 7, 11756 (2016).

    PubMed  PubMed Central  Google Scholar 

  16. Camunas-Soler, J. et al. Patch-seq links single-cell transcriptomes to human islet dysfunction in diabetes. Cell Metab. 31, 1017–1031.e4 (2020).

    CAS  PubMed  PubMed Central  Google Scholar 

  17. Dai, X.-Q. et al. Heterogenous impairment of α cell function in type 2 diabetes is linked to cell maturation state. Cell Metab. 34, 256–268.e5 (2022).

    CAS  PubMed  PubMed Central  Google Scholar 

  18. Kahn, S. E., Cooper, M. E. & Del Prato, S. Pathophysiology and treatment of type 2 diabetes: perspectives on the past, present, and future. Lancet 383, 1068–1083 (2014).

    CAS  PubMed  Google Scholar 

  19. 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).

    PubMed  PubMed Central  Google Scholar 

  20. Shapira, S. N., Naji, A., Atkinson, M. A., Powers, A. C. & Kaestner, K. H. Understanding islet dysfunction in type 2 diabetes through multidimensional pancreatic phenotyping: the Human Pancreas Analysis Program. Cell Metab. 34, 1906–1913 (2022).

    CAS  PubMed  Google Scholar 

  21. Chen, T. & Guestrin, C. XGBoost: a scalable tree boosting system. In Proc. 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining 785–794 (Association for Computing Machinery, 2016).

  22. Sander, M. et al. Genetic analysis reveals that PAX6 is required for normal transcription of pancreatic hormone genes and islet development. Genes Dev. 11, 1662–1673 (1997).

    CAS  PubMed  Google Scholar 

  23. Xin, Y. et al. Pseudotime ordering of single human β-cells reveals states of insulin production and unfolded protein response. Diabetes 67, 1783–1794 (2018).

    CAS  PubMed  Google Scholar 

  24. Efron, B. & Tibshirani, R. On the testing of significance of sets of genes. Ann. Appl. Stat. 1, 107–129 (2007).

    Google Scholar 

  25. Cahan, P. et al. CellNet: network biology applied to stem cell engineering. Cell 158, 903–915 (2014).

    CAS  PubMed  PubMed Central  Google Scholar 

  26. Wang, G. et al. A tumorigenic index for quantitative analysis of liver cancer initiation and progression. Proc. Natl Acad. Sci. USA 116, 26873–26880 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  27. Sansbury, F. H. et al. SLC2A2 mutations can cause neonatal diabetes, suggesting GLUT2 may have a role in human insulin secretion. Diabetologia 55, 2381–2385 (2012).

    CAS  PubMed  Google Scholar 

  28. Vlacich, G., Nawijn, M. C., Webb, G. C. & Steiner, D. F. Pim3 negatively regulates glucose-stimulated insulin secretion. Islets 2, 308–317 (2010).

    PubMed  PubMed Central  Google Scholar 

  29. Stancill, J. S. et al. Chronic β-cell depolarization impairs β-cell identity by disrupting a network of Ca2+-regulated genes. Diabetes 66, 2175–2187 (2017).

    CAS  PubMed  PubMed Central  Google Scholar 

  30. Ye, R. et al. Inositol 1,4,5-trisphosphate receptor 1 mutation perturbs glucose homeostasis and enhances susceptibility to diet-induced diabetes. J. Endocrinol. 210, 209–217 (2011).

    CAS  PubMed  PubMed Central  Google Scholar 

  31. Martina, J. A., Diab, H. I., Brady, O. A. & Puertollano, R. TFEB and TFE3 are novel components of the integrated stress response. EMBO J. 35, 479–495 (2016).

    CAS  PubMed  PubMed Central  Google Scholar 

  32. Ohta, Y. et al. Clock gene dysregulation induced by chronic ER stress disrupts β-cell function. eBioMedicine 18, 146–156 (2017).

    PubMed  PubMed Central  Google Scholar 

  33. Eizirik, D. L., Pasquali, L. & Cnop, M. Pancreatic β-cells in type 1 and type 2 diabetes mellitus: different pathways to failure. Nat. Rev. Endocrinol. 16, 349–362 (2020).

    CAS  PubMed  Google Scholar 

  34. Lytrivi, M., Castell, A.-L., Poitout, V. & Cnop, M. Recent insights into mechanisms of β-cell lipo- and glucolipotoxicity in type 2 diabetes. J. Mol. Biol. 432, 1514–1534 (2020).

    CAS  PubMed  Google Scholar 

  35. Pratt, E. P. S., Harvey, K. E., Salyer, A. E. & Hockerman, G. H. Regulation of cAMP accumulation and activity by distinct phosphodiesterase subtypes in INS-1 cells and human pancreatic β-cells. PLoS ONE 14, e0215188 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  36. Bryan, J. et al. ABCC8 and ABCC9: ABC transporters that regulate K+ channels. Pflugers Arch. 453, 703–718 (2007).

    CAS  PubMed  Google Scholar 

  37. Yang, Y. et al. The phosphatidylserine flippase β-subunit Tmem30a is essential for normal insulin maturation and secretion. Mol. Ther. 29, 2854–2872 (2021).

    CAS  PubMed  PubMed Central  Google Scholar 

  38. Palu, R. A. S. & Chow, C. Y. Baldspot/ELOVL6 is a conserved modifier of disease and the ER stress response. PLoS Genet. 14, e1007557 (2018).

    PubMed  PubMed Central  Google Scholar 

  39. Tang, N. et al. Ablation of Elovl6 protects pancreatic islets from high-fat diet-induced impairment of insulin secretion. Biochem. Biophys. Res. Commun. 450, 318–323 (2014).

    CAS  PubMed  Google Scholar 

  40. Gaulton, K. J. Mechanisms of type 2 diabetes risk loci. Curr. Diab. Rep. 17, 72 (2017).

    PubMed  Google Scholar 

  41. Nkonge, K. M., Nkonge, D. K. & Nkonge, T. N. The epidemiology, molecular pathogenesis, diagnosis, and treatment of maturity-onset diabetes of the young (MODY). Clin. Diabetes Endocrinol. 6, 20 (2020).

    PubMed  PubMed Central  Google Scholar 

  42. Alonso, L. et al. TIGER: the gene expression regulatory variation landscape of human pancreatic islets. Cell Rep. 37, 109807 (2021).

    CAS  PubMed  PubMed Central  Google Scholar 

  43. Kirkpatrick, C. L. et al. Hepatic nuclear factor 1α (HNF1α) dysfunction down-regulates X-box-binding protein 1 (XBP1) and sensitizes β-cells to endoplasmic reticulum stress. J. Biol. Chem. 286, 32300–32312 (2011).

    CAS  PubMed  PubMed Central  Google Scholar 

  44. Szabat, M. et al. Reduced insulin production relieves endoplasmic reticulum stress and induces β cell proliferation. Cell Metab. 23, 179–193 (2016).

    CAS  PubMed  Google Scholar 

  45. Li, H. et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078–2079 (2009).

    PubMed  PubMed Central  Google Scholar 

  46. Lareau, C. A., Ma, S., Duarte, F. M. & Buenrostro, J. D. Inference and effects of barcode multiplets in droplet-based single-cell assays. Nat. Commun. 11, 866 (2020).

    CAS  PubMed  PubMed Central  Google Scholar 

  47. Wolf, F. A., Angerer, P. & Theis, F. J. SCANPY: large-scale single-cell gene expression data analysis. Genome Biol. 19, 15 (2018).

    PubMed  PubMed Central  Google Scholar 

  48. Amemiya, H. M., Kundaje, A. & Boyle, A. P. The ENCODE blacklist: identification of problematic regions of the genome. Sci. Rep. 9, 9354 (2019).

    PubMed  PubMed Central  Google Scholar 

  49. The ENCODE Project Consortium. An integrated encyclopedia of DNA elements in the human genome. Nature 489, 57–74 (2012).

    PubMed Central  Google Scholar 

  50. Traag, V. A., Waltman, L. & van Eck, N. J. From Louvain to Leiden: guaranteeing well-connected communities. Sci. Rep. 9, 5233 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  51. Korsunsky, I. et al. Fast, sensitive and accurate integration of single-cell data with Harmony. Nat. Methods 16, 1289–1296 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  52. Satpathy, A. T. et al. Massively parallel single-cell chromatin landscapes of human immune cell development and intratumoral T cell exhaustion. Nat. Biotechnol. 37, 925–936 (2019).

    CAS  PubMed  PubMed Central  Google Scholar 

  53. Pliner, H. A. et al. Cicero predicts cis-regulatory DNA interactions from single-cell chromatin accessibility data. Mol. Cell 71, 858–871.e8 (2018).

    CAS  PubMed  PubMed Central  Google Scholar 

  54. Ji, X., Li, W., Song, J., Wei, L. & Liu, X. S. CEAS: cis-regulatory element annotation system. Nucleic Acids Res. 34, W551–W554 (2006).

    CAS  PubMed  PubMed Central  Google Scholar 

  55. Schep, A. N., Wu, B., Buenrostro, J. D. & Greenleaf, W. J. chromVAR: inferring transcription-factor-associated accessibility from single-cell epigenomic data. Nat. Methods 14, 975–978 (2017).

    CAS  PubMed  PubMed Central  Google Scholar 

  56. Fornes, O. et al. JASPAR 2020: update of the open-access database of transcription factor binding profiles. Nucleic Acids Res. 48, D87–D92 (2020).

    CAS  PubMed  Google Scholar 

  57. 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).

    CAS  PubMed  PubMed Central  Google Scholar 

  58. Kuleshov, M. V. et al. Enrichr: a comprehensive gene set enrichment analysis web server 2016 update. Nucleic Acids Res. 44, W90–W97 (2016).

    CAS  PubMed  PubMed Central  Google Scholar 

  59. Taliun, D. et al. Sequencing of 53,831 diverse genomes from the NHLBI TOPMed Program. Nature 590, 290–299 (2021).

    CAS  PubMed  PubMed Central  Google Scholar 

  60. Das, S. et al. Next-generation genotype imputation service and methods. Nat. Genet. 48, 1284–1287 (2016).

    CAS  PubMed  PubMed Central  Google Scholar 

  61. Hinrichs, A. S. et al. The UCSC Genome Browser Database: update 2006. Nucleic Acids Res. 34, D590–D598 (2006).

    CAS  PubMed  Google Scholar 

  62. van de Geijn, B., McVicker, G., Gilad, Y. & Pritchard, J. K. WASP: allele-specific software for robust molecular quantitative trait locus discovery. Nat. Methods 12, 1061–1063 (2015).

    PubMed  PubMed Central  Google Scholar 

  63. Mahajan, A. et al. Multi-ancestry genetic study of type 2 diabetes highlights the power of diverse populations for discovery and translation. Nat. Genet. 54, 560–572 (2022).

    CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

This publication includes data generated at the UCSD IGM Genomics Center using an Illumina NovaSeq 6000 that was purchased with funding from a National Institutes of Health (NIH) Shared Instrument Grant (#S10 OD026929). This work was supported by NIH U01DK105541 and R01DK122607 to M.S. and K.J.G.; R01DK114650 to K.J.G.; R01DK068471 to M.S.; U01DK120447 to P.E.M.; U01DK123716 to S.K.K. and P.E.M.; and UC4-DK112217, UC4-DK112232 and P30 DK116074 to S.K.K. Work at the UCSD Center for Epigenomics was supported by the UCSD School of Medicine. The funders had no role in study design, data collection and analysis, decision to publish or preparation of the manuscript. We thank H. Gao, Y. Shi, K. A. Frazer, B. Ren, and members of Sander lab for scientific discussions and input on the project. We also thank the organ donors and their families for their contributions to make this study possible.

Author information

Author notes

  1. These authors contributed equally: Gaowei Wang, Joshua Chiou, Chun Zeng.

  2. These authors jointly supervised this work: Sebastian Preissl, Kyle J. Gaulton, Maike Sander.

Authors and Affiliations

  1. Department of Pediatrics, University of California San Diego, La Jolla, CA, USA

    Gaowei Wang, Joshua Chiou, Chun Zeng, Ileana Matta, Nikita Kadakia, Mei-Lin Okino, Elisha Beebe, Medhavi Mallick, Kyle J. Gaulton & Maike Sander

  2. Pediatric Diabetes Research Center, University of California San Diego, La Jolla, CA, USA

    Gaowei Wang, Joshua Chiou, Chun Zeng, Ileana Matta, Nikita Kadakia, Mei-Lin Okino, Elisha Beebe, Medhavi Mallick, Kyle J. Gaulton & Maike Sander

  3. Biomedical Graduate Studies Program, University of California San Diego, La Jolla, CA, USA

    Joshua Chiou

  4. Center for Epigenomics, University of California San Diego, La Jolla, CA, USA

    Michael Miller, Jee Yun Han & Sebastian Preissl

  5. Department of Bioengineering, Stanford University, Stanford, CA, USA

    Joan Camunas-Soler

  6. Department of Pharmacology, University of Alberta, Edmonton, Alberta, Canada

    Theodore dos Santos, Xiao-Qing Dai, Cara Ellis & Patrick E. MacDonald

  7. Alberta Diabetes Institute, University of Alberta, Edmonton, Alberta, Canada

    Theodore dos Santos, Xiao-Qing Dai, Cara Ellis & Patrick E. MacDonald

  8. Department of Developmental Biology, Stanford University School of Medicine, Stanford, CA, USA

    Yan Hang & Seung K. Kim

  9. Departments of Medicine and of Pediatrics, Stanford University School of Medicine, Stanford, CA, USA

    Yan Hang & Seung K. Kim

  10. Stanford Diabetes Research Center, Stanford University School of Medicine, Stanford, CA, USA

    Seung K. Kim

  11. Department of Clinical Diabetes, Endocrinology & Metabolism, City of Hope, Duarte, CA, USA

    Fouad R. Kandeel

  12. Institute of Experimental and Clinical Pharmacology and Toxicology, Faculty of Medicine, University of Freiburg, Freiburg, Germany

    Sebastian Preissl

  13. Institute for Genomic Medicine, University of California San Diego, La Jolla, CA, USA

    Kyle J. Gaulton & Maike Sander

  14. Department of Cellular and Molecular Medicine, University of California San Diego, La Jolla, CA, USA

    Maike Sander

  15. Max Delbrück Center for Molecular Medicine in the Helmholtz Association, Berlin, Germany

    Maike Sander

Authors
  1. Gaowei Wang
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Joshua Chiou
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Chun Zeng
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Michael Miller
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Ileana Matta
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Jee Yun Han
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Nikita Kadakia
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Mei-Lin Okino
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. Elisha Beebe
    View author publications

    You can also search for this author in PubMed Google Scholar

  10. Medhavi Mallick
    View author publications

    You can also search for this author in PubMed Google Scholar

  11. Joan Camunas-Soler
    View author publications

    You can also search for this author in PubMed Google Scholar

  12. Theodore dos Santos
    View author publications

    You can also search for this author in PubMed Google Scholar

  13. Xiao-Qing Dai
    View author publications

    You can also search for this author in PubMed Google Scholar

  14. Cara Ellis
    View author publications

    You can also search for this author in PubMed Google Scholar

  15. Yan Hang
    View author publications

    You can also search for this author in PubMed Google Scholar

  16. Seung K. Kim
    View author publications

    You can also search for this author in PubMed Google Scholar

  17. Patrick E. MacDonald
    View author publications

    You can also search for this author in PubMed Google Scholar

  18. Fouad R. Kandeel
    View author publications

    You can also search for this author in PubMed Google Scholar

  19. Sebastian Preissl
    View author publications

    You can also search for this author in PubMed Google Scholar

  20. Kyle J. Gaulton
    View author publications

    You can also search for this author in PubMed Google Scholar

  21. Maike Sander
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

M.S., K.J.G. and S.P. conceived and supervised the research in the study; M.S., K.J.G., G.W. and J.C. wrote the manuscript; G.W. and J.C. performed analyses of single-cell and genetic data; C.Z., I.M., N.K., J.Y.H. and M.-L.O. performed experiments; M. Miller performed 10x single-cell assays; E.B. and M. Mallick contributed to data analyses; F.R.K. provided human islets; J.C.-S., T.d.S., X.-Q.D., C.E., Y.H., S.K.K. and P.E.M. provided Patch-seq data.

Corresponding authors

Correspondence to Sebastian Preissl, Kyle J. Gaulton or Maike Sander.

Ethics declarations

Competing interests

K.J.G. holds stock in Vertex Pharmaceuticals and Neurocrine Biosciences. J.C. is now employed by and holds stock in Pfizer Inc. The other authors declare no competing interests.

Peer review

Peer review information

Nature Genetics thanks Vivek Swarup, Judith Zaugg, and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.

Additional information

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Extended data

Extended Data Fig. 1 Quality control of snATAC-seq data.

(a) Steps for snATAC-seq data processing and quality control. (b) Representative quality control metrics for each donor. Log10 total reads, fraction of reads overlapping promoters, fraction of reads overlapping peaks, and fraction of reads overlapping mitochondria DNA distribution of cells from library JYH809 as example. Blue vertical lines denote thresholds of 1000 minimal fragment number, 15% fragments overlapping promoters, 30% fragments overlapping peaks, and 10% fraction of reads overlapping mitochondria DNA, respectively. Red vertical lines denote thresholds to identify top 1% barcodes with extremely high total fragment number and fraction of reads overlapping promoters and peaks, respectively. (c) Representative cell clustering from library JYH809. (d) Promoter chromatin accessibility in a 5 kb window around TSS for endocrine marker genes in individual nuclei library JYH809. Total counts normalization and log-transformation were applied. (e) Cell clustering of chromatin accessibility profiles from all donors. (f) Representative low-quality cluster and subcluster. Cells in cluster 14 (top, highlighted in red) have significantly lower unique fragment than cells in other clusters (p = 2.3e-9, n = 255,598 cells). Cells in subcluster 1 (bottom, highlighted in red) have significantly lower fraction of reads overlapping peaks than cells in other clusters (p = 4.8e-5, n = 16,296 cells). Data are shown as mean ± S.E.M., ANOVA test with sex, age, BMI, disease status as covariates. (g) Log10 total reads, fraction of reads overlapping peaks and fraction of reads in promoters of cells from each cluster in Fig. 1b. Data are shown as mean ± S.E.M. (h) Promoter chromatin accessibility in a 5 kb window around TSS for selected endocrine and non-endocrine marker genes for each profiled cell (alpha: GCG, beta: INS-IGF2, delta: SST, gamma: PPY, acinar: REG1A, ductal: CFTR, stellate: PDGFRB, endothelial: CLEC14A, immune: CCL3). The UMAP projection is the same as in the main Fig. 1b. (i) Genome browser tracks showing aggregate read density (scaled to uniform 1 × 106 read depth) for cells within each cell type cluster at hormone gene loci for endocrine islet cell types. The gene body of each gene is highlighted.

Source data

Extended Data Fig. 2 Identification of factors explaining donor variability in snATAC-seq data.

(a,d,g,j) Absolute Spearman correlation coefficient between the first 6 principle components (PCs) and each biological or technical variable in beta (a), alpha (d), delta (g), and gamma (j) cells. An absolute Spearman correlation threshold of 0.4 was used to identify factors having a high correlation with each PC. (b,e,h,k) Principal component analysis (PCA) based on cCREs in beta (b), alpha (e), delta (h), and gamma (k) cells from individual non-diabetic (n = 11), pre-T2D (n = 8), and T2D (n = 15) donors which are color-coded by disease status. Each donor in the space is defined by the first two principal components (left) and the two principal components (right) that show highest correlation with disease status. (c,f,i,l) Pairwise Spearman correlation coefficients between biological or technical variables in beta (c), alpha (f), delta (i), and gamma (l) cells.

Source data

Extended Data Fig. 3 Validation of beta cell T2D-differential cCREs in snATAC-seq data from an independent cohort of donor islets.

(a) Clustering of chromatin accessibility profiles from HPAP human islet snATAC-seq data (non-diabetic, n = 13; pre-T2D, n = 2; T2D, n = 5). Nuclei are plotted using the first two UMAP components. Clusters are assigned cell type identities based on promoter accessibility of known marker genes (see Extended Data Fig. 3b). The number of nuclei for each cell type cluster is shown in parentheses. (b) Promoter chromatin accessibility in a 5 kb window around TSS for selected endocrine and non-endocrine marker genes for each profiled nucleus (alpha: GCG, beta: INS-IGF2, delta: SST, acinar: REG1A). (c) Heatmap showing chromatin accessibility at differential cCREs identified in Fig. 1e in HPAP snATAC-seq data. Columns represent beta cells from each donor (non-diabetic, n = 13; T2D, n = 5) and all non-diabetic and T2D donors with accessibility of peaks normalized by CPM (counts per million).

Source data

Extended Data Fig. 4 T2D affects chromatin accessibility more profoundly in beta cells than in other endocrine cell types.

(a) Volcano plot showing differential cCREs in beta cells between type 2 diabetic (T2D) and non-diabetic donors. Panels show all beta cells (left), beta cells down-sampled to 15,000 (middle), and 5,000 cells (right). Each dot represents a cCRE. cCREs with FDR < .1 after Benjamini-Hochberg correction (red dots) were considered differentially accessible. (b) Volcano plot showing differential cCREs in alpha cells between T2D and non-diabetic donors. Panels show all alpha cells (left), alpha cells down-sampled to 15,000 (middle), and 5,000 cells (right). Each dot represents a chromatin accessible cCRE. cCREs with FDR < .1 after Benjamini-Hochberg correction (red dots) were considered differentially accessible. (c) Volcano plot showing differential cCREs in delta cells between T2D and non-diabetic donors. Panels show all delta cells (left) and delta cells down-sampled to 5,000 cells (right). Each dot represents a chromatin accessible cCRE. cCREs with FDR < .1 after Benjamini-Hochberg correction (red dots) were considered differentially accessible.

Source data

Extended Data Fig. 5 Machine learning identifies two beta cell subtypes.

(a) Clustering of chromatin accessibility profiles from 92,780 beta cells from non-diabetic, pre-T2D and T2D donor islets using Scanpy (resolution=0.5). Nuclei are plotted using the first two UMAP components. (b) Position of beta cells from representative non-diabetic (MM80), pre-T2D (MM55), and T2D (MM54) donors on the UMPA. (c) Illustration of the strategy for distinguishing gradual from subtype changes in beta cells using machine learning. Possible scenarios for cell changes during T2D progression and expected disease state prediction accuracies for each scenario. In the case of no T2D-associated changes, the prediction accuracy for each disease state would be random (scenario 1), gradual cell state changes would be reflected by high prediction accuracy in each disease state (scenario 2), and subtype changes would be reflected by median prediction accuracies (scenario 3, here shown for two cell subtypes). (d, f, h) Relative abundance of predicted disease state among beta (d), alpha (f), and delta (h) cells from each donor using XGBOOST. Each column represents nuclei from one donor. (e, g, i) Relative abundance of predicted disease state among beta (e), alpha (g), and delta (i) cells in non-diabetic, pre-T2D and T2D donor islets. Data are shown as mean ± S.E.M. (n = 11 non-diabetic, n = 8 pre-T2D, n = 15 T2D donors), dots denote data points from individual donors. (j) Illustration of the strategy for identifying a classifier capable of distinguishing the two beta cell subtypes.

Source data

Extended Data Fig. 6 Validation of beta cell subtypes using independent data and computational methods.

(a) Relative abundance of beta-1 and beta-2 cells in male and female donor islets. Data are shown as mean ± S.E.M. (n = 9 females, n = 25 males), dots denote data points from individual donors. ANOVA test with age, disease, BMI, and islet index as covariates. (b-d) Pearson correlation between relative abundance of beta-2 cells and BMI (b), islet index (c), age (d) across donors (n = 34 donors). The bands around the linear regression line represent the range of the 95% confidence interval, two-sided Pearson test. (e) Relative abundance of beta-1 and beta-2 cells in islet snATAC-seq data from an independent cohort (n = 13 non-diabetic, n = 5 T2D donors). Each column represents cells from one donor. (f) Relative abundance of each beta cell subtype in non-diabetic and T2D donor islets. Data are shown as mean ± S.E.M (n = 13 non-diabetic, n = 5 T2D donors). **P < .01; ANOVA test with age, sex, and BMI as covariates. (g) Clustering of chromatin accessibility profiles from 92,780 beta cells from non-diabetic, pre-T2D and T2D donors using beta cell differential cCREs between non-diabetic and T2D donors from Fig. 1e. Cells are plotted using the first two UMAP components. (h) Relative abundance of each beta cell cluster based on UMAP annotation. Each column represents cells from one donor. (i) Position of beta cells from non-diabetic, pre-T2D and T2D donors on the UMAP. (j) Position of beta cells from representative non-diabetic (MM80), pre-T2D (MM55) and T2D (MM54) donors on the UMPA. (k) Relative abundance of each beta cell cluster in non-diabetic, pre-T2D and T2D donor islets. Data are shown as mean ± S.E.M. (n = 11 non-diabetic, n = 8 pre-T2D, n = 15 T2D donors). **P < .01, *P < .05; ANOVA test with age, sex, BMI, and islet index as covariates. (l) Overlap between beta cell subtypes identified using machine learning and beta cell clusters from UMPA. The overlap is 76.6% between cluster 1 and beta-1 and 74.3% between cluster 2 and beta-2. P = 2.2e-16 (Two-sided Binominal test).

Source data

Extended Data Fig. 7 Validation and characterization of beta cell subtypes using multiome data.

(a, b) Clustering of chromatin accessibility profiles of nuclei from multiome data (a). (b) Clustering of gene expression profiles of cells from multiome data (b). Nuclei are plotted using the first two UMAP components. Clusters are assigned cell type identities based on expression levels or promoter chromatin accessibility of known marker genes (alpha: GCG, beta: INS, delta: SST, gamma: PPY). The number of nuclei for each cell type cluster is shown in parentheses. n = 20 donors (c) Clustering of gene expression profiles of beta cells from multiome data using genes linked to differential proximal (within ± 5 kb of a TSS in GENCODE V19) and distal (based on potential distal cCRE-promoter connections inferred from Cicero, see Methods) cCREs between non-diabetic and T2D beta cells from Fig. 1e. Nuclei are plotted using the first two UMAP components. (d) Plots of beta cell subtypes predicted from chromatin accessibility profiles of beta cells from multiome data by machine learning. (e-f) Correlation between changes in proximal cCRE (within ± 5 kb of a TSS in GENCODE V19) accessibility and gene expression differences between beta-1 and beta-2 cells for differentially expressed genes from Fig. 3b. There are 544 proximal cCREs and target gene pairs in total, of which 511 have consistent changes between proximal cCRE accessibility and gene expression (e). Correlation between changes in distal cCRE (potential distal cCRE-promoter connections inferred from Cicero, see Methods) accessibility and gene expression differences between beta-1 and beta-2 cells for differentially expressed genes from Fig. 3b. There are 85 distal cCREs and target gene pairs in total, of which 72 have consistent changes between distal cCRE accessibility and gene expression (f). Two-sided Pearson test. (g) Enriched gene ontology terms among genes (see Fig. 3b) with higher (proximal or distal) cCRE accessibility and expression in beta-1 compared to beta-2 cells (left) and higher (proximal or distal) cCRE accessibility and expression in beta-2 compared to beta-1 cells (right).

Source data

Extended Data Fig. 8 Beta-1 and beta-2 cell classification in scRNA-seq data from independent cohorts.

(a, d, g, j) Clustering of gene expression profiles of beta cells from cohort 15, cohort 212, cohort 322, and Patch-seq cohort using differentially expressed genes between beta-1 and beta-2 from Fig. 3b. Cells are plotted using the first two UMAP components. The number of donors for each cohort and cells for each cell cluster is shown in parentheses. (b, e, h, k) Heatmap showing pseudo-bulk expression levels of differentially expressed genes between beta-1 and beta-2 (see Fig. 3b) in beta cells from cluster 1 and cluster 2 of cohort 15, cohort 212, cohort 322, and Patch-seq cohort. Expression levels of genes are normalized by TPM (transcripts per million). (c, f, i, l) Relative abundance of each beta cell subtype in non-diabetic and T2D donor islets in cohort 15 (n = 5 non-diabetic, n = 4 T2D) cohort 212 (n = 12 non-diabetic, n = 6 T2D), cohort 322 (n = 10 non-diabetic), and Patch-seq cohort (n = 15 non-diabetic, n = 16 pre-T2D, n = 14 T2D). Data are shown as mean ± S.E.M., dots denote data points from individual donors. **P < .01, ***P < .001; ANOVA test with age, sex, and BMI as covariates.

Source data

Extended Data Fig. 9 Transcriptional programs distinguishing the two beta cell subtypes.

(a) Genome browser tracks showing aggregate RNA and ATAC read density at representative genes positively regulated by HNF1A, HNF4A or HNF4G. Beta cell cCREs with binding sites for HNF1A, HNF4A and HNF4G are shown. (b) Genome browser tracks showing aggregate RNA and ATAC read density at representative genes positively regulated by NEUROD1, NFIA or TCF4. Beta cell cCREs with binding sites for NEUROD1, NFIA and TCF4 are shown.(a, b) All tracks are scaled to uniform 1 × 106 read depth, differential expressed genes between beta-1 and beta-2 are indicated by grey shaded boxes. (c) Box plots showing accessibility at HNF1A, HNF4A and HNF4G proximal cCREs in beta-1 and beta-2 cells. Genomic coordinates of promotor regions are shown in parentheses. (d) Bar plots showing expression of HNF1A, HNF4A and HNF4G in beta-1 and beta-2 cells. (e) Bar plots showing accessibility at NEUROD1, NFIA and TCF4 proximal cCREs in beta-1 and beta-2 cells. Proximal region of genes were shown in parentheses. (f) Bar plots showing expression of NEUROD1, NFIA, and TCF4 in beta-1 and beta-2. (c-f) Accessibility of peaks is normalized by CPM, gene expression is normalized by TPM. Data are shown as mean ± S.E.M., n = 20 donors, Two-sided paired t-test. (g) Genome browser tracks showing aggregate RNA and ATAC read density at HNF1A, HNF4A and HNF4G in beta-1 and beta-2 cells. Beta cell cCREs with binding sites for HNF1A, HNF4A and HNF4G are shown. (h) Genome browser tracks showing aggregate RNA and ATAC read density at NEUROD1, NFIA and TCF4 in beta-1 and beta-2 cells. Beta cell cCREs with binding sites for NEUROD1, NFIA and TCF4 are shown. (i) Genome browser tracks showing aggregate RNA and ATAC read density at HNF1A and TCF4 in beta-1 and beta-2 cells. Beta cell cCREs with binding sites for HNF1A and TCF4 are shown. All tracks are scaled to uniform 1 × 106 read depth.

Source data

Extended Data Fig. 10 Transcriptional programs changed in both beta cell subtypes in T2D.

(a) Heatmap showing chromatin accessibility at cCREs with differential accessibility in beta-1 cells from non-diabetic and T2D donors. Columns represent beta cells from each donor (non-diabetic, n = 11; pre-T2D, n = 8; T2D, n = 15) with accessibility of peaks normalized by CPM (counts per million). (b) Heatmap showing chromatin accessibility at cCREs with differential accessibility in beta-2 cells from non-diabetic and T2D donors. Columns represent beta cells from each donor (non-diabetic, n = 11; pre-T2D, n = 8; T2D, n = 15) with accessibility of peaks normalized by CPM. (c) Heatmap showing expression of genes negatively regulated by TFs (green) with higher activity in non-diabetic compared to T2D beta-1 cells (see Methods) and TFs (red) with lower activity in non-diabetic compared to T2D beta-1 cells (n = 6 non-diabetic, n = 8 pre-T2D, n = 6 T2D donors). Representative target genes of individual TFs are highlighted and classified by biological processes. Gene expression is normalized by TPM (transcripts per million). # denotes TFs with decreased or increased expression in T2D in both beta-1 and beta-2 cells. (d) Heatmap showing expression of genes negatively regulated by TFs (green) with higher activity in non-diabetic compared to T2D beta-2 cells (see Methods) and TFs (red) with lower activity in non-diabetic compared to T2D beta-2 cells (n = 6 non-diabetic, n = 8 pre-T2D, n = 6 T2D donors). Representative target genes of individual TFs are highlighted and classified by biological processes. Gene expression is normalized by TPM (transcripts per million). # denotes TFs with decreased or increased expression in T2D in both beta-1 and beta-2 cells. (e-g) Pearson correlation of expression levels between indicated TFs across pseudo-bulk RNA profiles from each beta cell subtype (40 dots in total: 20 donors including n = 6 non-diabetic, n = 8 pre-T2D, n = 6 T2D). The bands around the linear regression line represent the range of the 95% confidence interval. Two-sided Pearson test.

Source data

Supplementary information

Reporting Summary

Peer Review File

Supplementary Table

Combined workbook including Supplementary Tables 1–14.

Source data

Source Data Figs. 1–6 and Extended Data Figs. 1–10

Statistical source data for Figs.1–6 and Extended Data Figs. 1–10. Source data for individual figures are labeled.

Rights and permissions

Springer Nature or its licensor (e.g. a society or other partner) holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Wang, G., Chiou, J., Zeng, C. et al. Integrating genetics with single-cell multiomic measurements across disease states identifies mechanisms of beta cell dysfunction in type 2 diabetes. Nat Genet 55, 984–994 (2023). https://doi.org/10.1038/s41588-023-01397-9

Download citation

  • Received:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1038/s41588-023-01397-9

This article is cited by

Search

Advanced search

Quick links

Nature Briefing

Sign up for the Nature Briefing newsletter — what matters in science, free to your inbox daily.

Get the most important science stories of the day, free in your inbox. Sign up for Nature Briefing