Multi-panel mass cytometry reveals anti-PD1 therapy-mediated B and T cell compartment remodeling in tumor-draining lymph nodes

Anti-PD1 therapy has become an immunotherapeutic backbone for treating many cancer types. While many studies have aimed to characterize the immune response to anti-PD1 therapy in the tumor and in the peripheral blood, relatively less is known about the changes in the tumor draining lymph nodes (TDLNs). TDLNs are primary sites of tumor antigen exposure that are critical to both regulation and cross-priming of the antitumor immune response. We employed multi-panel mass cytometry to obtain a high-parameter proteomic (39 total unique markers) immune profile of the TDLN in a well-studied PD1-responsive immunocompetent mouse model. Based on combined hierarchal gating and unsupervised clustering analyses, we found that anti-PD1 therapy enhances remodeling of both B and T cell compartments toward memory phenotypes. Functionally, expression of checkpoint markers was increased in conjunction with production of IFN, TNF, or IL2 in key cell types, including B and T cell subtypes and rarer subsets such as Tregs and NKT cells. A deeper profiling of the immunologic changes that occur in the TDLN milieu during effective anti-PD1 therapy may lead to the discovery of novel biomarkers for monitoring response and provide key insights toward developing combination immunotherapeutic strategies.


Introduction
Tumor draining lymph nodes (TDLN) are one of the primary sites of tumor antigen exposure and antigen-specific immune activation, a barrier between local and systemic immune responses. TDLNs are, however, also major sites of immune regulation, contributing to the expansion and suppressive action of immunosuppressive populations such as regulatory T cells (Tregs) and myeloid derived suppressor cells (MDSCs) (1)(2)(3)(4). Previous studies have shown that the tumor antigen experienced TDLNs exhibit a disproportionately increased B cell compartment relative to the T cell compartment (5)(6)(7)(8). Furthermore, regulatory B cells (Bregs) accumulate locally in the TDLNs and promote tumor growth (9,10). These factors recognize TDLNs as one of the primary targets of immune modulation for improving antitumor immunotherapeutic strategies through a better understanding of their role in tumor progression (11,12).
Immunotherapy agents targeting the programmed cell death protein 1 (PD 1) pathway has led to unprecedented responses in multiple cancer types (13). By antagonizing PD1, an inhibitory checkpoint molecule expressed on the surface of exhausted T cells, the immune system is reinvigorated to attack cancer cells (14). In contrast to the described immunoregulatory effects of TDLNs, a recent study has also demonstrated that TDLNs are a critical component of the antitumor immune response to anti-PD1 therapy (15). While there have been significant efforts to characterize the changes in the immune tumor microenvironment and peripheral blood with anti-PD1 therapy, the impact of anti-PD1 therapy on the immunologic makeup of TDLNs remains largely unexplored To study TDLNs, we established an immunocompetent syngeneic mouse model by injecting MC38 colon carcinoma cells subcutaneously in the right hind limb of C57Bl/6 mice, permitting direct interrogation of the TDLN at the right inguinal lymph node. The immunologic phenotype of this model is also well known for its responsiveness to anti-PD1 therapy and the changes occurring in the tumor-infiltrating immune cells have previously been characterized (16).
To study in depth both T and B cell subtypes and their functional status within the TDLNs, we employed multi-panel mass cytometry, "cytometry by time-of-flight" (CyTOF) (17), and performed single-cell proteomic profiling of the TDLN using a total of 39 unique markers. To enable a more robust downstream analysis, we also utilized a multiplexed staining approach based on CD45barcoding. By both unsupervised clustering and hierarchal gating analysis of the debarcoded CyTOF datasets, we found that anti-PD1 therapy preferentially modulates the B cell compartment, with skewing toward mature and memory phenotypes favored in both B and T cells. Analysis of functional marker expression in these B and T cell subsets demonstrated that CD69, checkpoints, and specific inflammatory cytokine expressions are enhanced by anti-PD1 therapy in memory B and T cells, and that these changes occur in association with an increase in a subset of Tregs that produce TNF, and NKT cells that produce IL2.

