The cellular composition of the human immune system is shaped by age and cohabitation

Abstract

Detailed population-level description of the human immune system has recently become achievable. We used a 'systems-level' approach to establish a resource of cellular immune profiles of 670 healthy individuals. We report a high level of interindividual variation, with low longitudinal variation, at the level of cellular subset composition of the immune system. Despite the profound effects of antigen exposure on individual antigen-specific clones, the cellular subset structure proved highly elastic, with transient vaccination-induced changes followed by a return to the individual's unique baseline. Notably, the largest influence on immunological variation identified was cohabitation, with 50% less immunological variation between individuals who share an environment (as parents) than between people in the wider population. These results identify local environmental conditions as a key factor in shaping the human immune system.

Figure 1

Data-driven analysis of immunological variation reveals biologically meaningful co-correlations between individual immune parameters. (a) On the complete dataset of 638 individuals (most recent sample only), a dendrogram of immune parameters was generated by hierarchical clustering on Euclidean distances of Spearman correlations between each parameter (left). Correlation plots using pairwise Spearman correlation coefficients between each two immunological parameters are shown (right). Coefficients are shown by the angle of eclipse (left-leaning, negative; right-leaning, positive) and colour (blue, negative; red, positive). Manually annotated thematic groups of immune parameters are shown by the colour bar next to the dendrogram. (b) Non-metric multidimensional scaling of pairwise Spearman correlations. The dataset is reduced from 54 immune parameters (54 dimensions) to a 2 dimensional representation, an exploratory approach to investigate the presence or absence of inter-relatedness between immune parameters. Each immune parameter is a point and the thematic groups are highlighted.

Figure 2

The human immune system is robustly maintained in multiple stable equilibriums. 177 individuals were sampled at least twice, allowing a dissection of inter- versus intra-individual variation. (a) Linear models were made for each immune parameter based on the multiple samples from each individual (ANOVA model: immune parameter ~ subject identifier + visit number). Open circles represent models built using all individuals (n=638) with multiple visits (177 individuals with up to 3 repeat visits; 921 visits in total). The filled circles represent models using only individuals who were continuously healthy between visits (152 individuals), filled squares represent models using only individuals who experienced acute gastroenteritis between their samples (24 individuals). For each cohort, R² values and (b) -log₁₀ of Bonferroni adjusted P values are shown for the linear models. (c) Proportion of the R² values from all volunteers attributable to either inter-individual differences or intra-individual differences. (d) Multidimensional scaling of the pre- and post-travel study visits. Each individual is represented twice; their first and second visits depicted with a dot or diamond respectively and linked by a grey line indicating immunological distance (50 individuals; 100 visits). Continuously healthy individuals (n=26) are shown in aqua, individuals with intervening acute gastroenteritis (n=24) are shown in orange. (e) Quantification of the immunological distance between the first and second visits for continuously healthy individuals (n=26) versus individuals with intervening acute gastroenteritis (n=24). A two-tailed Mann-Whitney test was used to compare the immunological distances.

Figure 3

Immunologial equilibria demonstrate elasticity following influenza vaccination. In a parallel cohort of 32 English individuals, volunteers were sampled prior to vaccination, with follow-up samples at day 7 and day 42 post-vaccination. (a) Samples were phenotyped, normalised to the day 0 value and assessed for change using paired t-tests. Unchanged variables are shown in grey, significantly modified variables are shown in red, with (b) boxplots for each significant immune parameter. Each parameter is labeled on the graph with an (uncorrected, two tailed) paired t-test P value. Boxes show median and interquartile ranges (IQRs), whiskers extend to 1.5 x IQR. (c) A linear model was made for each immune parameter based on the multiple samples from each individual (ANOVA model: immune parameter ~ subject identifier + visit number). For each parameter the proportion of the R² values from all volunteers attributable to either inter-individual differences or intra-individual differences was assessed for day 0 and 7 or (d) day 0 and 42. (e) Multidimensional scaling of the vaccination time-points. Each individual is represented at day 0 (green), 7 (red) and 42 (blue) and linked by a grey line indicating immunological distance (32 individuals; 96 visits). (f) Quantification of the immunological distance between day 0 and 7, 7 and 42, and 0 and 42 for each volunteer (n=32), with paired t-test. (g) For each volunteer, a Z-score was calculated for each parameter at day 0 and 42, indicating standard deviations from the mean value. Correlation analysis indicates the line of best fit. *, p<0.001; **,p<0.0001.

Figure 4

Age is a major determinant of immunological equilibria. (a) Each immune parameter was correlated with age and each other immune parameter using pairwise Spearman r values. The arrangement of the immune parameters is determined by their correlation with age. Inset panels show Spearman’s r (upper) and -log₁₀ Bonferroni corrected P values (lower) plotted against age. (b) R² for each immune parameter for models incorporating gender (filled circles), age (grey squares) or both gender and age (open circles) as the independent variable(s), with (c) accompanying -log₁₀ of Bonferroni corrected P values. (d) to (m) Individual scatterplots for each of the immune parameters with significant association with age. The percentage values of the flow parameters, or log₁₀ cytokine concentrations, are plotted against age, with women in red, men in blue and linear regression lines for the whole cohort in black. Data shown for (d) CD4⁺ RTE (p=3x10^-18), (e) transitional B cells (p=8x10^-7), (f) CD8⁺ RTE (p=1x10^-8), (g) T_H1 (p=3x10^-20), (h) CD4⁺ IL-2⁺ T cells (p=5x10^-23), (i) Tc1 (p=5x10^-31), (j) CD8⁺ IL-2⁺ T cells (p=8x10^-33), (k) CD8⁺ T cells (p=8x10^-11) (l) iNKT cells (p=6x10^-10) and (m) serum IL-6 (p=7x10^-11).

Figure 5

Immunoprofile is not influenced by BMI, depression, or anxiety. (a) R² for each immune parameter for models incorporating BMI (filled circles; 213 individuals), age (filled squares; 367 individuals) or both BMI and age (open circles; 213 individuals) as the independent variable(s), with (b) accompanying -log₁₀ of Bonferroni corrected P values. Analysis excluded children (<18 years). (c) Relationship between BMI and age in the analysed cohort. (d) The relative R² contributions of age and BMI to immune parameters that were significant in a model including both BMI and age (adjusted P<0.05). (e) R² and -log₁₀ adjusted P values for each immune parameter for a model incorporating HADS anxiety score and HADS depression score (235 individuals).

Figure 6

Parenthood shapes the immune system towards a shared equilibria. (a) 140 individuals were identified as adult (18-65 years) biological parents with a child still living at home. The immune profile (54 parameters) was compressed using multidimensional scaling (k=2) of the correlation matrix between individuals, visualising pairwise Spearman’s correlation coefficients between each individual. The immunological distance between each male:female pair is indicated by the connecting gray line. (b) The immunological distance, as measured by multidimensional scaling, between parental pairs versus random male:female pairs. To generate the random distribution, each male was computationally paired in a random fashion with 5 females from the parental dataset. Distributions were compared using a two-tailed Mann-Whitney test. *, p=8x10^-11.