Gene expression alterations in salivary gland epithelia of Sjögren’s syndrome patients are associated with clinical and histopathological manifestations

Sjögren’s syndrome (SS) is a complex autoimmune disease associated with lymphocytic infiltration and secretory dysfunction of salivary and lacrimal glands. Although the etiology of SS remains unclear, evidence suggests that epithelial damage of the glands elicits immune and fibrotic responses in SS. To define molecular changes underlying epithelial tissue damage in SS, we laser capture microdissected (LCM) labial salivary gland epithelia from 8 SS and 8 non-SS controls for analysis by RNA sequencing (RNAseq). Computational interrogation of gene expression signatures revealed that, in addition to a division of SS and non-SS samples, there was a potential intermediate state overlapping clustering of SS and non-SS samples. Differential expression analysis uncovered signaling events likely associated with distinct SS pathogenesis. Notable signals included the enrichment of IFN-γ and JAK/STAT-regulated genes, and the induction of genes encoding secreted factors, such as LTF, BMP3, and MMP7, implicated in immune responses, matrix remodeling and tissue destruction. Identification of gene expression signatures of salivary epithelia associated with mixed clinical and histopathological characteristics suggests that SS pathology may be defined by distinct molecular subtypes. We conclude that gene expression changes arising in the damaged salivary epithelia may offer novel insights into the signals contributing to SS development and progression.


Scientific Reports
| (2021) 11:11154 | https://doi.org/10.1038/s41598-021-90569-w www.nature.com/scientificreports/ expression microarray studies on labial salivary glands and peripheral blood detected dysregulation of type 1 interferon (IFN)-inducible genes 14 . Additional upregulated genes involved in antigen presentation, lymphocyte development and activation, as well as interferon-induced chemokines were also identified. Interestingly, IFNγ-inducible genes, such as Signal transducer and activator of transcription (STAT) family members STAT1 and STAT3, were also shown to be highly expressed in primary SS (pSS) 15 . Furthermore, IFN-mediated innate immune mechanisms have been implicated in the pathogenesis of pSS, where specifically 23 genes known to play a role in IFN signaling were identified 16 . Two GWAS studies have been performed in pSS and among the SS-associated non-HLA genes they identified STAT4 and IRF5 encoding transcription factors, BLK coding for B cell kinase as well as genes encoding the IL-12A cytokine. Genes involved in NF-κB signaling and CXCR5 chemokine production were also shown to be upregulated 9,17 . Collectively, these studies provided a strong support for the role of innate immunity, as well as adaptive immune mechanisms, in the pathogenesis of pSS. Despite these advances, treatment strategies for SS disease are almost inexistent, and molecular mechanisms that drive disease onset are poorly understood. One consistent feature of human autoimmunity is the finding that there is an initial phase in which patients have many of the serological abnormalities of autoimmune disease before they develop clinical pathology 18 . In this respect, the immunobiology of the glandular epithelial cells within the target tissue likely plays a critical role in initiating events in SS. Indeed, accumulating evidence suggests that impaired glandular secretion develops independent of lymphocytic infiltration, and defective secretion precedes the immune response during SS onset, suggesting a role for the epithelium in the etiology of SS 2 . For example, aberrantly high expression of BMP growth factors, such as BMP6, has been implicated in the induction of salivary and lacrimal gland dysfunction independent of the autoantibodies and immune infiltration associated with SS 19 . Recent studies, including our own, have also shown that the structural integrity of the epithelium of SS patients is lost, and that these defects are associated with aberrant apical-basal polarity [20][21][22][23][24][25] . Aberrant polarity can arise in various ways, including environmental toxins, viral infection, physical damage, and unbalanced hormone levels, and therefore it may be an early event that triggers disease states, including SS. Since structural changes of the epithelium of salivary glands may drive inflammation, the aim of this study was to elucidate molecular details of how epithelial structural dysfunction contributes to SS. We used laser capture microdissection of glandular epithelial regions (acini and ducts) from resected labial salivary glands of patients fulfilling or not fulfilling SS criteria (8 SS and 8 non-SS controls) and processed them for global transcriptomes using bulk RNAseq analysis. Our studies provide evidence that damage to salivary gland epithelia is associated with changes in gene expression signatures that provide insights into deregulated pathways and cellular processes which may underlie predisposition to, and early onset of, SS. Collectively, these findings provide support for the theory that alterations in gene expression arising in the salivary epithelia contribute to the etiology of SS. We postulate that the molecular sub-groups reflect varying stages of disease predisposition and that they may offer novel insights into the signals contributing to the progression of SS.

