Abstract

   Celiac disease (CD) is an autoimmune disease in which intestinal
   inflammation is induced by dietary gluten. The means through which
   gluten-specific CD4^+ T cell activation culminates in intraepithelial T
   cell (T-IEL) mediated intestinal damage remain unclear. Here, we
   performed multiplexed-single cell analysis of intestinal and
   gluten-induced peripheral blood T cells from patients in different
   celiac disease states and healthy controls. Untreated, active CD (ACD)
   and potential CD (PCD) were associated with an enrichment of activated
   intestinal T cell populations including CD4^+ follicular T-helper
   (T[FH]) cells, regulatory T cells (Tregs), and natural CD8^+ αβ and γδ
   T-IELs. Natural CD8^+ αβ and γδ T-IELs expressing activating natural
   killer cell receptors (NKRs) exhibited a distinct TCR repertoire in CD,
   which persisted in patients on a gluten-free diet (GFD) without
   intestinal inflammation. Our data further show that NKR-expressing
   cytotoxic cells, which appear to mediate intestinal damage in CD, arise
   from a distinct NKR-expressing memory population of T-IELs. Following
   gluten ingestion, both αβ and γδ T cell clones from this memory
   population of T-IELs circulated systemically along with gluten-specific
   CD4^+ T cells and assumed a cytotoxic and activating NKR-expressing
   phenotype. Collectively, these findings suggest that cytotoxic T cells
   in CD are rapidly mobilized in parallel with gluten-specific CD4^+ T
   cells following gluten ingestion.

One Sentence Summary:

   Gluten ingestion induces cytotoxic T-IELs from an NKR-expressing
   natural memory population in parallel with gluten-specific CD4^+ T
   cells.

INTRODUCTION

   Celiac disease (CD) is a common autoimmune disease characterized by
   dietary gluten-induced small intestinal mucosal injury in HLA-DQ2 or
   HLA-DQ8 individuals ([64]1). Beyond a lifelong gluten-free diet (GFD),
   current treatment options for CD are limited. In addition to HLA
   association, CD shares many hallmarks with other autoimmune diseases,
   including non-MHC genetics, immune-mediated tissue destruction, and
   auto-antibody production ([65]2). However, CD is one of the only
   autoimmune diseases for which the triggering antigen is well-defined
   and can be controlled, and thus provides an opportunity to study
   autoimmune T cell responses in different disease states.

   How gluten-dependent CD4^+ T cell activation culminates in intestinal
   epithelial cell (IEC) destruction remains incompletely understood
   ([66]3–[67]5). Gluten-specific CD4^+ T cell immunity is necessary for
   disease but is insufficient for IEC destruction in humans or in mouse
   models ([68]3, [69]6, [70]7). Instead, intestinal intraepithelial T
   cells (T-IEL) are thought to be the end effector cells that drive
   tissue damage in CD ([71]5). In patients with potential CD (PCD),
   characterized by positive serum celiac antibodies in the absence of
   mucosal damage, CD4^+ T cell responses are intact but appear to be
   uncoupled from T-IEL-mediated tissue damage ([72]6, [73]7). A
   substantial portion of PCD patients eventually develop CD, while others
   remain disease-free for years despite continued gluten consumption
   ([74]7).

   An increase in T-IELs is a pathologic hallmark of CD, and these cells
   persist in celiac mucosa despite GFD ([75]8). T-IELs in CD have not
   been demonstrated to recognize gluten and their cytolytic activity can
   be induced irrespective of antigen specificity ([76]8). They also
   express activating natural killer cell receptors (NKRs), which are
   considered critical for IEC destruction ([77]4, [78]9–[79]12). These
   findings have led to the prevailing model, which proposes that T-IELs
   within tissues are activated by stress and inflammation, presumably
   initiated by gluten-specific CD4^+ T cells ([80]4, [81]8).

   Re-introducing dietary gluten to celiac patients on a GFD rapidly
   induces gluten-specific CD4^+ T cells in the peripheral blood that
   express activation markers as well as the intestinal-trafficking β7
   integrin receptor ([82]13). This phenomenon likely represents the
   induction of a memory immune response to gluten, as seen in the
   peripheral blood en route from gut-associated lymphoid tissues to the
   intestine. In addition to gluten-specific CD4^+ T cells, oral gluten
   challenge leads to a parallel induction of activated, gut-homing CD8^+
   αβ and γδ T cells in the peripheral blood ([83]14). These cells
   phenotypically resemble IELs, express CD103, and share a TCR repertoire
   with intestinal T-IELs ([84]14, [85]15). However, the significance of
   these findings is unclear as the functional significance of
   gluten-induced circulating CD8^+ αβ and γδ T cells has not been
   established.

   Here, we leveraged the specific opportunities afforded to us with
   individuals who are in different states of CD and can be re-exposed to
   gluten to dissect mechanisms of autoimmunity using single-cell
   multiomics. We analyzed T cell responses and TCR repertoire at the site
   of T cell-driven autoimmunity in the presence or absence of the gluten
   trigger. Through oral gluten challenge, we were able to study the
   induction of CD through tracking intestinal T cell clones in the
   peripheral blood.

RESULTS

Single-cell profiling of duodenal T cells in CD

   We obtained duodenal biopsies from 11 patients with newly diagnosed or
   active CD (ACD), 19 celiac patients treated with a gluten-free diet
   (GFD), 7 PCD patients, and 17 healthy controls ([86]Fig. 1A). Marsh
   scores, which range from 0 (normal) to 3 (most severe), are used to
   grade disease severity in celiac biopsies through assessment of villous
   height and immune infiltration. All ACD patients had modified Marsh
   3a-3c scores, indicative of IEL infiltration and villous atrophy, while
   GFD/PCD patients had Marsh 0 (17/19 GFD, 6/7 PCD patients) or Marsh 1
   (2/19 GFD, 1/7 PCD patients) scores ([87]table S1) ([88]16). PCD
   patients were HLA-DQ2/8 with positive celiac serology
   (anti-tissue-transglutaminase 2 (TG2) IgA) in the absence of intestinal
   damage by histology ([89]table S1) ([90]6, [91]7). Live CD3^+ cells
   were flow sorted, and scRNA-seq, scTCR-seq, and 204 panel CITE-seq
   (cellular indexing of transcriptomes and epitopes by sequencing) were
   performed using the droplet-based 10X Genomics platform ([92]Fig. 1A,
   [93]fig. S1A, [94]tables S1 and [95]S2). Following quality control
   preprocessing to remove unidentifiable cells, multiplets, and
   low-quality cells, 41,285 intestinal T cells remained for analysis.