Results
Tumor draining lymph nodes are increased in size in response to anti-PD1 therapy. To study the changes in tumor draining lymph nodes in response to anti-PD1 therapy, we studied a wellestablished mouse model of PD1-responsiveness, MC38 colon cancer cells, in syngeneic C57Bl/6J background. Mice were injected subcutaneously in the right hind limb with MC38 cells and were treated with isotype or anti-PD1 antibody on days 14 and 18 when tumors reached approximately 500mm 3 in volume; and tumor volumes and TDLNs were assessed on day 21 ( Figure 1A). As expected, tumor volumes were smaller in the anti-PD1 group (Supplementary Figure 1A). Based on weight ( Figure 1B) or cellularity ( Figure 1C), TDLNs were significantly larger than the right inguinal lymph nodes from non-tumor bearing normal mice, and PD1-treated TDLNs were significantly larger than isotype-treated TDLNs. There were no significant differences in the TDLN size whether the mice were untreated or treated with the isotype antibody (Supplementary Figure 1B, C).
Increase in TDLN size related to anti-PD1 therapy is due to disproportionate expansion of the B cell compartment. To analyze the immune cell constituents of TDLNs, three groups were crosscompared: anti-PD1-treated TDLNs, isotype-treated control TDLNs, and naïve lymph nodes from non-tumor bearing normal mice as an additional control comparator. Lymph nodes were dissociated into single cells and subjected to PMA/ionomycin stimulation to simultaneously determine the immune cell capacity for cytokine production along with their subtyping markers. Samples belonging to each of the three groups were barcoded by staining with a CD45 antibody conjugated to a unique metal tag. CD45-labeled cells from each group were then combined into three-plex batches. The batches were subsequently aliquoted for multiplexed staining with either T or B cell-oriented CyTOF panels as shown in Figure 1D and Supplementary Tables 1 (for T  cells) and Table 2 (for B cells). Hierarchical gating on biaxial plots to identify T and B cell compartments identified the following subsets: naïve and memory cytotoxic T cells, naïve and memory helper T cells, regulatory T cells, T1-and T2-type transitional B cells, mature phenotype B cells, memory B cells, and regulatory B cells (gating strategies shown in Supplementary Figures  2 and 3). These analyses revealed that anti-PD1 therapy significantly expanded both B and T cell populations, but with more pronounced effects on B cells (Figure 2A). Compared to non-tumor bearing naïve lymph nodes, isotype-treated TDLNs exhibited 2.2-fold and 11.7-fold expansions in the T and B cell compartments, respectively, whereas anti-PD1-treated TDLNs showed 4.7fold T cell and 28.0-fold B cell compartment expansions. Further gating was done to profile subsets of both B and T cell compartments within the TDLNs (Figures 2B and Supplementary Figure 4). No significant differences in the TDLNs were attributable to the isotype antibody itself (Supplementary Figure 1D). In parallel, we also performed unsupervised clustering analysis using the FlowSOM algorithm. Using the T cell-oriented and B cell-oriented CyTOF panels, we identified 20 and 25 metaclusters that were then annotated into 7 T cell subtypes ( Figure 3A) and 10 B cell subtypes ( Figure 4A), respectively. All cell type annotations are listed in Supplementary Table 3. In terms of numbers of cells of each cell type in each lymph node, most cell types were significantly increased by anti-PD1 therapy compared to isotype controls ( Figure 3B, 4B). This was especially true for memory B and T cells and regulatory B and T cells ( Figure 3C, 4C, Supplementary Figure 5). Thus, our data suggest that anti-PD1 therapy stimulates T cell expansion and even greater B cell expansion, and differentiation leading to significant increases in the presence of both memory and regulatory subtypes.
Anti-PD1 therapy augments differentiation of B and T cell types away from naïve and toward memory and mature subtypes, especially for B cells. To further evaluate how each of the B and T cell compartments are remodeled by anti-PD1 therapy, we utilized the same CyTOF datasets to then analyze the subtypes as percentages within each of the T or B cell compartments. As visualized by a dimensionality reduction algorithm UMAP ( Figures 5A, 6A), the anti-PD1 therapy group demonstrated the highest percentage of memory T cells, both CD4 and CD8, and the lowest percentage of naïve T cells in the T cell compartment ( Figure 5B). Similarly, the highest percentage of memory B cell types and the lowest percentage of immature/transitional B cell types were generally noted in the anti-PD1 therapy group ( Figure 6B). Among all of the B and T cell types analyzed ( Supplementary Figures 6 and 7), statistically significant differences between isotype and anti-PD1 therapy were particularly noted for memory "B1" (CD21 + CD23 + CD38 + B220 hi IgD hi IgM hi , 1.55-fold vs. 1.37-fold change relative to non-tumor bearing LN) and memory "B3" (CD21 + CD23 + CD38 + B220 hi IgD hi IgM lo , 1.41-fold vs. 0.9-fold change relative to non-tumor bearing LN) subtypes.
Next, we performed a repeat run of the experiment with a different CyTOF panel to further validate these observations. This CyTOF panel was modified to (1) combine T and B subtyping markers into a single, simplified panel, (2) replace CD3+CD4+CD25+CD127-with CD3+CD4+Foxp3 as the canonical combination of markers to identify Tregs, and (3) exclude cytokine assaying to profiling the lymph nodes without a stimulation step (Supplementary Table 4). Based on the metal intensities from 18 subtyping markers, we performed FlowSOM and UMAP algorithms again for this dataset to identify and annotate the 25 metaclusters into seven T cell clusters, four B cell clusters, and one NK cell cluster (Supplementary Figure 9 and Table 5). We also performed, in parallel, hierarchal gating of the cell subtypes as performed previously. Based on the results from both analytical approaches, the anti-PD1 effects on the T cell compartment were consistent, again showing significant increases in memory and regulatory phenotypes ( Supplementary Figures 10-13). Moreover, the dramatic increases in B cells within the TDLNs was confirmed by this repeat dataset, i.e. the B cell compartment was the predominant component enhanced within the TDLNs. While the effects of anti-PD1 therapy on the B cell compartment remodeling within the TDLNs were less robust in this repeat set, mild trends could be observed for the memory subtype and regulatory B cells. The tumor itself appears to affect the B cell compartment within TDLNs and this may be enhanced with anti-PD-1 therapy.
Lymph node remodeling is specific to successful response to anti-PD1 therapy in tumor-bearing mice but not in normal mice without tumors. Since the observed changes in the B and T cell compartments could just represent non-specific changes occurring with PD1 inhibition, we also assessed by both CyTOF and flow cytometry the impact of anti-PD1 therapy in the lymph nodes from normal, non-tumor bearing mice and also contralateral (non-tumor draining; NTDLNs) lymph nodes from MC38-bearing mice. Except for non-significant minimal trends, treatment with anti-PD1 therapy did not have any notable impact in the lymph nodes from normal mice (Supplementary Figure 13). However, in the contralateral site in MC38-bearing mice, there were increased number of B and T cell subtypes in the tumor-bearing non-tumor draining lymph nodes (Supplementary Figure 9-12). Furthermore, when we assessed the TDLNs from mice bearing B16-F10 tumors, which are classically known to be PD1-resistant (Supplementary Figure 14A), no significant differences were detected between the isotype and anti-PD1 therapy groups within the TDLNs across all B and T cell types analyzed (Supplementary Figure 14B). These findings together suggest that the TDLN is an active site of response to PD1 therapy and that the drastic remodeling within the lymph nodes is unique to successful PD1 therapy in tumor-exposed immune responses. In addition to the checkpoint and costimulatory markers, production of three key inflammatory cytokines, interferon gamma (IFN), tumor necrosis factor alpha (TNF), and interleukin 2 (IL2) were also evaluated. IFN levels were the highest in memory CD8 T cells among the T cell subtypes studied and were significantly increased with anti-PD1 therapy ( Figure 5H

Discussion
Anti-PD1 therapy has transformed the management of many human cancers, but only a subset of patients responds to anti-PD1 therapy, and responses are rarely curative. Understanding the effects of anti-PD1 therapy on the composition and functional state of tumor-associated immune subsets can identify mechanisms of response and resistance to anti-PD1 therapy and provide insights into rational combination strategies to improve upon response rates to monotherapy. While the effects of anti-PD1 therapy within the tumor immune microenvironment in preclinical and human subjects has been previously reported (13), to our knowledge this is the first study to broadly characterize the effects of anti-PD1 therapy on TDLN. Studying the TDLN is particularly valuable since it is one of the primary sites at which the tumor antigen exposure occurs and impacts the immune response. Furthermore, to interrogate the changes in the B cell compartment, TDLN offers an opportune site since within the tumors of many syngeneic immunocompetent mouse models, B cell presence is very rare (18), whereas in the peripheral blood, the findings would be even more removed from fully capturing direct responses to tumor. Our results regarding the immune cell type composition, checkpoint expression, and cytokine expression within TDLN suggest that there are major remodeling events within the TDLN mediated by anti-PD1 therapy that are critical for conferring effective antitumor immune responses. We demonstrate that B and T cell compartments and other cell subsets present at lower frequencies, such as NKT cells, undergo significant activation, expansion, and maturation in response to anti-PD1 therapy, to shape antitumor programming in TDLNs. Importantly, our study provides a methodological demonstration that simultaneous use of complementary mass cytometry antibody panels enables deep interrogation of the immune profile and that performing unsupervised clustering and gating analyses in parallel empowers a more reliable and informative investigation.
Our results have important implications regarding the tumor immune response to anti-PD1 therapy. First, assaying immune cell activation states in the TDLNs may determine therapeutic responses and guide therapy development, including rational immunotherapy combinations. For example, we specifically identify a marked expansion of Bregs (PD-L1+ B cells and CD5+ B cells) in the TDLN in response to anti-PD1 therapy (19). Bregs have been implicated in tumor progression because they may facilitate an immunosuppressive environment and attenuate antitumor immune responses. This finding provides initial evidence that Bregs may be important for re-establishing tolerance in the setting of anti-PD1 therapy, suggesting that therapies that inhibit Bregs may improve upon response to anti-PD1 therapy. Secondly, our results show that anti-PD1 therapy broadly reprograms the B cell component in TDLN, highlighting the importance of exploring treatment strategies that target B cells to further enhance responses to PD1-targeted therapies. This idea is consistent with the results of a previous study looking at non-small cell lung carcinoma demonstrating that B cell-mediated antigen presentation to CD4 T cells is a critical feature of effective antitumor immune response (20). We recognize, however, that targeting B cells successfully is not a simple goal given how B cells are involved in both pro-tumor and antitumor programs as extensively reviewed elsewhere (21). Thirdly, our results provide new evidence of cytokine production changes within the TDLN milieu in response to PD1 therapy that likely represent augmented antigen processing and presentation to activate both B and T cell subsets. For example, we find a significant increase in IL2-producing memory B cells within the TDLN. Given the high number of B cells within the TDLN, our results indicate that B cells may be an important source of IL-2 that support the anti-tumor T cell response. These IL-2 producing B cells also express higher levels of CD69, consistent with increased activation.
Interestingly, our high-parameter analysis also revealed a Treg subset producing increased TNF, and NKT cells producing increased IL2 in response to anti-PD1 therapy. The role of TNF as an antitumor molecule and a stimulator of adaptive immunity is well established (22)(23)(24)(25). Recent studies, however, have demonstrated that Tregs are maintained by TNF via TNFR2 (26) and themselves produce TNF in response to activation, providing an autocrine immune-suppressive effect (27). Prior studies have also shown that TNF augments MDSC accumulation (28). In line with these reports, emerging evidence suggests blocking TNF may enhance anti-PD1 therapy (29,30). Thus, increased TNF in response to PD1 inhibition observed in our data may be represent both antitumor signaling and an immune-regulatory negative feedback loop. Prior studies have shown that NKT cells are able to stimulate the antigen-presenting functions of B cells (31) and also secrete IL2 upon activation to aid with reinvigoration of exhausted CD8 T cells (32). Thus, our data suggests that anti-PD1 therapy is associated with activation of NKT cells that in turn further stimulate the antitumor immunological processes mediated by both B and T cells. Furthermore, these results provide new evidence that IL2 and TNF accomplish this in TDLNs during anti-PD1 therapy (33).
The study has the following conceptual limitations. First, we have not investigated further into the T cell receptor and B cell receptor immune repertoires and/or utilized a tumor antigen-specific model, e.g. ovalbumin-expressing MC38, which would elucidate the antigen specificity of the TDLN response. Second, the study only evaluates a single time point, limiting our understanding of how the TDLN remodeling occurs over time and whether a particular time point would be critical to assess. Third, since our study focuses on TDLNs, in which the vast majority of cells are lymphoid cells, it does not yield any insight into the interplay between lymphoid and myeloid immune compartments, which is also very important in shaping the antitumor immune response. Thus, it would be worthwhile to simultaneously assess the changes that occur within the tumor infiltrating immune cells to understand the changes in the myeloid populations and also delineate the unique features of the TDLN. Fourth, we have not demonstrated any causality of these findings with regards to the efficacy of anti-PD1 therapy. This study thus prompts further investigations into the definitive roles of specific immune cell subtypes within the lymph nodes in the therapeutic outcome. Fifth, we demonstrate results from just one PD1-responsive (MC38) and one PD1-resistant (B16) model, restricting the generalized applicability of our conclusions. In light of these limitations, our findings are exploratory, warranting the need for more confirmatory and in-depth studies on changes that occur within the lymph nodes during anti-PD1 therapy by incorporating additional animal models and treatment regimens.
In summary, we broadly characterized the composition and functional state of immune cell subsets within TDLN in response to anti-PD1 therapy in a PD-1 therapy responsive tumor model. These data identified multiple immunological changes that may contribute to the activity of such agents, as well as potential mechanisms that may re-establish tolerance during anti-PD1 therapy. A causal interrogation of each of the observed changes is needed to uncover the mechanisms underlying the observed phenomena. Moreover, our findings would be strengthened by extending these investigations to studies utilizing clinical specimens. Our results suggest that the retrieval and analysis of tumor draining lymph nodes may enhance the evaluation of immunotherapy responses to this class of agents and should be integrated into neoadjuvant immunotherapy clinical trial designs. Supplementary Tables 1 and 2. Conjugation of primary antibodies (34) was performed using Maxpar Conjugation Kits according to the manufacturer's protocol. Briefly, purified antibodies were run through a buffer exchange protocol using 50kDa ultra filtration columns (Amicon) and then partially reduced with 4mM TCEP (Thermo Scientific). Polymers were loaded individually with isotopically enriched metals, 113In (Trace Sciences), 115In (Sigma), 163Dy (Fluidigm), and 175Lu (Fluidigm). Metal-loaded polymers were then conjugated to their respective antibodies and washed (see Supplementary Tables 1  and 2). 194Pt (Fluidigm) was directly conjugated to the reduced antibody. Antibody concentrations in the wash buffer were quantified using Nanodrop. The final antibody concentrates were then diluted in a stabilization buffer (Candor) containing 0.3% sodium azide. Each antibody was titrated by testing a range of 3-4 serial dilutions clear positive controls, e.g. stimulated lymphocytes for cytokine staining, and identifying the concentration that permits discrimination while minimizing spillover signals.

Antibodies: A list of mass cytometry antibodies, isotopes, concentrations used for a T-cell subtyping panel and a B-cell subtyping panel is listed in
Animal model: 7-week female C57Bl/6J mice were purchased from Jackson and housed in the Johns Hopkins animal facilities. Mice were allowed to acclimate for at least 1 week prior to experimentation. MC38 (Kerafast) cells were maintained in DMEM-based media containing 10% FBS, 1% L-glutamine, 100U/mL penicillin/streptomycin, 10mM HEPES, 1mM sodium pyruvate, and 0.1mM non-essential amino acids (Life Technologies) in 5% CO2 at 37C. B16-F10 (ATCC) cells were maintained in RPMI 1640 with glutamine containing 10% FBS, 100U/mL penicillin/streptomycin, 1mM sodium pyruvate, and 0.1mM non-essential amino acids. On day 0, 2.5x10 5 MC38 or 5x10 5 B16-F10 cells are injected subcutaneously in the right hind limb (16). Tumors were allowed to be established and grow for 2 weeks. Mice were then administered Isotype or PD1 antibody (clone RMP1-14, BioXCell) 10mg/kg intraperitoneally on days 14 and 18. Tumor draining right inguinal lymph nodes were identified and harvested for analysis on day 21. For assessment of non-tumor draining lymph nodes, left inguinal lymph nodes were harvested. Tumor volume was also calculated on day 21 based on caliper-measured major (D) and minor (d) diameters using the following formula: V = 1/2 * D * d 2 . Right inguinal lymph nodes from agematched non-tumor bearing normal mice were used as an additional control. Harvested lymph nodes were weighed and then processed into single cell suspension by manual dissociation through a 70uM filter in PBS with 4mM EDTA. Total cell numbers from individual lymph nodes were counted using a hemocytometer and recorded. 2x10 6 cells were plated into each well of a 96-well microplate for staining.
Mass cytometry staining and data acquisition: Cells were stimulated in RPMI with 10% FBS containing 1X PMA/ionomycin/brefeldin A cocktail (Biolegend) for 3 hours in 5% CO2 at 37C. In the final 15 minutes, cells were incubated in media containing 5mM EDTA (35). After stimulation, cells were washed with PBS with 4mM EDTA. Viability was marked by incubation in palladium chloride (Sigma) dissolved in DMSO and diluted in PBS to 500nM for 5 minutes at room temperature, after which complete media was used to quench the remaining palladium (36). TDLN cells from one mouse from each of the three groups were stained with a unique metal barcode tagged to CD45 for 25 minutes at room temperature ( Figure 1D). Cells from these 3 animals were then combined into a single tube for downstream processing. Each tube was then divided into two equal aliquots and blocked with 1ug anti-mouse Fc block (BD Biosciences) for 10 minutes at room temperature followed by either a cocktail of antibodies for the T cell or the B cell panel ( Figure 1D, Supplementary Tables 1 and 2). Surface marker staining was performed in cell staining buffers (Fluidigm) at room temperature for 30 minutes. Intracellular cytokine staining was performed using Cytofix/cytoperm kit (BD Biosciences) per manufacturer's protocol. Upon completion of staining, cells were stored in fresh 1% methanol-free formaldehyde in PBS (Thermo Scientific) until the day of data collection. Just before data collection, all cells were labeled with rhodium (Fluidigm) at 1:1000 for 45 minutes at room temperature. All events were acquired on a Helios™ mass cytometer (Fluidigm). Mass cytometry data were acquired at the University of Maryland School of Medicine Center for Innovative Biomedical Resources (CIBR) Flow Cytometry and Mass Cytometry Core Facility, Baltimore, Maryland.
Mass cytometry data preprocessing: Randomization, bead normalization, and bead removal of data collected were performed on CyTOF software (Fluidigm) v6.7. Using FlowJo (BD) v10.5, single cell events were identified by gating a tight population based on cell length and rhodium signal. Dead cells were then eliminated by manually gating out cells positive for 106Pd and 108Pd on a biaxial plot. Debarcoding was carried out by manual gating to select for events that are singly positive for the barcode of interest and doubly negative for the remaining two. Each individual preprocessed sample was exported as a separate fcs file for analysis. Briefly, single cell suspensions from lymph nodes were Fc blocked for 10 minutes at 4C after which the samples were incubated in the antibody cocktail in FACS buffer for 30 minutes on 4C. Cells were washed twice with FACS buffer and data was collected on CytoFLEX (Beckman). All supervised gating analyses were performed on Cytobank.
Identification and visualization of clusters from mass cytometry datasets: Both unsupervised and hierarchal gating analyses were performed in parallel. For unsupervised analysis, a modified analysis pipeline from Nowicka et al. (37), was utilized using R v3.5. FlowSOM algorithm (38) was used to define 20 metaclusters for T cell subtyping analysis and 25 metaclusters for B cell subtyping analysis. Resulting metaclusters from the respective datasets were annotated into 7 final T cell subtypes and 10 final B cell subtypes. These clusters were then visualized using a two-dimensional uniform manifold approximation and projection (UMAP) dimensional reduction algorithm (39) separately for T cell and B cell subtype datasets. 2000 cells per sample were used for visualization. Cluster proportions from each lymph node were also converted into absolute cell numbers of each cell type by multiplying the quantified total number of cells from the corresponding lymph node. Tables in csv format for cluster proportions were outputted for statistical analyses (Supplementary Tables 3-4). Expression of functional markers, CD69, CD40, CTLA4, PD1, PDL1, BTLA, IFN, TNF, and IL2, were compared as median intensities across clusters and groups (Supplementary Tables 5-6). For hierarchal gating analysis, preprocessed fcs files were imported into Cytobank. All gating hierarchies for T cell subtypes and B cell subtypes are shown in Supplementary Figures 2 and 3, respectively. Functional markers were analyzed by gating for positive populations based on biaxial plot visualizations for every T cell subtype and B cell subtype.
Statistical analysis: To compare the means of the three groups, as proportions, numbers, or median expression values, while taking into account any potential differences related to the multiplex batching, we performed repeated measures one-way ANOVA followed by pairwise testing between the groups. Adjusted p-values (FDR) of less than 0.05 were considered statistically significant.
Study approval: Experiments were performed in accordance with protocol approved by the institutional Animal Care and Use Committee.

Author contributions:
WJH performed the experiments, developed experimental and analytical methods, analyzed data, and wrote the manuscript. MY assisted with data analysis and edited the manuscript. SC performed the experiments, assisted with data analysis, and reviewed the manuscript. RM assisted with developing analytical methods and reviewed the manuscript. MBS assisted with experimental methods and reviewed the manuscript. LD developed analytical methods, and reviewed the manuscript. EJF supervised analytical methods development, and reviewed the manuscript. EMJ supervised the design, analysis of all experiments, and reviewed and edited the manuscript.        T   T  T  T  T  T  T