Results
Clinical characteristics of patients. Sixteen samples of labial salivary glands were obtained from female patients who presented with sicca symptoms to the Department of Otolayngology/Head and Neck Surgery at Haukeland University Hospital, Bergen. Patients were assessed for primary Sjögren's Syndrome (pSS) based on the American-European Classification Criteria 26 , where diagnosis of pSS required focal lymphocytic sialadenitis score (FS) of ≥ 1 or presence of anti-SSA/Ro or anti-SSB/La autoantibodies, as well as ocular symptoms of dry eyes, oral symptoms of dry mouth, positive Schirmer's test, and reduced unstimulated salivary flow ≤ 1.5 mL in 15 min (at least four positive tests including symptoms, or three objective tests). We categorized samples into two clinically distinct groups, patients diagnosed with SS (named SS, n = 8) or without SS (named non-SS, n = 8). Clinical characteristics of patients are shown in Table 1. The patient mean age was 53 ± 12.4 years, with three patients in the SS group presenting with anti-SSA/Ro and ANA antibodies, of which also two were positive for anti-SSB/La antibodies. All patients reported dry eyes and the majority reported dry mouth. Three SS patients had FS of 1, three had a FS of 2, with one patient presenting with FS of 3 and one with FS of 4. Interestingly, FS of either 1 or 4 were associated with positivity for all anti-SSA/Ro, anti-SSB/La and ANA, suggesting that these may be independent variables in our sample group. Likewise, there was no association between germinal centers (GCs) and autoantibody positivity. Furthermore, two of the non-SS patients had a negligible salivary flow rate, and there was little correlation between salivary flow rate and auto-antibody positivity in SS patients.
Histopathological evaluation of H&E sections of biopsies from non-SS and SS patients revealed that SS samples were characterized by acinar atrophy, ductal dilation, detectable fat infiltration, and lymphocytic aggregates with FS ≥ 1, whereas non-SS tissues were generally normal in appearance with some acinar atrophy, ductal dilation, and scattered inflammation (Fig. 1). From analysis of clinical presentation, we observed that there was no correlation between SS and non-SS groups for the parameters of positive Schirmer's test, saliva secretion, dry mouth or dry eyes. However, the presence of anti-SSA/Ro, anti-SSB/La, ANA autoantibodies, FS ≥ 1, and the presence of CD21+ GCs were only observed in the SS group. When assessed for age, the SS group was significantly older than the non-SS group (60.2 ± 10.7 years vs 45.6 ± 9.6 years respectively; p = 0.012).
Subclustering of epithelial-derived RNA-seq data reveals potential SS subgroups. To gain insight into the molecular characteristics of two groups with sicca symptoms, eight non-SS and eight SS, all samples were subjected to laser capture microdissection (LCM) of epithelial regions followed by RNAseq analysis to identify biomarkers restricted to the epithelial compartment of the labial salivary glands. We evaluated the entire set by principal component analysis (PCA) selected by the expression level of all differentially expressed (DE) genes. It was evident that the SS and non-SS samples did not entirely cluster together and rather three distinct transcriptomic subclusters were observed: a group containing non-SS samples (named Group A, n = 3), a group , and a group containing samples from both groups (named Group AB, n = 10) (Fig. 2a,b). The position of these groupings demonstrates different features of RNA expression in the epithelial tissues, suggesting that transcriptomic profiling of patient samples may characterize important molecular events in SS pathogenesis compared to grouping by clinical characteristics alone.
When assessed for age, Group B (69.3 ± 12.7 years) was significantly older than both the Group AB (52.2 ± 7.4 years; p = 0.020) and Group A (39.0 ± 6.6 years respectively; p = 0.002), whereas there was no significant age difference between Group A and Group AB. Clinical findings for the RNA-Seq groupings revealed that all samples in Group B had GCs with FS of ≥ 2. One of the patients in this group developed lymphoma during the course of this study. In the mixed Group AB, two of the five SS patients also had GCs and lower mean FS. There was no correlation observed for RNAseq grouping for the presence of anti-SSA/Ro, anti-SSB/La and ANA autoantibodies, salivary flow and Schirmer test. Likewise, no specific pattern regarding Group A patients was detected, except that all presented with positive Schirmer tests and a subjective feeling of dry eyes. Differential expression testing of SS subgroups. Genome-wide differential expression testing was performed using DESeq2 for multiple group comparisons. First, we tested the clinical grouping of SS against non-SS. Second, we compared the RNAseq grouping of Group B against Group A samples to determine if there were differences in gene expression between clinical diagnosis of Sjögren's and principal component analysis (PCA) for clustering. Third, the RNAseq Group A and Group B groups were compared to the Group AB.
RNAseq analysis of the epithelial regions for the groupings clinically diagnosed with SS showed no significantly altered differentially expressed (DE) genes compared to the non-SS group out of 25,444 total genes analyzed (Fig. 2c). However, when analyzed according to PCA clustering, there were 1247 DE genes between Group B and Group A at a FDR cut-off of 0.05 (Fig. 2c), further demonstrating a benefit of clustering the samples according to RNAseq data. Table 2 summarizes the top 10 DE upregulated and downregulated genes for both clinical Sjögren's (SS vs non-SS) and RNAseq (Group B vs Group A vs Group AB) groupings.
Inter-group analysis with the Group AB subcluster revealed the following comparisons: 767 DE genes compared to the Group A (FDR cut-off = 0.05) and 1480 DE genes compared to Group B (FDR cut-off = 0.05). While strictly not statistically meaningful, the FDR-value ranking was used to define signatures to be annotated by pathway enrichment analysis.
Pathway analysis identifies epithelial signals associated with SS. Gene set enrichment analysis (GSEA) was conducted for each of the groups, and hallmark gene sets were evaluated (Fig. 3a). Comparison of SS vs non-SS groupings revealed 8 signaling pathways that were statistically upregulated and 27 pathways downregulated in the SS group (Supplementary Table S1). Notably, immunological and autoimmune disease pathways were enriched in the SS group while cholesterol and protein-related metabolic pathways were downregulated (Fig. 3b,c). Interferon response genes were upregulated in both the clinical SS and RNAseq (Group B) groupings compared to non-SS and Group A samples, respectively, suggesting conservation of relevant processes reflected in RNA profiling. Further hallmark GSEA of comparisons between the subcluster groupings revealed 9 upregu- Table 1. Clinical characteristics of patients whose salivary gland biopsies were analyzed by RNAseq. ANA antinuclear antibodies; Saliva secretion (unstimulated flow); normal flow > 1.5 ml/15 min; Schirmer test positive when reduced tear secretion (objective test); Dry mouth and dry eyes, subjective symptoms from the patient, FS Focus Score, NT not determined. Germinal center-like structures may result from the deregulated proliferation and infiltration of lymphocytes in glandular tissue.   (Fig. 3a). Compared to Group AB, both Group A and Group B showed upregulated pathways in genes encoding components of the apical junction complex and apical surface of epithelial cells, epithelial-mesenchymal transition and coagulation, all indicative of changes in epithelial cell polarity. Significantly, the Group AB showed intermediary expression of genes involved in protein secretion and those activated through the mTORC1 complex as enrichment of genes in these pathways was low in Group B and high in Group A.
While the RNAseq data generated from the epithelia of non-SS and SS salivary glands showed complexity in gene signatures among patient samples, a closer comparison revealed genes that encode secreted proteins that may contribute to the disease. These included genes encoding lactotransferrin (LTF), matrix metalloproteinase 7 (MMP7), bone morphogenetic protein 3 (BMP3), interferon-inducible proteoglycan testican-2 (SPOCK2) as well as cytokines IL-15, IL-19 and BAFF (TNFSF13B). Both BAFF and IL-15 have been shown to be associated with B cell development, physiology and malignancy 13,27 . Furthermore, SPOCK2 was also identified as a downstream target of lncRNA-dependent miRNAs involved in the modulation of ECM genes in SS 13 .
RNAscope analysis reveals increased LTF, MMP7, and BMP3 expression in Sjögren's labial salivary gland epithelia. Among the top genes whose transcriptional expression was increased in the RNAseq analysis of the SS grouping were LTF and in Group B were MMP7 and BMP3. LTF is a glycoprotein found in exocrine fluids, such as milk, saliva, tears, and nasal secretions, as well as in secondary granules of polymorphonuclear leukocytes (PMNs) and was previously shown to be upregulated on the protein level in saliva of SS patients 28,29 . MMP7 belongs to a family of ubiquitous proteolytic enzymes, the MMPs, and has been shown to be www.nature.com/scientificreports/ elevated in corneal specimens of pSS patients 30 . BMP3 is a member of the TGF-β superfamily of proteins, associated with inflammation in rheumatoid arthritis patients with concomitant inhibition of bone repair 31 .
In order to validate the increased LTF, MMP7, and BMP3 expression identified in our RNAseq data, we performed independent RNAscope analyses for expression of these genes. This in situ hybridization approach is based on the amplification of target-specific signals but not non-specific background hybridization. Samples were counterstained with pan-cytokeratin in order to restrict analysis of RNA to epithelial regions and analyzed using the CellProfiler image analysis software 32-34 (http:// www. cellp rofil er. org). To increase our sample size, we included an additional group of labial salivary gland samples obtained from the University of Oslo. Figure 4 illustrates combined RNAscope analysis of LTF, MMP7 and BMP3 expression from the Bergen samples described in Table 1 and Oslo samples described in Table 3. We note that the age of SS patients in the Oslo cohort was not significantly different than the non-SS group (51.0 ± 19.1 years vs 44.5 ± 13.1 years respectively) although the new Oslo cohort included SS samples from two male patients and the SS groups were mostly positive for autoantibodies ( Table 3). Quantitation of the RNA showed significant increases in the epithelial percentage of LTF, MMP3, and BMP3 expression in the SS group compared to the non-SS group (Fig. 4a,b). The observed increases in LTF, MMP3 and BMP3 were not driven by age as there was no significant correlation between patient age and expression of either gene (Fig. 4, "Methods", "Statistical analysis"). These results correlate with our RNAseq findings of increased epithelial expression of LTF, MMP7, and BMP3 in SS.

Discussion
Sjögren's syndrome (SS) is a debilitating autoimmune disease with poorly understood etiology and limited biomarkers available for early detection and prognosis. In this study, we generated and annotated gene expression signatures of labial salivary gland epithelia from human SS and non-SS patients using RNAseq and computational analyses with a goal of obtaining new insights into molecular changes associated with this disease. Our studies www.nature.com/scientificreports/ provide evidence that damage to salivary gland epithelia is associated with changes in gene expression signatures that provide insights into deregulated pathways and cellular processes which may underlie predisposition to and early onset of SS. We postulate that the molecular sub-groups identified in our studies reflect varying stages of disease predisposition and that they may offer novel insights into the signals contributing to the development of SS. PCA analysis showed a clear segregation of gene signatures between some non-SS and SS patients (Fig. 2). In addition, it revealed a novel subgroup of gene signatures derived from both non-SS and SS tissues that may reflect a mixed patient population. Interestingly, however, analyses of clinicopathological features of salivary epithelia from SS and non-SS patient cohorts did not reveal any consistent association with GCs and focus scores; indeed, there were no substantial differences among these patient groups in autoantibody positivity, saliva secretion, tear secretion, and subjective feeling of dry eyes or dry mouth. Our identification of a mixed gene expression www.nature.com/scientificreports/ signature subgroup (Group AB) sheds novel insights into molecular changes in labial salivary glands that may reflect early disease in individuals whose symptoms do not meet the criteria for SS diagnosis. Better understanding of the molecular features that underlie different subsets of patients is essential for disease diagnostics, prognosis, and effective treatment interventions. The focal infiltration of the salivary gland provides the basis for the focus score used to describe the pathological changes of the salivary glands in the currently used classification criteria for pSS. While the European-American Consensus Group criteria 26 include symptoms of dry mouth and dry eyes, the ACR-criteria include only objective criteria of reduced secretion of saliva and tears, ocular staining, positive focus score, and presence of anti-SSA/Ro in serum 35 . Both classification criteria for pSS require that the patients have anti-SSA/Ro antibodies and/or have a positive focus score, but do not include germinal centers. The epithelial gene expression profiles we describe show some association with clinical/pathological features of SS, but such molecular similarities are also evident in non-SS patients. This may be a consequence of the stochastic nature of clinical/pathological characteristics currently used to stratify patients into non-SS and SS groups.
Changes in SS epithelial and stromal tissue biology have been increasingly shown to play a role in SS pathobiology 13,36 . Previous studies revealed alterations in subpopulations of epithelial cells in NOD/ShiLtJ mouse model of SS, revealing a decrease in the tissue area of secretory cells in favor of ductal epithelial cells and expansion of keratin 5-positive basal cells 37 . Additional studies highlighted the key roles of the innate and adaptive immune systems, as well as stromal cells, in promoting the pro-inflammatory milieu 38,39 . A more recent report using an LCM approach coupled with RNAseq analyses described identification of cell-type and disease-specific markers in SS patients 40 . However, the latter study did not strictly focus on salivary gland epithelia, and thus precluded identification of gene signatures associated only with the epithelial compartment. Additionally, recent analyses of PBMCs from SS patients provided evidence for the involvement of the perturbed long non-conding RNAs (lncRNAs)-microRNAs axes in the expression of genes involved in the pathogenesis of SS, including those with roles in B cell function and extracellular matrix remodeling, among others 13 . Here, our results identifying IL-15 and BAFF cytokines with important functions in B cell physiology, as well as SPOCK2 proteoglycan   Table 3. Additional patient samples used for RNAScope (Oslo cohort). ANA antinuclear antibodies; Saliva secretion (unstimulated flow); Schirmer test positive when reduced tear secretion (objective test); Dry mouth and dry eyes, subjective symptoms from the patient, FS Focus Score.

Germinal centers (CD21 + )
Non-SS www.nature.com/scientificreports/ contributing to the increasingly recognized roles of the stromal compartment in the pathogenesis of SS 36 , provide important new information into their epithelial origin in SS.
In addition to signals known to be enriched in SS, including IFN-γ and JAK/STAT-regulated genes, and genes encoding secreted factors implicated in immune responses, our studies identified LTF, also known as lactoferrin (LF), to be upregulated in expression in SS epithelia compared to non-SS tissues. Increased levels of salivary LTF with some relationship to the labial salivary gland histopathology grade were reported previously 28 , while no significant changes in LTF expression were detected in tear secretions in SS 41 . LTF is an iron-binding multifunctional glycoprotein known to function in the prevention of infections, late onset sepsis, and necrotizing enterocolitis 29 . LTF has been shown to stimulate a balanced T-helper-1/T-helper-2 cytokine immune response and to be a part of the innate defense, mainly at mucosa 42 . In our study, the observed increase in LTF gene expression in SS ductal epithelia suggests that it may be a relevant biomarker predictive of early disease with a potential causal effect on the induction of immune response.
Our data also show that MMP7 and BMP3 transcripts are upregulated in labial salivary glands from SS patients. MMP7 is a member of ubiquitous MMPs with important functions in tissue remodeling during development, but whose inappropriate expression has been linked to tissue destruction in autoimmune rheumatoid arthritis 43 and with corneal melting in patients with primary SS 30 . This suggests that upregulation of MMP7 in labial tissues from SS patients is likely to contribute to salivary gland destruction. In contrast, relatively less is known about the role of BMP3 in pSS. Previous studies documented increased expression of BMP6 in salivary and lacrimal glands from human SS patients independent of autoantibodies 19 . To date, increased expression of BMP3 has been linked to immune infiltration in rheumatoid arthritis and inhibition of bone repair 31 . Thus, our RNAseq-based identification of BMP3 in salivary epithelia from SS patients may provide additional mechanistic insights into how inflammation induces BMP3 expression to drive epithelial dysfunction in pSS.
Accumulating evidence suggests a role for the epithelium in the etiology of SS. Impaired epithelial secretion has been reported to develop independently of lymphocytic infiltration and frequently, defective secretion precedes the immune response during SS onset 2 . Epithelial defects are often accompanied by pro-fibrotic changes in the stroma indicating remodeling of the extracellular matrix and fibroblast activation 36 . Here, our findings align alterations in gene expression signatures in the salivary epithelia with different degrees of compromised salivary and lacrimal glands' function. Collectively, these studies offer a molecular explanation for exocrine tissue dysfunction as it progresses to SS and suggest that transcriptional signatures may serve as predictive biomarkers for early and advanced disease.

Methods
Patient selection and sample collection. The study protocol was approved by the Norwegian Regional Committee for Medical and Health Research Ethics of West (REK 2009/686) and South-East Norway (REK 2010/1292-1). The study was performed in compliance with the tenets of the Declaration of Helsinki. Prior to participation in the study, written informed consent was obtained from all participants. The data was deidentifed prior to analysis.
Labial salivary glands were obtained from a total of 27 patients with sicca symptoms. Based on the American-European Consensus Criteria 26 , these patients were sub-divided into 13 pSS (focus score (FS) 1-4) and 14 non-SS sicca controls (FS < 1). Inclusion criteria were: SS-patients fulfilling the AECG classification for pSS; non-SS controls were patients evaluated for SS, but not fulfilling the AECG classification. Importantly, none of the non-SS controls had a positive focus score or positive anti-SSA/Ro or anti-SSB/La autoantibodies. All participants provided informed consent prior to the initiation of this study. Clinical characteristics of patients whose minor salivary gland (SG) biopsies were analyzed for RNAseq are summarized in Table 1 (Bergen cohort) and for RNAscope analyses in Table 3 (Oslo cohort). These included assessment of objective tests of antinuclear antibodies (ANA, anti-SSA/Ro, anti-SSB/La), saliva secretion (unstimulated flow of whole saliva), tear secretion by Schirmer test, FS, and germinal center (GC)-like structures, as well as subjective symptoms noted by the patients (dry mouth and dry eyes).
Histopathology. Paraffin-embedded, formalin-fixed human salivary gland tissue sections, 5-µm-thick, were stained with hematoxylin and eosin (H&E) and focal lymphocytic infiltration areas in H&E slides were captured and assessed using Aperio digital technology (Leica Biosystem). Minor SG FS of SS patients was assessed based on the currently accepted criteria 4 . The sections were evaluated by one oral pathologist (KS) in order to determine their FS. The tissue sections were further assessed for the presence of GC-like structures both on H&E sections and with an additional staining with CD21 (a marker of follicular dendritic cells) to improve the reliability of GC identification 44 . Fresh frozen biopsies from the Bergen cohort characterized in Table 1, were used for RNAseq analyses. Additional minor SG biopsies from another patient cohort (Oslo cohort, Table 3) were used for immunofluorescence labeling and RNAscope analyses.
Laser capture microdissection (LCM) analysis. Fresh surgical specimens were preserved in the Tissue-TEK OCT compound and frozen in isopentane, followed by storage at − 80 °C. Tissue sections (5-μm-thick) were stained with an ethanol-based cresyl violet and eosin stain, followed by examination under the microscope and annotation of epithelium and inflammatory infiltrates on digital images. Acinar and ductal epithelium was subsequently captured by LCM using a Leica LMD6 microdissection system with regions containing inflammatory cells being excluded.
RNAseq analyses. RNA  Processing RNA-sequencing data. Paired-end RNA-sequencing reads for each sample were first trimmed for  Illumina TruSeq adapters (AGA TCG GAA GAG CAC ACG TCT GAA CT CCA GTC AC and AGA TCG GAA GAG  CGT CGT GTA GGG AAA GAG TGT AGA TCG GAA GA GCG TCG TGT AGG GAA AGA GTG TAG ATC TCG GTG  GTC GCC GTA TCATT) using Cutadapt tool 45 , setting the minimum read quality threshold to 20 for the start and end of the reads, and 15 for the minimum read length. Reads were then aligned to a reference hg19 genome using Tophat2 46 , with the number of allowed read mismatches, gap length and edit distance set to 1. Reads were then mapped to gene features using Ensembl hg19 (human genome release 19) genome annotations, and quantified as raw counts per feature using featureCount 47 using default parameters.
RNA differential expression testing. Gene feature raw RNA expression counts were normalized and tested for differential expression between different sample groupings using DESeq2 48 . In our analyses, we did not include age as a variable to be adjusted for when performing differential expression testing, although we are aware that age can be a confounding factor. Pairwise comparisons for the non-SS and SS RNAseq data provided sufficient coverage to estimate the differential expression of moderately to highly expressed genes.
Gene set enrichment analysis. For testing enrichment of pathways, genes were first ranked (genome-wide) based on their signed average log-fold change between group comparisons (e.g. Sjogren's vs control), using the pre-ranked GSEA tool v2.1.1 given a specific differential expression ranked list (based on the 't-statistic' of differential expression). Data were then tested for enrichment of annotated hallmark gene sets from the molecular signature database (MSigDB) 49 using the pre-ranked Gene Set Enrichment Analysis (GSEA) tool 50 . Pathwaybased enrichment testing allowed the identification of coordinated sets of genes, and GSEA computationally determined if a defined set of genes was of statistical significance. We note that we did not use any fold-change cut-off here in addition to an FDR threshold.
RNAscope analyses. Formalin-fixed, paraffin-embedded non-SS and SS patient minor SG biopsies from the Bergen and Oslo patient cohorts were cut into 4-μm-thick sections and mounted onto Superfrost Plus slides. Samples were deparafinized and in situ hybridization was conducted using the RNAscope 2.5 HD Assay Reagent Kit (ACDBio), probing for lactotransferrin Hs-LTF (LTF) RNA, matrix metalloproteinase 7 Hs-MMP7 (MMP7) RNA, and bone morphogenetic protein 3 (BMP3) RNA, Positive Control Probe Hs-PPIB, and Negative Control probe Dap-B according to the manufacturer's instructions (ACDBio). After the final amplification, fast-red chromogenic detection was performed. The epithelial localizations of LTF, MMP7, and BMP3 were determined by counterstaining the tissue sections with pan-cytokeratin protein. Sections were blocked with 5% donkey serum-TBST (DS-TBST) and incubated in primary antibody 1:100 in 5% DS-TBST for 60 min, followed by incubation in 1:500 secondary goat anti-mouse conjugated with Alexa 647, counterstained with DAPI, and mounted in ProLong Gold Antifade and coverslipped. Image processing of the LTF/BMP3/MMP7 and pan-cytokeratin (pCK) signals was carried out on the free, open-source software CellProlifer v4.0.6 across all 2048X2048px images 32-34 (http:// www. cellp rofil er. org). We utilized the ability to mask nuclei and epithelial cells, followed by merging the two masks to define a full epithelial cell. Nuclei outside the epithelial mask was discarded from further analysis. Following the estimation of epithelial cells, we measured the area of pixels of RNAscope probes within the of area covered by the epithelial mask to derive the percentage of area covered by a probe signal over all the epithelial areas (Supplementary Figure 1). This analysis was carried out on all patients on a minimum of 3 images per patient.