Fig 1. Distinct CD4^+ Treg and T[FH] cell phenotypes enriched in untreated CD
mucosa.

   Fig 1.
   [96]Open in a new tab

   (A) Schematic representation of experimental design and sorting
   strategy. (B-D) UMAP plot of intestinal T cells colored by (B) patient
   condition, (C) annotated cell type, and (D) identified clusters with
   T[FH] and Treg subclusters highlighted (encircled). (E) Quantification
   of CD4^+ T cells within T[FH] and Treg subclusters in control (n = 17),
   PCD (n = 7), ACD (n = 11), and GFD (n = 19) mucosa. (F) Heatmap of mean
   RNA (left) and ADT (right; ADT-stained samples only) expression of
   selected genes and cell surface markers identified through differential
   expression across T[FH] and Treg subclusters and remaining CD4^+ T
   cells. (G) Heatmap of mean RNA expression of select genes, identified
   through differential expression, within T[FH] subcluster cells across
   states of CD. (H) Scatterplot comparing percentage of T[FH] and Treg
   subcluster cells in untreated CD (n = 18) mucosa. Shaded region
   indicates 95% confidence interval following linear regression.
   Two-sided Wilcoxon rank sum test (E). *P < 0.05, **P < 0.01, ***P <
   0.001, and ****P < 0.001. Pearsons rank correlation coefficient (H).

   To resolve CD3^+ T cell subtypes, we used scTCR-seq, CITE-seq antibody
   derived tags (ADTs), and reference-based scRNA-seq cell type
   recognition to annotate CD4^+ TCRαβ (CD4^+), CD8^+ TCRαβ (CD8^+), and
   γδ T cells corroborating relative percentages determined during FACS
   sorting ([97]fig. S1, [98]A to [99]C) ([100]17–[101]19).

   Transcriptomic data from all patient samples was integrated using
   canonical correlation analysis, and dimensionality reduction was
   performed to visualize unsupervised clustering of single cells as
   Uniform Manifold Approximation and Projection (UMAP) embeddings
   ([102]Fig. 1, [103]B and [104]C, [105]fig. S1, [106]D to [107]F). Eight
   clusters were identified and annotated by means of distinct
   transcriptomic signatures and cell surface marker expression ([108]Fig.
   1D, [109]fig. S1G, [110]table S3). CD4^+ T cell populations were
   identified as T-helper 1 (Th1), T-helper 17 (Th17), and a shared
   cluster of follicular T-helper (T[FH]) and regulatory T cells (Tregs)
   (T[FH]/Treg cluster) ([111]Fig. 1, [112]C and [113]D). CD8^+ T cell
   populations were classified as effector memory (T[EM]), resident memory
   IELs (T[RM]), and natural IELs (natural-memory; T[NM],
   natural-effector; T[NE]) ([114]20), which shared phenotype with γδ T
   cells ([115]Fig. 1, [116]C and [117]D). We also observed a small
   cluster of naïve T cells, containing both CD4^+ and CD8^+ T cells
   ([118]Fig. 1, [119]C and [120]D).

Induction of activated CD4^+ T cell populations in CD

   The study of CD4^+ T cells in CD has largely focused on gluten-specific
   T cells, which were reported to have a distinctive T[FH] phenotype
   ([121]21). Several groups have also reported an increase in FOXP3^+ T
   cells in active celiac mucosa ([122]22, [123]23).

   We found that a distinct population of activated CD4^+ T cells
   (T[FH]/Treg cluster), comprised of both T[FH]-phenotype and
   FOXP3^+CD25^+CD4^+ Tregs, was highly enriched within CD mucosa relative
   to controls ([124]Fig. 1, [125]B to [126]D, [127]fig. S2A, [128]table
   S4). The T[FH]/Treg cluster contained activated T cells with elevated
   expression of CD39, TIGIT, CCR4, and ICOS and low expression of CD49a,
   IL7R, and CD73 ([129]fig. S2B). We subclustered the T[FH]/Treg cluster,
   which yielded a subcluster of Tregs and a subcluster of T[FH]-phenotype
   T cells ([130]Fig. 1D). Both activated T[FH] and Treg-phenotype cells
   were enriched in ACD and PCD relative to GFD and controls ([131]Fig.
   1E).

   T[FH]-subcluster CD4^+ T cells matched the previously reported,
   pro-inflammatory phenotype of gluten tetramer-reactive cells ([132]21),
   expressing CXCL13, CD200, CD84, IFNG, IL-21, PD-1, CD39, and
   β7-integrin ([133]Fig. 1F). We also observed an upregulation of
   transcription factors (TFs) BCL6 and LEF1, both shown to be critical
   for T[FH] differentiation ([134]Fig. 1F) ([135]24, [136]25). Notably,
   T[FH]-subcluster CD4^+ T cells in inflamed ACD mucosa had elevated
   expression of IFNG relative to T[FH]-subcluster cells in non-inflamed
   PCD mucosa, which highly expressed IL-21 ([137]Fig. 1G).

   While FOXP3-expressing CD4^+ T cells have now been well-described in
   CD, their functional significance is unclear ([138]26, [139]27).
   Consistent with a regulatory role, cells within the Treg subcluster
   exhibited an activated and suppressive phenotype, expressing CD39,
   GITR, and CD27 in addition to IL10, CD25, CTLA4, and PD-L1 ([140]Fig.
   1F) ([141]28, [142]29). This CD39^+TIGIT^+IL10^+ phenotype matches that
   of Tregs recently described in Crohn’s disease ([143]30).

   To investigate the regulatory networks characterizing the distinct
   T[FH]and Treg phenotypes, we performed gene co-expression analysis
   followed by the identification of TF cis-regulatory motifs ([144]31,
   [145]32). 310 active TF modules were found and scored individually for
   each cell. T[FH] and Treg subclusters shared several active TF modules
   including LEF-1, BATF, ETV7 and c-MAF ([146]fig. S2C). In contrast,
   other CD4^+ T cell populations used RORγt, RORα, KLF3, and RUNX2 TF
   networks ([147]fig. S2C).

   It has been reported that circulating gluten-induced CD4^+ T cells in
   CD co-expressed FOXP3 and CD39, suggesting that both effector and
   regulatory CD4^+ T cell responses in CD may be gluten-specific
   ([148]26). We investigated the TCR repertoire across T[FH] and Treg
   subclusters. In CD mucosa, T[FH]-phenotype cells exhibited elevated
   usage of V-gene biases previously identified in gluten-specific CD4^+ T
   cells ([149]fig. S2D) ([150]33–[151]36). We next applied the GLIPH2
   algorithm (grouping of lymphocyte interaction by paratope hotspots) to
   cluster TCRs by potential epitope specificity ([152]37). Within the
   T[FH] cluster, GAET, LRGG, R(X)PSTDT, and QNTG motifs were identified
   and shared across multiple patients ([153]table S5). Within the LRGG
   motif, we observed TCRs utilizing arginine at CDR3β position 5,
   previously implicated in binding the immunodominant gliadin α-II
   epitope ([154]35, [155]38) ([156]table S5). Collectively, these results
   corroborate the presence of gluten-specific CD4^+ T cells within the
   T[FH] population.

   Within the Treg subcluster, different CDR3β motifs, including S(X)YRE,
   GALS, and SLGLAGA(X)E, were enriched, and there was negligible sharing
   between T[FH] and Treg subclusters, suggesting they are selected by
   different antigens and arguing against the gluten specificity of Tregs
   ([157]table S5). Despite distinct TCR repertoires, we observed a strong
   correlation between the percentage of CD4^+ T cells exhibiting T[FH]
   and Treg phenotypes within untreated CD patients, suggesting both
   populations are induced together in CD ([158]Fig. 1H). Accordingly, in
   PCD, where gluten-induced CD4^+ T cell responses are intact, both T[FH]
   and Treg populations were enriched ([159]Fig. 1, [160]E and [161]H).

Natural CD8^+ αβ and γδ T-IELs exhibit an NK-polarized phenotype

   We next examined populations of intestinal CD8^+ and γδ T cells
   ([162]Fig. 1, [163]C and [164]D). We observed the presence of effector
   memory CD8^+ T cells (T[EM] cluster; CD28^+CD44^+) in addition to
   several populations of CD103^+ T-IELs, which we termed resident memory
   (T[RM], CD73^+IL7R^+), natural-memory (T[NM,] CD38^+CD39^+NKR^+IL7R^+)
   and natural-effector (T[NE,] CD38^+CD39^+NKR^+IL7R^−) ([165]Fig. 1D,
   [166]Fig. 2, [167]A to [168]C, [169]fig. S1G). Notably, T[NM] and T[NE]
   populations contained both CD8^+ and “naturally-occurring” γδ T-IELs
   ([170]Fig. 1C) ([171]20).

