Research was conducted under a protocol approved by the University of Virginia (UVA) institutional review board (study #18519), and informed consent was obtained prior to conducting any study-related procedures. The first 150 patients seen at UVA rheumatology meeting classification criteria for RA as defined by the American College of Rheumatology (ACR) or European League Against Rheumatism (EULAR) (2, 20) were selected. Exclusion criteria included prisoner status, non-english speaking, anemia with hemoglobin <10g/dL, and pregnancy. The 150 participant number was targeted to maintain adequate sample size post processing based on prior power calculations.

Medical records were reviewed using UVA’s EPIC EMR system. Outside records were requested in CareEverywhere through the UVA EPIC system. Demographics represent the most recent data available at the time of data mining.

Date-dependent clinical parameters were collected on or as close to the date of research sample collection as possible. In cases where laboratory results or other data for the same test existed in approximately the same timeframe both before and after research sample collection, results obtained on the date preceding sample collection were chosen. Raw values and positivity status were collected for each parameter, with positivity evaluation based on laboratory clinical standards as noted for each patient.

Clinical data acquired from multiple physicians, laboratories, and hospitals were harmonized into a uniform format. Either Disease Activity Score in 28 joints (DAS28) or Clinical Disease Activity Index (CDAI; whichever was available) was recorded for each patient, and values were converted to standard scores of remission, low, moderate, and high disease activity (</= 2.8, 2.9-10, 10.1-22, >22 for CDAI; <2.6, 2.6-3.0, 3.1-5.0, >5.0 for DAS28). Double positivity refers to patients with available data for both RF and historical anti-CCP and both were positive.

Complete blood counts, erythrocyte sedimentation rate (ESR), and C-reactive protein (CRP) values are routine tests for ongoing RA management, so most were on or near the date of research sample collection. RF and anti-CCP testing are used diagnostically and thus collected less frequently, so were often found at dates farther from sample collection. As RF and anti-CCP values tend to be less labile, this increased gap between clinical test and research sample for RF and anti-CCP was deemed allowable. Of note, these values mined from patient charts are denoted as historical anti-CCP. Current anti-CCP values were derived from sera collected at the same timepoint as that used for ddPCR analysis and determined with the QUANTA Lite^® CCP3.1 IgG/IgA ELISA kit (positivity cutoff ≥20 U/mL).

Fresh peripheral blood mononuclear cells (PBMCs) from 150 RA patients and 9 HCs (Research Blood Components) were isolated using Ficoll-Paque gradient separation per manufacturer’s instructions (Cytiva). Then fresh CD8+ T-cells were isolated using RosetteSep reagent according to manufacturer’s instructions (StemCell Technologies). Genomic DNA (gDNA) was extracted from isolated CD8+ cells using the Anaprep Cultured Cell DNA Extraction Kit (Biochain) and stored at -20°C.

gDNA from isolated CD8+ T-cells was not sufficient to allow replicates of multiple assays; therefore, gDNA was pre-amplified as a single amplicon suitable as input for all downstream assays. Samples were tested with a quantification assay to determine starting copy number, and 40,000 genomic copies were used as input for pre-amplification. After amplification, the D661 wild-type (WT) assay was performed to quantify pre-amplification outputs. 40,000 amplified copies were then used as input for the D661Y and Y640F mutation-specific assays and the S614-G618 dropout assay designed to detect any mutation in that region. ddPCR was performed using the QX200 and digital droplet generator, ddPCR Supermix for probes (No dUTP), and FAM/HEX probes under manufacturer’s recommendations (BioRad) ( Supplementary Figure 1 ). The three samples with the highest frequency mutations were confirmed by next-generation sequencing as previously described (21).

At least two separate runs using 40,000 amplified copies as input were performed for each mutation-specific assay. A third confirmatory assay was performed using 20,000 amplified copies as input. Biorad’s Quantasoft software provided the variant allele frequencies (VAFs) based on a Poisson distribution, and thresholding was placed manually based on the overlay of all samples on each plate with a threshold cutoff of at least 3 positive droplets per well as per manufacturer recommendations. Sample VAFs were considered replicated and used for statistical analysis if they had at least two runs with VAFs in the same bin (Bins: <threshold, 0.008-0.05%, >0.05-5%, >5%). Assuming heterozygosity, the maximum bin cutoff is equivalent to 1 mutation in every 10 cells, and the minimum threshold is equivalent to 1 mutation in every 6,250 cells. Following all pre-amplification and ddPCR assay steps, 98 patients had sufficient sample volume for analysis.

Statistical testing was performed in R (ver 4.3.0). Chi-squared tests with p-values based on Monte Carlo simulation were used for categorical variables, comparing different VAF bins, WT to STAT3 mutant patients, and Y640F patients to D661Y. ANOVA compared clinical features between WT, D661Y, and Y640F groups based on mutation status. Welch’s ANOVA or Kruskal-Wallis was used depending on assumptions. T-tests compared quantitative variables between WT and STAT3 mutant samples. Either Welch’s T or Wilcoxon Rank Sum tests were used depending on assumptions. All tests and p-values are in Supplementary Tables 1 - 3 . The cutoff for p-value significance was <0.05. Sample sizes varied based on clinical data availability. Adjustment of p-values for multiple comparisons was done for each test type using the Holm method ( Supplementary Tables 1 - 3 ). Sensitivity analyses investigating clinical differences between patients with or without STAT3 mutation were performed using either a cutoff VAF≥0.008% or VAF>0.05%, with the testing from VAF>0.05% reported in all tables and figures except Supplementary Figure 2 which uses ≥0.008% as the cutoff.

We modeled mutation presence (VAF≥0.008%) as a function of RA status and age using logistic regression. The model incorporated data from our ddPCR analysis as well as from a published dataset of 99 HCs that utilized deep amplicon sequencing of the two main exons (20 and 21) of the STAT3 SH2 domain, with 2x300bp reads optimized for high sequencing depth (>25,000x). Their resulting VAF range (0.007%-1.2%) was comparable to the 0.008%-5.93% we observed using our targeted ddPCR assays (22).

Due to the putative role of cytotoxic T-cells in RA, symptomatic overlap, and the frequent co-occurrence of LGL leukemia with RA (6, 23), we hypothesized that CD8+ T-cells from RA patients may exhibit a higher prevalence of STAT3 mutations than healthy individuals. To address this question, we isolated CD8+ cells from the whole blood of RA patients, extracted, and then amplified DNA from these cells. The DNA was assayed by ddPCR to detect the following STAT3 mutations: Y640F, D661Y, or any mutations in the region of S614-G618 of STAT3 ( Figure 1 ). Each assay showed distinct regions of positive droplets that were used to determine VAF ( Figure 2 ). Example ddPCR outputs of zero, near threshold, and high VAF are shown in Supplementary Figure 3 . Y640F and D661Y ddPCR assays utilized highly sensitive allelic discrimination assays as previously reported (12). The S614-G618 assay was designed to detect any mutation in that region. As a dropout assay, less sensitivity is expected as separation of mutant from WT droplets is not as distinct ( Figure 2 ).

STAT3 ddPCR VAFs were determined for 98 RA patients and 9 HCs ( Table 1 ). Overall, we observed that 52/98 RA patients had either Y640F or D661Y mutations (VAF ≥0.008%), significantly more than 1/9 in HCs (p=0.029, Chi-Squared test). The resulting VAFs were divided into four different ranges, <threshold, 0.008-0.05%, >0.05%-5%, and >5%. Of note, no mutations were detected between 1% and 5% VAF. The upper bin was set at >5% to represent a minimum 1:10 cells with a mutation assuming heterozygosity.

The majority of HCs (8/9 Y640F assay, 9/9 D661Y assay) had VAFs below the 3-droplet threshold of the assay. Of the RA patients sampled, 3.1% (3/98) had VAFs above 5% for any of the three mutations analyzed, and 10/98 Y640F and 8/98 D661Y mutant samples had STAT3 mutations above 0.05% VAF. Y640F and D661Y mutations were frequently detected in the RA samples in the 0.008-0.05% VAF range (Y640F 31/98 and D661Y 17/98), but not at all in the HCs (Y640F 0/9 and D661Y 0/9). Several RA patients had detectable levels of both mutations (14/98). Interestingly, 7/8 of the highest VAF (>0.05%) D661Y patients also had some level of Y640F mutation. Of note, the two mutations detected with the S614-G618 dropout assay were confirmed via rhAmpSeq, and the variants were S614R and G618R, both of which have been previously reported in LGL leukemia (12, 24).

To validate this difference in mutation presence between RA samples and HCs, we utilized a publicly available dataset from Valori, et al. (22) This study utilized deep amplicon sequencing of the two main exons of the STAT3 SH2 domain to determine the prevalence of STAT3 mutations in healthy donor CD8+ T-cell DNA with nearly identical sensitivity to our ddPCR assay ( Figure 3A ). In their screen for mutations across the STAT3 SH2 domain, Valori et al. observed that 24/99 controls exhibited any mutation. However, only 7/99 (VAF 0.07-0.18%) had Y640F or D661Y mutations, with only 1/99 having a Y640F substitution. The prevalence of either Y640F or D661Y mutations in this larger control population did not differ significantly from our control dataset (p=0.17 and p=1, respectively).

Comparison of the Y640F and D661Y STAT3 distribution of this new cohort of 108 healthy donor samples (99 individuals plus 9 from our study), to the 98 RA patients in our study showed a statistically different VAF distribution of STAT3 mutations in RA patients compared to HCs (Y640F mutation p=0.0005 and D661Y mutation p=0.0005, Pearson’s Chi-Squared test), with the majority of the difference stemming from mutations in the 0.008%-0.05% VAF range ( Figure 3A ).

To account for the potential confounding factor of patient age differences between cohorts, logistic regression was used to explain the presence of any detectable mutation (≥0.008% VAF) using RA status and age. Even after taking age into account, there was a significantly higher probability of STAT3 mutation in RA patients compared to HCs across both mutation types (Either Y640F or D661Y p=1.74x10^-6, Y640F alone p=4.57x10^-5, D661Y alone p=0.00465) ( Figures 3B–D ).

In LGL leukemia, STAT3 mutations are correlated to specific disease phenotypes and response to treatment (11, 13, 15, 25). Therefore, we investigated the association of STAT3 mutations with clinical and demographic characteristics of RA patients to determine whether STAT3 mutation presence was related to disease manifestation. Many of the mutations detected in this study are in a very small proportion of CD8+ T-cells. Initially we looked at the clinical differences between WT patients and those with any detectable STAT3 mutation (VAF≥0.008%) and observed no strong differences between groups (p-values in Supplementary Tables 1 - 3 ). With the assumption that there may be a minimal level of mutation necessary to impact clinical features, we used a higher cutoff (VAF >0.05%, equivalent to 1:1000 heterozygous mutant cells) that captured approximately the upper quartile of all mutations. All samples with a VAF >0.05% for any of the three mutation assays were considered STAT3 mutant (n=16), and the rest were considered WT (n=82). Four main categories were examined: patient characteristics, measures of disease state, hematologic parameters, and treatment ( Table 2 ).

The overall cohort exhibited similar demographics to what had previously been reported for RA. Clinical features were similar between mutant and WT groups for most parameters with exceptions such as time since diagnosis graphed in Supplementary Figure 4 . While 59% (57/97) of patients exhibited erosive disease, the majority of patients with available data had well-controlled disease, as indicated by CDAI/DAS28 values in the remission and mild disease categories for 67% (62/92) of individuals. CRP serum levels differed across WT, Y640F, and D661Y, but the average values (1.67 mg/L, 0.29 mg/L, 1.12 mg/L respectively, p=0.04 Welch’s ANOVA) in each category were well below the positivity cutoff of 10 mg/L ( Supplementary Figure 4 ). Blood counts were typically in the normal range, and most patients were exposed to methotrexate (86%, 84/98), while <50% were exposed to TNF-α inhibitors.

Utilizing a threshold cutoff of any detectable mutation (VAF ≥0.008%) for comparison across WT and STAT3 mutant groups, we observed a correlation between the presence of STAT3 mutation and seropositivity for either RF or anti-CCP ( Supplementary Figure 2 ). The majority of patient samples in the >0.05% VAF bin were above the threshold for seropositivity ( Supplementary Figure 5 ). Using this refined 0.05% cutoff, we observed that seropositivity rates were significantly higher in patients with STAT3 mutations vs WT. Specifically, RF positivity was significantly more prevalent (p=0.047) in patients with mutant STAT3 (87.5%) than WT STAT3 (58.5%) ( Figure 4A ), with the median value of RF in patients with mutant STAT3 (Median: 152.5 U/mL, Q1: 64.6 U/mL, Q3: 234.3 U/mL) trending higher compared to WT STAT3 (Median: 59.6 U/mL, Q1: 20 U/mL, Q3: 176.8 U/mL, p=0.099) ( Figure 4B ). Historical anti-CCP positivity was also significantly more common (p=0.031) in patients with mutant STAT3 (92.3%) compared to WT STAT3 (57.5%) ( Figure 4C ). Anti-CCP levels trended higher (p=0.339) in the STAT3 mutant group (Median: 173.5 U/mL, Q1: 30.4 U/mL, Q3: 250 U/mL) compared to WT (Median: 38.5 U/mL Q1: 0.5 U/mL, Q3: 300 U/mL), but due to the lack of titration in clinical testing, the full potential range of serum levels was not represented in the data available ( Figure 4D ). Current anti-CCP levels measured from the same time point as sequencing were completed in-house to standardize and determine the full range of serum values in patients. Both the current anti-CCP positivity (STAT3 mutant 87.5% vs. WT 69.5%, p=0.216) and levels (STAT3 mutant Median: 305.7 U/mL, Q1: 47 U/mL, Q3: 649 U/mL vs. WT Median: 94.3 U/mL, Q1: 9.6 U/mL, Q3: 456.3 U/mL p=0.177) exhibited similar trends to the historical values but lacked significance ( Figures 4E, F ). However, double positivity for RF and anti-CCP was significantly more common (p=0.016) in mutant STAT3 (84.6%) vs. WT RA patients (48.8%) ( Figure 4G ).