Differential immune- and apoptosis-related gene signatures in pancreatic alpha and beta cells contribute to their fate in type 1 diabetes

Differential immune- and apoptosis-related gene signatures in pancreatic alpha and beta cells contribute to their fate in type 1 diabetes

Collection of raw scRNA-seq data

For the human islets datasets, raw scRNA-seq FASTQ files (10X Genomics) were downloaded from the Human Pancreas Analysis Program (HPAP) portal (https://hpap.pmacs.upenn.edu) (Table S1). HPAP-young and HPAP-older donors represent matched T1D and T2D controls curated by the HPAP consortium based on similar age and body mass index (BMI) characteristics. We analyzed 15 HPAP-young donors (HPAP-022, HPAP-026, HPAP-027, HPAP-034, HPAP-035, HPAP-036, HPAP-037, HPAP-039, HPAP-040, HPAP-044, HPAP-047, HPAP-056, HPAP-082, HPAP-099 and HPAP-104), 13 HPAP-older donors (HPAP-052, HPAP-053, HPAP-054, HPAP-059, HPAP-063, HPAP-074, HPAP-075, HPAP-077, HPAP-080, HPAP-093, HPAP-101, HPAP-103 and HPAP-105), 9 AAB1+ donors (HPAP-024, HPAP-029, HPAP-038, HPAP-045, HPAP-049, HPAP-050, HPAP-072, HPAP-092 and HPAP-114), 2 AAB2+ donors (HPAP-043 and HPAP-107), 11 T1D donors (HPAP-020, HPAP-021, HPAP-023, HPAP-032, HPAP-055, HPAP-064, HPAP-071, HPAP-084, HPAP-087, HPAP-113 and HPAP-123) and 10 T2D donors (HPAP-051, HPAP-057, HPAP-058, HPAP-061, HPAP-065, HPAP-070, HPAP-079, HPAP-083, HPAP-085 and HPAP-088). The clinical characteristics of these donors are available at the HPAP portal (https://hpap.pmacs.upenn.edu/explore/download?donor) by referencing their respective donor IDs. The raw scRNA-seq FASTQ files (10X Genomics) of hiPSC-derived and human embryonic stem cells (hESC)-derived islet-like cells were retrieved from the Gene Expression Omnibus (GEO) database under the accession numbers: GSE190726 (Chandra et al., hESC stage 6), GSE203384 (Szymczak et al., iPSC stage 7), GSE20083 (hESC-derived cells matured retrieved in vivo after transplantation into humanized mouse) and GSE20084 (hESC-derived cells differentiated in vitro), as summarized in Table S1.

scRNA-seq data processing

Raw scRNA-seq data (10x Genomics) were processed using Cell Ranger (v6.1.2) [18] with default parameters. The cellranger count pipeline was used to perform read alignment to the human reference genome (GRCh38), barcode and UMI quantification, and generation of the gene-by-cell expression matrix. Ambient mRNA contamination was removed using SoupX v1.6.1 [19], leveraging marker genes (INS, GCG, SST, TTR, IAPP, PYY, KRT19, and TPH1) representative of the major cell types identified in the initial clustering of each sample. The corrected gene count matrix was then imported into Seurat v4.3.0 [20] for downstream quality control. Cells were retained based on the following criteria: (1) genes expressed in at least three cells, and each cell expressing at least 200 genes; (2) potential doublets identified using scDblFinder v1.12.0 [21] were removed; (3) cells with fewer than 9000 detected genes and fewer than 10 000 total UMI counts were kept; and (4) cells with >5% mitochondrial gene expression were excluded. Additional information is provided in Supplementary Methods.

Bulk RNA-seq analysis processing

Raw FASTQ sequencing reads for the FACS-sorted alpha-cell and beta-cell samples from human islets were downloaded from the Human Pancreas Analysis Program (HPAP) data portal (https://hpap.pmacs.upenn.edu/) (Table S1). All samples underwent quality control using fastp v0.19.6 [22] and were subsequently analysed as described in Supplementary Methods.

Gene set enrichment analysis (GSEA)

Gene set enrichment analysis (GSEA) was performed using the fgsea v1.20.0 package [23] in R. The ranked gene list was generated based on the log2 fold change from differential expression analysis of both single-cell and bulk RNA-seq data. Gene sets of Reactome were downloaded from the Molecular Signatures Database (MSigDB v7.2, https://data.broadinstitute.org/gsea-msigdb/msigdb/release/7.2/) in GMT format. For each analysis, enrichment scores and associated statistics were calculated using the fgsea() function with the following parameters: minSize = 15, maxSize = 500, and nperm =  10 000. P values were adjusted using the Benjamini-Hochberg method. Pathways with adjusted P < 0.05 were considered significantly enriched. A custom R script was used to visualize the top 15 significantly up- and downregulated pathways in each analysis, with pathway ranking based on normalized enrichment score (NES) and plotting implemented using ggplot2. If less than 15 pathways were significantly enriched (adjusted P-value < 0.05) in either direction, all significant pathways were displayed.

Rank-rank hypergeometric overlap (RRHO) analysis

To compare the global transcriptomic similarity between HPAP-young (type 1 diabetes control) and HPAP-older (type 2 diabetes control) in the context of alpha-cell versus beta-cell differences, we performed Rank-Rank Hypergeometric Overlap (RRHO) analysis. We applied an enhanced RRHO pipeline, RedRibbon [24], which offers improved speed and accuracy compared to the original RRHO method [25], enabling efficient identification of shared up- or downregulated genes between two independent datasets [24]. Detailed information on the procedure is provided in Supplementary Methods.

EndoC-βH1 cell line culture and transduction

The human pancreatic beta cell line EndoC-βH1 [16] (RRID:CVCL_L909) was kindly provided by R Scharfmann, Institut Cochin, University of Paris, France, and is available commercially (https://www.humancelldesign.com/endoc-bh1/). According to the vendor: “This high-quality human pancreatic beta cell line displays a homogeneity that exceeds 99% and an insulin content of 0.5 to 1 µg/million cells. Functionality of EndoC-βH1 cells has been validated using glucose stimulated insulin secretion (GSIS) assay”. The presence of mycoplasma infection is regularly controlled using the MycoAlert Mycoplasma Detection kit (Lonza). Cells are periodically refreshed from frozen stocks to maintain a low passage number. EndoC-βH1 cells were cultured in Dulbecco’s modified Eagle medium with 5.6 mM glucose (Gibco, Thermo Fisher Scientific) as previously described [26]. All the results shown for EndoC-βH1 cells refer to independent biological samples, i.e., different passages. MEG3 gene expression was silenced using 3 different siRNAs [Thermo Fisher Scientific]; siMEG3 #1 (siRNA n272552): [5ʹ-GUGUUCACCUGCUAGCAAAtt-3ʹ] and siMEG3 #2 (siRNA n272559): [5ʹ-UCUUAUUUAUUCUCCAACAtt-3], targeting 15 different transcript variants at exon 3 (siMEG3 #1 at position 460; siMEG3 #2 at position 883), and siMEG3 #3 (siRNA n503814): [5′-GAACUGCGGAUGGAAGCUGtt-3′], targeting 15 transcript variants at exons 4 to 7 (positions ranging between 978 and 1342). AllStars Negative Control siRNA (siCTRL) (Qiagen) was used as a negative control; the siRNA control does not interfere with beta-cell gene expression, function or viability [27]. Cells were transfected using Lipofectamine RNAiMAX lipid reagent (Invitrogen, Life Technologies) in Opti-MEM (Gibco, Thermo Fisher Scientific) with 30 nM of each siRNA, overnight as previously described [10]. After transfection, cells were kept in culture for a 24 h recovery period and subsequently exposed or not to human IFN-α (2000 U/ml; PeproTech), IFN-γ (1000 U/ml; PeproTech) alone or together with IL-1β (50 U/ml; PeproTech) for 24 h, 48 h or 6 days. In some experiments, cells were also exposed to TNF-α (1000 U/ml; PeproTech). Cells were also exposed to the ER stressor Thapsigargin (1 μM) for 48 h. These conditions are based on our previously published time- and dose-response experiments [7, 10, 28].

Human islets microtissues (hIsMTs) culture and transduction

Human islet microtissues (hIsMTs) (InSphero, MT-04-002-01-60) were generated after dispersion and reaggregation of primary human islets from a single donor (UNOS ID #ALIR128; male; 41 y.o.; Hispanic; BMI 27.77; positive serologies: CMV IgG, EBV IgG; HbA1c: 5.3%; cause of death: head trauma). hIsMTs were transduced with AAVs encoding specific shRNA targeting MEG3 (siMEG3 #1 sequence) from Vector Biolabs during aggregation for 5 days in Akura hanging drop plates (InSphero, CS-PF24) and then released into Akura 96-well plates (InSphero, CS-09-001-00). hIsMTs were cultured in standard culture medium (InSphero, CS-07-005-01) or standard culture medium containing cytokines: TNF-α (25 ng/mL, R&D Systems, 10291-TA-050); IL-1β (5 ng/mL, Sigma, H6291-10UG); IFN-γ (25 ng/mL, R&D Systems, 10067-IF-100) and IFN-α (10 ng/mL, Peprotech, 300-02AA) starting at day 9 after reaggregation, for 6 days. Medium was exchanged and cytokines were re-dosed every 2–3 days throughout the experiment.

Cell viability assessments

EndoC-βH1 cell viability was determined by fluorescence microscopy using the nuclear dyes propidium iodide (10 μg/ml; Sigma-Aldrich) and Hoechst 33342 (10 μg/ml; Sigma-Aldrich), as previously described [27, 29]. A minimum of 500 cells was counted per condition. Viability was evaluated by two independent observers, one of them being unaware of the sample identity. The agreement between the two observers was >90%. The results are expressed as percentage of apoptosis, calculated as the number of apoptotic cells/total number of cells and presented as fold change compared to siCTRL cytokines-treated cells. In EndoC-βH1 and hIsMTs, Caspase 3/7 activity was assessed using the Caspase-Glo 3/7 Assay System (Promega, G8091). Caspase reagent was added to wells containing cells or just culture medium (used as background control) and total Caspase 3/7 activity was evaluated by luminescence following manufacturer’s instructions.

Quantitative real-time PCR, protein extraction and Western blot analysis

EndoC-βH1 were washed with PBS, detached and collected in lysis buffer. Poly(A) + RNA was isolated using the Dynabeads mRNA DIRECT Kit (Invitrogen) and reverse-transcribed using the Reverse Transcriptase Core Kit (Eurogentec), following the manufacturer’s protocol. hIsMTs were pooled (6 per well), washed twice in PBS and transferred to PCR-clean 2 mL Eppendorf tubes. QIAzol lysis buffer (Qiagen, 79306) was added and hIsMTs were vortexed and sonicated to aid thorough lysis. Additional details and primers composition are provided in Supplementary Methods.

Total protein was extracted using Laemmli buffer supplemented with phosphatase and protease inhibitors (Roche) and separated on 10% SDS–PAGE. The nitrocellulose membranes were probed using specific primary antibodies: Phospho-Stat1 (Tyr701) (58D6) Rabbit mAb #9167 and Stat1 (9H2) Mouse mAb #9176 Cell Signaling Technology, or Human/Mouse/Rat GAPDH Antibody #2275-PC-100 R&D Systems diluted 1:1000 in TBST (TBS, 0.1% Tween 20) with 5% BSA (GAPDH was diluted 1:10 000). The densitometric values were quantified by ImageLab software version 6.1 (Bio-Rad Laboratories, RRID:SCR_014210) and normalized to GAPDH or the respective total protein form, after background subtraction. Full length uncropped original western blots used are provided in the Supplementary Material.

Functional analysis and measurements of Insulin / CXCL10 secretion by ELISA

Glucose-stimulated insulin secretion was performed in hIsMTs sequentially in KRHB containing 2.8 mM and 16.7 mM glucose for 2 h. ELISA (Human Insulin, Mercodia, Uppsala, Sweden for EndoC-βH1 cells, STELLUX® Chemi Human Insulin ELISA (Alpco, 80-INSHU-CH01) for hIsMTs or Human CXCL10/IP-10 Immunoassay, Quantikine ELISA kit, R&D Systems for both cell models) was used to quantify insulin or CXCL10 secretion to the culture medium (by 50 000 EndoC-βH1 cells / 6 hIsMTs / 200 μl of culture medium) as indicated in the figure legends. The cell culture supernatants were collected at the end of the experiments and immediately stored at −80 °C until samples’ processing. For insulin secretion, all samples were diluted ten-fold. For CXCL10 secretion, samples under non-treated conditions were diluted two-fold and samples exposed to cytokines were diluted 70-fold. The assay procedures and calculation of the results were conducted following the manufacturer’s recommendation.

3D staining and imaging

hIsMTs were washed twice in PBS (with Mg2+Ca2+), fixed for 15 min in 4% PFA, washed twice more with PBS and kept in PBS with 0.05% sodium azide until staining and imaging as described in Supplementary Methods.

Quantification and statistical analysis

Data are presented as mean ± Standard Error of the Mean (SEM) and analysed by paired Student’s t test or one-way ANOVA followed by Bonferroni post hoc correction for multiple comparisons, as indicated in the figure legends. P values less than or equal to 0.05 were considered statistically significant. Normality was evaluated by Shapiro-Wilk test. Sample size is displayed in each graph and it was calculated based on the variability observed in preliminary experiments and power analysis using GPower 3.1 software. Grubbs method was used to identify data outliers. Statistical analyses were performed using GraphPad Prism software version 8 (GraphPad Software, La Jolla, USA). Statistics analysis for single-cell and bulk RNA-seq data are described in the corresponding sections above.