Fig 2. Natural αβ and γδ T-IELs exhibit NK-reprogramming in CD.

   Fig 2.
   [172]Open in a new tab

   (A) Dotplot heatmap showing mean expression of select ADT markers
   across identified clusters (ADT-stained samples only). (B) ADT
   expression profile of selected cluster-defining cell surface markers
   (ADT-stained samples only). (C) RNA expression profile of selected
   cluster-defining marker genes. (D) UMAP plot of intestinal T cells
   colored by AUC scores of individual cells for selected canonical
   pathway activities. (E) Quantification of CD8^+ and γδ T cells across
   identified T[RM], T[EM], and natural clusters in control (n = 17), PCD
   (n = 7), ACD (n = 11), and GFD (n = 19) mucosa. Two-sided Wilcoxon rank
   sum test (E). *P < 0.05, **P < 0.01, and ***P < 0.001.

   T[RM] cluster T-IELs expressed CD73, CD127, and CD103, in combination
   with high expression of KLRC1 (NKG2A) ([173]Fig. 2, [174]A to [175]C).
   Consistent with the capacity of NKG2A to restrict T cell activation
   ([176]39), T[RM] T-IELs scored highly for gene ontology (GO) pathways
   involved in negatively regulating both T cell immunity and cytotoxicity
   ([177]Fig. 2D). These cells comprised the largest fraction of CD8^+ T
   cells from healthy controls, and as such, represent the “baseline”
   CD8^+ T-IEL population within healthy duodenum ([178]Fig. 2E).

   T[EM] cluster cells exhibited a highly cytotoxic phenotype, expressing
   GZMK, EOMES, NKG2D, PD-1, and KLRG1 with limited CD103 expression,
   suggestive of their presence within the lamina propria (LP) as opposed
   to the epithelium ([179]Fig. 2, [180]A to [181]C, [182]fig. S3A)
   ([183]20). Cytotoxic CD8^+ T cells expressing EOMES and KLRG1 have been
   implicated in mediating tissue damage in inflammatory bowel disease
   (IBD) ([184]32, [185]40). However, we did not observe an enrichment of
   T[EM] cells across CD disease states suggesting they do not similarly
   drive tissue damage in CD ([186]Fig. 2E, [187]fig. S3B).

   T-IELs that we termed “natural” highly expressed activating NKRs,
   including KLRC2 (NKG2C), KLRC3 (NKG2E), and KLRC4 (NKG2F), in addition
   to killer immunoglobulin-like receptors (KIRs), including KIR2DL4,
   analogous to previously described “natural” murine T-IELs ([188]20,
   [189]39) ([190]Fig. 2C, [191]fig. S3, [192]C and [193]D). Natural T-IEL
   populations exhibited GO signatures consistent with NK polarization
   ([194]3, [195]4, [196]9–[197]12), including NK cell-mediated immunity,
   cytokine production, and NKR-driven interactions with non-classical MHC
   class I molecules ([198]39) ([199]Fig. 2D, [200]fig. S3A). Notably, we
   observed an increased frequency of “natural” CD8^+ T-IELs in all celiac
   disease states, including PCD ([201]Fig. 2E).

Natural T-IEL populations define stages of CD

   Given the importance of NKR-expressing T-IEL populations in CD ([202]3,
   [203]4, [204]9–[205]12), we re-clustered 28,380 CD8^+ and γδ T-IELs,
   and identified ten distinct clusters of CD8^+ and γδ T-IELs ([206]Fig.
   3, [207]A and [208]B, [209]fig. S4, [210]A to [211]C). Subsequent
   analysis was performed on T-IELs, which were distinguished through
   expression of CD103 ([212]20).

Fig 3. Distinct Natural T-IEL phenotypes define CD progression.

   Fig 3.
   [213]Open in a new tab

   (A-B) UMAP plot of intestinal CD8^+ and γδ T-IELs colored by (A)
   patient condition and (B) identified clusters. (C) Heatmap of mean RNA
   expression of selected genes identified through differential expression
   across natural T-IEL clusters and remaining CD8^+ and γδ T-IELs. (D)
   Quantification of CD8^+ and γδ T-IELs across natural T-IEL clusters in
   control (n = 17), PCD (n = 7), ACD (n = 11), and GFD (n = 19) mucosa.
   (E) Heatmap of mean RNA expression of select genes, identified through
   differential expression, within T[NE] T-IELs across states of CD. (F)
   Volcano plot of differentially expressed genes (log2FC expression)
   between γδ (left) and CD8^+ (right) T[NE] T-IELs from (n = 37) CD
   patients. (G) Scatterplot comparing percentage of CD8^+ and γδ T-IELs
   in CD (n = 37) mucosa across select natural clusters. Shaded region
   indicates 95% confidence interval following linear regression.
   Two-sided Wilcoxon rank sum test (D). *P < 0.05, **P < 0.01, and ***P <
   0.001. Pearsons rank correlation coefficient (G).

   We observed substantial changes within NKR-expressing populations
   across different stages of CD. Within the previously described
   natural-memory population, we observed three distinct clusters, which
   we termed Innate-like (Innate), natural-memory 1 (T[NM]1), and
   natural-memory 2 (T[NM]2) ([214]Fig. 3B). The most abundant population
   of γδ T-IELs in healthy mucosa (T[NM]1 cluster) expressed TGFB1
   (TGF-β), AREG (Amphiregulin), the TFs ID2 and ID3, and 4–1BB,
   suggestive of immuno-regulatory and tissue repair-promoting roles
   ([215]Fig. 3, [216]B to [217]D, [218]fig. S4D) ([219]41). Therefore, as
   T[RM] (described above) represents the baseline CD8^+ T cell phenotype
   ([220]Fig. 2, [221]A to [222]E, [223]fig. S4C), T[NM]1 similarly
   represents the baseline γδ T cell phenotype in healthy duodenum.

   Despite being the baseline γδ T cell phenotype, there were clearly
   CD8^+ T cells present within the T[NM]1 cluster ([224]Fig. 3D,
   [225]fig. S4A). In contrast to γδ T cells, there were relatively few
   T[NM]1-phenotype CD8^+ T cells in healthy control mucosa compared to
   celiac mucosa ([226]Fig. 3D). However, both T[NM]1-phenotype CD8^+ and
   γδ T cells were less abundant within ACD mucosa compared to
   non-inflamed PCD and GFD mucosa ([227]Fig. 3D). This population, which
   has not been well characterized previously, represents a distinct
   memory T-IEL population, as the TF ID3 is critical for long-lived
   memory CD8^+ T cell generation ([228]42).

   In healthy mucosa, we also observed the presence of γδ T-IELs (Innate
   cluster) highly expressing TYROBP (DAP12), CD247, and SH2D1B in
   addition to elevated expression of CD122 and natural cytotoxicity
   receptors (NCRs) NKp46 and NKp44 ([229]Fig. 3, [230]A to [231]D,
   [232]fig. S4D). These innate-like γδ T-IELs were diminished in CD
   mucosa and highly utilized TRGV4, matching a previously described
   population of Vδ1^+ γδ T-IELs lost in CD mucosa with disease onset
   ([233]Fig. 3D, [234]fig. S4E) ([235]43).

   It was previously reported that a distinct Vδ1^+ γδ T-IEL population
   persists in celiac patients despite GFD, consistent with “permanent
   reshaping” or “scarring” of T-IELs in CD ([236]43). Consistent with
   this report, T-IELs exhibiting a similar NKR-expressing phenotype
   (T[NM]2 cluster), were significantly increased in both ACD and GFD,
   irrespective of gluten exposure or intestinal damage ([237]Fig. 3,
   [238]C and [239]D). In addition to Vδ1^+ γδ T-IELs, we clearly observed
   both CD8^+ and non-Vδ1^+ (~52% of T[NM]2 γδ T cells) T-IELs within this
   population ([240]Fig. 3D, [241]fig. S4, [242]A and [243]F). Both αβ and
   γδ T-IELs within this T[NM]2 population exhibited similar enrichment in
   CD mucosa relative to controls ([244]Fig. 3D), and highly expressed
   activating NKRs, KLRC2, KLRC3, and KLRC4, and markers indicative of
   cytotoxic potential including HCST (DAP10) and KLRK1 (NKG2D) ([245]Fig.
   3, [246]C and [247]D). However, despite their cytotoxic potential, this
   population was present and enriched within GFD mucosa which lacked any
   discernable tissue damage, suggesting that they were not actively
   mediating cytotoxicity ([248]Fig. 3D).

   Within inflamed ACD mucosa, we observed a population of
   pro-inflammatory, IFNG-expressing CD8^+ and γδ T-IELs, which we termed
   natural-effector (T[NE]) ([249]Fig 3, [250]A to [251]D). T[NE] T-IELs
   exhibited a highly cytotoxic and activated phenotype through expression
   of activating NKR KLRC2 in addition to PRF1, NKG7, GZMB, CD71, CD52,
   and CD161 ([252]Fig. 3C, [253]fig. S4D). Their highly cytotoxic
   phenotype and enrichment in ACD mucosa strongly suggest that T[NE]
   T-IELs are responsible for IEC destruction in CD ([254]Fig. 3, [255]A
   to [256]D). Interestingly, T[NE] T-IELs were also enriched in PCD
   relative to GFD or control mucosa, in the absence of histologic tissue
   damage ([257]Fig. 3D). Compared to ACD, T[NE] T-IELs in PCD were
   present at lower frequency and exhibited reduced expression of PRF1,
   IFNG, KLRC2, KLRC3, and KLRC4 relative to ACD T[NE] T-IELs ([258]Fig.
   3, [259]D and [260]E). These findings imply additional mechanisms of
   T[NE] activation may differentiate ACD from PCD ([261]3–[262]5,
   [263]11, [264]44).

   Aside from their TCR genes, CD8^+ and γδ T-IELs within natural T-IEL
   populations were phenotypically indistinguishable, suggestive of
   redundant roles ([265]Fig. 3F, [266]fig. S4G, [267]table S6). Within
   T[NE], T[NM]1, and T[NM]2 T-IEL clusters, we observed a moderate to
   strong correlation between the percentage of CD8^+ and γδ T-IELs
   exhibiting each phenotype, implying these NKR-expressing CD8^+ αβ and
   γδ populations are induced together in CD ([268]Fig. 3G).

