GZMK+CD8+ T cells Target a Specific Acinar Cell Type in Sjögren’s Disease

doi:10.21203/rs.3.rs-3601404/v1

Download PDF

Article

GZMK+CD8+ T cells Target a Specific Acinar Cell Type in Sjögren’s Disease

https://doi.org/10.21203/rs.3.rs-3601404/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this older preprint version

Read the latest preprint version →

Sjögren's Disease (SjD) is a systemic autoimmune disease without a clear etiology or effective therapy. Utilizing unbiased single-cell and spatial transcriptomics to analyze human minor salivary glands in health and disease we developed a comprehensive understanding of the cellular landscape of healthy salivary glands and how that landscape changes in SjD patients. We identified novel seromucous acinar cell types and identified a population of PRR4+CST3+WFDC2- seromucous acinar cells that are particularly targeted in SjD. Notably, GZMK+CD8 T cells, enriched in SjD, exhibited a cytotoxic phenotype and were physically associated with immune-engaged epithelial cells in disease. These findings shed light on the immune response's impact on transitioning acinar cells with high levels of secretion and explain the loss of this specific cell population in SjD. This study explores the complex interplay of varied cell types in the salivary glands and their role in the pathology of Sjögren's Disease.

Health sciences/Diseases/Immunological disorders/Autoimmune diseases

Health sciences/Medical research/Translational research

Sjögren's Disease (SjD) is the second most common systemic autoimmune disease (0.1-4% prevalence) with a striking 9:1 female to male preponderance and a dramatically increased lifetime risk of lymphomagenesis¹^,². Despite important advances in understanding the pathogenesis of SjD, neither a unifying etiology nor an efficacious medication that meaningfully improves quality of life or preserves organ function has been established. The inherent clinical and biological heterogeneity of SjD is cited as the main factor that has limited the clinical success of identified targets and developed agents³. SjD is characterized by chronic lymphocytic infiltration of the exocrine organs and reduced secretory function. The associated clinical sequelae include xerostomia (i.e., dry mouth) and xerophthalmia (i.e., dry eyes). Additional organ systems may be involved as salivary-type exocrine organs support hydration, lubrication, and antimicrobial activity for multiple organ systems and tissues (e.g., upper and lower airway [nasal cavity, trachea, bronchus], genitourinary tract [Bartholin’s, endocervical glands], digestion [pancreas]) exemplifying the systemic nature of the disease. SjD patients also experience debilitating extra-glandular systemic effects including fatigue, widespread musculoskeletal pain, and polyarthritis⁴^,⁵. There are no effective and approved therapies for SjD and management remains suboptimal but includes palliative and supportive care for symptoms and immunosuppressants for more severe and systemic symptoms.

Our current understanding of the etiological basis of SjD suggests a complex interaction between genetics, environmental insults, and hormonal factors. Multiple biological pathways, including both the innate and adaptive immune systems, are reported to be involved in SjD. Genome-wide association studies (GWAS)⁶^,⁷ have now identified 22 single nucleotide polymorphisms (SNPs) associated with SjD. Using polygenic risk scores, the genetic risk of SjD divides the SjD population into anti-SSA autoantibody positive patients, with genetic association dominated by human leukocyte antigens (HLA) SNPs, and anti-SSA autoantibody negative patients, whose SNPs are limited in number and have a more muted effect size⁶. How these risk haplotypes translate to cellular function and autoimmune pathogenesis in targeted organs is not clear.

The salivary glands are composed primarily of acinar cells and ducts, with associated stromal tissues and tissue-resident inflammatory cells. Fluid and the majority of protein secretion in the salivary glands initiates in the acinar cells⁸. In humans, there are three sets of major salivary glands: parotid, composed purely of serous acinar cells; submandibular gland, a mixture of serous/seromucous cells and mucous acinar cells; and sublingual, predominantly mucous cells. In addition to these three pairs of major glands, hundreds to thousands of minor (submucosal) glands are distributed throughout the oral cavity and upper aerodigestive tract. Minor salivary gland (salivary glands) biopsies are used to evaluate the histopathological component of SjD classification criteria and serve as biological surrogates for major gland involvement^9,10. The SjD-affected salivary glands exhibits a chronic, mixed, predominantly lymphocytic infiltrate that surrounds ducts and interdigitating with acini in areas not exhibiting non-specific damage (i.e., interstitial fibrosis, ductal dilation, and acinar atrophy)^11–13. The immune cells contributing to this infiltrate are T lymphocytes (CD4 > CD8), B lymphocytes, and plasma cells; other immune cells are also found, although in fewer numbers including natural killer (NK) cells, macrophages, and dendritic cells^10,14,15. In some patients, tertiary lymphoid-like structures with germinal center formation can be found¹⁶. The immunopathology of SjD is modulated by both genetics^17,18 and clinical features of disease (e.g., anti-SSA+¹⁹, Age¹). Advances in genetic, epigenetic, transcriptomic, and proteomic technologies have advanced our understanding of the pathogenesis of SjD. These studies have shown that multiple inflammatory cytokines are involved including: Type I and Type II interferons (IFNs: IFNɑ, IFN𝛃, and IFNγ)^20,21, interleukins (ILs: IL6, IL-7, IL-12, IL-21, and IL-23)^20,22, and cytokines and chemokines (CXCL10, -11, -12)^23–26. Emerging technologies (e.g., single cell and spatial transcriptomics) are unraveling disease-specific cellular networks and signaling pathways, and spotlight effector cells and therapeutic pathways.

To date, most attention in SjD has been towards understanding the B cell and CD4 T cell involvement in SjD. The direct involvement of CD8 T cells in SjD requires further study. CD8 T cells, also known as cytotoxic T cells, are a subset of T cells involved in the immune response to viral infections, cancer, and autoimmune diseases. In Sjögren's disease, CD8 T cells are thought to target the acinar and/or ductal epithelial cells in the salivary gland, leading to tissue damage and glandular dysfunction. However, direct evidence of CD8 T cells exerting cytotoxic effects on epithelial cells is lacking, their phenotypic profile is unknown, and the specific targeted epithelial cells (e.g., a specific type of acinar or ductal cell, or another cell type) or the antigens initiating or driving disease.

Most of what is understood about the pathogenesis of SjD in the salivary glands is based on unbiased bulk approaches (e.g., microarray, RNAseq) or on low-throughput methods using focused panels of target genes or proteins. There is a critical need for an unbiased understanding of the molecular changes in the salivary glands at the cell level within the tissue context. Herein, we leveraged the strengths of scRNA-seq and spatial transcriptomics, ex vivo cellular assays, and high-throughput microscopy methods to achieve a comprehensive molecular mapping of the cellular composition and cell-state changes in human salivary glands in SjD patients. We then show the disease-specific cellular interactions involved in the immune insult in the gland and novel pathways driving the loss of gland function and parenchymal maintenance. This evidence further supports the importance of T cells in the pathogenesis of SjD including colocalization with, and direct effects upon, seromucous acinar cells of the exocrine epithelia.

Single cell and spatial transcriptomics of the minor salivary glands in Sjögren’s Disease.

We collected research biopsies for single cell RNA sequencing from 11 Sjögren’s and 14 nonSjD patients (Supplemental Table 1). After quality control, 94,000 single cells were collected and sequenced (Fig. 1). Six thousand spots of spatial RNA-sequencing were collected from 38 patients: 19 with a diagnosis of Sjögren's (SjD), 14 without a diagnosis of Sjögren's (non-SjD) and 5 healthy volunteers (HV) (Supplemental Table 2). These patients included clinical presentations with variation in age, sex, severity of disease and positivity for autoantibodies (Fig. 1b). A principal component analysis (PCA) of clinical data from patients for which single-cell data was generated clearly discriminated between anti-SSA autoantibody positive (SSA+) and anti-SSA autoantibody negative (SSA-) patients along the second principal component (Fig. 1c). Using the single-cell data to define cell types, 31 unique clusters were identified by Leiden clustering followed by manual annotation based on defined cell populations (Fig. 1e). These included seven unique types of seromucous acinar cells, two sets of mucous cells, and five groups of T cells (Fig. 1e, Figure S1a,b). Some cell types were enriched in disease or in patients that were SSA+ (e.g., CD8 + Exhausted T cells, CD4 + T cells, B cells, and Natural Killer (NK) cells; Figure S1a). These unique and complementary datasets were helpful to establish an unbiased portfolio of alterations in the SjD-affected salivary glands and a foundation for discovering and testing important cell types, cell states, and cell interactions in SjD.

Spatial transcriptomics identifies cell interactions and architectural alterations in SjD.

Important disease and tissue-specific changes in cellular composition and gene expression state changes can be gleaned from single cell approaches in autoimmune diseases. Although high resolution at the single cell level, this comes at the expense of losing important cell interactions in the tissue context which can only be inferred bioinformatically. To determine cell interactions, spatially-encoded RNA-sequencing (spRNAseq) was performed (Fig. 2a). Within-spot cell proportions were estimated using Cell2Location and expression patterns from scRNAseq (Fig. 2b). Structural characteristics of the gland, at a resolution of 55 µm, were compared between disease and health and recapitulate typical architecture of the glands (Fig. 2a). Similar cell types cluster together in both health and diseased glands (Fig. 2a,c). The differences in within-spot cell proportions reveal important disease-specific changes in cellular interactions such as the absence of Tregs with collocated ductal progenitor cells, increased association between immune cell types within lymphocytic foci, periductal fibrosis, and heightened contact between seromucous acinar cell types (Fig. 2c, Figure S2a,b). Scoring cooccurrence for all cell types, and comparing SjD to nonSjD supported these differences, highlighting a notable increase in within-spot cooccurrence of CD8 + exhausted T cells and PRR4 + CST3 + WFDC2- seromucous acinar cells (Fig. 2d).

Cells grouped by index expression suggest two paths to dysfunction.

