Genetic and Environmental interactions contribute to immune variation in rewilded mice

Abstract

The relative and synergistic contributions of genetics and environment to inter-individual immune response variation remain unclear, despite its implications for understanding both evolutionary biology and medicine. Here, we quantify interactive effects of genotype and environment on immune traits by investigating three inbred mouse strains rewilded in an outdoor enclosure and infected with the parasite, Trichuris muris. Whereas cytokine response heterogeneity was primarily driven by genotype, cellular composition heterogeneity was shaped by interactions between genotype and environment. Notably, genetic differences under laboratory conditions can be decreased following rewilding, and variation in T cell markers are more driven by genetics, whereas B cell markers are driven more by environment. Importantly, variation in worm burden is associated with measures of immune variation, as well as genetics and environment. These results indicate that nonheritable influences interact with genetic factors to shape immune variation, with synergistic impacts on the deployment and evolution of defense mechanisms.

Fig. 1:. Interactions between Genotype and Environment contribute to variation in immune composition in murine PBMCs.

(A) 51 C57BL/6, 41 129S1, 43 PWK/PhJ different mice strains (total =135) were used in these experiments. Some were kept in conventional vivarium (n = 63), and some were housed in the wild enclosure (Rewilded), (n = 72) for total of 5 weeks. In addition, some were exposed to 200 T. muris L3 eggs (n =61) while the others were left unexposed (n = 74). Experiments was repeated twice in July, Block 1 (n = 61) and August, Block 2 (n = 74). Blood was collected for Flow cytometry, Plasma cytokine assessment and CBC analysis. MLN cells were also collected for Single cell RNA sequencing (ScRNA-seq), Flow cytometry with two different panels - lymphocyte and myeloid panel as well as cytokine profiling of supernatants from MLN stimulated cells. For Flow cytometry visualization, only Block 2 was used for downstream analysis. (B) Bar plots showing the pseudo R² measure of effect size of predictor variables and interactions as calculated by multivariate distance matrix regression analysis (MDMR) (C) PCA of immune cell clusters identified by unsupervised clustering in the blood and the density of each population along the PCA and (D) Bar plot showing variance on PC1 and PC2 axis of PCA plots in (B). (E) Barplot showing GMFI of CD44 on blood CD4+ T cells, (F) % of Tbet+ CD4 T cells of Live, CD45+ T cells. Bar plots showing percentage of Neutrophils, Lymphocytes, Monocytes, Eosinophils, and Basophils out of total at (G) 2 weeks post re-wilding and (H) 5 weeks post re-wilding based on assessment by CBC with differentials. Statistical significance was determined by a based on MDMR analysis with R package (B) or based on one-way ANOVA test between different groups with Graph-Pad Software (D), (E) and (F). For (E) and (F) direct comparison was done between groups of interest with one-way ANOVA test. Data are displayed as mean ± SEM. ns p>0.05; * p<0.05; ** p<0.01; *** p<0.001; **** p<0.0001

Fig. 2:. Interactions between Genotype, Environment and Infection contribute to variation in immune composition in murine MLNs.

(A) Bar plots showing the pseudo R² measure of effect size of predictor variables and interactions as calculated by multivariate distance matrix regression analysis (MDMR). (B) PCA of immune cell clusters identified by unsupervised clustering in the MLN with the lymphoid panel and the loading factor of each population along the PCA. (C) Bar plot showing variance on PC1 and PC2 axis of PCA plots in (A). (D) MLN cell count from each mouse group and (E) pseudo R² measure of effect size of predictor variables and interactions as calculated by MDMR analysis based on MLN cell count (F) Bar plots depicting GMFI of CD44 on MLN CD4+ T cells. n = 5–15 mice per group, Block = 2. (G) Representative Histogram from Block 1 and 2 showing concatenated files from T. muris infected and rewilded mice of each mice strain (H) Representative Histogram showing concatenated files from different groups of mice in Block 2 with corresponding (I) Bar plots depicting proportion of B cells expressing CD44 on MLN cells. n = 6–18 mice per group, Block = 2. Statistical significance was determined by a based on MDMR analysis with R package for (A) and (E) or based on one-way ANOVA test between different groups with Graph-Pad Software for (C) and (D). For (F) and (I), one-way ANOVA test was used to test statistical significance between the different groups of interest. Data are displayed as mean ± SEM. ns p>0.05; * p<0.05; ** p<0.01; *** p<0.001; **** p<0.0001.

