Detailed annotation of networks based on interpretation of gene set parameters and leading edge gene patterns
As described in the Methods section, the number of LE genes shared between GSs determined their position, interconnectivity and cluster membership in correspondence with all other GSs of a network. Thus, GSs in close proximity to each other share distinct similarities in their LE gene patterns. To resolve redundancies, commonly caused by GSs representing complex pathways, highly interconnected network areas require additional interpretation. The same accounts for large GSs annotated with terms too general to reflect the true theme shared by their LE members.
The basis for this curated annotation, written in italic type in Figures 3, 5A, 6A and 7A, is formulated upon analysis of (1) the LE members of each MCL cluster displayed as LE gene clouds generated by using a vector graphics–capable adaptation of the WordCloud Cytoscape plugin  (Figures 4, 5B, 6B and 7B), (2) each GS’s LE genes (Additional file 2) and (3) current literature–based interactome maps (Additional file 1: Figure S8). Additional file 3 comprises the networks displayed in Figures 3-7 as infinitely scalable and electronically searchable vector graphics, thereby allowing the visualization of network detail.
Transcriptional changes underlying themes being enriched during progression from pre- to subclinical Sjögren’s syndrome–like disease
In Figure 4, the LE gene cloud for Cluster_01-01, in combination with the percentage of each pathway covered by its LE members (Additional file 2; TAGS), points toward two pattern recognition receptors, namely, Toll-like receptor 3 (TLR3) and IFN-induced helicase C domain–containing protein 1 (IFIH1), also known as MDA5. Both these receptors are key molecules upstream of IFN regulatory factor 3 (Irf3) and signal transducer and activator of transcription 1 (STAT1) (Figure 4, Cluster_01-01). Also delineated by this cluster are upregulation of Tlr4 and its coreceptors Cd14 and lymphocyte antigen 96 (Ly96). These may, via their upregulated signaling cascade, deliver the strongest trigger for the observed canonical activation of nuclear factor κ-light-chain-enhancer of activated B cells (NF-κB) and mitogen-activated protein kinase 8 (MAPK8) [4, 25]. In addition, the gene nerve growth factor (Ngf), which encodes another important inducer of NF-κB, was upregulated, even though NGF’s crucial receptor p75 neurotrophin receptor (p75NTR)  was absent from the list of LE genes for RE_P75NTR SIGNALS VIA NF-KB; TAGS = 73% (Additional file 2).
Determining the effect of the interconnecting NC_TRAIL SIGNALING PATHWAY; TAGS = 63% and BI_MET PATHWAY; TAGS = 74% GSs is more difficult because of incomplete coverage of the different arms of the TNF-related apoptosis-inducing ligand (TRAIL) pathway  and the hepatocyte growth factor receptor (MET) pathway  by their respective LE genes (Additional file 2). However, their involvement in determining cell fate and proliferation is reflected by their central position in Figure 3.
The GSs and LE genes associated with Cluster_02 (Figures 3 and 4) suggest that INSR and IGFR1, via their shared downstream signaling cascade involving PI3K and AKT, upregulate mammalian target of rapamycin complex 1 (mTorC1) and mTorC2. The mTor system in turn is pivotal in determining cell fate . Increased autophagy is inferred by the presence of PI3K pathway members and several autophagy-related protein (Atg) encoding genes (Figure 4, Cluster_02).
In close proximity, Cluster_03-01 and Cluster_03-02 delineate cell-matrix adhesion complexes that transmit regulatory signals and mechanical forces (Figures 3 and 4) [30, 31]. Cluster_03-01, including GO_0031581 (HEMIDESMOSOME ASSEMBLY); TAGS = 73%, defines the hemidesmosome-mediated, laminin-5-dependent anchorage of epithelial cells’ intermediate filaments to the basal lamina of the ECM. Cluster_03-02 represents, in large part, signaling pathways that are activated by alterations in a cell’s immediate surroundings and are transmitted via actin cytoskeleton–anchored FAs, such as NC_SIGNALING EVENTS MEDIATED BY FOCAL ADHESION KINASE; TAGS = 68% . Thus, 31 of 77 genes annotated in the CC GO_0005925 (FOCAL ADHESION) GS were located in its LE (TAGS = 40%).
Matching the integrin genes (that is, Itgav, Itgb1, Itgb4 and Itgb5) with the dominant growth factor receptor genes (that is, Met, Insr, Igfr1, fibroblast growth factor receptor 1 (Fgfr1) and Tgfbr1) in the LE profiles displayed in Figure 4 suggests that integrin αvβ5, upstream of the enriched integrin-linked kinase signaling-associated GSs, provides a basis for IGFR1-integrin cross-talk . Similarly, αvβ5 and αvβ1 may allow for TGFBR1 signaling by collaborating with integrin pathways (Figure 3) . Tgfb1, Tgfb2 and Tgfb3, together with the TF Smad family member 2 (Smad2) and Smad4 downstream of Tgfbr1 and the negative feedback–associated Smad7, are all present in the LE of Cluster_04 (Figure 4). The presence of osteopontin (Spp1), another ligand of integrins αvβ1 and αvβ5 in the LE of Cluster_03-02 (Figure 4), indicates that FA maturation may also occur in relation to innate immune cells.
Supporting a critical role of FA remodeling during this transition from pre- to subclinical SS-like disease, calpain 1 (Capn1) and Capn2 (Figure 4), which regulate the dynamics of FA assembly and disassembly, are at the center of the two calpain-specific GSs (Figure 3). In addition, all other effector phases of non-muscle-cell movement are represented by GSs and LE genes of Cluster_03-02, Cluster_03-03 and the intercalated section of Cluster_04, respectively (Figures 3 and 4) .
The fourth biological theme shares 14.3% of its LE genes with Cluster_03-02 described above. This is due to molecular similarities between CC GO_0030055 (CELL-SUBSTRATE JUNCTION); TAGS = 41% and CC GO_0005913 (CELL-CELL ADHERENS JUNCTION); TAGS = 46%. Multiple LE genes belonging to the claudin and the occludin gene families further indicate increased formation of tight junctions at CC GO_0016327 (APICOLATERAL PLASMA MEMBRANE); TAGS = 51% . These two types of cell–cell junction complexes depend critically on CDH1 expressed by epithelial cells . Correspondingly, NC_E-CADHERIN SIGNALING IN THE NASCENT ADHERENS JUNCTION; TAGS = 60% and CDH1 anchorage-related GS GO_0017166 (VINCULIN BINDING); TAGS = 60% were significantly enriched and are mapped at the center of Figure 3. Furthermore, enrichment of GS GO_0030057 (DESMOSOME); TAGS 40%  delineates a third class of intercellular junction complexes associated with genes that are upregulated approximately 8 weeks prior to the onset of SS-like disease in C57BL/6.NOD-Aec1Aec2 mice.
Transcriptional changes underlying themes being depleted during progression from pre- to subclinical Sjögren’s syndrome–like disease
In Figure 5A, CC term GO_0031012 (EXTRACELLULAR MATRIX) is located at the center of the first major biological theme, becoming depleted during this time period (Figure 5A). It was the largest GS that yielded significance in this study, with 159 of its 320 members contributing to its significance (TAGS = 50%). The LE genes grouped in Cluster_01-01 delineate broad downregulation of genes encoding collagens of the ECM (Figure 5B) . Contributed by GS BI_INTRINSIC PATHWAY; TAGS = 73% and suggesting endothelial cell activation, this cluster also includes coagulation factor–encoding genes.
The LE gene cloud of Cluster_01-02 in Figure 5B lists genes associated with all categories of specialized ECM proteins . These include laminin (LAM) encoding subunits, such as Lama4 and Lamb1; proteoglycans, such as versican (Vcan); and glycoproteins, such as fibrillin 1 (Fbn1) and Fbn2. Genes coding for all matrix metalloproteinases (MMPs) capable of degrading collagens, as well as distinct members of the disintegrins and metalloproteinases with thrombospondin motif (ADAMTS) family, also contributed to the significance of the GSs grouped in Cluster_01-02 (Figure 5A). ADAMTS peptidases catalyze procollagens (for example, Adamts3) and inhibit angiogenesis (for example, Adamts5, Adamts8, Adamts9 and Adamts20) . Genes annotated as inducers of wingless-type mouse mammary tumor virus (MMTV) integration site family members (Wnt) (for example, Norrie disease (Ndp)), several Wnt genes (for example, Wnt1) and all Wnt1-inducible signaling pathway proteins (Wisp1, Wisp2 and Wisp3)  completed the LE of Cluster_01-02. These changes complement the marked and broad deceleration of ECM turnover as a potential consequence of the ongoing innate immune response and/or delayed conclusion of developmental processes in the SGs of C57BL/6.NOD-Aec1Aec2 mice.
The second theme delineates downregulation of genes associated with GSs annotating gap junction core proteins (for example, PF_00029 (CONNEXIN); TAGS = 67%) (Figures 5A and 5B) and thus represents the only class of cell-cell junctions not enriched at 8 weeks of age.
The third theme is dominated by genes coding for ligand-gated ion channels essential for neurotransmission (Figure 5B, Cluster_03) . The largest part of these genes encodes subunits of anionic Cys loop receptors (GABAA 12/12, GABAA-ρ 2/3 and GlyR 5/5), cationic Cys loop receptor subunits (serotonin-gated 5-HT3A and 5-HT3B and nicotinic ACh receptor 14/16 subunits), 18 of 20 ionotropic glutamate receptor subunits and ATP-gated channels P2X purinoceptors P2X1, P2X3, P2X5 and P2X6, as well as subsets of voltage-gated and acid-sensing potassium channels (for example, amiloride-sensitive cation channels 1 to 3 (ACCN1 to ACCN3) and ACCN5). The remaining clusters of this gene cloud represent mainly metabotropic receptors involved in sensory perception, whereas the LE genes associated with GO_0051971 (POSITIVE REGULATION OF TRANSMISSION OF NERVE IMPULSE; TAGS = 62%) also include inflammatory mediators such as IFNγ, tumor necrosis factor (TNF) and interleukin 6 (IL-6), all of which are known to decrease the threshold for nerve impulse generation (Additional file 2) .
Transcriptional changes underlying stabilization of subclinical disease between 8 and 12 weeks of age
LE genes associated with the continued depletion of GO_0005581 (COLLAGEN); TAGS = 46% encode all peptide chains for collagen type I, the most abundant collagen of the ECM, and collagen type III (Additional file 1: Figure S2B and Additional file 2). Collagen type IV, which has coverage of 67%, is associated with basal membranes .
Transcriptional changes underlying themes being enriched during transition from subclinical to overt Sjögren’s syndrome–like disease
In Cluster_01 of Figure 6A, GO_0009897 (EXTERNAL SIDE OF PLASMA MEMBRANE); TAGS = 33% interconnects the EM-associated GSs that delineate the adaptive effector immune response. The LE pattern of integrins (Cluster_01; Figure 6B) suggests an increase of CDH1 adhesive integrin αEβ7-expressing intraepithelial T cells, whereas CD11c, encoded by Itgax and Itgb2, points toward antigen-presenting cells (APCs) of myeloid origin . The latter represent the most probable source for the concomitant increase in transcription of various INFα-encoding genes (Infα-1, Infα-5, Infα-9 and Ifitm1) in the SGs of C57BL/6.NOD-Aec1Aec2 mice (Figures 6A and 6B).
The establishment of a NK cell population in the targeted tissues is supported by several distinct LE members (Cluster_01; Figure 6B). Cytotoxicity-triggering receptors NKG2-D type II integral membrane protein (Klrk1), Cd244 and its ligand encoded by UL16-binding protein 1 (Ulbp1) represent three key components of NK cells’ effector pathway. In contrast, the NK cell receptor complex encoded by killer cell lectinlike receptor subfamily D member 1 (Klrd1) and G member 1 (Klrg1) exert a regulatory anticytotoxic effect [38, 39]. The LE genes Cd244 and Cd48, in conjunction with Cd2 and intercellular adhesion molecule 2 (Icam2), may further suggest regulation of CD8+ T cells by NK cells. Expression of major histocompatibility complex (MHC) and MHC-related genes, however, were skewed toward upregulation of MHC class II (MHCII) and MHCII invariant chain (Cd74) expression (Figure 6B).
The chemokine receptor-ligand profile characterizes emigration of multiple APC and lymphocyte populations (Cxcr4:Cxcl12), as well as reinforced recruitment of T-helper type 1 (TH1) cells, NK cells and plasmacytoid dendritic cells (Cxcl3:Cxcl9/Cxcl10 and Ccr7:Ccl19) (Figure 6B) . Immune cell homing may also be facilitated by increased expression of LE genes that encode mucosal vascular addressin cell adhesion molecule 1 (MAdCAM-1) and lymphocyte function–associated antigen 1 (Lfa-1), encoded by Itgal and Itgb2 and Icam1.
ICAM1 and LFA-1 ligation is also critical for Cd28-dependent T-cell activation . The pattern of LE genes encoding costimulatory molecules assigns importance to both the activating Cd28-dependent pathway and the inhibitory cytotoxic T-lymphocyte antigen 4 (Ctla4)–dependent pathway (Figure 6B; Cluster_01). Regarding the T-cell-associated central component of the immunological synapse, T-cell receptor (TCR) accessory proteins (for example, Cd3), TCR coreceptors Cd4 and Cd8 and TCR-associated molecules (for example, Cd45 (Ptprc)) are also covered by the LE gene cloud of Cluster_01. The concomitant upregulation of Il2, Il2rb and Il2rg, as well as the presence of Cd69, represent effects downstream of T-cell activation . Regulating activation of T-cell effector lineages at this stage may thereby involve the two LE gene B7 family members B and T lymphocyte attenuator (Btla) and Tnf receptor superfamily 18 (Tnfrsf18) (Figure 6B).
With respect to late costimulatory signals, Cd40:Cd40lg and inducible T-cell costimulator (Icos):IcosL are the receptor-ligand pairs present in the LE of Cluster_01 (Figure 6B). Both these systems, together with LE-gene Il4, are critical for mounting effective TH2 responses .
B-cell-specific genes (for example, immunoglobulin heavy constant μ (Igh-6), Cd79a, Cd79b, Cd19 and Cd22) are highly represented in the LE gene cloud of Cluster_01 (Figure 6B). Increased transcription of Tnfrsf13C (that is, Baffr), Tnfrsf17 (that is, Bcma) and Tnfrsf13B (that is, Taci), together with their common ligand Tnfsf13b (that is, Baff), as well as the activation-induced cytidine deaminase gene (Aicda) (Figure 6B and Additional file 2; LE gene list for KE_04672 (INTESTINAL IMMUNE NETWORK FOR IGA PRODUCTION); TAGS = 60%), further indicates strong signaling for survival, proliferation and differentiation of B cells in SGs marked by overt disease [38, 40].
The second biological theme enriched by 16 weeks of age pertains to neurotransmission and marks a partial reversal of changes that occurred earlier in the disease course. Of the 78 LE genes defining enrichment at this later stage (Figure 6B; Cluster_02), 40 were previously associated with depletion of GSs concerning neurotransmission at 8 weeks of age (Figure 5B; Cluster_03). Reinitiating transcription are mainly ACh-, GABA- and Gly-gated ionotropic receptors coding genes (Additional file 2). Overlaps were also found for genes encoding receptors for dopamine (that is, Drd1a) and substance P (that is, Tacr1). Unique to enrichment at 16 weeks of age were genes encoding for metabotropic receptors specific for ACh (that is, Chrm4) and somatostatin (that is, Sstr1, Sstr2, Sstr4, Sstr5) (Figure 6B and Additional file 2) .
Transcriptional changes underlying themes being depleted during transition from subclinical to overt Sjögren’s syndrome–like disease
GSs (Figure 7A) and their LE genes (Figure 7B) depleted and downregulated, respectively, during this time period, predominantly signify the reversal of previous enrichments in FAs and cell-cell junction–associated GSs observed at 8 weeks of age (Figures 3 and 4). Pairwise comparison of the overlapping GSs revealed that, on average, 54% of LE genes contributing to depletion at 16 weeks of age also contributed to these GSs’ prior enrichment at 8 weeks of age. The highest percentage of LE members following this pattern was identified for NC_SIGNALING EVENTS MEDIATED BY FOCAL ADHESION KINASE, with 79%, and the lowest percentage was found for GO_0043296 (APICAL JUNCTION COMPLEX), with 43% (Additional file 2). The LE genes not included in these LE overlaps did not define additional biological themes, but instead contributed predominantly to the increased average coverage of the EM-related GSs at 16 weeks (TAGS 45%) compared to 8 weeks of age (TAGS 38%) (Additional file 2).