Distinct natural T-IEL TCR αβ and γδ repertoires in CD

   We next analyzed scTCR-seq data to investigate both the CD8^+ αβ and γδ
   TCR repertoires across different CD states. Across identified
   populations, we observed substantial variations in the degree of clonal
   expansion in both CD8^+ and γδ T cells ([269]Fig. 4, [270]A and [271]B,
   [272]fig. S5, [273]A to [274]D). Interestingly, in healthy mucosa,
   “baseline” populations (CD8^+, T[RM] 1–4; γδ, T[NM]1), contained
   elevated fractions of highly expanded TCR clones (33.1% αβ, 46.4% γδ)
   ([275]Fig. 4, [276]A and [277]B). Conversely, natural (T[NM]1, T[NM]2,
   and T[NE]) T-IELs, the predominant αβ and γδ T-IEL populations in CD
   mucosa, were clonally diverse in CD, affirming previously reported
   findings of increased clonal diversity in CD ([278]Fig. 4, [279]A and
   [280]B) ([281]45, [282]46).

Fig 4. Distinct TCR repertoire of Natural αβ and γδ T-IELs in CD.

   Fig 4.
   [283]Open in a new tab

   (A-B) (A) CDR3β and (B) CDR3δ clonal expansion quartiles across and
   select clusters within healthy (HC) and CD mucosa. (C-D) UpSet plot
   displaying shared (C) αβ or (D) γδ clonotypes across T-IEL clusters
   within untreated CD mucosa. Barplot on left indicates count of unique
   (C) CDR3β or (D) CDR3δ within clusters. Shared clonotype pairing is
   visualized as a black circle joined by black lines to indicate sharing
   sets, with number of unique shared (C) CDR3β or (D) CDR3δ between
   clusters (independent of expansion) annotated above the barplot.
   Highlighted bars indicate sharing between natural clusters. (E-F) Chord
   diagram displaying TRGV-TRDV pairing usage across natural T-IEL
   clusters within (E) control and (F) CD mucosa.

   We examined TCR sharing across clusters within untreated CD mucosa,
   where we observed ~15% of unique clonotypes were present across
   multiple clusters (570/3,737 shared/total αβ clonotypes; 151/971
   shared/total γδ clonotypes) ([284]Fig. 4, [285]C and [286]D). Broadly,
   within untreated CD mucosa containing all natural, T[RM], and T[EM]
   populations, natural T-IELs showed extensive TCR repertoire sharing
   with other natural populations, as opposed to T[RM] and T[EM]
   populations ([287]Fig. 4, [288]C and [289]D). For both CD8^+ and γδ
   T-IELs, we observed frequent TCR sharing between T[NM]1 and T[NM]2
   clusters (Most frequent αβ, 140 clonotypes, Most frequent γδ, 81
   clonotypes) and T[NM]2 and T[NE] populations (2^nd most frequent αβ, 50
   clonotypes; 2^nd most frequent γδ, 18 clonotypes) ([290]Fig. 4, [291]C
   and [292]D). TCRγδ sharing between innate and natural (T[NM]1, T[NM]2,
   and T[NE]) populations was very limited, suggesting these populations
   arise independently ([293]Fig. 4D). Thus, our TCR data show that in CD,
   natural T-IEL populations, including T[NM]1, T[NM]2 and T[NE],
   extensively shared TCR repertoire with one another, but not with T[RM],
   T[EM], or Innate populations. Accordingly, we frequently observed
   clones present within each of the T[NM]1, T[NM]2 and T[NE] clusters for
   both αβ (27 clonotypes) and γδ (11 clonotypes) T-IELs ([294]Fig. 4,
   [295]C and [296]D).

   Among γδ T-IELs, healthy controls preferentially utilized TRDV3 (64.5%
   Control; 32.2% CD) compared to CD patients, who preferentially utilized
   TRDV1 (25.1% Control; 50.9% CD) ([297]fig. S5, [298]E to [299]G). In
   both healthy and CD mucosa, innate T-IELs exhibited a semi-invariant
   usage of both TRGV4/TRDV1 and TRGV4/TRDV3, and were found at levels
   consistent with those of previously described butyrophilin-reactive
   Vδ1^+ γδ T-IELs ([300]Fig. 3D, [301]Fig. 4, [302]E and [303]F)
   ([304]43). In contrast, within natural clusters, CD patients exhibited
   highly distinct TRGV/TRDV pairings when compared to healthy controls,
   with increased diversity of pairings observed in CD mucosa ([305]Fig.
   4, [306]E and [307]F). As such, natural T-IELs in CD exhibited greater
   αβ and γδ clonal diversity when compared to healthy controls ([308]Fig.
   4, [309]A and [310]B). Collectively, these results suggest TCR-mediated
   selection of disease-associated γδ T-IELs in CD.