To understand disease- and cell-type-specific pathway utilization, we hypothesized that scRNAseq data could be clustered based on indices of gene expression. Expression data was converted from log-normalized counts to Z-scores and summed across annotated indices to form S-scores. Subsequently, new clusters were identified, each a composite of previously named cell types (Fig. 3a-c). Given these clusters encompassed multiple cell types and were characterized by a regulatory program or shared function — exhibiting contractility, or undergoing cell division — they remained unnamed (Fig. 3c). Like the clusters based on expression (Fig. 1e-g, Figure S1a), clusters based on index scores were enriched for clinical features of interest (Fig. 3c-f, Figure S3a,b). Cells from nonSjD and SSA- patients’ cells shared expression of indices, indicating their lack of dissimilarity —sensory taste perception and pathological glycosylation. Conversely, cells from SSA + patients were associated with many enriched pathways involving effector immune cell function, antigen processing and presentation, and IFN signaling (Figure S3a,b). To illustrate this more clearly, a 21-gene Type I interferon-stimulated genes score was constructed and found to be higher in every annotated cell type in SSA + patients, irrespective of SjD status (Figure S3c,d), and shown to be highest in endothelial cells, macrophages, dendritic cells, and ducts. Calculating a ratio of cells in each clinical feature, two index clusters have notable enrichments with disease or SSA+ (Fig. 3e,f). These cells were almost exclusively identified previously as either T-lymphocytes or antigen-presenting cells (i.e., B-lymphocytes, dendritic cells, or macrophages; Fig. 1e, Fig. 3b-e, Fig. 4a). Comparing the indices directly, indices associated with T cell activation scored higher in the SSA + cluster; no indices were significantly enriched in the Sjögren's-diagnosed cluster (Figure S3b).

CD8-GZMK T cells are enriched in SjD and exhibit effector functionality.

The data above highlight important cell types in SjD pathogenesis. The vast majority of these cells are found within dense collections of inflammatory cells termed “foci” (Fig. 1b, Fig. 4a). SjD immune foci are characterized by predominantly lymphocytic infiltrates composed of T- and B-lymphocytes. To visualize and quantify the expression, measure the function of T cells infiltrating SjD salivary glands, and to understand how these cells interact and organize across the gland, we analyzed spRNAseq, scRNAseq, and flow cytometry data. T cell subclusters are annotated based on their expression profiles into CD4 + T cells (CD40L), Regulatory T cells (FOXP3), Progenitor T cells (MKI67, CD25), Exhausted CD8-T cells (Tex; CD8A, GZMK, PDCD1), Cytotoxic CD8-T cells (Tc; CD8A, GZMB), and NK cells (NKG7, GNLY) (Fig. 4b). The most disease-enriched immune cell type in SjD salivary glands are Tex (4.75-fold enriched Tex, q = 0.0009; Fig. 4b,c; Fig. 1e; Figure S1a). Focusing on this population we wanted to understand their role in salivary gland pathology in SjD. Clustering by disease feature shows enrichment of T cell subtypes in SjD and anti-SSA + specific clustering along the periphery of the individual subclusters (Fig. 4c). Looking at gene expression in T cells, GZMA and GZMK mRNA were higher in both SSA- SjD and SSA + SjD, while GNLY and GZMB exhibited lower in expression in the same patients, indicative of the shift towards more GZMK + Tex in the salivary glands (Fig. 4d). Tc and Tex displayed similar pathway activation, with greater activation in SSA + patients (Fig. 4d,e). Although GZMB + Tc cells were increased in number in SjD, they are fewer than Tex and not significantly increased; however, they do display similar pathway enrichment (as Tex).

To directly test the potential cytotoxicity of T cells infiltrating the salivary glands, we isolated infiltrating immune cells from patients’ salivary glands (nonSjD: n = 5 patients; SjD: n = 6 patients) and performed an ex vivo T-lymphocyte cytotoxicity assay to measure intracellular cytokine expression and surface expression of CD107a. Expectedly, SjD salivary glands show marginally increased proportions of T cells in the glands (p = 0.08; two-sided t-test), with increased CD4 + T cells (1.6-fold, p < 0.05; two-sided t-test) and decreased CD8 + T cells (0.8-fold, p < 0.01; two-sided t-test) ratios (Fig. 4f). A high proportion of CD8 + T cells from salivary glands express IFN𝛄; however, SjD glands show increased proportions of CD8 + T cells positive for surface CD107a in response to PMA/ionomycin stimulation, indicative of enhanced degranulation (e.g., cytotoxic potential; 2.1-fold, p < 0.001; two-sided t-test; Fig. 4f). While it was not possible to separate Tc and Tex subsets from human salivary glands for this assay, it supports that both Tc and Tex exhibit enhanced effector function in the glands of SjD patients. Based on this and our data above, this suggests that Tex, like Tc, can exert effector function in the SjD salivary gland. This raised important questions including regarding the arrangement of GZMK + CD8-Tex cells in the gland and their role in SjD pathogenesis?

Spatial phenotyping demonstrates Sjögren’s-Specific Cellular Interactions.

To answer this question, we used multiplex spatial phenotyping to map the positions and interactions of cells and to provide in the salivary glands of SjD and healthy patients. Representative salivary gland FFPE sections containing 2–5 glands from SjD (n = 5) and nonSjD (n = 5) patients were analyzed using 35-plex spatial phenotyping platform to identify 16 cell types: Th, Treg, Tc, NK, DC/mac, LEC, VEC, progenitor VEC, myofibroblasts, fibroblasts, myoepithelial, ducts, acini (Figure S4a, Supplemental Methods), although this initial panel was unable to discriminate all cell types resolved in scRNAseq (Fig. 1e; e.g., clear separation of seromucous from mucous acinar cells). Spatial phenotyping corroborates the general scRNAseq-based cellular proportion changes including: decreased numbers of acini and myoepithelial cells, greater proportions of ducts and immune cells including statistically more Tc (2-fold; t-test, p < 0.05) and Th (2-fold; two-sided t-test, p < 0.05), and marginally more Treg (5-fold; two-sided t-test, p = 0.052) (Figure S4b). Th and Tc cells in SjD express higher CD45, IFNγ, HLA-DR, HLA-A, CD107a, and reduced expression of GZMB and immune checkpoint markers, PD-1 and IDO1 (Figure S4c). SjD ducts, acini, and myoepithelial cells show greater HLA-A and HLA-DR expression. Importantly, SjD acini show altered expression of acinar markers: increased CD66 (greater mucous than seromucous expression), lower Galectin-3 (greater seromucous than mucous expression), and decreased Ki67 (Figure S4c). To directly visualize the arrangement of immune cells in the salivary glands of SjD, we constructed spatial plots of the identified cell types. Consistent with above, SjD glands exhibit reduced acinar density and dense immune cell clusters (i.e., foci), circumscribing ducts and inter-acinar spaces (Fig. 4g). Th, Tc, B cells congregate around epithelial structures with increased expression of HLA-A at the immune interface (Fig. 4g). Using neighborhood analysis, SjD increases the interactions between multiple immune cell types (e.g., B, DC, Treg, Tc; Fig. 4g, Figure S4d). At the immune-epithelial interface, acini show enhanced interactions with Treg, Th, DC, and Tc; ducts exhibit increases in B, Treg, Th, DC, and Tc; and myoepithelial cells show increased interactions with DC, TC, Th, but loss of interactions with other cell types (e.g., acini, ducts, stromal, and endothelial cells) (Fig. 4g, Figure S4d) illustrating the complexity of cell interaction changes in SjD.

Concentrating on T cell interactions in SjD, we performed receptor-ligand analysis of our scRNAseq data (Fig. 1c) using CellChat. Examining interactions that were found in anti-SSA + SjD but absent in anti-SSA- SjD show “MHC-I signaling” as the most dominant inferred interaction in SSA + SjD salivary glands (Fig. 4h). Additional SSA + SjD-specific receptor-ligand interactions included “CD22” and “BAFF” (Figure S4e), as well as “GRN”, “GAS”, “SELL”, “KIT”, “ICAM”, and “ADGRE5”. These results support that autoantibody-positive disease exhibits signaling centered on T cell/B cell signaling pathways (costimulation and maturation of B cells [CD22, BAFF], the Tc immunological synapse [GRN], and gamma-associated sequences [GAS], as well as, pathways involved in signaling infiltration of lymphocytes from the blood to peripheral organs [SELL, ICAM] (Fig. 4g, Figure S4e,f). These results highlight the importance of effector (Tc) interaction with target cells (ducts and acini) and further support the biological distinction between SSA + from SSA- SjD in the glandular pathology. What remains to be answered is what is the potential effector function of these SjD-enriched GZMK-Tex cells in the salivary glands?

GZMK Drives Type I IFN Signaling in Salivary Epithelial Cells.

To answer this question, we narrowed our analyses on those cells that were most discriminant in SjD infiltrates from health, the presence of GZMK + CD8-Tex. Yet, very little is known about how these cells exert an effector function or contribute to SjD pathology. Effector engagement with target cells occurs at the immunological synapse whereby a regulated exocytosis of cytotoxic granules, containing perforin-1 and granzymes, induce target cell effects including apoptosis. Unlike Tc cells that secrete cytotoxic granules containing GZMB, other granzymes (e.g., GZMA, GZMK) have unique targets and downstream functions. Specifically, GZMK does not directly induce caspase-dependent apoptosis²⁷. Alternately, cytosolic GZMK has been shown to cleave BID in the mitochondrial membrane to drive the generation of reactive oxygen species (ROS) and potentially, caspase-independent cell death²⁷. We hypothesize that in addition to the potential generation of ROS, cytosolic GZMK may affect mitochondrial membrane permeability promoting the relocalization of mitochondrial nucleic acids to the cytosol and the initiation of the pathogen recognition response (PRR)^28,29. In the context of autoimmunity, epithelial cells expressing autoantigenic peptides on their MHC-I (Fig. 4, Figure S4) may initiate effector cell engagement of Tex cells to drive subacute effects on inflammatory signaling thereby explaining 1) minimally increased apoptosis in the glands (i.e., the immunopathology does not seem to be directly Tc), and 2) the canonical increased Type I IFN signature in SjD-affected ducts and acini.

