Human islet donors
Information about the donors included in the study is presented in Supplementary Data 1. Human islets for dataset 1 were obtained from the Nordic Network for Islet Transplantation (www.nordicislets.org). Written informed consent was obtained from the donors or their relatives, and the ethics committees at Uppsala and Lund Universities approved all procedures. Human islets from dataset 2 were purchased from Prodo Laboratories Inc. (Irvine, CA, USA), providing pancreatic islets isolated from donors obtained with research consent from Organ Procurement Organizations (OPOs)75. The use and storage of human islets and tissue samples were performed in compliance with the Declaration of Helsinki, ICH/Good Clinical Practice, and were approved by the Regional Ethics Committee (Gothenburg, Sweden).
Sample processing and single-cell RNA sequencing
Human islets were stored in 30 ml HEPES-buffered medium. Islets were transferred to a 50 ml tube and centrifuged at 150 g for 2 min at RT. Supernatant was removed, and cells were washed once with 5 ml of Accutase. Another 5 ml of Accutase was added for incubation at 37 °C for 8–10 min. Thereafter, 5 ml of cold islet media was added, and cells were suspended by pipetting. Cells were passed through a 40-μm cell strainer and centrifuged at 500 × g for 5 min and washed twice with PBS at 500 × g for 5 min. Finally, cells were resuspended in HBSS and kept on ice before being FACS-sorted as single cells into 384-well plates. cDNA libraries were generated using the Smart-seq2 protocol27 and sequenced on a HiSeq 2500 sequencer.
Filtering and cell type clustering
For dataset 1. STAR 2.4.2a was used for alignments with the reference genome hg38 using 2-pass alignment for improved performance of de novo splice junction reads, filtered for uniquely mapping reads only, and saved in .bam files. The count matrix showed the individual counts aligning to each gene per cell. Expression values were computed as reads per kilobase of the gene model and million mappable reads (RPKMs). Sequenced cells with mRNA reads <100,000, percent of uniquely mapping reads <50%, or percent of uniquely exonic reads <40% were removed, obtaining 3645 transcriptomes.
For dataset 2, raw sequencing data of human pancreatic islets were demultiplexed and converted into fastq files using Illumina bcl2fastq with default settings. Reads from human data were aligned to the human genome hg38 using the STAR aligner. Uniquely aligned reads mapping to the RefSeq gene annotations were used for gene expression estimation at reads per kilobase transcript and million mapped reads (RPKMs) using rpkmforgenes. Low-quality cells were excluded from downstream analysis when they failed to meet the following criteria for retaining cells: (1) ≥50,000 sequence reads; (2) ≥40% of reads uniquely aligned to the genome; (3) ≥40% of these reads mapping to RefSeq annotated exons; and (4) ≥1000 genes with RPKM > = 1. In addition, doublets detected by Scrublet were further removed, resulting in a total of 4866 single cells for downstream analysis. For cell type assignment, the top 1000 most variable genes were used for PCA dimension reduction, followed by clustering using affinity propagation, resulting in seven non-endocrine cell clusters (acinar, ductal, PSC, MHC II, mast, stem-like, endothelial cells) and one endocrine cell cluster. To further annotate endocrine cell types, we applied the same approach to endocrine cells only and identified five clusters representing alpha cells, beta cells, delta cells, gamma cells, and epsilon cells, respectively. For the number of cells in each cluster and dataset by disease status, see Supplementary Data 1.
Integration of datasets
For datasets 1 and 2, the expression of transcript variants with the same gene name was summed. The datasets were integrated and clustered using Conos, and cell types were assigned using label propagation based on previous clustering analysis and known markers. We used the basicP2proc function from Pagoda to preprocess each dataset and facilitate integration. A conos object was created with ‘Conos::new()‘ to analyze the joint dataset. Clustering and dimensional reduction were performed using the buildGraph function from pagoda. A shared nearest neighbor (SNN) graph was built using principal component analysis (PCA) as the dimensionality reduction method, retaining the top 40 components. Global clusters were detected using the Leiden clustering algorithm implemented in the findCommunities function from Pagoda. The graph was embedded into a two-dimensional space using UMAP for visualization, and markers were selected based on a combination of differential expression and biological expertise.
Identification of cell type-specific gene expression
For Fig. 1C, we plotted canonical markers for each cell type in the integrated dataset (1 + 2). For Fig. 1D, we ran “getDifferentialGenes” function in Conos and plotted the top 5 Z-scores ordered by the specificity of expression.
Differential Gene Coordination Network Analysis (dGCNA)
Cell type-specific gene-gene correlation matrices for each cell type with more than 30 cells/donor were constructed separately for T2D and non-T2D donors using Pearson correlation between transcriptional levels in log2(RPKM + 1) of each gene pair. We filtered non-expressed genes (RPKM + 1 < 2) and genes detected (RPKM + 1 > 2) in less than 5 cells/donor. The correlations were calculated using a linear mixed-effects model including random effects to correct for donor-specific expression differences. To build cell type-specific differential networks between disease and healthy states, we calculated the difference in pairwise Pearson correlations for each gene pair between T2D and non-T2D donors. To remove false positive weights in the differential network, we used a bootstrapping approach, generating random differential networks by selecting 512 pairs of random groups of donors. We filtered all correlations with a relative frequency lower than 0.975, keeping only links that were significantly stronger than random. This approach provides an internally scaled thresholding that takes the quality of the underlying dataset into account. We then performed clustering of the differential networks (one per cell type) using TOMdist (WGCNA package) and hierarchical clustering (hclust) to reveal gene communities altered between the conditions. Networks of differentially coordinated genes (NDCGs), i.e., gene modules, were determined using the cutreeDynamic function (WGCNA package) at different thresholds. Next, we evaluated whether these T2D-related NDCGs were associated with known biological functions using a hypergeometric test (Bonferroni correction, p < 0.05) with Gene Ontology terms (Cellular component, Biological function, and Molecular function). Protein-coding genes detected in our dataset for each cell type were used as the background gene-set. Network scores were calculated with the function eigen_centrality (igraph R package) for the positive (hyper-coordinated) subnetwork (network with just positive weights) and for the negative (de-coordinated) subnetwork (network with just negative weights) separately.
Differential expression analysis
We used the DESeq2 method to calculate the differential expression of genes. Expression of all single-cell counts of the same cell type and the same donor was summed before calculating differential expression. We removed unwanted variation from the RNA-seq counts with the package RUVSeq using control genes. Thresholds of the volcano plots were established to a fold change of 2 and an adjusted p-value of 0.05. For bulk RNA-Seq data, the expression count matrix was generated with Salmon v1.5.0, and differential expression analysis was performed with DESeq2 v1.30.1.
Evaluating the replicability of the dCGNA analysis
We evaluated similarity between the picked modules in datasets 1 and 2 by calculating the pairwise Jaccard index, considering the overlap in significant GO-terms.
Enrichment analysis of published datasets (Patch-Seq, GWAS, and whole genome CRISPR screen)
We used a hypergeometric test for enrichment of overlap between the dCGNA modules in T2D- and islet-related gene sets. To decrease the loss of power due to multiple comparison correction, we required a minimum overlap size of 2 or more genes for a test to be performed. All expressed genes in the cell type corresponding to the module were considered as background. We computed the adjusted p-value within each gene-set for all modules using the FDR method.
Selection of genes for functional validation
We selected genes for functional follow-up based on the following criteria: high eigencentrality in the RDN from different modules. We then applied filters: (1) no previous publication with search terms “beta cell” or “islet” on PubMed, (2) robust expression in INS-1 832/13 cells (assessed with bulkRNAseq). One exception was made for DBI, which had two older papers in a hamster cell line and rat islets76,77. TMEM176A/B and CEBPG were selected for further in vivo studies due to the availability of KO mouse models.
INS-1 832/13 cell culture and siRNA-mediated gene silencing
INS-1 832/13 cells47 were cultured at 5% CO2 and 37 °C in RPMI1640 medium (Sigma Aldrich, St Louis, MO) containing 2 g/L D-glucose, supplemented with 10% fetal bovine serum, 10 mM HEPES, 1 mM sodium pyruvate, and 50 μM β-mercaptoethanol (Sigma Aldrich). Gene silencing in INS-1 832/13 cells was performed using Lipofectamine RNAiMAX (Life Technologies, Waltham, MA) and 60 nM siRNA targeting Atraid mRNA (s151757 and s235878), Cebpg mRNA (s129415 and s129417), Dbi mRNA (s128869 and s128871), Dstn mRNA (282408 and 282409), Gabarap mRNA (s133040 and s133041), Spire1 (s158071 and s158073), Tmem176a mRNA (s150574) and Tmem176b mRNA (s140062). The sequences for scrambled siRNA were sense: 5´-GAGACCCUAUCCGUGAUUAtt-3´ and antisense: 5´-UAAUCACGGAUAGGGUCUCtt-3´ (Silencer Select, rat negative control #1; scrambled siRNA and all siRNAs were from Ambion, Life Technologies). Transfection complexes were prepared in accordance with the instructions provided by the manufacturer and added to 2 × 105 cells seeded per well in 24-well plates.
Human islets were transfected with 60 nM siRNA targeting human CEBPG (s2901, Silencer Select Pre-designed siRNA, Ambion, Life Technologies). The sequence for CEBPG siRNA was sense: 5’-AGAGCCGGUUGAAAAGCAAtt-3´ and antisense: 5’-UUGCUUUUCAACCGGCUCUtt-3´.
RNA extraction
Total RNA was isolated from INS-1 832/13 cells 48 h after transfection using the NucleoSpin extraction kit (Macherey-Nagel, Düren, Germany). The amount of isolated RNA was measured using the NanoDrop system. For bulk RNA sequencing, samples were also analyzed using TapeStation. Total RNA from human islets was extracted using Nucleospin RNA XS (Macherey-Nagel) three days post-transfection.
cDNA synthesis and real-time qPCR
For the experiment, 1 μg of isolated total RNA was reverse-transcribed to cDNA using the RevertAid First Strand cDNA synthesis kit (Life Technologies). qPCR was performed with 25 ng cDNA template using TaqMan Expression Master Mix (Life Technologies) according to the instructions provided by the manufacturer. Tbp and Hprt1 (Rn01527840_m1 and Rn01455646_m1, respectively) were used as housekeeping genes. Expression levels were calculated using the 2-ΔΔCt-method. TaqMan assays used were: Atf4 (Rn00824644_g1), Atf6 (Rn01490844_m1), Atraid (Rn01468747_g1), Cebpg (Rn01764319_m1), Dbi (Rn00821402_g1), Dstn (Rn01415640_g1), Gabarap (Rn00490680_g1), Ins1 (Rn0212433_g1), Ins2 (Rn01774648_g1), Spire1 (Rn01491339_m1), Tmem176a (Rn01451723_m1), Tmem176b (Rn00508100_m1), Xbp1s (Rn03464499), Cdc42 (Rn00696671), Nkx2.2 (Rn04244749), Nkx6.1 (Rn01450076), NeuroD1 (Rn00824571), Pdx1 (Rn00755591), Rab3a (Rn07311159), Rab27a (Rn00568995), Snap25 (Rn00578534), Syt11 (Rn01218287), and Xbp1 (Rn01443523_m1). For human islets, 300 ng of RNA was reverse-transcribed to cDNA using RevertAid First Strand cDNA synthesis kit. RT-qPCR for target genes and two endogenous controls (HPRT1 and TBP) was performed using 2.5 ng cDNA. TaqMans for human islet were: ATF4 (Hs00909569_g1), ATF6 (Hs00232586_m1), CEBPG (Hs01922818_s1), DDIT3 (Hs99999172_m1), HPRT1 (Hs4326321_m1), INS (Hs02741908_m1), TBP (Hs00427620_m1), and XBP1 (Hs00231936_m1). All TaqMan assays were from Life Technologies, and qPCR reactions were run on the Viia7 real-time PCR system (Applied Biosystems, Foster City, CA).
Insulin secretion experiments
INS-1 832/13 cells were seeded in 24-well plates (2 × 105 cells/well) and allowed to adhere for approximately 5 h before transfection. On the day of insulin secretion, cells were washed in PBS, incubated with 2.8 mM glucose for 2 h, and then incubated for 60 min with 2.8 mM glucose, 16.7 mM glucose, 16.7 mM glucose +0.1 mM IBMX, 2.8 mM glucose +35 mM K+, 2.8 mM glucose +10 mM α-ketoisocaproic acid, 16.7 mM glucose +10 mM arginine, or 16.7 mM glucose +1 mM palmitate. For the experiments, glucose was dissolved in the secretion assay buffer (SAB; 114 mM NaCl, 4.7 mM KCl, 1.2 mM KH2PO4, 1.16 mM MgSO4, 1.16 mM CaCl2, 20 mM HEPES, 25.5 mM NaHCO3, and 0.2% BSA). The cells were then placed on ice, and the culture media were collected for insulin and protein determinations. Insulin concentration was determined by ELISA from Mercodia (Uppsala, Sweden), and protein concentration was determined using Protein Assay Dye Reagent (BioRad, Hercules, CA).
Thapsigargin treatment
INS-1 832/13 cells were seeded in 24-well plates (2 × 105 cells/well). Cells were incubated with 100 nM thapsigargin (Sigma Aldrich) for 24 h. In a second experiment, cells were seeded in 24-well plates and Cebpg mRNA was silenced using siRNA (s129415) for 24 h. 100 nM thapsigargin was then added to these cells for 24 h. Total RNA was extracted from the cells from both above-mentioned experiments.
Immunohistochemistry
Primary antibodies: anti-rabbit ATF6 (code NBP1-76675, dilution 1:200, Novus Biologicals, Centennial, CO), anti-rabbit CEBPG (code HPA012024, dilution 1:20, Sigma Aldrich), anti-rabbit DSTN (code ab186754, dilution 1:100, Abcam), anti-rabbit GABARAP (code HPA-78365, 1:500 dilution, Sigma Aldrich), anti-guinea pig insulin (code A0564, 1:1000 dilution, Dako, Glostrup, Denmark), anti-rabbit KIAA1324 (code PA5-67123, dilution 1:800, Life Technologies), anti-rabbit SPIRE1 (ab130403, 1:1000 dilution, Abcam), anti-rabbit TMEM176A (code NBP1-83283, 1:20 dilution, Novus Biologicals), anti-rabbit TMEM176B (code CSB-PA023758LA01HU, 1:50 dilution, CusaBio, Houston, TX), anti-rabbit DBI (code DBI-AK3, 1:800 dilution (kind gift from Prof. Efendik, Karoliska Institute, Stockholm, Sweden)), anti-rabbit ATRAID (code bs-9806R, 1:200 dilution, Thermo Scientific), and anti-rabbit RPS3 (code ab128995, 1:200 dilution, Abcam).
All antibodies were diluted in PBS with 0.25% BSA and 0.25% Triton X-100. Slides were incubated with primary antibodies overnight at 4 °C, then washed twice in PBS with 0.25% Triton X-100. Slides were then incubated for 1 h with secondary antibodies at RT. Secondary antibodies used: Donkey anti-guinea pig AlexaFluor594 (for insulin), donkey anti-rabbit Cy2 (for ATF6, CEBPG, DSTN, GABARAP, KIAA1324, SPIRE1, TMEM176A, TMEM176B, DBI, ATRAID, RPS3), and goat anti-guinea pig AlexaFluor405 (for insulin). For phalloidin staining, Phalloidin-iFluor488 (Abcam) was diluted in PBS (+0.25% Triton X-100 and 0.25% BSA; 1:1000) and incubated for 60 min at RT. Slides were then washed twice in PBS with 0.25% Triton X-100 and mounted. Immunofluorescence was examined in an epifluorescence microscope (Olympus BX60, Olympus, Tokyo, Japan) or a Zeiss LSM800 confocal microscope (Oberkochen, Germany). Images were acquired with a digital camera (Olympus DP74, Olympus) using the CellSens Dimensions software (Olympus) or Zen System 2.6 (Zeiss). Pixel intensity was measured using ImageJ software, and colocalization (expressed as Pearson´s colocalization coefficient) analysis was performed using the CellSens Dimensions software in a blinded manner.
Animals
Tmem176a/b double knockout (DKO) mice have been described in detail elsewhere53. For this study, the mice were kept in Macrolon cages in a temperature-controlled environment (21 °C) with a relative humidity of 40–60%, on a 12 h light/dark cycle and with free access to diet (3 weeks of standard rodent chow or high-fat diet (60% fat, 20% protein, 20% carbohydrates, Special Diet Services, Witham, UK)) and tap water. All experiments were performed in accordance with institutional guidelines for animal care and use for the University of Nantes. The Cebpg KO mice48 were kept in Macrolon cages in a temperature-controlled environment (21 °C) with a relative humidity of 40–60%, on a 12 h light/dark cycle, and with free access to standard rodent chow and tap water. All animal experiments were approved by the animal ethics committee in Malmö and Lund, Sweden.
Beta cell mass measurements in Cebpg KO mice and Tmem176a/b DKO mice
Male Tmem176a/b DKO and WT (n = 8 per group, age: 12–14 weeks, weight: 28.5 + 3) and Cebpg KO (n = 6) and WT (n = 16) (50% males, age: 12 weeks) mice were analyzed. Beta cell mass, average islet size, and total islet number measurements in the KO mouse models were performed as detailed78 at three depths of the pancreas (100 μm difference between a series of three sections). Insulin was visualized using an anti-insulin antibody (code A0564, 1:1000 dilution, Dako) and anti-guinea pig AlexaFluor594 (1:400 dilution), and the signal was used to determine the total beta cell area, which was divided by the total pancreas area to generate a measure of beta cell mass. Immunofluorescence was examined under an epifluorescence microscope (Olympus BX60, Olympus, Tokyo, Japan) with a digital camera (Olympus DP74, Olympus). Measurements were made at 40× magnification using the CellSens Dimensions software (Olympus) in a blinded manner.
Intraperitoneal glucose tolerance test (IpGTT)
Male Tmem176a/b DKO mice fed HFD (n = 10, age: 15 weeks, weight: 28 + 2.7 g) or a control diet (n = 10, age 12 weeks, weight: 25 + 2.9) and male WT mice fed HFD (n = 10, age 14 weeks, weight: 27 + 0.9) or a control diet (n = 10, age 13 weeks, weight: 24.6 + 2.7) were used, as previous studies show that males develop more pronounced insulin resistance and beta cell stress compared with females79. After a 4-h fasting period, basal blood samples were drawn via retro-orbital puncture. The mice were then injected intraperitoneally with glucose (2 g/kg bodyweight). Blood samples were collected after 10, 20, 30, 60, and 90 min into chilled EDTA-tubes. Samples were maintained on ice and centrifuged at 500 × g for 5 min, and plasma was stored at −80 °C until analysis.
Insulin secretion and glucose measurement
Plasma insulin from the Tmem176a/b DKO mice subjected to IpGTT and fasting plasma insulin from the Cebpg KO mice were determined by ELISA (Mercodia) according to the manufacturer’s instructions. Glucose from mouse plasma was analyzed using the Infinity Glucose (Ox) kit from Thermo Scientific according to the instructions provided by the manufacturer. Acute insulin response (AIR) in the Tmem176a/b DKO mice was calculated as the difference between 10-min insulin levels and basal insulin levels. Integrated area under the curve (iAUC) for the glucose levels in the Tmem176a/b DKO mice subjected to IpGTT was calculated using GraphPad Prism 8 (GraphPad Prism, San Diego, CA).
Mitochondrial oxygen consumption measurements
Oxygen consumption rates (OCRs) were measured by Seahorse XFe24 Extracellular Flux Analyzer (Agilent Technologies, USA). INS-1 832/13 cells (120,000 cells/well) in poly-L-lysine-coated Seahorse XF cell culture plates were preincubated in a CO2-free incubator for 2 h at 37 °C in assay buffer (114 mM NaCl, 4.7 mM KCl, 1.2 mM KH2PO4, 1.16 mM MgSO4, 20 mM HEPES, and 2.5 mM CaCl2, pH 7.2) with 2.8 mM glucose. To determine changes in mitochondrial respiratory response, stabilized cellular respiratory rate was measured before and after injection of 16.7 mM glucose and subsequently after sequential injection of 4 µg/ml oligomycin, 2 µM FCCP, and 1 µM rotenone. Aside from injection timings, experiments were performed as previously described80. Primary data analysis was performed using the Seahorse Analytics web tool (https://seahorseanalytics.agilent.com/; version 1.0.0-699).
Metabolomic analysis
Metabolites were extracted and profiled using gas chromatography/mass spectroscopy, as previously described81. Metabolite data were mean-centered and scaled, followed by analysis by principal component analysis (PCA; prcomp) and linear models (lm; Anova, car R package). Plots were produced in ggplot (ggplot2 R package). The analysis provided 623 features, out of which 110 were putatively identified by mass spectrum and retention index using the MS-DIAL associated libraries and in-house libraries.
Ca2+ imaging
INS-1 832/13 cells were cultured in imaging-specific plates coated with poly-L-lysine (1μg/ml working solution) from Nunc (Lab-Tek Chambered #1.0 Borosilicate Coverglass System) in regular RPMI1640 media (Sigma Aldrich). Ca2+ measurements were performed as detailed82 72 h after transfection. Briefly, cells were washed in PBS and then preincubated in 2.8 mM glucose for 2 h at 37 °C. 1.5 h into the incubation, cells were loaded with 2.5 μM, 0.5 mM sulfinpyrazone and incubated for the remaining period. Prior to imaging, the pre-incubation buffer was discarded, and the cells were washed with imaging buffer. Upon imaging, a 2-min baseline (2.8 mM glucose) was recorded before the addition of a highglucose (16.7 mM) buffer. Images were analyzed using ImageJ software.
Western blotting for Grp78/BiP
INS-1 832/13 cells were lysed in RIPA buffer 72 h after transfection. Protein content was measured using the BCA Protein assay kit from Pierce. 20 μg of protein was loaded onto Mini PROTEAN TGX stain-free gels (BioRad). Precision Plus Protein Kaleidoscope Prestained Protein Standards (#1610375, BioRad) were used. Gels were activated using a ChemiDoc MP CCD camera from BioRad. Blotting was performed using a Transblot Turbo Transfer System (BioRad). Membranes were activated and then blocked for 1 h using TBS-T with 5% non-fat dry milk. Then, membranes were incubated overnight at 4 °C with Grp78 primary antibodies (#3183S, Cell Signalling, Danvers, MA), diluted 1:1000 in TBS-T with 5% non-fat dry milk. Membranes were then washed three times in TBS-T, incubated with secondary antibody diluted in TBS-T with non-fat dry milk for 60 min, washed three times in TBS-T, developed using SuperSignal West Pico Maximum Sensitivity Substrate (Thermo Scientific), and imaged using a ChemiDoc MP CCD camera. Proteins were quantified based on total protein normalization with stain-free gels in Image Lab version 6.0.1 (BioRad).
Protein aggregation assay
Thioflavin T (ThT) is a benzothiazole dye that exhibits enhanced fluorescence when binding to ER stress-induced protein aggregates83. ThT was dissolved in 70% ethanol (60 mM) and stored at −20 °C. Immediately before use, ThT was diluted in RPMI1640 to create a 5 mM stock solution. The stock solution was further diluted in assay media to a final concentration of 5 μM. INS-1 832/13 cells were seeded into 96-well plates at a total density of 1 × 104 cells/well for Cebpg KD. Then, 72 h post-transfection, cells were treated with thapsigargin (500 nM) and ThT (5 μM), and fluorescence (excitation 437 nm, emission 485 nm) was read after a 1-h incubation at 37 °C.
Phalloidin staining in INS-1 832-13 cells
Cells were transfected with Tmem176a and Tmem176b siRNAs in 8-well chamber slides (Sarstedt). 72 h after transfection, cells were incubated with anti-guinea pig insulin antibodies (Dako) overnight at 4 °C, washed twice in PBS (with 0.25% BSA and 0.25% Triton X-100), incubated with donkey anti-guinea pig AlexaFluor594 (Jackson ImmunoResearch) for 60 min at room temperature, again washed twice with PBS (with 0.25% BSA and 0.25% Triton X-100). Phalloidin (iFluor488 from Abcam) was diluted 1:1000, and cells were stained for 60 min at RT. Phalloidin intensity was quantified using ImageJ software in a blinded manner.
Analysis and statistics
All experimental data were analyzed using 2-way ANOVA, followed by Tukey’s test post hoc, or by two-tailed Student’s t-test. Differences with p < 0.05 were considered significant. All replicates mentioned were biological replicates.
Reporting summary
Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.