CITE-seq clustering of intestinal T cells

   To further profile cell surface marker repertoire of previously
   described populations, protein level data from a 204 antibody CITE-seq
   panel across ADT-stained patient samples (n = 26,066 IEL and LP T
   cells) was integrated using reciprocal principal component analysis to
   cluster intestinal T cells by cell surface marker phenotype ([311]fig.
   S6, [312]A and [313]B and [314]table S1). CITE-seq largely recaptured
   scRNAseq-based populations of intestinal T cells ([315]fig. S6C). We
   again observed that activated T[FH] and Treg populations (CITE-seq
   Tfh/Treg cluster) were highly enriched in ACD and clustered distinctly
   from other CD4^+ T cell populations ([316]fig. S6, [317]B to [318]D).

   Confirming our transcriptomic analysis, healthy mucosa was enriched
   with KLRC1^+CD73^+CD103^+CD8^+ T[RM] cells, while cytotoxic
   KLRG1-expressing T[EM] cells did not show significant changes across
   disease states ([319]fig. S6, [320]A to [321]G). Increased levels of
   activated CD8^+ T cells expressing CD38, CD39, CD71, IFNG, and PRF1
   (CITE-seq T[NE] cluster) were again observed in ACD mucosa ([322]fig.
   S6, [323]A to [324]G). While surface protein level clustering was
   unable to fully resolve scRNA-seq-defined T[NM]1 and T[NM]2 T-IEL
   clusters, we observed an enrichment of “natural” CD38^+CD39^+NKR^+CD8^+
   T-IELs (CITE-seq T[NM] cluster) in GFD mucosa, providing further
   evidence for the persistence of natural, NKR-expressing CD8^+ T cells
   in CD mucosa ([325]fig. S6, [326]B to [327]H).

Gluten challenge induces rapid reprogramming of natural CD8^+ and γδ T-IELs

   Following oral gluten challenge in CD patients on a GFD, CD8^+ and γδ T
   cells expressing the activation marker CD38 and gut-homing integrins
   CD103 and β7 integrin appear in peripheral blood ([328]14). Their
   functional significance, however, remains unclear ([329]47). We
   analyzed 3,924 CD38^+CD103^+ CD8^+ and γδ T cells from five patients at
   day 6 following three-day oral gluten challenge ([330]Fig. 5A,
   [331]fig. S7, [332]A and [333]B). Four of the five patients underwent
   intestinal biopsy prior to gluten challenge, enabling us to track the
   phenotype of shared T cell clones in the intestine and peripheral blood
   ([334]Fig. 5A, [335]fig. S7C, [336]table S1). The relative frequency of
   gluten-induced CD38^+CD103^+CD8^+ (0.16–1.3%) and CD38^+CD103^+γδ
   (0.08–3.73%) T cells in the peripheral blood varied across patients,
   consistent with prior reports ([337]fig. S7B) ([338]14).

Fig 5. Gluten challenge induces circulating T[NE]-phenotype cells from
intestinal T[NM]1 clones.

   Fig 5.
   [339]Open in a new tab

   (A) Schematic representation of experimental design and sorting
   strategy of peripheral blood, gluten-induced (CD38^+CD103^+) CD8^+ and
   γδ T cells. (B-C) Projection of transcriptomic profile of
   gluten-induced T cells onto UMAP plot of intestinal T-IELs colored by
   (B) gluten challenge sample and (C) predicted cluster using unimodal
   UMAP projection. (D) Quantification of predicted cluster of
   gluten-induced circulating T cells from (n = 5) gluten challenge
   patients. (E) RNA expression profile of selected gluten-induced T
   cell-enriched genes identified through differential expression within
   gluten-induced circulating and intestinal CD8^+ and γδ T cells. (F)
   UMAP plot colored by shared TCR clones of pre-challenge intestinal and
   post-challenge peripheral blood T cells, split by predicted phenotype
   of gluten-induced T cells. (G) Volcano plot of differentially expressed
   genes (log2FC expression) between shared TCR clones of pre-challenge
   intestinal (left) and post challenge peripheral blood (right) T[NE]
   T-IELs from (n = 4) CD patients. (H) Chord diagram displaying TRGV-TRDV
   pairing usage within gluten-induced γδ T cells. Two-sided Wilcoxon rank
   sum test (D), **P < 0.01.

   To investigate the phenotype of these gluten-induced peripheral blood T
   cells, we projected their transcriptomic profile atop our previously
   described intestinal CD8^+ and γδ T-IEL populations using unimodal UMAP
   projection ([340]Fig. 5B). Gluten-induced CD38^+CD103^+ CD8^+ and γδ T
   cells exhibited a highly cytotoxic phenotype consistent with T[NE]
   intestinal T cells, highly enriched in untreated CD ([341]Fig. 5,
   [342]B to [343]E, [344]fig. S7D). These cells similarly exhibited a
   pathogenic, NK-reprogrammed phenotype, highly upregulating IFNG, PRF1,
   CD52, KLRC2, HCST, KLRK1, in addition to MKI67 (KI67) and cell surface
   expression of CD94 ([345]Fig. 5E, [346]fig. S7, [347]E and [348]F).
   Confirming a prior report, gluten-induced CD8^+ and γδ T cells also
   highly expressed activation markers observed among T[FH] cells
   including HLA-DR and CD71 ([349]Fig. 1F, [350]fig. S7F) ([351]47).
   Thus, gluten induces circulating CD8^+ and γδ T cells in peripheral
   blood exhibiting a NK-reprogrammed cytotoxic T[NE] phenotype, similar
   to T-IELs enriched in PCD and ACD mucosa.

   For patients who underwent both biopsy and gluten challenge, we
   observed significant sharing of αβ and γδ TCRs across the intestine and
   peripheral blood ([352]15). Interestingly, we found that these common
   clones had distinct phenotypes in intestine vs. peripheral blood,
   suggestive of gluten-induced reprogramming ([353]Fig. 5, [354]F and
   [355]G, [356]table S7). Intestinal αβ and γδ T-IEL clones from the
   T[NM]1 cluster, the baseline γδ T cell population in healthy controls,
   were later seen in the blood exhibiting the cytotoxic T[NE] phenotype
   ([357]Fig. 3, [358]B to [359]D, [360]Fig. 5F, [361]table S8). While we
   also observed a limited number of intestinal T[EM] and T[RM] clones
   circulating in peripheral blood after gluten challenge, in contrast to
   natural T-IELs, these clones largely retained their intestinal
   phenotype in the blood ([362]Fig. 5F, [363]table S8). Therefore, as
   suggested by TCR sharing between T[NM] and T[NE] T-IELs within
   untreated CD mucosa ([364]Fig. 4, [365]C and [366]D), these results
   from gluten challenge affirm the common ancestry of natural-memory and
   natural-effector αβ and γδ populations. These data further show that
   intestinal αβ and γδ T cell clones present in the T[NM]1 cluster
   undergo rapid reprogramming upon gluten ingestion, assuming the
   cytotoxic and pro-inflammatory T[NE] phenotype in the peripheral blood.

   It was previously reported that gluten-induced circulating CD8^+ and γδ
   T cell populations show TCR similarity and clones are recalled upon
   repeat gluten challenge, suggesting antigen-driven selection ([367]14).
   To assess the possibility that these gluten-induced peripheral blood T
   cells were selected by specific antigens, we used the GLIPH2 algorithm
   to cluster TCRs to identify motifs critical for antigen-binding
   ([368]37). We observed convergence within the CDR3β region of
   gluten-induced CD8^+ T cells and identified numerous motifs shared
   across patients ([369]table S9). Further, like natural celiac γδ
   T-IELs, gluten-induced circulating γδ T cells did not show the biased
   TRGV usage seen in natural T-IELs from healthy mucosa ([370]Fig. 4,
   [371]E and [372]F, [373]Fig. 5H, and [374]fig. S5G). Collectively these
   results suggest that natural αβ and γδ T-IELs in CD are selected by
   antigens distinct from natural T-IELs in healthy mucosa. Although
   gluten recognition by T-IELs has not been demonstrated, these
   circulating αβ and γδ T[NE] T-IELs are induced by gluten with the same
   kinetics as gluten-specific CD4^+ T cells ([375]8, [376]14).