To answer this, we then confirmed the location and proportion of GZMK + and GZMB + CD8-T cells in the peri-epithelial immune infiltrates using immunofluorescence microscopy. GZMB + CD8-T cells were fewer in number in the infiltrates and more distant from epithelial structures in SjD (Fig. 4h, Fig. 5a); GZMK + CD8-T cells were more numerous and reside in close proximity to, or within, involved ducts and acini (Fig. 5a). To model the effector activity of cytotoxic granules containing GZMK from Tex on targeted epithelia in SjD, we transfected GZMB, GZMK into the cytosol of human acinar cells (Figure S5a). Transfection of GZMK did not induce apoptosis by staining with Annexin V/PI flow cytometry (Figure S5a); transfection of GZMB increased the amount of dead/floating cells in the supernatant, but not significant differences by flow cytometry. Western immunoblot analysis of GZMK shows the dose-dependent differences in GZMK protein levels in HSG cells (Figure SbB), and shows the time-dependent presence of cytotoslic GMZK and GZMB protein in NS-SV-TTAC cells out to 48 hours (Figure S5c). Then, using immunofluorescence microscopy, GZMK transfection, but not GZMB, results in mtDNA relocalization (loss of collocation with MitoTracker) to the cytosol (2.6-fold, p < 0.0001; ANOVA; Fig. 5b) in acinar cells. This resulted in an increase in pIRF3/IRF3 ratio suggesting IRF Pathway activation and Type I IFN signaling (Figure S5d). The ability for GZMK to drive cytosolic mtDNA in reporter cells was confirmed using fluorescent microscopy of THP1-macrophages (Figure S5e). To show that GZMK-induced cytosolic mtDNA drives Type I IFN through PRR^28,29, we used reporter cell lines to show IRF Pathway activation and Type I IFN responses, (1.6-fold, p < 0.0001; ANOVA; Figure S5e) and Type I IFN response (1.3-fold, p = 0.023; ANOVA) (Figure S5f), respectively. To further confirm this effect, we transfected GZMK into the cytosol of human primary salivary gland epithelial cells (pSGEC). GZMK induced cytosolic mtDNA and nuclear translocation of pIRF3 (1.2-fold, p < 0.0001; two-sided t-test; Figure S5g) in pSGECs.

GZMK + Tex-like cells have been reported to be involved in multiple target organs (e.g., kidneys, synovium) in autoimmune diseases^27,30. Contemporaneous to our analyses, GZMK + Tex-like cells were found to be gland-enriched and clonally-expanded in SjD patients’ glands, possibly suggesting a salivary-specific autoantigen³¹. Adding upon this foundation, our results explain a unique effector function of these cells by demonstrating the potential for autoreactive GZMK + Tex cells to secrete GZMK on target cells to drive Type I IFN immune signaling in the salivary gland epithelia. This may partially explain Type I IFN signaling characteristic of SjD including in duct and acinar cells and implicate these cells in modulating immune signaling in foci with sustained antigen engagement. In our data, this effect is enhanced in the SSA + SjD subjects. However, the real-time direct observation of these events in situ has not yet been possible.

Sjögren's is defined by loss of a seromucous acinar cell population

The principal target cells and antigens have yet to be elucidated in SjD. To study the SjD-specific effects on epithelial cells, we used our single cell analyses and spatial transcriptomics (Fig. 1, Fig. 6a) to study the composition and architectural changes in the context of SjD. Seven distinct seromucous acinar cell (SMAC) populations were identified using Leiden clustering (Fig. 1e, Fig. 6b). Instead of subdividing these populations manually, two overlapping populations were collapsed into a single “PRR4- Transitioning SMACs” population. The most “stem-like” cells identified were a subpopulation of cells clustering within the ductal cell cluster subsequently named ‘Ductal Progenitors’ (KRT5)³². Exploring the SMAC cells across this trajectory score, genes of ductal differentiation fall in expression (WFDC2) and genes associated with seromucous acinar function increase in expression (ZG16B, PRR4; Fig. 6b-d). MUC7, a seromucous acinar protein, is highest in the CST3 + SMACs. The cells seem to follow two trajectories terminating in PRR4 high and ZG16B high cell clusters (Fig. 6d). When examining the differences in proportions of cells in the salivary glands of SjD patients, the PRR4 + CST3 + acinar cell cluster, which constitute the largest single population of SMACs in nonSjD glands, are significantly reduced in SjD (Figure S6b) while the other acinar cell populations remain statistically unaffected (Figure S1b).

scRNAseq reveals striking and underappreciated heterogeneity of the acinar cell composition of the salivary glands and that in SjD seems to result in preferential loss of SMACs (Figure S1b; Figure S6a,b)³³. To confirm changes in heterogeneity from scRNAseq among the salivary glands cell populations, we used hybrid fluorescent 12-plex ISH/3-plex immunofluorescent spatial phenotyping. We applied markers that discriminate SMAC clusters (MUC7, MUC5B, ZG16B, WFDC2, CST3, and PRR4); and applied canonical markers of salivary differentiation (BPIFB2, PIP, AMY1A), markers of response to Type I IFN (ISG15, B2M), and protein markers (pan-cytokeratin [pCTK], CD45) (Fig. 6a, Figure S6a,b). Following recursive fluorescent whole slide imaging and registration, cells are segmented and phenotyped. Then, dimensionality reduction methods identify clusters of cells based on mRNA expression (Fig. 6e, Figure S6a) and spatial plots of cell phenotypes were constructed (Fig. 6f). Visualizing the expression of MUC7 and ZG16B together (green) separates SMACs from mucous acinar cells (MUC5B, magenta) (Fig. 6f). Cells expressing elevated pCTK but negative for MUC7, ZG16B, and MUC5B are identified as ductal cells (Fig. 6a, 6f). Fluorescent signals for ZG16B (in green), MUC5B (in magenta), PRR4 (in white), WFDC2 (in red), and CST3 (in yellow) showcased a higher diversity in SMAC types than previously reported (Fig. 6f). Examining specific areas of the gland, highlighted in yellow (nonSjD) or orange (SjD) dashes, it became evident that striated ductal cells exhibit positivity for WFDC2, but not other SMAC markers (Fig. 6e,6f). In contrast, intercalated ducts express WFDC2 and CST3, and lower, but not absent, levels of other SMAC markers (e.g., MUC7, ZG16B) (Fig. 6f). In the seromucous acinar cells (enlarged images in white), some cells displayed elevated levels of PRR4 and ZG16B, with variable levels of CST3 (Fig. 6e, 6f). Notably, PRR4 was expressed in isolated seromucous cells and separately within the “serous demilunes” of mucous acini indicating a unique cell type. Mucous cells (MUC5B+) expressed CST3 at lower levels (Fig. 6a,f). These observations demonstrate distinct populations representing pure seromucous acini which are distinct from “serous demilune” cells associated on the periphery of predominantly mucous acini.

To visualize annotated cell types in their tissue context, segmentation and phenotyping were performed (6f, spatial plots) on identified cells positive for MUC7 or MUC5B and the log-normalized expression of six target genes: MUC7, MUC5B, ZG16B, WFDC2, CST3, and PRR4 was calculated, followed by unsupervised clustering. This analysis identifies 15 cell clusters: eight clusters are MUC5B+/MUC7- (mucous acinar cells), and seven are MUC5B-/MUC7+ (SMAC cells). ZG16B expression mirrors MUC7, while WFDC2 and CST3 are expressed in 11 clusters. Alignment with the phenotypes with our reference scRNAseq (Fig. 1e; Fig. 6a; Figure S6a) to identify six SMAC types and at least two distinct mucous acinar cell types confirming the observed microscopic and single cell transcriptomic heterogeneity.

Considering the spatial context, we observe differences between nonSJD and SjD salivary glands (Fig. 6f). In nonSjD glands, a full complement of scRNAseq cell types are identified illustrating a non-pathological state. However, in SjD, some SMACs subtypes are reduced or absent. Closer inspection of the glands and representative spatial plots in health and SjD revealed fewer pure seromucous acini and analogous structures between health and SjD show reductions in PRR4 + expression and increased WFDC2 expression. Mucous acini seem to be increased in proportion and thus the appearance of more “serous demilune”-type SMACs is apparent. Logistic regression analysis of the proportions of cells in the glands show significant changes: increased proportions of terminal mucous acini (+ 3.8-fold, p < 0.001), and decreased proportions of PRR4 + CST3 + WFDC2 + SMAC (-0.7-fold, p < 0.000), and High-ZG16B + SMACs (-0.4-fold, p < 0.000) in SjD glands (Fig. 6f; Figure S6c). PRR4 + CST3 + WFDC2 + SMACs were found in demilunes of seromucous acini, and pure seromucous acini primarily had PRR4-negative cells (Fig. 6f). PRR4 + CST3 + WFDC2-SMACs cells are concentrated at gland lobule peripheries. These spatial transcriptomics data roughly align with scRNAseq data. One minor difference from scRNAseq data includes significant changes in PRR4 + CST3 + WFDC2 + and ZG16B + SMACs, but not PRR4 + CST3 + WFDC2- SMACs which may depend on gene expression patterns in this transitional state. Generally, these results suggest that SMAC populations are more sensitive to SjD-mediated immune insult and can be illustrated by a proportionate increase in mucous cells, loss of seromucous acinar cells, and increased expression of MHC-I and MHC-II (Fig. 4, Figures S4) in putative target cells. Furthermore, these data provide a novel high-resolution map of acinar cell heterogeneity in health and highlight architectural disturbances in SjD glands beyond simply compositional changes emphasized by scRNAseq. These results underscore the complex interplay of different epithelial cell populations in salivary gland pathophysiology, but questions about what changes occur in targeted cells are raised.

To understand the changes in these putative target cells, we analyzed the per-cell changes in pathway utilization in the most SjD-lost clusters from SjD patients compared to the same cell type in nonSjD patients (i.e., PRR4 + CST3 + WFDC2-). The members of this cluster found in SjD individuals does have a distinct gene expression pattern, with cells from nonSjD individuals exhibiting higher levels of expression of indices associated with differentiation and development in a range of tissues led by an index called post-anal tail morphogenesis (PATM, among several indices involved in patterning and differentiation [Figure S7a]). These developmental indices are populated primarily by transcription factors (TFs). The distribution of expression of these indices in nonSjD fits within the lower end of the distribution in SjD, suggesting the cells with high levels of expression are being targeted in disease (Figure S7b). The remaining cells in this cluster in SjD also have a higher chance of colocalizing with immune cells characteristic of the disease, including IgG plasma cells and GZMK + T cells (Fig. 6e; Fig. 5b; Fig. 2). Furthermore, spatial phenotyping illustrates a relationship between proximity to inflammatory cells and expression of MHC-I protein (Fig. 4), indicative of CD8-T cell/MHC-I signaling.