Fig. 3:. Genotype has a bigger effect on cytokine responses as a fixed predictor variable.

(A) Bar plots showing the pseudo R² measure of effect size of predictor variables (left) and interactions (right) as calculated by multivariate distance matrix regression analysis (MDMR) from plasma cytokine data, Block 1 and Block 2 (B) PCA showing circulating plasma cytokine levels from individual mice and (C) their loading factors (D) Heatmap depicting levels of circulating cytokines in plasma (E) Bar plots showing the pseudo R² measure of effect size of predictor variables and interactions as calculated by multivariate distance matrix regression analysis (MDMR) from cytokine supernatant data of MLN cells stimulated with either CD3/CD28 beads, LPS, Candida albicans, Clostridium perfringes, Bacteroides vulgatus or T. muris antigens. (F) Bar plots showing the pseudo R² measure of effect size of predictor variables and interactions as calculated by multivariate distance matrix regression analysis (MDMR) from cytokine supernatant data of MLN cells stimulated with T. muris antigen. Barplot showing transformed IFN-γ cytokine levels in the supernatant for (G) controls as well as (H) following stimulation with T. muris antigen. Statistical significance was determined by a based on MDMR analysis with R package for (A), (E) and (F). Plasma cytokine data included samples from both Block 1 and Block 2 while MLN Cytokine supernatant data included samples from Block 2 alone due to technical problems with stimulation assays from Block 1. Data are displayed as mean ± SEM. ns p>0.05; * p<0.05; ** p<0.01; *** p<0.001; **** p<0.0001

Figure 4:. Single Cell Sequencing Analysis for assessing immune variation in cellular composition and cytokine profiles.

(A) UMAP visualization of scRNASeq data identifying 23 major immune cell subset (B) Bar plot showing the pseudo R² measure of effect size of predictor variables and interactions as calculated by multivariate distance matrix regression analysis (MDMR) based on cellular composition of cells identified in (A) in each mouse (C) PCA of MLN cellular compositional data as determined by scRNAseq analysis Bar plots showing percentage of B Follicular cells (D) and CD4 T cells (E) based on the scRNA-Seq identified in (A). (F) Cytokine-expressing cell clusters. (G) Proportion of cells expressing cytokine-related genes of those identified in (F). (H) Bar plots showing the pseudo R² measure of effect size of predictor variables and interactions as calculated by multivariate distance matrix regression analysis (MDMR) based on data from proportion of cytokines expressing cells identified in (F). (I) Bar plot showing cytokine expressing cells identified in the (F). Statistical significance was determined by a based on MDMR analysis with R package for (B) and (H). Data are displayed as mean ± SEM. ns p>0.05; * p<0.05; ** p<0.01; *** p<0.001; **** p<0.0001.

Figure 5:. Genetics and Environmental factors predicts outcomes during exposure with T. muris parasite.

Significant variation in worm burden among exposed mice, 3 weeks after inoculation with 200 eggs of Trichuris muris per host. (A) Worm burden (number of nematodes remaining in the caecum at that timepoint) followed a negative binomial distribution. (B) Worm burden depicted as Number of worms per mouse and Percentage of mice (Prevalence) still infected by worms. Each was predicted by a combination of Genetic and Environmental effects, including Gen*Env for worm burden (see text). (C) When we used PC2 from the scRNAseq data (Fig 3B) as an index of immune variation among individuals in our statistical models, we found that Gen*Env was no longer significant. Instead, the best model included main effects of host strain (C57BL/6 vs 129SL vs PWK/PhJ), environment (Lab vs RW), and PC2. The figure depicts 1000 model-estimated values for the effect of each predictor on worm burden. (D) Loading Factors for PC2 of the scRNAseq dataset.