DISCUSSION

   Here, we leveraged opportunities afforded in CD to comprehensively
   study T cell phenotype and TCR repertoire across different disease
   states in CD. As CD shares traits with other autoimmune diseases, these
   findings should prove informative beyond CD.

   Intestinal gluten-reactive CD4^+ T cells were previously ascribed a
   T[FH]-phenotype through mass cytometry and bulk RNAseq ([377]21).
   Through our analysis, we identified a similar population of CD4^+ T
   cells, which were enriched in untreated CD mucosa, but largely absent
   in controls. Consistent with gluten-reactivity, we observed motifs seen
   in gluten-reactive TCRs and biased variable gene usage within this
   population ([378]33–[379]36, [380]38). The T[FH] phenotype of these
   cells are consistent with the idea that gluten-specific T cells
   interact with antigen-specific B cells (anti-gliadin, anti-TG2) serving
   as antigen presenting cells, reciprocally amplifying humoral and
   cell-mediated immune responses in CD ([381]48). These T[FH] cells
   phenotypically resembled those reported in rheumatoid arthritis and
   Crohn’s disease, and thus may have a similar role in these diseases
   ([382]30, [383]49).

   Tregs have been shown to be expanded within CD mucosa ([384]22,
   [385]23), but their functional significance is unclear ([386]27). Our
   data show CD39^+CD25^+FOXP3^+ Tregs in CD shared activation markers and
   TF networks with T[FH]-phenotype cells, but exhibited a suppressive
   phenotype (IL10, CD25, CTLA4, PD-L1). We further observed a positive
   correlation between the relative frequencies of T[FH] cells and Tregs
   in ACD and PCD mucosa, suggesting they are induced together. It was
   reported that gluten challenge induces circulating gluten-specific
   CD39^+FOXP3^+ Tregs ([387]26), although a more recent study
   demonstrated that FOXP3-expressing gluten-specific T cells do not
   co-express CD25 or phenotypically resemble Tregs ([388]21). Our TCR
   data show established gluten-reactive motifs present in T[FH] but not
   Treg cells, and negligible repertoire overlap between the two
   populations. Therefore, our data are consistent with the idea that
   inflammatory (T[FH]) and regulatory CD4^+ T cells are induced together
   in CD through gluten exposure, but are selected by different antigens.

   The importance of activating NKR-expressing T-IELs in celiac
   pathogenesis has now been well established ([389]4, [390]8–[391]12).
   The current prevailing model proposes that these T-IELs, which are
   responsible for tissue damage, are activated downstream of CD4^+ T
   cell-driven inflammation ([392]5, [393]8). Our data provide several
   novel insights regarding how NKR-expressing T-IEL responses are
   modulated in CD.

   An increase in T-IELs is a pathologic hallmark of CD, and known to
   persist even with dietary gluten exclusion. We describe two populations
   of NKR-expressing, natural memory T-IEL, T[NM]1 and T[NM]2, that are
   modulated in CD. While celiac-associated T-IEL populations were
   previously reported to consist of Vδ1^+ γδ T cells ([394]43), we found
   both CD8^+ αβ and non-Vδ1^+ γδ T-IELs present across all CD-modulated
   natural populations. The T[NM]1 population, the most abundant
   population of γδ T-IELs in healthy mucosa, expressed TGFB1 and AREG
   suggestive of immuno-regulatory and tissue repair-promoting roles. The
   T[NM]2 phenotype exhibited an NK-reprogrammed phenotype and highly
   upregulated HCST and KLRK1, suggestive of cytotoxic potential. Despite
   their cytotoxic potential, T[NM]2 T-IELs were relatively enriched in
   non-inflamed CD mucosa, suggesting they were not actively mediating
   intestinal damage.

   Our analysis instead implicates cells within the natural-effector
   (T[NE]) cluster as the cells responsible for tissue damage in CD. This
   cluster contained activated T-IELs with a highly cytotoxic and
   pro-inflammatory phenotype, expressing PRF1, GZMB, and IFNG, and was
   significantly enriched in damaged ACD mucosa. While this population was
   predominantly comprised of CD8^+ T-IELs, we also observed the presence
   of phenotypically indistinguishable γδ T[NE] T-IELs, suggesting
   redundant roles of natural CD8^+ and γδ T-IELs in CD. These
   natural-effector cells extensively shared TCR repertoire with
   CD-enriched natural-memory T-IELs, consistent with the prevailing model
   whereby exposure to IL-15 and CD4^+ T cell-mediated inflammation
   induces activation of NKR-expressing, natural-memory T-IELs to assume
   the fully cytotoxic natural-effector phenotype ([395]2–[396]5, [397]8,
   [398]44).

   Notably, we also observed an enrichment of natural-effector T-IELs in
   noninflamed PCD mucosa when compared to both control and GFD patients.
   When compared to ACD mucosa, T[NE] T-IELs were present at lower levels
   in PCD mucosa and expressed lower levels of pro-inflammatory mediators.
   While it is possible that a threshold level of T[NE] T-IELs is
   necessary for tissue damage and not reached in PCD, we speculate that
   mucosal expression of IL-15 and markers of epithelial stress, which
   distinguish PCD from ACD, may be required to unleash the inflammatory
   potential of T[NE] T-IELs ([399]3–[400]5, [401]11, [402]44).

   Our TCR data show that in CD, natural-memory and natural-effector cells
   show evidence of TCR-mediated selection and shared TCR repertoire,
   establishing that these populations share ancestry. One caveat to our
   study is that we were unable to obtain biopsies from the same patients
   in different disease states, limiting our ability to infer
   directionality of differentiation of T cell populations that share TCR
   repertoire. However, we were able to obtain gluten-challenged
   peripheral blood from four patients that previously underwent
   intestinal biopsy. These data establish directionality and timing of
   differentiation, showing that intestinal natural-memory clones,
   predominantly from the T[NM]1 phenotype, are reprogrammed to assume the
   natural-effector phenotype rapidly, within 6 days, upon gluten
   ingestion. This reprogramming occurs in parallel with gluten-specific
   CD4^+ T cell activation, before significant intestinal inflammation is
   seen histologically, suggesting it occurs deliberately and not
   downstream of CD4^+ T cell-induced inflammation ([403]14). Further
   evidence of the gluten-induced transition from T[NM]1 to
   natural-effector phenotype is seen in the diminished frequency of
   T[NM]1 cells in ACD compared to GFD mucosa. Beyond showing TCR
   repertoire sharing, our data do not show directionality of
   differentiation between T[NM]2 and T[NM]1/T[NE] populations. We
   speculate that once reprogrammed T[NE] cells are recruited to the
   intestine through gluten exposure, clones are retained as T[NM]2
   T-IELs, which can be re-induced to mediate cytotoxicity upon gluten
   re-encounter.

   Although T-IELs in CD are not thought to recognize gluten themselves
   ([404]8), their activation may still be dependent on gluten-specific
   CD4^+ T cells. Multiple aspects of effector CD8^+ T cell responses are
   known to depend upon CD4^+ T cell help ([405]50–[406]52). We propose
   that gluten-specific CD4^+ T cells similarly license the activation and
   recruitment of cytotoxic natural-effector T-IELs in CD. This licensing
   may occur via IL-21, implicated in CD8^+ T cell migration and effector
   function ([407]53), which was highly produced by gluten-reactive CD4^+
   T cells ([408]8, [409]21). Interestingly, IL-21R was highly expressed
   in both T[NM]1 T-IELs and Tregs in CD. If gluten-specific CD4^+ T cells
   do indeed interact with T[NM]1 CD8^+ T cells, this interaction would
   likely occur in gut-associated lymphoid tissue, as these two cell types
   reside in different locations within the intestine (epithelium vs.
   lamina propria). As our study was entirely performed on dissociated
   tissues, spatial context is lacking and we relied on established
   markers to infer the presence of cells in the intestinal epithelium
   ([410]20).

   The question of how gluten induces parallel inflammatory responses by
   gluten-specific CD4^+ T cells and T-IELs that have not been
   demonstrated to recognize gluten is an important one, and likely
   applicable in other HLA class II-associated autoimmune diseases
   ([411]2). A similar phenomenon of parallel activation of circulating
   CD4^+, CD8^+ and γδ T cells was reported in a murine model of CD4^+ T
   cell-initiated autoimmunity ([412]54). In this model, circulating CD8^+
   T cells functioned as regulatory CD8^+ T cells, suppressing
   proliferation of autoreactive CD4^+ T cells ([413]54). In contrast, our
   data strongly suggest that gluten-induced circulating natural CD8^+ and
   γδ T cells in CD have a pro-inflammatory, rather than regulatory, role.
   We note, however, that the pro-inflammatory role we ascribe to natural
   CD8^+ and γδ T cells in CD and the regulatory role seen in the mouse
   model are not mutually exclusive. One possibility that reconciles these
   two models is that natural CD8^+ and γδ T-IEL populations are
   self-reactive. Indeed, T-IEL populations in humans and mice have been
   shown to display self-reactive properties ([414]20, [415]55, [416]56).
   Their self-reactivity would enable both rapid mobilization to kill
   tissue in response to infection as well as feedback inhibition of
   pro-inflammatory CD4^+ T cells.