Building on our prior work,^34,35 in this study RNA-seq data from single cells and spatial RNA-sequencing data were collected from patients with and without a diagnosis of Sjögren's Disease. The single-cell data was annotated to identify unique clusters, including seven unique seromucous acinar cell types not previously annotated and five different T cell types. These clusters were used to quantify the percent of cell type in each spot in the spatial data, providing evidence of cell-cell colocalization changes in disease and health. Spatial phenotyping and ex vivo cellular assays further confirm that SjD-enriched GZMK + CD8 Tex cells are enriched in the glands, colocate with ducts and acini, and exhibit cytotoxic potential. SSA + and Sjögren’s diagnosis shared a group of associated indices, including NK- and T-cell-mediated cytotoxicity and antigen presentation, and a list of indices lower in expression. However, the manifold built on index S-scores revealed two distinct paths to T cell dysfunction in Sjögren’s Disease, one associated with autoantibody positivity, and one negative. Changes characteristic of SSA + included increased T cell activation, antigen presentation accompanied by changes in MHC expression and an upregulation across the gland of genes inducible by interferon signaling; these changes were absent in SSA- disease confirming the distinct biological heterogeneity between antibody positive and antibody negative disease.

For the first time we show both direct and indirect evidence that GZMK + CD8 Tex cells may exert potent and a specific immune response in target cells that express markers of acinar differentiation but lie transcriptionally nearest to a progenitor duct-like cell type origin. This effect ultimately reduces the proportions of terminally differentiated acinar cells, with appreciably less effect upon the mucous acinar cells of the glands. This can be appreciated microscopically whereby mucous acinar cells and ducts dominate the typical histopathological profile in Sjögren’s disease MSG biopsies.

What was definitional of disease in both SSA + and SSA- Sjögren’s patients was a loss of PRR4 + CST3 + WFDC2- seromucous acinar cells characterized by the highest levels of expression of MUC7, which fell from the largest single cell population in nonSjD patients to representing less than half that level in SjD. The proportions of other SMAC populations remained statistically unchanged by disease; yet exhibited lower proportions. This can be explained by gland atrophy which is a feature of more longstanding disease. Moreover, it is not possible to quantify the cells that were lost prior to biopsy. One hypothesis for this vulnerability comes from animal models of autoimmune thyroid disease where CD8 T cells remove hypersecreting mutants from glandular epithelia in a theory called “autoimmune hypersecreting mutants (AHSM)”³⁶. In Sjögren’s Disease, we surmise that autoreactive GZMK + CD8 + T cells may selectively remove or restrain the differentiation of transitioning seromucous acinar cells with high levels of secretion. Unsurprisingly for this autoimmune disease, invasive lymphocytes replaced PRR4 + CST3 + WFDC2- seromucous acinar cells nearly one-to-one. Surviving members of this seromucous acinar cluster in Sjögren’s Disease individuals have a distinct gene expression pattern that offers clues to the vulnerability of this population, with cells from non-Sjögren's individuals having higher levels of expression of indices associated with development and proliferation which may have led to sensitivity in the Sjögren’s immune microenvironment. Looking at cell-cell type correlations within spots in the spRNAseq, negative correlations between unrelated cell types in the non-Sjögren's individuals gave way to zero or low positive correlations in Sjögren's Disease, suggesting a higher prevalence of interactions between previously distinct cell types. Differential colocalization scores revealed a pathogenic adaptive immune infiltrate adjacent to vulnerable seromucous acinar cells suggesting effector cell targeting of putatively antigen-expressing epithelial cells.

In conclusion, Sjögren’s Disease is a chronic systemic autoimmune disease characterized by loss of exocrine gland function. Using multiple new transcriptomic and proteomic datasets, we identify both the most sensitive target cells in the salivary epithelia of Sjögren’s Disease patients and an important population of GZMK + CD8 + T cells defined by a novel effector function of this cell type: altering mitochondrial membrane permeability and eliciting cellular innate immune signaling. These results provide a potential rationale for non-viral Type I IFN signaling from the affected epithelia in Sjögren’s Disease.

By taking a sample at a single time point this study represents a snapshot of dysfunctional glands and does not explore how Sjögren's disease progresses. The disease duration in each patient can only be roughly estimated based on symptom prevalence. Some of our data suggest that pre-clinical disease (those with SSA + without satisfying classification criteria) exhibit altered pathway utilization similar to disease-affected patients and may provide clues about “pre-clinical” disease. Moreover, undiagnosed comorbidities, changes in the microbiome, changes in hormone signaling, changes in metabolism and so on represent blind spots of variables that were not measured and cannot be easily controlled. Gene expression changes result from changes in histone modification, DNA accessibility, and regulatory element methylation, none of which were measured as part of this study. Protein quantification of a limited set of gene products was performed, mostly for immune cells and state markers. Despite sampling more than 90,000 unique cells with scRNAseq across our patients, certain cell types known to be present in the salivary gland, most notably nervous tissue, are not recovered by this method of single-cell library preparation. To overcome limitations associated with modest numbers of cases, multiple orthologous (e.g., spatial proteomics, spatial transcriptomics) and experimental approaches (in vitro assays using multiple cell lines including patient-derived primary cells) were employed. PhenoCycler-Fusion spatial phenotyping, in this dataset, exhibited limited ability to segment in the densest immune infiltrates, attenuating the true number of inflammatory cells present in the glands and the interactions amongst those cells. Thus, direct comparisons between spatial phenotyping and scRNAseq-based shifts in cellular proportions are not possible. However, the directions of the shifts are maintained, thus spatial phenotyping differences may be of lower magnitude albeit with maintained tissue context representing a potential trade-off between resolution and in situ location. While not a limitation, per se, our GZMK transfection experiments, using multiple cell lines (AC, THP-1 macrophages, and HSG, we were not able to reliably induce detectable ROS or cell death (early or late apoptosis) up to the maximum of what we could transfect into cells. This may reflect the cell lines tested or that the specific effects on targeted cells may be different (e.g., fibroblasts or salivary epithelial cells). This highlights opportunities to further explore the role of granzymes (e.g., Granzyme B, Granzyme A, Granzyme K) on intercellular innate immune signaling in a cell-type specific contexts.

Potential Conflicts of Interest: BMW has Cooperative Research Award and Development Agreements [CRADA] from Pfizer, Inc., and Mitobridge, Inc. (A subsidiary of Astellas Pharma, NIAMS has CRADAs with Astra Zeneca and Bristol Myers Squibb. In the past 3 years, SAT has been a member of scientific advisory boards of Element Biosciences, Qiagen and ForeSite Labs and is a co-founder and equity holder of Transition Bio and EnsoCell. These agreements did not financially support the experimental results presented herein. SAT, KMB, QTE, BFM, and BMW, are all active members of the Human Cell Atlas.

Funding/Support: This research was principally supported through research awards to BMW from the Division of Intramural Research (DIR) Program of the National Institute of Dental and Craniofacial Research of the National Institutes of Health (NIH/NIDCR ZIA: DE000704) with additional support through (Z01-DE00695, PI-Chiorini), and the DIR Program of the National Institute of Arthritis and Musculoskeletal and Skin Diseases of the NIH (NIH/NIAMS ZIA: AR041199). TFJP is supported through the Graduate Partnership Program. This work also received support from the VCU Massey Comprehensive Cancer Center Bioinformatics Shared Resource, funded, in part, by NIH-NCI Cancer Center Support Grant P30 CA016059.

ACKNOWLEDGMENT: We owe our gratitude to the research participants who provide their support of our clinical studies through volunteerism. We thank the generous support of the NIH/UMD Graduate Partnership Program for supporting TFJP. We would also like to thank the support and the staff of the NIDCR Office of the Clinical Director, NIDCR Dental Clinic, the National Eye Institute Ophthalmology Clinic, and the NCI Center for Cancer Research Anatomic Pathology. This work utilized the computational resources of the NIH HPC Biowulf cluster. (http://hpc.nih.gov). We would also like to thank the generous support of the NIH/UMD Graduate Partnership Program for TFJP. We would like to thank the following NIH core facilities which enabled this study: NIDCR/NICHD Computational Biology and Genomics Core, the NIDCR Microscopy Core,

DATA AVAILABILITY. Data were collected based on existing resources, and no additional software was generated for this manuscript. Analysis notebooks will be available at the time of publication before at the time of upload to a pre-print server (e.g., medrxiv), whichever comes first. Single cell, spatial transcriptomics, and spatial phenotyping data, as well as the associated metadata, will be uploaded to dbGAP or GEO, consistent with consent language of the parent clinical protocols. All other data is available upon reasonable request.

Clinical Evaluation of Patients.

Research participants were consented to NIH-approved IRB protocols (15-D-0051, NCT00001390) prior to any study procedure. All participants were evaluated and classified according to 2016 American College of Rheumatology (ACR) and the European League Against Rheumatism (EULAR) classification criteria⁹. Comparator tissues included subjects (nonSjD) who were otherwise healthy, but did not meet 2016 ACR-EULAR criteria. All subjects were screened for evidence of systemic autoimmunity and received comprehensive oral, salivary, rheumatological, and ophthalmological investigations³⁶. Clinical investigations were conducted in accordance with the Declaration of Helsinki principles.

In vitro immune cell stimulation and flow cytometry.

Blood was collected in EDTA tubes and PBMCs were isolated by Ficoll-Paque density centrifugation (GE Life Sciences). Minor salivary glands were minced using a GentleMACS Tissue Dissociator (Miltenyi Biotec) and were enzymatically digested in a shaker incubator at 37°C for 40 min in RPMI medium containing 1 mg/ml Collagenase D, 1 mg/ml hyaluronidase and 50 U/ml DNase I (Sigma Aldrich). Cell suspensions were then passed through a 100 μm cell strainer and enriched for lymphocytes using a 40% Percoll density gradient centrifugation (GE Life Sciences). Cells were stimulated with phorbol myristate acetate and ionomycin (ThermoFisher) at 37°C for 4 hours in the presence of anti-CD107a (BioLegend), brefeldin A, and monensin (ThermoFisher) in Complete RPMI media containing 10% fetal calf serum (FCS), 1% Sodium pyruvate, 25 mM Hepes, and 2 mM l-glutamate (ThermoFisher). Cells were stained with various combinations of fluorochrome-labeled antibodies and Fixable Viability Dye eFluor 780 (purchased from BioLegend or ThermoFisher) in PBS containing 1% FCS and 1X Brilliant Stain Buffer (BD Biosciences) for 20 min at 4°C. After fixation and permeabilization using Foxp3/Transcription Factor Staining Kit (ThermoFisher), intracellular cytokine staining was performed in PBS containing 1% FCS, 1X Brilliant Stain Buffer, and 1X Permeabilization buffer for 1 hour at 4°C. Samples were acquired using a FACSymphony cytometer (BD Biosciences) and analyzed using FlowJo software (FlowJo, LLC).

Cell type analysis and further annotation.

Analysis of the single-cell RNA-sequencing was done in Python 3.9.4 using numpy 1.21.2³⁷, scipy 1.6.3³⁸, statsmodels 0.12.2³⁹, scikit-learn 0.24.2⁴⁰, pandas 1.2.4⁴¹, anndata 0.8.0 and scanpy 1.7.2⁴². Figures were made using matplotlib 3.4.2⁴³ and seaborn 0.11.1⁴⁴. Cells were gated on 15% mitochondrial reads, 50% ribosomal reads, 5% hemoglobin gene reads and a minimum of 10 genes and 100 counts. Nearest neighbors were identified using the BBKNN algorithm⁴⁵ followed by UMAP manifold representation⁴⁶ and Leiden clustering⁴⁷. Leiden clusters were named by gene expression using known markers as well as scanpy’s rank_genes_groups tool. In multiple cases a single Leiden cluster corresponded to multiple known cell types - when that occurred, clusters were divided based on gene expression and UMAP embeddings to generate new annotation classes. When multiple clusters corresponded to a single previously-known cell type annotation, individual pairwise calls to rank_genes_groups were used to identify genes discriminating between these clusters. Differentially expressed genes were identified using a rank_genes_groups call between two clinical features.

Proportion changes.

To test changes in cell type proportions between clinical features the crosstab function from pandas was used to get a per-patient proportion for each cell annotation. T-tests between clinical features produced q-values after Benjamini-Hochberg multiple test correction, and q-values below 0.05 were deemed statistically significant.

Calculating Index S-scores.

A comprehensive set of gene indices including KEGG terms^48,49, GO terms^50,51, Reactome Pathways⁵², mSigDB Hallmark Pathways^53–55, four pathways that represent senescence-associated secretory phenotypes⁵⁶, a previously-published pathway of ECM factors involved in Sjögren's fibrosis³⁴ and a set of genes regulated by interferon was used to transform from a gene expression space to an index S-score space. Expression data from each gene were converted to Z-scores and this Z-scores matrix was multiplied by an encoding matrix representing the relationships between genes and indices. The new matrix has a S-score for each index which is a sum of the Z-scores of the genes that make up that index. These S-scores were then treated as if new gene expression data and a new UMAP was generated and Leiden clustering was performed. The relationship between the new and old clustering methods was quantified using the pandas crosstab function.

Trajectory analysis and differentiation.

Trajectory analysis in the single-cell was performed using the Wishbone algorithm on the subset of epithelial tissues. The initial cell was set as the cell with the highest MKI67 expression, a marker of stemness, which happened to be in the Ductal Progenitor population. All of the epithelial cells were arranged by Wishbone distance from the initial cell and the expression of key markers of epithelial identity were plotted after normalization with a 1000-cell cosine rolling window.

Spatial RNA-sequencing analysis.

The spatial RNA-sequencing data was analyzed in Python 3.9.7 using numpy 1.21.5, scipy 1.8.0, statsmodels 0.13.2, pandas 1.4.1, networkx 2.7.1⁵⁷ and scanpy 1.8.2 and figures were made using matplotlib 3.5.1 and seaborn 0.11.2. Cells were gated identically to the single-cell data. A cell2location model⁵⁸ was trained using the previous scRNAseq data and applied to the spRNAseq data followed by normalization to proportions. The geojson package (2.5.0) was used to import manual annotations. Within- and next-spot correlations were calculated as well as a within-spot co-occurrence score C, which was quantified as:

$$C\left({\kappa }_{1},{\kappa }_{2}\right)=\frac{\left({\kappa }_{1}+{\kappa }_{2}\right){e}^{-3\left|{\kappa }_{1}-{\kappa }_{2}\right|}}{2}$$

where ${\kappa }_{1}$ is the proportion of the first cell type and ${\kappa }_{2}$ the proportion of the second cell type within a spot. This score increases as the proportions increase and is higher if the proportions are close to the ${\kappa }_{1}$ = ${\kappa }_{2}$ axis.

Multiplex proteomics (Phenocycler Fusion).

The multiplex was performed using 5 µm FFPE sections mounted on SuperFrost plus (ThermoFisher MA, USA), which underwent deparaffinization and rehydration. Slides were soaked in a Coplin jar with 1:20 AR9 buffer (Akoya Biosciences, MA – USA). The jar was seated in a pressure cook for 15 minutes at low pressure. Samples went to tabletop cooling for 30 min, followed by 30 seconds in deionized water and 100% EtOH for 3min. All slides went through pre-staining procedures by immersing the slides in the Hydration Buffer for 2 minutes and staining buffer for 20 minutes (Akoya Biosciences, MA – USA). The antibody cocktail for primary incubation was prepared using 4 blockers (G, S, J, and N), including 9.5 microliters of each in 362 µL of staining buffers. For each slide, we aliquoted 150ml and included 1 microliter of each of the antibodies (list below). The slides were removed from the staining buffer and accommodated in a humidity chamber (StainStray, Sigma-Aldrich, MO – USA), and primary incubation was done overnight on 4°C refrigerator. After incubation, the slides were followed by a post-staining fixing solution for 10 minutes. After fixation, sequential 1-minute PBS washes were performed, followed by an immersion in ice-cold MeOH for 5 minutes. The sections were treated with a 200 µL of final fixative solution for 20 minutes. Additional washes were performed to remove the final fixative solution. Slides were dried and mounted with the Akoya flow using a press that seals the flow cell/coverslip on the slides for 30 seconds.

The slides were removed from the press and soaked in a 1X PCF buffer (Akoya Biosciences, MA – USA). To prepare the PCF reporter wells, a 15 mL Falcon tube was initially covered with aluminum foil. Subsequently, into this Falcon tube, we introduced 6.1 mL of nuclease-free water, along with 675 µL of 10X PCF buffer, 450 µL of PCF assay reagent, and 4.5 µL of a concentrated DAPI solution that we had prepared in-house. This DAPI addition resulted in a final DAPI concentration of 1:1000. The reporter stock solution was then pipetted into 18 amber vials, for each slide, with each vial containing 235 µL of the solution. To each vial, we added 5 µL of reporter per cycle. The total volume per vial was either 245 µL for a cycle involving 2 reporters or 250 µL for a cycle with 3 reporters (see Supplemental Methods Table 2.).

A specific criterion was established for each cycle, allowing for a maximum of 3 reporters. These reporters were chosen from a pool of Atto550, AlexaFluor 647, and AlexaFluor 750, as deemed appropriate based on the experimental requirements (specific reporter FluorChannels and antibody barcodes in the list below).

Distinct pipette tips were employed to transfer the contents of each amber vial into a 96-well plate. Vials containing DAPI were pipetted into wells within the H-row, while vials containing reporters were distributed into wells in other rows of the 96-well plate.

After completing the well-filling process, the wells were securely sealed using adhesive aluminum foil (Akoya Biosciences, MA – USA). Subsequent imaging procedures were conducted using a PhenoImager Fusion system, which was connected to a PhenoCycler (specifically, the PhenoCycler Fusion system from Akoya BioSciences). The imaging was performed with a 20X objective lens from Olympus.

The necessary solutions for instrument operation included ACS-grade DMSO from Fisher Chemical, nuclease-free water, and 1X PCF buffer with the addition of a buffer additive. This latter solution was prepared by mixing 100 mL of 10X PCF buffer and 100 mL of buffer additive with 800 mL of nuclease-free water.

Image acquisition and Segmentation.

A multi-step procedure was employed to analyze the protein multiplex images following image acquisition using Qupath⁶⁰. The whole slide images were analyzed using the 16-bit image directly provided as the raw file from PhenoImager. Individual slides were evaluated by a trained pathologist (BM) for each of the 40 antibody markers to set manual thresholds and use as a comparison point for the computational assignment described in the next session. The images underwent segmentation using a pre-trained model based on Cellpose 2.0⁶¹. Segmentation was refined iteratively, utilizing 80 fluorescent images sourced from post-mortem biopsies. The model underwent multiple training iterations until precise delineation of cell expansion was achieved in both the acinar cells and ducts. Each of the minor salivary gland biopsies were segmented separately. Following nuclei-based cell segmentation, raw data was exported, and the number of pixel value units per cell and channel was quantified based on the threshold.

Multiplex Proteomics Image Analysis: Preparing cell type signature matrix.

The matrix has cell type markers in rows and cell types of interest in columns. An element of this matrix has a value 1 for every cell type marker considered to be highly expressed in the given cell type (i.e., column) and 0 otherwise. Markers distinguishing deeper levels of cell granularity (e.g., CD4 for Helper T cells which are a subclass of T cells, or CD8 for Cytotoxic T cells which are another subclass of T cells) assume values corresponding to their cell granularity level in the signature matrix (e.g., 2 for hierarchy level 2).

Cell type annotation using TACIT.

The TACIT algorithm takes two main inputs: a matrix consisting of z-normalized marker intensities for each cell, and the cell type signature matrix described above. To start, the two matrices are multiplied to derive the Cell Type Relevance scores (CTR) for all cell types in all cells. The CTR score becomes a linear combination of the intensities of markers defining the specific cell type in the cell. The higher the CTR score, the stronger the evidence that the cell is representing a given cell type. As next, TACIT examines the distribution of the CTR scores across all cells and conducts primary cell type assignments in multiple steps: (1) Using the z-normalized marker intensities, cells are first grouped into a large number of seed clusters, each of which averages 0.5% of total cells and represents a small but tight community of cells with highly similar overall marker profiles; (2) For each cell type, the distribution of the cluster median CTR scores across all seed clusters is plotted in ascending order; (3) Piecewise regression is then employed to identify breakpoints that separate data (here, clusters) into different linear trends over different data regions (here, cluster ranks). The piecewise regression model is fit to the data allowing for up to three breakpoints, and the model that best represents the data is chosen based on the Akaike Information Criterion (AIC) score; (4) Once breakpoints are determined, clusters are categorized into the “low relevance group (LRG)” and “high relevance group (HRG)”, with the former encompassing all clusters ranking below the lowest breakpoint, and the latter group encompassing clusters ranking above the highest breakpoint; (5) CTR value that minimizes classification error between LRG and HRG is determined and used as a cutoff to assign primary cell type. Steps (1) to (5) are performed on each considered cell type individually. A binary cell type matrix (CTM) is formed, where cells (in rows) are assigned 1 for all plausible cell types (in columns) whose CTR scores pass the cutoffs determined in step (5).

Mixed cell type deconvolution.

CTM columns representing lower granularity cell types (columns for T-helper cells or cytotoxic T cells) are updated in accordance with the individual marker positivity (CD4 or CD8, respectively, using the cited example) and an AND logical gate. The individual marker positivity (0 if below or 1 is above the cutoff) is derived as it was the case of CTR score described in step 5 above, except that this time CTR score is a single marker (not a linear combination of multiple markers). The updated binary CTM matrix will contain a fraction of cells uniquely defined (rows for which only once a value of 1 has been registered across the columns), with the remaining cells exhibiting multiple cell type assignments (multiple columns have a value of 1 in a row, indicating the assignment of various cell types to the cell in that row). To resolve the conflict, clean cells from each cell type involved in the mixture group together with the ambiguous cells will be projected onto the feature space composed of relevant markers only (i.e., those that define considered cell types). A cell with mixed identities is reassigned to the cell type exhibited by the majority of its k nearest clean cell neighbors using KNN classification approach. The cells with mixed identities are then reassigned to the most similar cell type of their closest neighbors using KNN classification approach. Uniform Manifold Approximation Projection (UMAP) is used for visualization and verification of the cell type assignment results.

Cell-cell interactions and neighborhood analysis.

PhenoCycler data from each individual tissue was processed with Squidpy (version 1.3.1)⁶² that describes cellular interactions as graphs with nodes representing individual cells and edges potential cellular interactions as determined by Delaunay triangulation. A 99th percentile distance threshold for each tissue was set to remove edges representing improbably long cell-to-cell distances. Non-deconvoluted cells (classified under “Unknown” cell types) were excluded from the analysis before performing Delaunay triangulation. An interaction matrix was constructed, with an element ai,j representing the number of edges shared between the cell type i and cell type j. Healthy interaction matrix is an aggregate of connections in healthy tissues, while disease interaction matrix is derived from all diseased tissues. The matrices were column-normalized and used for comparative analysis of cellular interactions between healthy and diseased tissues. The elements of the diseased interaction matrix was divided by the elements of the healthy interaction matrix, and the resultant fold-change interaction matrix was examined for differences in cell type interactions between the two conditions, with FC values above 1 indicating more pronounced interactions in diseased tissues, and FC values below 1 suggesting more interactions in healthy tissues. To visually represent these differences, a hierarchically clustered heatmap with Euclidean distance was used.

Cell Lines and Cell Culture.

THP1-Lucia ISG and IFN-α/β Reporter HEK 293 Cells were purchased (Invivogen, San Diego, CA, USA); the NS-SV-TTAC Acinar Cells (aka: AC cells) were provided by Dr. Jay Chiorini. THP1 and 293 cells were cultured as previously reported in RPMI (Corning Cellgro, Corning, NY, USA, #15-040-CV) and DMEM (Gibco, New York, USA, #10313-021), respectively, supplemented with 10% FBS (Gibco, #A47668-01) and penicillin/streptomycin 100 U/mL, 2 mM L-glutamine, 100 ug/mL Normocin. The NS-SV-TTAC cell line was cultured in Defined K-SFM (ThermoFisher, #10744019), supplemented with Defined Keratinocyte-SFM Growth Supplement. All cell lines were maintained in a humidified 37° C incubator with 5% CO₂, and cell lines were tested and confirmed negative for mycoplasma, validated with MycoStrip (InvivoGen, San Diego, CA, USA, #rep-mys-50).

IRF Pathway Activation.

THP1 ISG Lucia cells were cultured according to manufacturer’s guidelines and used to measure the activation of the IRF Pathway. Briefly, cells were differentiated into macrophages using 150 nM Phorbol 12-myristate 13-acetate (PMA; Invivogen, San Diego, CA, USA) for 24 hrs. After differentiation, PMA-media was removed and replaced with basal growth media. After 72 hrs, the media was replaced with OPTIMEM (Gibco, New York, USA, #31985-070) prior to transfection. Cells were transfected with XFect Protein Transfection Reagent (TakaRa, Shiga, Japan, #631324). Cells were pre-treated, for 30 min with the 10 uM of caspase inhibitor, Q-VD-OPH (Sigma-Aldrich, Burlington, MA, USA, #SML0063), as indicated and transfected with either recombinant human Granzyme K (Sino Biological, #19732-H08H), recombinant human Granzyme B (Sino Biological, #10345-H08H) at indicated concentrations in OPTIMEM media for 3 hours. OPTIMEM was replaced with growth media and luciferin (a direct readout of IRF pathway) was collected from the cellular supernatant and measured at 24 and 48 hours post transfection using QUANTI-Luc: Luciferase Detection Reagent (Invivogen, San Diego, CA, USA, #rep-qlc4lg1).

IFN-α/β Pathway Activation.

IFNα/β Reporter HEK 293 cells were cultured according to the manufacturer's guidelines to detect Type I interferons. Cells were transfected with XFect Protein Transfection Reagent (TakaRa, Shiga, Japan, #631324) with a caspase inhibitor, 10 uM Q-VD-OPH (Sigma-Aldrich, Burlington, MA, USA, #SML0063). Cells were transfected with either recombinant human Granzyme K (Sino Biological, #19732-H08H) or Granzyme B (Sino Biological, #10345-H08H) in OPTIMEM. STAT1/2-inducible secreted embryonic alkaline phosphatase (SEAP), was quantified utilizing the manufacturer's guidelines at 24 and 48 hours transfection using QUANTI-Blue Solution (Invivogen, San Diego, CA, USA, #rep-qbs).

Immunofluorescence microscopy and analysis.

THP1 Macrophages and NS-SV-TTAC Acinar Cells were plated on an 8-well Lab-Tek II Chamber Slide (Thermo Scientific, 154534). After indicated treatment, cells were loaded with Quant-iT PicoGreen Reagent (1:50) (Invitrogen, #P7589A) and MitoTracker Red CMXRos (Invitrogen, #M7512) (1:10000) and incubated for 45 minutes at 37°C. The cells were fixed for 20 minutes at room temperature in 4% PFA. The cells were blocked in 3% BSA diluted in 1X PBS for 1 hour. Cells were washed with 1X PBS and then stained with primary antibodies overnight and then stained by secondary antibodies for 1 hour at room temperature (Supplemental Table 2). Images were acquired on Nikon A1 HD (Nikon) confocal microscope and processed with CellProfiler in ImageJ (Broad Institute).

Primary salivary gland epithelial cell (pSGEC) generation and culture.

pSGEC were generated according to methods published by Jang et al (2015)⁶³. The pSGEC cell lines were not authenticated. pSGECs were plated on collagen-coated 8-well chamber slides and transfected as indicated. Cells were washed with 1X PBS and immediately fixed with 4% PFA for 30 minutes at 37°C. The reaction was quenched via the addition of glycine (100 mM) for 20 minutes at room temperature. Cells were washed with 1X PBS and permeabilized with 1% TritonX-100 for 12 minutes at room temperature. The cells were blocked in 5% nonfat dry milk in PBS-T (blocking buffer) for 1 hour at room temperature. Cells were washed with PBS-T and then incubated with primary antibodies in blocking buffer and washed three times with PBS-T (10 min each) before incubated with the secondary antibodies (in blocking buffer) for 1 hour each at room temperature (Supplemental Table 2). After washing, cells were mounted with Fluoro-Gel II with DAPI (EMS) and images were acquired on an Olympus (FV1200 MP) confocal microscope and processed with CellProfiler in ImageJ.

AnnexinV Staining.

The potential for apoptosis induced by Granzyme K or B was measured using the flow-cytometry-based Dead Cell Apoptosis Kit with Annexin V Alexa Fluor 488 and Propidium Iodide (Invitrogen) following the manufacturer’s protocol on a FACS Symphony flow cytometer.

Western Immunoblot:

Extracted total protein was resolved using sodium dodecyl sulfate-polyacrylamide gel electrophoresis (SDS PAGE) (4% stacking, 12% resolving, Invitrogen) and transferred onto Polyvinylidene fluoride (PVDF) membranes. Membranes were incubated with primary and secondary antibodies (online Supplemental Methods Table 2). The signal was detected by ChemiDocMP Imaging System (BIO-RAD), and the density of the bands was analyzed using Fiji.

Supplemental Table 2: Reagents and Dilutions

Reagents		Manufacturer	Dilution
Immunofluorescence Confocal Microscopy
Figure 5
	Quant-iT PicoGreen Reagent	Invitrogen	1:50
	MitoTracker Red CMXRos	Invitrogen	1:10000
	6-His Tag Antibody (A190-113A) Goat Polyclonal	Bethyl	1:200
Figure 5	Phospho-IRF-3 (Ser396) (4D4G) Rabbit mAb	Cell Signaling Technology	1:100
Figure 5	Granzyme K Antibody (PA5-50980) Rabbit Polyclonal	Invitrogen	1:100
Western Blots
HRP-conjugated anti β-actin Mouse mAb		Invitrogen	1:10000
Phospho-IRF-3 (Ser396) (4D4G) Rabbit mAb		Cell Signaling Technology	1:300
IRF-3 (D83B9) Rabbit mAb		Cell Signaling Technology	1:300
Granzyme K Antibody (PA5-50980) Rabbit Polyclonal		Invitrogen	1:300
Granzyme B Antibody (MA1-80734) Mouse mAb		Invitrogen	1:300
RDye® 800CW Donkey anti Rabbit IgG (H + L)		LI-COR Biosciences	1:10000
IRDye® 680RD Donkey anti-Mouse IgG (H + L)		LI-COR Biosciences	1:10000
SuperSignalTM West Pico PLUS Chemiluminescent Substrate		Thermo Scientific	-

HiPlex RNAscope In Situ Hybridization.

12-plex RNAscope fluorescent in situ was performed on 5 µm-thick FFPE sections, using the RNAscope HiPlex12 Reagents Kit (Advanced Cell Diagnostics, Newark, CA, USA), according to the manufacturer’s protocols. The protocol was modified to incorporate antibody-based immunofluorescence. The following modifications were introduced: the slides were baked in a dry air oven for 30 min at 60°C. After rehydration the slides were processed for target retrieval using 1X Co-Detection Target Retrieval, during 20 min. The slides were quickly rinsed in water and then washed in PBS-T (PBS 1X, 0.1% Tween-20, pH 7.2). A hydrophobic barrier was created using ImmEdge pen and primaries antibodies were added to detect pCTK (Novus Bio, NBP2-29429, dilution 1:200 in Co-detection Antibody Diluent from ACD), ACTA2 (Antibodies.com, #AB2445, dilution 1:300 in Co-detection Antibody Diluent from ACD) and CD45 proteins (NB100-77417SS, dilution 1:300 in Co-detection Antibody Diluent from ACD) and incubated at 4°C overnight. The next day primary antibodies were washed in PBS-T, three times for 2 min at room temperature. Tissues and primary antibodies were post-fixed in 10% Neutral Buffered Formalin (NBF) for 30 min at room temperature and washed in PBS-T thrice for 2 min. Tissues were then digested using RNAscope Protease III during 30 min at 40°C. Immediately continuing with the HiPlex assay according with standard protocol. The following RNAscope HiPlex probes were used: Hs-CST3(T1), Hs-PIP(T2), Hs-WFDC2 (T3), Hs-BPIFB2 (T4), Hs-CFTR(T5), Hs-B2M(T6), Hs-ZG16B(T7), Hs-ISG15(T8), Hs-PRR4(T9), Hs-AMY1A(T10), Hs-MUC7(T11), Hs-MUC5B(T12). All tissues were first run with positive and negative control probes (species-specific) to ensure tissue and RNA quality for RNAscope ISH.

Recursive Image Acquisition.

After each round of HiPlex RNAscope hybridization, fluorescent images were acquired on Axio Scan-Z1 (Zeiss) fluorescent slide scanner. Whole slide images were collected with a 40x/0.95 objective, with a 0.5 µm step size. Images for each individual in the study were collected from at least 2–4 salivary glands. After each round of imaging for HiPlex RNAscope, coverslips were removed by submerging the slides in SSC 5X buffer until coverslips fall off the slide. Subsequently, fluorophores were cleaved, autofluorescence quenched, and then new fluorophores were added (Advanced Cell Diagnostics, Newark, CA, USA). Immediately after imaging was repeated on the same tissue slices using the additional ISH probes. We performed this for three rounds of RNAscope, after final ISH imaging step, fluorophores were cleaved, autofluorescence quenched and the following secondary antibodies were added: donkey-antiMouseAF488, donkey-antiRabbitAF446 and donkey-antiRat647, Jackson Immunoresearch, all in dilution 1:300 in Co-detection Antibody Diluent from ACD. Finally, antibodies were washed in PBS-T and mounted in prolong-Diamond antifade with DAPI. A final round of imaging was performed to detect the signal of the secondary antibodies.

Image Processing: Registration, Segmentation, Quantification of RNAscope Puncta.

After collecting the images, the following process was used to register the individual fluorescent images: First the set of whole slide images to analyze was defined as a QuPath project. The project was then opened in Fiji, and the registration was performed using both elastix⁶⁴, and a semi-automated approach using BigWarp to interactively improve the results of the automated registration. The transformations were retrieved and applied in QuPath to produce on-the-fly transformed images that correspond to superimposition of all four rounds of images per slide. Once the overlayed images were produced in QuPath, new OME-tiff files were exported, the new 19-channel images were used to create a new Qupath project. Cells were segmented using CellPose developing our own model^65,66. Once the cells were segmented, we quantify the number of dots per cell using the subcellular detection in QuPath.

HiPlex ISH analysis.

The matrix of counts of dots per cell obtained from Qupath was cleaned and parsed in R. MUC7 mRNA positive cells were subset to identify SMACs subtypes. Log2-transformed gene expression of WFDC2, CST3, PRR4, ZG16B, MUC7, and MUC5B were used as inputs UMAPS were obtained using the uwot package in R, the clusters were classified based on the expression of WFDC2, CST3, PRR4, ZG16B, MUC7, and MUC5B, and their similarity to reference cell types identified in scRNAseq. With this set of cell type identification, we ran a co-occurrence test using SPIAT⁶⁷. After reclassification, cell type proportions were calculated for each of the patients in the study and the cell proportions were compared between nonSjD and SjD using a mixed model logistic regression on the package lme4⁶⁸.

General statistical methods.

All non-high throughput data (e.g., analysis of biological assays) were analyzed using Prism 9. Appropriate statistical tests were selected based on data type and are described throughout the text. Figures were constructed using Prism 9 when not otherwise described elsewhere.

Mavragani, C. P. & Moutsopoulos, H. M. Sjögren syndrome. CMAJ 186, (2014).
Mavragani, C. P. & Moutsopoulos, H. M. Sjögren’s syndrome: Old and new therapeutic targets. J. Autoimmun. 110, 102364 (2020).
Hall, J. C. et al. Molecular Subsetting of Interferon Pathways in Sjögren’s Syndrome. Arthritis Rheumatol. 67, 2437–2446 (2015).
Moutsopoulos, H. M. Sjögren’s syndrome: autoimmune epithelitis. Clin. Immunol. Immunopathol. 72, 162–5 (1994).
Mariette, X. & Criswell, L. A. Primary Sjögren’s Syndrome. N. Engl. J. Med. 378, 931–939 (2018).
Khatri, B. et al. Genome-wide association study identifies Sjögren’s risk loci with functional implications in immune and glandular cells. Nat. Commun. 13, 4287 (2022).
Taylor, K. E. et al. Genome‐Wide Association Analysis Reveals Genetic Heterogeneity of Sjögren’s Syndrome According to Ancestry. Arthritis Rheumatol. 69, 1294–1305 (2017).
Chibly, A. M., Aure, M. H., Patel, V. N. & Hoffman, M. P. Salivary gland function, development, and regeneration. Physiol. Rev. 102, 1495–1552 (2022).
Shiboski, C. H. et al. 2016 American College of Rheumatology/European League Against Rheumatism Classification Criteria for Primary Sjögren’s Syndrome: A Consensus and Data‐Driven Methodology Involving Three International Patient Cohorts. Arthritis Rheumatol. 69, 35–45 (2017).
Verstappen, G. M., Pringle, S., Bootsma, H. & Kroese, F. G. M. Epithelial–immune cell interplay in primary Sjögren syndrome salivary gland pathogenesis. Nat. Rev. Rheumatol. 17, 333–348 (2021).
Chisholm, D. M. & Mason, D. K. Labial salivary gland biopsy in Sjögren’s disease. J. Clin. Pathol. 21, 656–660 (1968).
Greenspan, J. S. The histopathology of Sjögren’ssyndrome in labial salivary gland biopsies. 37, (1974).
Tarpley, T. M., Anderson, L. G. & White, C. L. Minor salivary gland involvement in Sjögren’s syndrome. Oral Surg. Oral Med. Oral Pathol. 37, 64–74 (1974).
Maria I Christodoulou, Efstathia K Kapsogeorgou, & Haralampos M Moutsopoulos. Characteristics of the minor salivary gland infiltrates in Sjögren’s syndrome. J Immunol. 34, 400–7 (2010).
Gottenberg, J.-E. et al. Activation of IFN pathways and plasmacytoid dendritic cell recruitment in target organs of primary Sjögren’s syndrome. Proc. Natl. Acad. Sci. 103, 2770–2775 (2006).
Pipi, E. et al. Tertiary Lymphoid Structures: Autoimmunity Goes Local. Front. Immunol. 9, 1952 (2018).
N Gestermann et al. STAT4 is a confirmed genetic risk factor for Sjögren’s syndrome and could be involved in type 1 interferon pathway signaling. Genes Immunol. 11, 432–8 (2010).
Thompson, N., Isenberg, D. A., Jury, E. C. & Ciurtin, C. Exploring BAFF: its expression, receptors and contribution to the immunopathogenesis of Sjögren’s syndrome. Rheumatology 55, 1548–1555 (2016).
Yoshimi, R., Ueda, A., Ozato, K. & Ishigatsubo, Y. Clinical and Pathological Roles of Ro/SSA Autoantibody System. Clin. Dev. Immunol. 2012, 1–12 (2012).
Gandolfo, S. & Ciccia, F. JAK/STAT pathway targeting in primary Sjögren syndrome. Rheumatol. Immunol. Res. 3, 95–102 (2022).
Nezos, A. et al. Type I and II interferon signatures in Sjögren’s syndrome pathogenesis: Contributions in distinct clinical phenotypes and Sjögren’s related lymphomagenesis. J. Autoimmun. 63, 47–58 (2015).
Nguyen, C. Q., Hu, M. H., Li, Y., Stewart, C. & Peck, A. B. Salivary gland tissue expression of interleukin‐23 and interleukin‐17 in Sjögren’s syndrome: Findings in humans and mice. Arthritis Rheum. 58, 734–743 (2008).
Lee, Y. J. et al. Salivary chemokine levels in patients with primary Sjögren’s syndrome. Rheumatology 49, 1747–1752 (2010).
Hernandez-Molina, G., Michel-Peregrina, M., Hernandez-Ramirez, D. F., Sanchez-Guerrero, J. & Llorente, L. Chemokine saliva levels in patients with primary Sjögren’s syndrome, associated Sjögren’s syndrome, pre-clinical Sjögren’s syndrome and systemic autoimmune diseases. Rheumatology 50, 1288–1292 (2011).
McCoy, S. S. et al. Ruxolitinib inhibits IFNγ-stimulated Sjögren’s salivary gland MSC HLA-DR expression and chemokine-dependent T cell migration. Rheumatology 61, 4207–4218 (2022).
Barone, F. et al. CXCL13, CCL21, and CXCL12 Expression in Salivary Glands of Patients with Sjögren’s Syndrome and MALT Lymphoma: Association with Reactive and Malignant Areas of Lymphoid Organization. J. Immunol. 180, 5130–5140 (2008).
Jonsson, A. H. et al. Granzyme K ⁺ CD8 T cells form a core population in inflamed human tissue. Sci. Transl. Med. 14, eabo0686 (2022).
Bouwman, A. C., Van Daalen, K. R., Crnko, S., Ten Broeke, T. & Bovenschen, N. Intracellular and Extracellular Roles of Granzyme K. Front. Immunol. 12, 677707 (2021).
Riley, J. S. & Tait, S. W. Mitochondrial DNA in inflammation and immunity. EMBO Rep. 21, e49799 (2020).
Shi, M. et al. Single-cell RNA sequencing shows the immune cell landscape in the kidneys of patients with idiopathic membranous nephropathy. Front. Immunol. 14, 1203062 (2023).
Xu, T. et al. Single-cell profiling reveals pathogenic role and differentiation trajectory of granzyme K+CD8+ T cells in primary Sjögren’s syndrome. JCI Insight 8, e167490 (2023).
Costa-da-Silva, A. C. et al. Salivary ZG16B expression loss follows exocrine gland dysfunction related to oral chronic graft-versus-host disease. iScience 25, 103592 (2022).
Huang, N. et al. SARS-CoV-2 infection of the oral cavity and saliva. Nat. Med. 27, 892–903 (2021).
Yin, H. et al. Sclerosing Sialadenitis Is Associated With Salivary Gland Hypofunction and a Unique Gene Expression Profile in Sjögren’s Syndrome. Front. Immunol. 12, 699722 (2021).
Yael Korem Kohanim, Avichai Tendler, Avi Mayo, Nir Friedman, & Uri Alon. Endocrine Autoimmune Disease as a Fragility of Immune Surveillance against Hypersecreting Mutants. Immunity 52, 872–884 (2020).
Warner, B. M. et al. Sicca Syndrome Associated with Immune Checkpoint Inhibitor Therapy. The Oncologist 24, 1259–1269 (2019).
Harris, C. R. et al. Array programming with NumPy. Nature 585, 357–362 (2020).
Virtanen, P. et al. SciPy 1.0: fundamental algorithms for scientific computing in Python. Nat. Methods 17, 261–272 (2020).
Seabold, S. & Perktold, J. Statsmodels: Econometric and Statistical Modeling with Python. in 92–96 (2010). doi:10.25080/Majora-92bf1922-011.
Pedregosa, F. et al. Scikit-learn: Machine Learning in Python. Mach. Learn. PYTHON.
McKinney, W. Data Structures for Statistical Computing in Python. in 56–61 (2010). doi:10.25080/Majora-92bf1922-00a.
Wolf, F. A., Angerer, P. & Theis, F. J. SCANPY: large-scale single-cell gene expression data analysis. Genome Biol. 19, 15 (2018).
Hunter, J. D. Matplotlib: A 2D Graphics Environment. Comput. Sci. Eng. 9, 90–95 (2007).
Waskom, M. seaborn: statistical data visualization. J. Open Source Softw. 6, 3021 (2021).
Polański, K. et al. BBKNN: fast batch alignment of single cell transcriptomes. Bioinformatics 36, 964–965 (2020).
McInnes, L., Healy, J. & Melville, J. UMAP: Uniform Manifold Approximation and Projection for Dimension Reduction. Preprint at http://arxiv.org/abs/1802.03426 (2020).
V. A. Traag, L. Waltman, & N. J. van Eck. From Louvain to Leiden: guaranteeing well-connected communities. Sci. Rep. 9, (2019).
Kanehisa, M. & Goto, S. KEGG: Kyoto Encyclopedia of Genes and Genomes.
Kanehisa, M., Sato, Y., Kawashima, M., Furumichi, M. & Tanabe, M. KEGG as a reference resource for gene and protein annotation. Nucleic Acids Res. 44, D457–D462 (2016).
Ashburner, M. et al. Gene Ontology: tool for the unification of biology. Nat. Genet. 25, 25–29 (2000).
The Gene Ontology Consortium. The Gene Ontology knowledgebase in 2023. Genetics 224, (2023).
Jassal, B. et al. The reactome pathway knowledgebase. Nucleic Acids Res. gkz1031 (2019) doi:10.1093/nar/gkz1031.
Liberzon, A. et al. Molecular signatures database (MSigDB) 3.0. Bioinformatics 27, 1739–1740 (2011).
Subramanian, A. et al. Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles. Proc. Natl. Acad. Sci. 102, 15545–15550 (2005).
Liberzon, A. et al. The Molecular Signatures Database Hallmark Gene Set Collection. Cell Syst. 1, 417–425 (2015).
Basisty, N. et al. A proteomic atlas of senescence-associated secretomes for aging biomarker development. PLOS Biol. 18, e3000599 (2020).
Hagberg, A. A., Schult, D. A. & Swart, P. J. Exploring Network Structure, Dynamics, and Function using NetworkX. (2008).
Kleshchevnikov, V. et al. Cell2location maps fine-grained cell types in spatial transcriptomics. Nat. Biotechnol. 40, 661–671 (2022).
Bankhead, P. et al. QuPath: Open source software for digital pathology image analysis. Sci. Rep. 7, 16878 (2017).
Pachitariu, M. & Stringer, C. Cellpose 2.0: how to train your own model. Nat. Methods 19, 1634–1641 (2022).
Palla, G. et al. Squidpy: a scalable framework for spatial omics analysis. Nat. Methods 19, 171–178 (2022).
S I Jang et al. Establishment of functional acinar-like cultures from human salivary glands. JDR 94, 304–11 (2015).
Klein, S., Staring, M., Murphy, K., Viergever, M. A. & Pluim, J. elastix: A Toolbox for Intensity-Based Medical Image Registration. IEEE Trans. Med. Imaging 29, 196–205 (2010).
Stringer, C., Wang, T., Michaelos, M. & Pachitariu, M. Cellpose: a generalist algorithm for cellular segmentation. Nat. Methods 18, 100–106 (2021).
Cutler, K. J. et al. Omnipose: a high-precision morphology-independent solution for bacterial cell segmentation. Nat. Methods 19, 1438–1448 (2022).
Feng, Y. et al. Spatial analysis with SPIAT and spaSim to characterize and simulate tissue microenvironments. Nat. Commun. 14, 2697 (2023).
Bates, D., Mächler, M., Bolker, B. & Walker, S. Fitting Linear Mixed-Effects Models Using lme4. J. Stat. Softw. 67, (2015).

There is NO Competing Interest.

NMSupplementaryDataTable1.docx
Supplemental Data Table 1
NMSupplementaryDataTable2.docx
Supplemental Data Table 2
NMSupplementalDataTable3.docx
Supplemental Data Table 3
SF1.pdf
Figure S1. (a) Barplot of changes in cell proportion in disease and SSA positivity across all cell types. (b) Dotplot of genes found to be specific to Leiden clusters and used to annotate those clusters. (c) Dotplots of differentially expressed genes in Sjögren's and the cell types those genes are differentially expressed in. (d) Dotplots of differentially expressed genes in SSA+/- and the cell types those genes are differentially expressed in. (e) Dotplots of genes differentially expressed across Sjögren's and anti-SSA positivity and the cell types those genes are differentially expressed in.
SF2.pdf
Figure S2. The quantification of cell types per spot according to Cell2Location. (a) Shown are CD8+ exhausted T cells and PRR4+CST3+WFDC2- seromucous acinar cells (SMACs) along with manually annotated ducts and immune infiltrates. Top tissue is from a healthy volunteer and bottom tissue is from a patient with a diagnosis for Sjögren’s. (b) Cell cooccurrence scores most different in disease. Immune cell colocalization with the PRR4+CST3+WFDC2-seromucous population is highest in SSA+ Sjögren's patients.
SF3.pdf
Figure S3. Interferon response increased across cell types with SSA positivity. (a) Cells in the index UMAP colored by patient’s diagnosis and SSA positivity. (b) Changes in index expression with disease diagnosis and SSA positivity. (c) Feature plots of IFN score in SSA+ and SSA- patients. Except in M2 macrophages and endothelial populations, IFN response score is low across cell types in SSA negative individuals. (d) Barplot of IFN response S-score across cell types in SSA positive and negative patients. IFN response increases in every cell population.
SF4.pdf
Figure S4. PhenoCycler-Fusion 35-plex spatial proteomics confirms changes in glandular composition and cellular interactions. (a) Dot plot of multiple protein expression of cellular phenotypes (left-side of dot plot; EPCAM to CD20) and cellular state markers (right side of dot plot; HHP3 to IFNG). (b) Cellular proportion changes in SjD compared with nonSjD. (c) Per-cell log2 fold-changes in protein expression in SjD compared with non-SjD were calculated for all phenotype and segmented cells. (d) Cellular neighborhoods in SjD compared with non-SjD. (e) CellChat ligand-receptor analysis corroborating neighborhood analysis and the importance of B- and T-lymphocyte co-stimulatory signaling in SjD. (f) Representative PhenoCycler-Fusion immunofluorescent images from a lymphocytic focus showing infiltration with T- and B-lymphocytes around ducts and acini.
SF5copy.pdf
Figure S5. (a) Evaluation of the induction of apoptosis after transfection of recombinant GZMK, GZMB, β-Gal, or Naive (mock). (b) Optimization of protein transfection with GZMK in HSG cells. (c) Transfection of GZMK and GZMB into NS-SV-TTAC exhibits time-dependent degradation. (d) GZMK transfection induced dose-dependent increases in pIRF3 by Western immunoblot. (e) THP1-macrophages were transfected with proteins (i.e., GZMK, GZMB, GFP, and mock) and analyzed by fluorescence microscopy showing GZMK induces cytosolic relocalization of mtDNA to the cytosol. (f) Engineered cell lines were used to measure IRF Pathway and Type I IFN induction using THP1-ISH Lucia and HEK-IFNa/b SEAP cells, respectively.
SF62.pdf
Figure S6 (a) Expression of acinar cell type markers in the HiPlex ISH data. (b) Statistically significant changes in cell proportions between SjD and nonSjD with PRR4+CST3+WFDC2- SMACs presented with the greatest loss in SjD. (c) Cell proportion changes in the ISH clusters with in nonSjD and SjD. The proportion of terminal mucous acini increases as seromucous acinar cell types are diminished. (d) Neighborhood analysis from the ISH data. As demonstrated with spRNAseq, the neighborhood analysis demonstrates involution with acinar cell types associating more frequently in SjD (greater mucous acinar cell connections and decreased seromucous to mucous acinar cell connections).
SF7.pdf
Figure S7. Surviving members of the lost acinar cell population show reduced expression of TFs found in development and growth. (a) Differential index S-scores in the lost acinar cell population reveal indices associated with differentiation and growth specific to cells found in non-Sjögren's patients. (b) KDE plot of one index, “Post Anal Tail Morphogenesis” from the Gene Ontology database. Cells from Sjögren’s patients fit within the distribution of cells from individuals without Sjögren’s disease.

Download PDF

Version 1

posted

You are reading this older preprint version

Read the latest preprint version →

GZMK+CD8+ T cells Target a Specific Acinar Cell Type in Sjögren’s Disease

Status:

Version 1

Abstract

Figures

INTRODUCTION

RESULTS

Sjögren's is defined by loss of a seromucous acinar cell population

DISCUSSION

LIMITATIONS

DECLARATIONS

Online Methods

References

Additional Declarations

Supplementary Files

Status:

Version 1