MATERIALS AND METHODS

Study Design

   The goal of our study was to leverage unique traits of celiac disease
   to study the recruitment, activation, and phenotypic alterations of T
   cells across distinct states of disease, at the single-cell level. To
   this end, we generated a single-cell multiomics dataset from T cells
   isolated from duodenal biopsies of eleven ACD, nineteen GFD, seven PCD,
   and seventeen control patients in addition to gluten-induced peripheral
   blood T cells from five CD patients following oral gluten challenge
   (see [417]table S1).

Human Specimens

   Duodenal biopsies from patients undergoing endoscopy at Columbia
   University Irving Medical Center were collected with informed consent,
   as approved by the Columbia University IRB (AAAB2472). All celiac
   disease patients had a biopsy-confirmed diagnosis. Patients were
   categorized as potential celiac disease based on positive TTG serology
   and HLA-DQ2/8 in the absence of histologic findings of celiac disease
   ([418]7). Gluten-free diet patients were treated for at least 6 months
   and had confirmed mucosal healing based on pathology. Peripheral blood
   from volunteers who underwent oral gluten challenge was collected with
   informed consent, as approved by the Columbia University IRB
   (AAAR1608). For gluten challenge, volunteers with biopsy-proven CD,
   adherent to a strict GFD for at least 6 months, consumed four slices of
   white bread (~15g gluten) per day for 3 consecutive days (days 1–3),
   and donated blood on day 6, as described ([419]14).

Biopsy Sample Processing

   Fresh tissue was dissociated into single-cell suspensions using
   sequential enzymatic digestions to isolate both IELs (1x HBSS, 1mM
   EDTA, 1mM DTT, and 5% FBS: 30 min shaking) and LPLs (1x HBSS, 5% FBS,
   and 0.025mg/mL Liberase: 90 min shaking) or with Multi Tissue
   Dissociation Kit 1 (Miltenyi Biotec) ([420]Table S1). Single-cell
   suspensions were cryopreserved (40% FBS 10% DMSO) and stored in liquid
   nitrogen.

Peripheral Blood Mononuclear Cell Isolation

   Whole blood was collected in BD Vacutainer Collection tubes. Peripheral
   blood mononuclear cells were isolated by density gradient
   centrifugation (Ficoll-Paque, Fisher Scientific). Single-cell
   suspensions were cryopreserved (40% FBS 10% DMSO) and stored in liquid
   nitrogen.

CITE-Seq Antibody-Oligo Hyperconjugation

   Conjugation of purified antibodies for cell hashing was performed as
   described ([421]https://cite-seq.com/protocol/). Antibody clones and
   barcodes are available in [422]table S2.

CITE-seq

   Cryopreserved single-cell suspensions were thawed at 37°C, aliquots of
   IEL/LPL fraction pooled (when applicable), and rested for 8 hours in
   RPMI/10% FBS. Single-cell suspensions were incubated with Fc-receptor
   block for 10 minutes at 4°C (TruStain FcX, BioLegend) followed by
   hashing antibody (0.5μg per 2 million cells) for 30 minutes at 4°C.
   Following multiple washes, samples were stained with anti-CD3 (clone
   OKT3, BioLegend) FITC, anti-CD4 (clone SK3, BD-Horizon) BUV737,
   anti-CD8α (clone RPA-T8, BioLegend) BV785, anti-CD45 (clone HI30,
   BioLegend) BV711, and anti-TCRγδ (clone 5A6.E9, Thermo Fisher) PE.
   Gluten challenge samples were stained with anti-CD38 (clone HIT2,
   BioLegend) APC, anti-CD103 (clone Ber-ACT8, BioLegend) BV711, and
   anti-Integrin β7 (clone FIB504, BioLegend) PE/Dazzle^™ 594. Samples
   were incubated with pre-titrated TotalSeq-C CITE-seq antibody cocktail
   (BioLegend, PN900000115) for 30 minutes at 4°C. LIVE/DEAD Fixable Red
   (ThermoFisher) and DAPI were used to gate viable cells. T cells were
   sorted into RPMI by FACS AriaII or Propel Labs Bigfoot Cell Sorter by
   gating viable CD3^+ intestinal cells and CD3^+CD4^−CD38^+CD103^+ cells
   from gluten challenge samples. FACS analysis performed using FlowJo
   (v.10.7.1) software.

Droplet-based CITE-seq

   FACS-sorted T cells were loaded onto the Chromium 10x Genomics V1 or
   v1.1 platform (10X Genomics) per manufacturer’s instructions. Single
   cell 5’ transcript, TCRαβ, and CITE-seq libraries were prepared using
   Chromium Single Cell 5’ V(D)J with Feature Barcoding Reagent Kits
   ([423]CG000186) or Chromium Next GEM Single Cell V(D)J Reagent Kits
   v1.1 with Feature Barcode Technology for Cell Surface Protein
   ([424]CG000208). Cell hashtag library generation was performed as
   described ([425]https://cite-seq.com/protocol/). TCRγδ libraries were
   prepared as described ([426]17).

   cDNA libraries were quantified using High Sensitivity D5000 ScreenTape
   (Agilent) and sequenced using Illumina MiSeq, NextSeq 550, and NovaSeq
   6000 platforms.

Data Processing

   To process gene expression and CITE-seq sequencing reads, Cell Ranger
   v5.0.0 software (10x Genomics) was used. Raw base call files were
   demultiplexed using cellranger mkfastq and aligned to GRCh38 reference
   genome and CITE-seq barcode list for 5’ mRNA and CITE-seq respectively
   with cellranger count. TCRαβ, and TCRγδ sequencing reads were aligned
   to GRCh38 reference genome with cellranger vdj version 3.0.2 (10x
   Genomics).

   Raw cell hashtag reads were processed using CITE-seq-Count to generate
   cell hashing count matrix
   ([427]https://github.com/Hoohm/CITE-seq-Count). Hashes were normalized
   to total counts across cells (total UMI counts set at 100). Cell
   barcode identity was assigned to most highly expressed hashtag that
   encompassed a minimum of 2/3 the recovered UMIs. Shared CDR3 nucleotide
   sequence was used to rescue unidentified cells. Remaining unidentified
   cells were removed from subsequent analysis.

Multimodal Single-cell Analysis

   Preprocessing and clustering analysis steps were performed using the
   Seurat R package. The percentage of genes mapping to mitochondria,
   ribosomes, and total UMIs were used to remove “low-quality” cells (>15%
   mitochondrial genes OR < 350 genes detected). Global-scaling
   normalization and identification of 2000 most variable genes were
   performed prior to integration across runs using canonical correlation
   analysis available in the Seurat R package. A second round of filtering
   was performed to remove additional “low-quality”
   (mitochondrial-gene-expressing, ribosomal-gene-expressing,
   HSP-gene-expressing, cycling cells) clusters of cells. CITE-seq count
   matrix was normalized using centered log-ratio transformation following
   adjustment to UMI counts to account for isotype controls in an
   isotype-specific manner. CITE-seq clustered data was integrated using
   reciprocal PCA analysis as available in the Seurat R package. Cells
   were visualized as UMAP embeddings and clustered with resolutions of
   0.2 for all intestinal T cells, 0.55 for CD8^+ and γδ T-IELs, and 0.35
   for CITE-seq-clustered data. Subclustering was performed on UMAP
   derived all intestinal T cells on a shared cluster T[FH]/Treg/T-Naïve
   cluster to isolate a cluster of naïve T cells. Subclustering was
   performed on T-IEL UMAP T[RM] cluster to isolate additional
   subpopulations of T[RM] T-IELs.

   Gluten-induced peripheral blood T-IELs were additionally filtered using
   CITE-seq ADTs (GC 1–3) and scRNA-seq (GC 1–5) to remove cells that did
   not express CD38, CD103, and β7-Integrin. Unimodal UMAP projection was
   performed to project transcriptomic profile of gluten-induced T-IELs
   onto UMAP structure of intestinal CD8^+ and γδ T-IELs as available in
   the Seurat R package. For TCR sharing, shared T cell clones were
   assigned phenotype (T[NE]/T[NM]1/T[NM]2/T[EM]/T[RM]/MAIT/Innate) based
   on most abundant predicted cell phenotypes per clone.

Identification of Cell Types

   TCRαβ and TCRγδ cells were annotated based on scTCR-seq output.
   Remaining T cells without a detected TCR were annotated as αβ or γδ T
   cells through CITE-seq expression of TCRαβ and TCRγδ in addition to
   expression of TCRγδ genes (i.e. TRDC, TRGC1, TRDV1). αβ T cells were
   then annotated as CD4^+ or CD8^+ based on combination of CITE-seq CD4
   and CD8 antibodies followed by annotation with the SingleR package with
   Monaco Immune Dataset for unbiased reference-based scRNA-seq cell type
   annotation. T-IELs were annotated through transcriptomic and CITE-seq
   expression of CD103.

scRNA-seq Differential Expression Analysis

   Differential expression between populations was performed using the R
   package MAST via the Seurat wrapper function FindMarkers. Minimum
   log-fold threshold was set to default parameters (0.25) with exception
   of generation of volcano plots where minimum log-fold threshold was set
   to 0. Differentially expressed genes available in [428]Supplementary
   Materials.

TCR Analysis

   TCRαβ and TCRγδ sequences annotated as “productive” were kept for
   further analysis. TCR counts were based on CDR3β, CDR3δ, CDR3γ amino
   acid sequences in a patient-dependent manner. Expanded TCRs were split
   into quartiles based on clonal frequency. Variable gene pairings
   visualized using the chordDiagram function available in the circlize R
   package. Sharing of CDR3β and CDR3δ sequences across clusters
   visualized using the UpSetR R package.

   CDR3β, TRBV gene, TRBJ gene, and CDR3α sequence were input into the
   GLIPH V2 algorithm ([429]http://50.255.35.37:8080/) using cell type
   matched reference. Fisher-exact test performed to assess confidence of
   clusters.

Gene Ontology and Pathway Enrichment Analysis

   Individual cells were scored for pathway activities using
   area-under-the-curve (AUC) analysis with the AUCell R package. Gene
   expression rankings for each cell was created using the
   AUCell_buildRankings function with default parameters. The canonical
   pathway database was downloaded from the Broad Institute, with gene
   sets used to score each cell individually for each gene set using
   calculated AUC values with the function AUCell_calcAUC.

Gene Regulatory Network Analysis

   The R package SCENIC was used to examine transcription factor module
   activity across intestinal T cells. Log-transformed count matrix was
   filtered to exclude genes detected in fewer than 50 cells, minimum
   counts per gene less than 130 (~1 average UMI count across all cells),
   and genes not detected in RcisTarget databases. Positive and negative
   target associations between transcription factors and genes were
   identified using the runCorrelation function on the filtered count
   matrix. Gene regulatory network analysis using stochastic Gradient
   Boosting Machine regression with early-stopping regularization was
   performed using GRNBoost2, as available in the arboreto python program.
   Detected coexpression modules for transcription factors were scored
   using AUCell package as previously described in “[430]Gene Ontology and
   Pathway Enrichment Analysis”.

Statistical Analysis

   Statistical tests performed in [431]Fig. 1E, [432]1H, [433]2E, [434]3D,
   [435]3G, [436]5D, [437]figs. S2A, [438]S3B, [439]S4B, and [440]S6D,
   were performed in R using the ggpubr package. Differential expression
   was performed in R using Seurat and MAST packages. Statistical tests
   performed in [441]tables S5 and [442]S9 used GLIPH2 program. Unpaired
   Mann-Whitney tests were performed for non-parametric data (min. 10
   CD4^+, CD8^+, γδ T cells/patient). Pearson rank correlation coefficient
   in [443]Figs. 1H and [444]3G. [445]Table S4 contains data input into
   statistical calculations. Fisher’s exact test was used to determine
   enrichment of GLIPH2-identified TCR motifs.

Supplementary Material

   Table S3
   [446]NIHMS1924047-supplement-Table_S3.xlsx^ (10.1KB, xlsx)
   Table S2
   [447]NIHMS1924047-supplement-Table_S2.xlsx^ (21.6KB, xlsx)
   Table S4
   [448]NIHMS1924047-supplement-Table_S4.xlsx^ (81.1KB, xlsx)
   Table S5
   [449]NIHMS1924047-supplement-Table_S5.xlsx^ (17.7KB, xlsx)
   main supplementary
   [450]NIHMS1924047-supplement-main_supplementary.pdf^ (1.7MB, pdf)
   Table S8
   [451]NIHMS1924047-supplement-Table_S8.xlsx^ (9.8KB, xlsx)
   Table S9
   [452]NIHMS1924047-supplement-Table_S9.xlsx^ (47.4KB, xlsx)
   Table S6
   [453]NIHMS1924047-supplement-Table_S6.xlsx^ (1.5MB, xlsx)
   Table S7
   [454]NIHMS1924047-supplement-Table_S7.xlsx^ (611.8KB, xlsx)
   MDAR checklist
   [455]NIHMS1924047-supplement-MDAR_checklist.docx^ (68.2KB, docx)
   Table S1
   [456]NIHMS1924047-supplement-Table_S1.xlsx^ (14.1KB, xlsx)

Acknowledgements: