Despite the importance of the immune system in many diseases, there are currently no objective benchmarks of immunological health. In an effort to identifying such markers, we used influenza vaccination in 30 young (20–30 years) and 59 older subjects (60 to >89 years) as models for strong and weak immune responses, respectively, and assayed their serological responses to influenza strains as well as a wide variety of other parameters, including gene expression, antibodies to hemagglutinin peptides, serum cytokines, cell subset phenotypes and in vitro cytokine stimulation. Using machine learning, we identified nine variables that predict the antibody response with 84% accuracy. Two of these variables are involved in apoptosis, which positively associated with the response to vaccination and was confirmed to be a contributor to vaccine responsiveness in mice. The identification of these biomarkers provides new insights into what immune features may be most important for immune health.
A systems analysis of immune biomarkers in 89 young and older adults revealed age‐dependent and age‐independent features, including markers of apoptosis that correlated with antibody responses to a seasonal influenza vaccine.
Influenza hemagglutinin peptide arrays reveal age‐associated effects that correlate with both pre‐existing and vaccine‐induced antibody titers.
Age‐dependent and age‐independent baseline immune parameters correlate with and substantially predict the serological response to a seasonal influenza vaccine.
Soluble FasL and gene modules associated with apoptosis are predictors of the serological response to an influenza vaccine, which was abrogated in Fas‐deficient mice.
The role of the immune system in human health is widely recognized. This is particularly evident in aging, where the immune response is often attenuated as observed in the decreased protection following vaccination and high prevalence and severity of infectious diseases in the elderly. Of the current estimated 36 000 influenza deaths in the United States each year, 90% occur in people older than 65 years of age (Thompson et al, 2003), and in this age group influenza alone contributes to 64% of the estimated total economic annual burden of $10.4 billion (Molinari et al, 2007). Moreover, the most common causes of death in the older population, including coronary heart disease, stroke and cancer, may also be due, at least in part, to immune system dysfunction (Hansson, 2005; de Visser et al, 2006).
The changes that occur with age affect both the innate and the adaptive components of the immune system. These include major alterations in cell subset composition, partly due to thymic involution (Fagnoni et al, 2000), defects in antigen presentation, reduced cytotoxic function, a restricted T‐cell repertoire and defects in antibody responses (Pawelec and Larbi, 2008; Weiskopf et al, 2009).
Because of the complexity of the immune system, the majority of the studies to date have measured a relatively small number of variables. However, recent advances in technology enable many parameters to be measured quickly and at low cost, allowing much broader immune profiling. This ‘systems immunology’ approach to vaccination has been pioneered in the studies of Sékaly and colleagues (2008) and Pulendran and colleagues (2009) in analyses of yellow fever vaccine, and more recently in the study of Nakaya et al (2011) in influenza vaccination, who used gene expression data from post‐vaccination samples to identify differences between responders and non‐responders. Here we have used pre‐vaccination samples to determine whether it is feasible to find baseline biomarkers from blood that could predict the serological response to the trivalent, inactivated seasonal influenza vaccine (TIV) in advance of its administration. We comprehensively characterized baseline young and older individuals’ blood using whole‐genome DNA microarrays, peptide microarrays for detection of influenza‐specific antibodies to linear epitopes, multiplexed serum protein assays comprising 50 cytokines and chemokines, and high‐resolution immune phenotyping for 15 immune cell subsets. We also measured the antibody responses to TIV and the signaling responses to cytokines in various cell subsets from peripheral blood, using antibodies to phosphorylated STAT proteins, also known as phosphoFlow (pFlow).
A strong negative effect of pre‐existing anti‐influenza antibodies with hemagglutinin (HA) inhibition (HAI) activity on the fold‐change HAI titer following vaccination has been reported in several studies (Beyer et al, 1996; He et al, 2008). We characterized the reactivity of influenza‐specific antibodies recognizing HA peptides from baseline serum samples with the aim of (1) investigating the HA linear epitope specificities in pre‐vaccine antibodies and (2) finding the minimal set of pre‐existing antibodies that could potentially explain the negative effect of baseline HAI titers on the response to the influenza vaccine. We found a set of pre‐existing antibodies recognizing viral HA peptides that can distinguish individuals with high versus low pre‐vaccine HAI titers. We also identified antibodies against HA peptides that associate with the HAI response to the vaccine. Using a machine learning approach to identify other relevant predictors of vaccine response, we found age‐related and ‐independent variables, which correlated with and predicted the serological response with 84% accuracy. Two of these markers were involved in apoptosis, which was confirmed to be important in influenza vaccine responses in mice.
The immune predictors found here may also be important in the response to other vaccines and help to define metrics of immunological health (Davis, 2008).
Antibody responses to influenza vaccination in young versus older vaccine recipients
We carried out an influenza vaccine study in 91 young and older ambulatory subjects to better characterize immune parameters from peripheral blood that associate with vaccine responsiveness. Eighty‐nine individuals’ serum specimen sets (pre‐ and post‐vaccine) were complete for testing the vaccine response. The number of previous influenza immunizations, history of influenza infection, demographic data and other characteristics of the subjects are depicted in Table I. We measured a comprehensive set of immune parameters in pre‐vaccine samples (day 0) (Figure 1) and the antibody responses by the HAI assay at day 0 and 28±7 days after vaccination. We used the HAI assay, because it is a direct measure of antibodies produced in vivo against viral HA proteins in response to infection or vaccination and is generally accepted as the gold standard for gauging the effectiveness of influenza vaccines. Seroconversion (defined as an HAI titer post‐/pre‐vaccine ⩾4) was computed for each strain in the vaccine. The individuals were classified according to the breadth and magnitude of the response as poor responders (PR), if seroconversion was achieved for one strain or none, or good responders (GR), if seroconversion was achieved for two or all three strains in the vaccine (Figure 2A). By this definition, 80% of young individuals were GRs and only 38% were GRs in the older cohort (Figure 2B). These results are in agreement with several previous studies.
An important feature of the immune response to influenza in humans is the presence of pre‐existing HAI antibodies, the titer of which negatively correlates with responsiveness (fold increase) to the vaccine (Beyer et al, 1996; He et al, 2008). Pre‐existing flu‐specific memory CD4+ T cells seem to have an important role by activating (CD56dim)‐NK cells that are able to inhibit antigen presentation by dendritic cells, thereby suppressing the subsequent CD4+ T‐cell help to B cells (He et al, 2008). Therefore, we calculated for each subject the pre‐vaccine HAI geometric mean titer for all three strains in the vaccine (pre‐GMT) and compared pre‐vaccine titer in young individuals with older individuals. The older cohort had significantly less pre‐GMT than the young cohort (P=0.0017; Figure 2C), which suggests that pre‐GMT is a stronger determinant of the vaccine response in young than in older adults. To address this in detail, we compared the pre‐GMT levels in PRs versus GRs in young or older individuals. The differences between PRs and GRs were more pronounced for the young group (P=0.005) than for the older group (P=0.039; Supplementary Figure 1).
Pre‐existing antibodies to HA peptides and pre‐vaccine HAI titer
To characterize the pool of pre‐existing HAI antibodies with more precision, we developed peptide microarrays for profiling the antibody repertoire directed against viral HA peptides. These studies were restricted to 76/89 subjects because of limitations in the availability of serum from some subjects. Our peptide microarray technology is able to cover ∼40% of the HA protein sequence for each vaccine strain. Thus, to maximize the number of positive hits, we selected an HA region found in an H5N1 influenza strain containing epitopes for neutralizing antibodies against different influenza strains (broadly neutralizing antibodies; Khurana et al, 2009).
To determine whether the reactivity against these HA peptides at baseline was associated with the presence of pre‐existing HAI antibodies, we compared individuals with high versus low pre‐GMT (above versus below the group average, respectively). As anticipated, there were significantly more GR within the subjects with low pre‐GMT compared with those with high pre‐GMT (P<0.01 by Fisher's exact test; Figure 3A). We identified four peptides to which serum reactivity was significantly more prevalent in individuals with high compared with low pre‐GMT (P<0.05, false discovery rate (FDR) Q<0.1; Figure 3B). This result indicated that individuals with a high pre‐vaccine HAI titer, associated with poor responses here and elsewhere, are more likely to have antibodies directed against peptides H1_23, H3_5, H3_8 and BH14 corresponding to potential linear HA epitopes from each of the virus strains in the vaccine.
Pre‐existing antibodies to HA peptides and prediction of the HAI response
We next addressed the possibility that the pre‐existing antibodies to HA peptides were also associated with and could predict the HAI titer in response to the vaccine. To do so, we used the ‘elastic net’ (Friedman et al, 2010) with fivefold cross validation, a machine learning procedure that enables the discovery of relevant features from high‐dimensional data and estimates the performance of the resulting prediction model, such that the most informative factors are selected out of the pool of available parameters in an automated manner, to generate the model with the lowest error. In the fivefold cross validation, the samples are randomly partitioned into five sets, of which four sets are used to ‘train’ the algorithm and blindly predict the outcome of the fifth set (test set). This process is repeated iteratively five times until all sets are tested. Then a cross‐validated area under the ROC curve (cvAUROC) can be computed to assess the goodness of the model without the known problem of overfitting.
By this procedure, we built a first model (Model 1) that included age, and identified a total of 19 potential epitopes recognized by pre‐vaccine antibodies, resulting in an cvAUROC=0.754 (Figure 4A), which is significantly better than the model including age alone (Model 2), which yielded an cvAUROC=0.686 (P<0.0001 by likelihood‐ratio test). Not surprisingly, in Model 1 age had the largest regression coefficient (contribution to the model; Figure 4B). Of note, two of the peptide predictors identified in Model 1 (H1_23, NH2‐FALSRGFGSGIINSNAPMD‐COOH and H3_8, NH2‐SSGTLEFNNESFNWTGVTQ‐COOH) were significantly more prevalent in subjects with high pre‐GMT (Figure 3B) and neither of these peptide predictors were age‐associated (see below). Thus, we identified two HA peptides that contribute to distinguish individuals with a high versus low pre‐vaccine HAI titer and to which high antibody levels are negatively associated with the HAI activity in response to influenza vaccination.
We then tested the predictability and changes in the selection of predictors when age is removed from the model (Model 3) and only peptide variables are used to predict the HAI response. By performing this procedure, 26 peptides were selected, yielding a cvAUROC=0.711, only moderately better than the model where age is used as the only predictor (P<0.001, by likelihood‐ratio test). Intriguingly, pre‐vaccine reactivity to several sets of adjacent peptides was computationally selected for the prediction model when the contribution of age was removed (Figure 4C). The fact that clusters of adjacent peptides (both by sequence and within the three‐dimensional structure of HA) were selected in these models highlights potential HA epitopes important for inhibiting the binding of influenza virus to cells.
Age‐related pre‐existing antibodies against HA peptides
To investigate whether reactivities to HA peptides were different in young versus older individuals, we performed simple regression analysis on each peptide variable and age, and corrected for multiple testing using permutation tests (Benjamini and Hochberg, 1995). Older subjects showed differences in the pre‐existing reactivity to 11/72 peptides at an FDR<20% (P <0.05, FDR Q<0.2; Supplementary Table 1); three of which were selected as predictors of vaccine response in Model 1 and four in Model 3 (Figure 4C, asterisks).
These findings demonstrate that the levels of some pre‐existing antibodies targeting viral HA peptides correlate with pre‐vaccine HAI titers and with the magnitude of the HAI antibody response to each of the three strains in TIV. The fact that the antibody levels to some of these peptides are lower in older subjects might reflect differences in influenza strain exposure or changes in their antibody repertoire with age (Yager et al, 2008).
Baseline parameters and prediction of the HAI response: Models 1 and 2
We next tested the hypothesis that the baseline levels of immune parameters and subjects’ characteristics are indicative of an ability to respond to vaccination. To do so, we first generated gene modules from whole‐blood gene expression data and estimated their corresponding regulatory programs using a method similar to what has previously been described (Segal et al, 2003). Of a total of 48 771 gene probes in the microarray per sample, we derived 109 gene modules and their corresponding regulatory programs (Supplementary Table 2 and Methods). In addition to these gene modules, we also explored other immune parameters, including serum proteins, immune cell subsets, signaling events and the subjects’ characteristics, such as demographics and CMV and EBV serological status. A total of 241 variables were used for feature selection and prediction of the antibody response. As with the peptide microarray data, we performed cross validation using the elastic net algorithm (Friedman et al, 2010), which is well suited for scenarios having more features than samples.
We first generated two models as with the peptide array data: Model 1 (including all features) that selected nine relevant features, which correlated with the HAI response and resulted in a cvAUROC=0.844 (Figure 5A), and Model 2 in which the immune features are excluded and only age is used to predict the HAI response. The latter model with age alone predicted the HAI response with an accuracy of cvAUROC=0.694, which is significantly less accurate than our full model (P=0.000012 by likelihood‐ratio test).
Only one of the nine features identified in Model 1 (gene module 47) had a positive regression coefficient (Figure 5B). This gene module is enriched in genes involved in apoptosis (P<0.01; module ‘APO’), including (1) GSTP1, which suppresses the anti‐apoptotic protein JNK (Gate et al, 2004; Yu et al, 2004), (2) the STAT4 inhibitor PIAS4, which stimulates transcriptional activity of p53 and increases Rb‐dependent corepression through recruitment to E2F‐responsive promoter (Bischof et al, 2006; Zhang et al, 2004), (3) IL17D, which acts by downregulating T cell‐mediated immune hypereactivity and proinflammatory cytokine production (Hamano et al, 2003; Villarino et al, 2003; Owaki et al, 2008; Yoshimoto et al, 2008) and (4) ZNF‐148 (also known as ZBP‐89), which encodes for a protein that binds to p53 and p21 (waf1), enhancing their transcriptional activity (Bai and Merchant, 2001; Zhang et al, 2010). Our regulatory program analysis indicated that the gene module APO is likely repressed by IRF9, XBP‐1 and Bcl‐3, among other transcription factors (Supplementary Figure 2) that promote cell proliferation and survival in diverse cell types, and negatively associate with apoptosis (Takaoka et al, 2003; Tsuno et al, 2009; Thorpe and Schwarze, 2010). Of note, XBP‐1 also controls expression of IL‐6, which in turn, promotes plasma cell growth and survival (Iwakoshi et al, 2003).
Age and pre‐GMT were selected as negative predictors of the HAI response in Model 1, which is consistent with previous reports (Beyer et al, 1996, 1999; Beyer et al, 2004; He et al, 2008). Gene module 85 had a negative regression coefficient as well, and is composed of genes associated with the regulation of B‐cell proliferation (P<0.0001; module ‘PROL’), such as CDK6, BCS1L, PBX4, CD320 and CD40L. The regulatory program of the PROL module indicated that it is likely repressed by TSC‐22d1 (TGF‐stimulated clone 22), implicated in the TGF‐mediated apoptosis (Uchida et al, 2000) and the androgen receptor (AR) coregulator ARA70 (NCOA4), which is known to be repressed by Akt and has been found to act as a tumor suppressor by inducing AR‐dependent apoptosis (Ligr et al, 2010). Activators of this gene module include the following: (1) PHB2, an anti‐apoptotic factor upregulated by Myc (O'Connell et al, 2003) (2) ID3, a target of Blimp‐1 and essential for B‐cell proliferation (Pan et al, 1999) and (3) STAT4. Of note, PIAS4, a suppressor of STAT4, was found in the APO module (Supplementary Table 2), a positive predictor of vaccine response (see above).
This suggested that upregulation of genes in the APO module and the concomitant downregulation of genes found in the PROL module is beneficial for cell renewal and optimal antibody production in response to immunization, which is consistent with the notion that the limited niche availability restricts the number of cells to respond to antigenic challenges (Discussion and Franceschi et al, 2000). Indeed, interesting data from Xiang et al (2007) demonstrated that crosslinking the IgG receptor FcγRII induces apoptosis of bone marrow plasma cells in mice and of plasmablasts in humans likely contributing to cell turnover (Xiang et al, 2007). This prompted us to compare the expression of the FcγRIIB gene (FCGR2B) in PRs versus GRs. However, there was not sufficient variation in the expression levels of FCGR2B for its selection and gene module construction (Methods). Therefore, we focused on FcγRIIA (FCGR2A), as in mice FcγRIIB and FcγRIIIA exert similar pro‐apoptotic characteristics (Fang et al, 2012), and mouse FcγRIIIA is predicted to be the counterpart of human FcγRIIA in terms of functionality and sequence similarity in the extracellular portion.
Strikingly, FCGR2A clustered with a number of genes participating in cell apoptosis in module 54 (Supplementary Table 2), with differential expression between PRs and GRs (Supplementary Figure 3) and selected for prediction of the HAI response in Model 3 (see below). This suggests that diverse pathways implicated in apoptosis must function properly for a high cell turnover and the establishment of new antibody‐producing cells.
To determine whether apoptosis genes were also linked to the presence of pre‐existing antibodies to the potential linear epitopes identified using the peptide array technology, we divided the subjects according to their gene expression levels of the module APO into high or low, if they were above or below the group average, respectively (Figure 6A). As anticipated, the fraction of GRs was much higher in subjects with high expression of genes in the module APO (P=0.00007). Significant differences were also observed in many of the pre‐existing antibodies to linear determinants (P<0.05, FDR Q<0.1; Figure 6B). Of note, some of these antibodies (Figure 6B, asterisks) were identified to be predictive of the HAI response in our previous analyses (Figure 4).
For serum cytokines, we found that soluble Fas ligand (sFasL) and IL‐12p40 were negative predictors (Figure 5B). sFasL prevents membrane‐bound FasL from binding to Fas, thereby decreasing susceptibility to apoptosis (Knox et al, 2003). This strengthens our previous findings and indicates that the presence of soluble apoptosis inhibitors could affect the antibody response to vaccination.
Additional immune features that were selected in Model 1 were the frequency of central‐memory CD4+ and effector‐memory CD8+ T cells, and the baseline STAT1 levels in CD8+ T cells, which are elevated in the elderly and negatively correlate with the pSTAT response to diverse cytokine stimuli.
Baseline parameters and prediction of the HAI response: Model 3
We next tested the possibility that new immune features may emerge from the classification model when age is excluded from the list of candidate predictors (Model 3). Indeed, excluding age resulted in the identification of a new set of features in addition to those previously discovered (Figure 7A), and yielded a cvAUROC=0.700, which is similar to Model 2 (age alone). Model 3 included four positive and four negative additional associations (Figure 7A and B). The positive predictors were the CD8+ naive T‐cell frequency, gene module 34 (cell‐to‐cell signaling and interaction, P<0.01), gene module 43 (RNA post‐transcriptional modification, P<0.01) and gene module 54 (cell death by apoptosis, P<0.01). The negative predictors incorporated in Model 3 were the frequency of CD8+ CD28− cells and of NK cells, and modules 79 and 92 both enriched for genes associated with metabolism of carbohydrates (P<0.01). In addition to age and pre‐GMT, the frequency of CD8+ CD28− cells was the only feature identified in previous studies to associate with poor vaccine response in the elderly (Goronzy et al, 2001).
The fact that three apoptosis‐related features (module 47, module 54 and sFasL) correlated with a strong vaccine response motivated us to investigate these findings in a mouse model. Thus, we tested two strains of mice lacking CD95 (Fas); MRL/Mpj‐Faslpr/J) and B6.MRL‐Faslpr/J. In addition to the lack of Fas, plasma cells from these mice do not express FcγRIIB and, therefore, they are also protected from FcγRII‐mediated apoptosis (Xiang et al, 2007). After vaccination with TIV, mice lacking Fas on either genetic background showed a significant reduction in their ability to produce influenza‐specific antibodies (P<0.05 by t‐test) (Figure 8A and B), which suggests a critical role of the apoptosis system in the antibody response in vivo.
Age‐dependent and independent predictors of vaccine response
It has been well established that many older individuals fail to generate protective antibody titers following vaccination with TIV. Therefore, given our model for feature selection where features are not necessarily orthogonal, i.e., they can correlate with one another (Methods); some of the predictors could be correlated with age. To identify age‐dependent and ‐independent features in our data, for each predictor we performed a multiple regression analysis on all samples by age, gender and CMV status (Supplementary Table 3). We used an FDR lower than 5% (P<0.05; FDR Q<0.05) as a cut‐off value for significance. Of the 16 immune predictors selected in Model 3, only 2 were not associated with age (Supplementary Figure 3, asterisks). These are the levels of sFasL and expression of the PROL module (module 85). The elderly expressed low levels of positive predictors and high levels of negative predictors (Supplementary Table 3), with the exception of pre‐GMT, which were reduced in older compared with young subjects. This suggests that the effect of age on the antibody response overcomes the effect of pre‐existing immunity to the vaccine.
In this study, we comprehensively characterized the immune system of young and older subjects across multiple technological platforms and identified new candidate traits that associate with and partially predict the humoral response to vaccination in both a young and an old cohort.
We first characterized the reactivity of influenza‐specific antibodies recognizing HA peptides from baseline serum samples with the aims of (1) examining the HA‐linear epitope specificities in the pre‐vaccine antibodies and (2) finding the minimal set of pre‐existing antibodies that could potentially explain the negative effect of baseline HAI titers on the response to the influenza vaccine observed here and elsewhere (Beyer et al, 1996; He et al, 2008). We found a set of four peptides, the reactivities of which correlated with the pre‐vaccination HAI titer. In two of these, robust activity negatively correlated with the HAI response to influenza vaccination. Thus, our system can rapidly identify GRs from PRs based on their pre‐vaccination antibody repertoire against HA protein regions with known neutralization activity (Khurana et al, 2009), as conserved sequences are characteristic for neutralization epitopes on HA across different influenza virus strains (Yamashita et al, 2010). The negative effect of pre‐existing antibodies on the response to vaccination has been suggested to be due to pre‐existing flu‐specific memory CD4+ T cells that inhibit antigen‐presentation by dendritic cells, thereby suppressing the subsequent CD4+ T‐cell help to B cells (He et al, 2008). Therefore, it is possible that such memory CD4+ T cells are able to maintain a background T‐cell‐dependent antibody production represented by the reactivities against the linear epitopes we found here. Further studies are needed to clarify the mechanisms underlying these observations.
We identified both negative and positive peptide predictors of the HAI response, many of which were age‐associated. These may partially explain the overall poor response observed in aging. However, neither of the reactivities to peptides associated with pre‐vaccine HAI changed with age, although the levels of baseline HAI titers were lower in older subjects. This could be due to the limitations of the assay that (1) does not detect conformational epitopes and (2) does not allow for the identification of antibodies that depend on glycosylation for binding. Several previous studies have compared the pre‐GMT (baseline titers) with influenza in elderly versus young individuals with no general consensus (Clark et al, 2009; Vajo et al, 2010; Chen et al, 2011; Ehrlich et al, 2012). These discrepancies have been examined in detail by Sasaki et al (2011) who used plasmablast‐derived polyclonal antibodies and found greater influenza‐specific broad reactivity of antibodies from older versus young vaccine recipients. Thus, the pre‐GMT in the elderly is more dependent on their previous exposure history than that of young subjects. This is relevant for vaccinology, as it suggests differential strategies in young versus older subjects based on the historical findings of circulating strains.
Additional predictors of vaccine response included various immune features connected with apoptosis function. Consistent with this, apoptosis‐deficient (lpr) mice exhibited poor serologic responses to the influenza vaccine. It is known that the lpr mice have defects in T‐cell and B‐cell development. Thus, our results might not be only due to apoptosis defects in these mice (although this is suggested by the apoptosis gene modules and sFasL), but these other factors. However, there is no agreement as to whether the germinal center (GC) formation and B‐cell responses in these mice are normal (Smith et al, 1995; Takahashi et al, 2001). In particular, our results argue against a previous study that reported memory and antibody‐forming cell populations appear to be normal in lpr mice (Smith et al, 1995). However, more recent studies have shown that these mice have defects in clonal selection and the establishment of the memory B‐cell repertoire (Takahashi et al, 2001).
Apoptosis has a fundamental role in lymphocyte development and in the termination of the immune response. In the GC reaction, an optimal apoptotic machinery ensures the survival of only high‐affinity plasmablasts (Smith et al, 2000; Takahashi et al, 1999), which we found to decrease with age in a recent study (Sasaki et al, 2011). In addition, a functional Fas‐mediated pro‐apoptotic program is required for clearance of reactive T and B cells after an immune response takes place (Nagata, 1999), as well as for augmentation and maintenance of naive cells, which ultimately results in an optimal balance between the naive and memory cell pools (Zhou et al, 1995). When an immune response ends, the expanded clones of effector cells must be reduced in size. This allows the immune system to cope with new influenza virus challenges, such as those originating from antigenic drift and/or found in a new influenza vaccine preparation. Thus, the appropriate regulation of apoptosis and rapid removal of reactive memory cells may improve immune responses to newly encountered antigenic challenges as observed in the study by Haynes et al (2005), in which new functional CD4+ naive T cells developed after CD4+ T‐cell depletion in aged and young mice, restoring optimal responses to antigen (Haynes et al, 2005). This has a number of implications for vaccination and suggests possible ways to overcome, at least partially, the immune system defects observed here. This would also suggest competition and a restricted immunological niche largely occupied by memory cells that, unless cross‐reactive, cannot respond optimally to novel antigens.
We also find an age‐dependent reduction in the gene modules associated with apoptosis. However, at the level of the serum cytokine sFasL, there was no difference between young and older subjects. This indicates that aging results in the accumulation of factors contributing to less cell renewal and complements previous observations of decreased susceptibility to apoptosis in cells from older subjects (Zhou et al, 1995; Salvioli et al, 2001; Hsu and Mountz, 2010). Indeed, Hsu et al (2006) showed that T cells from nonagenarians have much higher levels of Fas and increased activation‐induced cell death than subjects 65 to 89 years old (Hsu et al, 2006), suggesting that apoptosis is beneficial in maintaining good immune health and attaining a longer life span. Although the role of apoptosis in the maintenance of immune tolerance has been well established, to the best of our knowledge, its association with a robust antibody response to a vaccine has not been shown previously.
The accumulation of memory cells in those older subjects who do not respond to vaccination could also be due to chronic antigenic stimulation. Thus, it has been postulated that the constant activation of clonal T‐cell populations may lead to increasing numbers of memory cells filling the immunological space and compromising responses to newly encountered antigens (De Martinis et al, 2005). Our results indicate that the increasing size of the memory compartment in older individuals that do not respond well could result from a loss of pro‐apoptotic activity and increased levels of cell survival factors. Indeed, recent investigations have shown that memory CD4+ T cells are more resistant to apoptosis than naive CD4+ T cells (Grayson et al, 2002; Jaleco et al, 2003) and, at the gene expression level, memory CD4+ T cells express lower levels of pro‐apoptotic‐related genes and higher levels of cell survival and proliferation‐related genes (Liu et al, 2001).
Our study also identifies a link between cell survival and proliferation, and defective antibody responses. We found that the CD40L gene clustered in a module enriched for proliferation genes (module PROL), which is a negative indicator of the vaccine response. Although molecules in the gene module PROL (e.g., CD40L, CD320) are critical for plasmablast growth and survival, results from different studies suggest that elevated proliferation signals from these molecules have a detrimental effect. For example, heightened CD40 signaling causes B cells to shunt into an extrafollicular plasma cell fate and this prevents the generation of long‐lived bone‐marrow plasma cells. This has consequences for the B‐cell response, including the premature termination of the humoral immunity and the disruption of GC formation in vivo (Erickson et al, 2002). Therefore, it is appealing to hypothesize that weaker responses to the vaccine in individuals with augmented CD40L expression result from the preferential generation of extrafollicullar plasma cells that compromise the accumulation of somatic mutations in the GC B cells, dampening late‐appearing high‐affinity antibodies, as demonstrated in transgenic CD40L mice (Kishi et al, 2010). An interesting finding in support of this hypothesis is the impairment in GC formation found in old mice with poor responses to vaccine (Eaton et al, 2004). Further studies in humans could potentially address this.
Our results also suggest the presence of other genes in the PROL module with similar functions to that of CD40L, which could contribute to the diminished antibody responses. For example, the gene CD320 has been shown to participate in GC differentiation by directing B cells to mature into plasma cells (Zhang et al, 2001).
The large variation in the responses to vaccines in humans creates a number of challenges. We have utilized a ‘systems biology’ approach that allows us to embrace this variation and to integrate diverse measurements in the same individuals to generate hypothesis on how a system functions. Systems biology approaches, which have been extensively applied in the study of metabolism networks and genetics, have only recently been utilized in vaccinology (Gaucher et al, 2008; Querec et al, 2009; Nakaya et al, 2011). For instance, a recent study by Nakaya et al (2011) has been successful in finding features of the innate immune response, 3 and 7 days after influenza vaccination, which partially predict the subsequent antibody response in young, healthy adults. One of the genes the authors found in their predictive signatures and confirmed in a mouse model, encoding for the Ca2+/calmodulin kinase IV (CaMKIV), was inversely correlated with the HAI response to the influenza vaccine. This kinase has been implicated in the regulation of CRE‐dependent transcription (Sato et al, 2006) and, consistent with our findings, a recent study has shown the requirement of CaMKIV for the accumulation of the anti‐apoptotic molecules Bcl‐2 and Bcl‐x and the survival of DCs (Illario et al, 2008). Despite the fact that we did not find CaMKIV in our predictive gene modules, we hypothesize that CaMKIV could also prevent apoptosis in other immune cells, resulting in an overall low cell renewal and poor immune response. The study of Nakaya et al (2011) is therefore complementary to the work described here, which identifies biomarkers before vaccination in different compartments of the human immune system, and across age groups. Being able to predict whether or not an individual will respond productively to a given vaccine is important clinically, because if the prognosis is poor, one might choose not to vaccinate or choose a more potent formulation. For example, in seasonal influenza vaccines, there is a 4 × formulation that has been reported to be more efficacious than the standard one used here (DiazGranados et al, 2013), and adjuvanted influenza vaccines, although not yet licensed in the United States, are available in many other countries.
In summary, we have used immune aging as a model for impaired immunity and identified biomarkers that point to what factors have the greatest role in the response to a seasonal influenza vaccine, now recommended annually for all persons older than 6 months of age with few exceptions. An important feature of our approach is that each individual is profiled for many facets of the immune system. This enables us to observe within single individuals or across groups of individuals the relationship between multiple factors as well as the effects these factors have on the immune response. Ultimately, this methodology, applied to further studies of influenza or other vaccines and infections, should enable us to gain a more complete understanding of immune system function and dysfunction, as well as identifying key variables and mechanisms of immunological health (Davis, 2008).
Materials and methods
Subjects and sample collection
Ninety‐one healthy donors (ages 20 to >89 years) were enrolled in an influenza‐vaccine study at the Stanford‐LPCH Vaccine Program during the fall of 2008 of which eighty‐nine completed the study. The protocol of this study was approved by the Institutional Review Board of the Research Compliance Office at Stanford University. Informed consent was obtained from all subjects in the study. All individuals were ambulatory and generally healthy as determined by clinical records. Females of childbearing potential were tested for pregnancy by a urine sample. Volunteers had no acute systemic or serious concurrent illness, no history of immunodeficiency, nor any known or suspected impairment of immunologic function, including clinically significant liver disease, diabetes mellitus treated with insulin, moderate to severe renal disease, blood pressure >150/95 at screening, chronic hepatitis B or C, and recent or current use of immunosuppressive medication. In addition, none of the volunteers were recipients or donors of blood or blood products within the past 6 months and 6 weeks, respectively, nor showed any signs of febrile illness on the day of enrollment and baseline blood draw. A total of ∼40 ml per visit peripheral blood were obtained at day 0 (pre‐vaccine) and 28±7 days after receiving a single intramuscular dose of trivalent seasonal influenza vaccine Fluzone (Sanofi Pasteur). Each dose of the vaccine contained 15 μg HA each of H1N1, H3N2 and B strains of the virus. Whole blood was used for gene expression analysis (below). Peripheral blood mononuclear cells (PBMCs) were obtained by density gradient centrifugation (Ficoll‐Paque) and frozen at −80°C for 24–48 h before transferring to LN2. Serum was separated by centrifugation of clotted blood and was stored at −80°C before use. Whole blood, PBMCs or serum from the first visit (baseline: day 0) were processed and used for determination of gene expression, leukocyte subset frequency, signaling responses to stimulation, serum cytokine and chemokine levels, and CMV and EBV serostatus by ELISA (Calbiotech, San Diego, CA). Serum samples from day 0 and day ∼28 were used for HAI titer determination.
HA microarray design and fabrication
We printed streptavidin‐surface glass slides (Arrayit, Sunnyvale, CA) with biotinylated 19‐mer peptides (Sigma Aldrich, St Louis, MO) at 0.5 mg/ml in PBS plus 2.5% glycerol (Sigma) in triplicate using a VersArray ChipWriter Compact robotic microarrayer and ChipWriter Pro software (BioRad). We then blocked printed peptide arrays with biotin at 0.1 mg/ml (Sigma) in PBS for 10 min at room temperature (RT) with 500 r.p.m. orbital agitation. After several washes with peptide‐binding buffer (PBB): 50 mM Tris pH 7.5, 150 mM NaCl, 0.05% NP‐40 plus 2.5% FCS, we incubated the arrays with the subject's sera diluted in PBB. Following primary incubation for 2.5 h at RT with light rocking agitation, we rinsed the arrays with PBB, then with PBST (PBS, pH 7.5 with 0.05% Tween‐20), before incubating with goat anti‐mouse IgG or goat anti‐human IgG conjugated to Cy3 (Jackson ImmunoResearch) diluted to 1.5 μg/ml in PBST plus 20% FCS for 30 min at RT. We then rinsed arrays several times in PBST, then with PBS and finally with diH20, and dried arrays by centrifugation at 300 g for 5 min at RT. We immediately scanned processed microarrays with an Axon digital scanner and analyzed scanned images with Genepix Pro 6.1 software (Molecular Devices, Sunnyvale, CA).
Whole‐blood microarray analysis of gene expression
Total RNA was extracted from PAXgene blood RNA tubes (PreAnalytiX, USA) using the QIAcube automation RNA extraction procedure according to the manufacturer's protocol (Qiagene, Valencia, CA). The amount of total RNA, and A260/A280‐ and A260/A230‐nm ratios were assessed using the NanoDrop 1000 (Thermo Fisher Scientific, Wilmington, DE). RNA integrity was assessed using the Agilent 2100 Bioanalyzer (Agilent Technologies, Santa Clara, CA). For each sample, 750 ng of total RNA were hybridized to Beadchips (HumanHT‐12v3 Expression Bead Chip, Illumina, San Diego, CA) that contains 48 771 probes for around 25 000 annotated genes. The hybridized Beadchips were scanned on an Illumina BeadScan confocal scanner and were analyzed by Illumina's GenomeStudio software, version 2.0. After checking the quality of each individual array, the Feature Extraction Files were imported into R Bioconductor and analyzed using the Beadarray package for probe filtering, quantile normalization, replicate probe summarization and log2 transformation. The original microarray probe‐level data files were entered into the GEO repository under accession number GSE41080.
Leukocyte subset frequency determination
PBMCs were thawed in warm media, washed twice and stained with three separate anti‐human antibody cocktails containing: (1) anti‐CD3 AmCyan, CD4 Pacific Blue, CD8 APCH7, CD28 APC; (2) CD3 AmCyan, CD4 Pacific Blue, CD8 APCH7, CD27 PE, CD45RA PE‐Cy5; and (3) CD3 AmCyan, CD19 Alexa Fluor700, CD56 PE, CD33 PE‐Cy7, TCR APC, all reagents from BD Biosciences. Incubation with antibodies was performed for 40 min at 4°C. Cells were washed with FACS buffer (PBS supplemented with 2% FBS and 0.1% Na Azide) and resuspended in 200 μl FACS buffer. Data were collected using DIVA software in an LRSII instrument (BD Biosciences). Data analysis was performed using FlowJo 8.8.6 by gating on live cells based on forward versus side‐scatter profiles, then using double gating for singlet discrimination, followed by cell subset‐specific gating.
Phosphorylation of intracellular proteins by pFlow
Cells were thawed in warm media and rested at 37°C in RPMI with 10% FBS. Cells were then distributed in 96 deep‐well blocks (2 ml) and stimulated with IFN‐γ, IL‐2, IL‐6, IL‐7, IL‐10 or IL‐21 at 50 ng/ml or with 104 U/ml IFN‐α for 15 min. After stimulation, cells were immediately fixed with 1.5% PFA for 10 min at RT, washed with an excess of plain PBS and permeabilized with 95% ice‐cold methanol for 20 min on ice. Different stimuli conditions were barcoded using a 3 × 3 matrix with Pacific Orange and Alexa Fluor 750 (Invitrogen Corp.) at 0.03 and 0.04 μg/ml for low and 0.2 and 0.3 μg/ml for high staining, respectively. Incubation with barcoding dyes was performed at 4°C for 30 min. After several washes with FACS buffer, stimulated and barcoded cells were pooled into single tubes and stained for 30 min at 4°C, with an antibody cocktail containing anti‐pSTAT1 Alexa Fluor 488, pSTAT3 Alexa Fluor 647, pSTAT5 PE, CD3 Pacific Blue, CD4 PerCP‐Cy5.5, CD20 PerCP‐Cy5.5 and CD33 PE‐Cy7 (all from BD Phosflow). After washing, cells were resuspended in FACS buffer and acquisition was performed on an LSRII instrument (BD Biosciences). Data were collected using DIVA software. Data analysis was performed using FlowJo 8.8.6. by gating on live cells, then using double gating for singlet discrimination, followed by cell subset‐specific gating. Phosphorylation of STAT1, 3 and 5 proteins in B cells, CD4− or CD4+ CD3+ T cells, or monocytes was analyzed by deconvolution of stimuli‐specific gating. Baseline levels and fold increase between stimulated and unstimulated conditions were calculated using the 90th percentile fluorescence intensity of the pSTAT1, 3 or 5 signals.
Phosphorylation of Akt and PLC‐γ was assessed in B cells by crosslinking of the B‐cell receptor. After resting PBMC samples at 37°C (as conducted for cytokine stimulations), cells were distributed in V‐bottom 96‐well plates at 0.5 × 106 cells per well, and incubated for 4 min at 37°C in CO2 incubator with anti‐IgM and anti‐IgG both at 10 μg/ml (BD Biosciences) and 3% H2O2 for phosphatase inhibition. Cells were then fixed with 1.5% PFA for 10 min at RT. After washing twice with plain PBS, cells were permeabilized by 20 min incubation in 95% ice‐cold methanol. Cells were then washed with FACS buffer and stained with an antibody cocktail containing: CD3 Pacific Blue, CD20 PerCp Cy5.5, CD27 PE Cy7, PLGγ2 (BD Biosciences) and pAkt‐S473 Alexa Fluor 488 (Cell Signaling Technologies, Danvers MA). After 30 min incubation at 4°C, cells were washed in FACS buffer and analyzed by flow cytometry (as for cytokine stimulation). Median fluorescence intensity was recorded and used for the calculation of baseline levels of phosphorylated proteins and fold increase after BCR stimulation.
Serum cytokine levels determination
Cytokines were measured by a Luminex system (Luminex Corp., Austin, TX). The 50‐plex kits were purchased from Millipore and were used according to manufacturer's recommendations, with modifications as described below. Briefly, serum samples were mixed with antibody‐linked polystyrene beads on 96‐well filter plates and incubated at RT for 2 h followed by overnight incubation at 4°C. Plates were then vacuum filtered and washed twice before 2 h incubation with biotinylated detection antibody. Samples were filtered as above, washed twice and incubated with streptavidin‐PE for 40 min, then filtered and washed twice again before resuspending in reading buffer. Each sample was measured in duplicate. Plates were read using a Luminex LabMap200 instrument with a lower bound of 100 beads per sample per measured cytokine. The Luminex LabMap200 outputs the fluorescence intensity of each bead measured for a given cytokine in a sample. For each well, we considered the median fluorescence intensity (MFI) of all beads measured for a given cytokine and averaged the MFI of the two replicates. Values were normalized to a control sample ran in each of the plates.
The HAI assay was performed on sera from day 0 and day 28 using a standard technique (Webster et al, 2002); serially diluted 25‐μl aliquots of serum samples in PBS were mixed with 25‐μl aliquots of virus, corresponding to 4 HA units, in V‐bottom 96‐well plates (Nunc, Rochester, NY, USA) and incubated for 30 min at RT. At the end of the incubation, 50 μl of 0.5% chicken red blood cells were added and incubated for a minimum of 45 min before reading for HAI activity. The HAI titer of a given sample was defined as the reciprocal of the last serum dilution with no HA activity. A titer of five was assigned to all samples in which the first dilution (1:10) was negative. Seroconversion has been defined as an increase of at least fourfold in the antibody titer from pre‐ to post‐vaccination. Categories of vaccine PRs, if seroconversion was achieved for 0 or 1 strain, and vaccine GRs, if seroconversion was achieved for 2 or all 3 strains were utilized to account for both the strength and the breadth of the response.
Vaccination of CD95−/− mice and ELISA for detection of influenza‐specific antibodies
Mutant mice MRL/Mpj‐Faslpr/J (Mpj/lpr) and B6.MRL‐Faslpr/J (B6/lpr), and MRL/MpJ (Mpj) and C57BL/6L (B6) control mice were obtained from Jackson Laboratories (Bar Harbor, ME). Young mice (5–7 weeks of age) were housed at the Stanford Animal Facility until the protocol was completed. All mice were handled in accordance with APLAC and Stanford University animal care guidelines. Four male mice per group in the Mpj background, six male B6/Lpr and three B6 mice were given a single intramuscular dose (1/50 dilution, ∼1 μg HA) of the trivalent 2008 seasonal influenza vaccine Fluzone (Sanofi Pasteur). Samples were obtained from tail blood draws before and ∼4 weeks after vaccination. Serum was separated by centrifugation and stored at −80°C until collection of all samples. Vaccine‐induced influenza‐specific antibodies levels were determined by ELISA as follows. Ninety‐six‐well Vinyl Microtiter Microplates were coated with 2008 TIV (Fluzone) at 1 mg/ml in PBS. Plates were incubated overnight at 4°C, washed with wash buffer (0.1% of Tween 20 in PBS) and blocked with blocking buffer (3% BSA in PBS) for 1 h at 37°C. Serum samples were serially diluted with blocking buffer, added to the wells of coated/blocked plates and incubated for 1 h at 37°C. Wells incubated with blocking buffer only served as negative control. The plates were then washed, incubated for 1 h at 37°C with peroxidase‐conjugated goat anti‐IgG (KPL) diluted 1:2000 with blocking buffer. After washing, the plates were developed with TMB substrate (KPL) and the OD450nm of each well determined with an ELISA plate reader.
pFlow normalization scheme
We analyzed 89 samples across 13 different plates. Fold change due to stimulation was computed as the ratio of the cell, cytokine stimulation, phospho‐protein measure to the raw, unnormalized, cell phospho‐protein matching baseline that was measured on the same plate. Fold‐change values were then normalized by the median fold‐change difference of a given cell cytokine stimulation phospho‐protein measure within a given plate. We tested each assay for plate‐dependent differences and no significant differences between plates were detected post‐day normalization.
Gene module construction
Of a total of 48 771 gene probes in the microarray per sample, we first selected 6234 (s.d. cutoff value of 0.24) and subsequently normalized their expression by centering and scaling the expression so that each gene's expression across all subjects had euclidean norm equal to 1 for purposes of clustering. We utilized hierarchical agglomerative clustering with average linkage, euclidean distance and a height cutoff value of 1.5 to derive 109 modules. For each gene module, we assigned a set of regulatory genes (regulatory program), based on regression analysis of genes in the modules onto expression of transcription factors using a Akaike Information Criterion (AIC;Akaike, 1974). To do so, we used a set of candidate regulators composed of known signaling and transcription factors, 394 in total (Supplementary Table 4), of which 188 met the s.d. cutoff value. Briefly, we performed linear regression with elastic net penalty of each module's expression onto the set of regulators using LARS‐EN algorithm with l2 penalty weighted by 0.01. The LARS‐EN algorithm provides fits of increasing number of predictors. To select the best model among the outputs of LARS‐EN, we assessed quality of the resulting models by AIC, with sample specific terms weighted by within‐module variance. The fit with the best AIC score was selected for each module. Gene modules and their regulatory programs can be accessed at http://www.cs.unc.edu/~vjojic/fluy2/.
Cross validation and feature selection for prediction of antibody response to vaccination
An integral part of our training algorithm is a procedure for fitting a logistic regression model with l1 and l2 penalties, the elastic net penalty, a regularization algorithm that uses cyclical coordinate descent in a pathwise fashion as described (Friedman et al, 2007, 2010). The optimization cost can be stated as: where n is the number of donors in the sample, p is the number of predictors, xt denotes a vector of predictor values for subject t and yt is the observed outcome (PRs versus GRs based on the seroconversion to 0–1 or 2–3 strains, respectively). We assume all of our predictors are standardized to mean 0 and s.d. 1. The result of our fitting procedure is the set of predictor weights β and intercept α for the logistic regression model. In practice, penalty weights λ and κ in Eq. 1 are set by a data‐driven procedure, such as cross validation. The minimum λ was chosen to yield the highest AUC with the minimum set of features. For prediction of vaccine response using peptide predictors, we applied sevenfold cross validation with the aim of including at least 10 subjects per set (these experiments were conducted in 76 subjects only). For prediction of HAI response using blood measurements (gene module expression levels, serum proteins, pFlow conditions, immune cell subsets), all 89 subjects were included. Hence, eightfold cross validation was performed to include at least 10 subjects per test set.
Gene enrichment analysis
Genes from each module were explored by using Ingenuity Pathway Analysis for function enrichments. Data was imported and Core Analysis was performed with the following setting: Data Source: Ingenuity Expert Findings; Confidence: Experimentally Observed, TarBase, Protein–protein Interactions, Additional interactions; Species: Human. Most significant function enrichment for each module were explored in selected modules and were used for quality control. As an example, regression of each module by age and gender (Supplementary Table 3) identified a gene module (106) highly correlated with gender (r=0.714, P‐value=7.9 × 10−14) of which 17 genes (42% of total genes in module) are Y chromosome‐linked genes.
Function enrichment with highest significance, cell death and cell cycle, is reported from predictive modules APO and PROL, respectively. In addition, all genes in these modules and module regulators with highest regression coefficient or with known function were manually curated using a variety of sources including PubMed, IPA and BIOBASE Knowledge Library.
We are grateful to all our volunteers for agreeing to participate in this study. We also thank Yueh‐hsiu Chien, Rob Tibshirani, Jerry Friedman, Karen Sachs, Harry Greenberg, Janet McElhaney (University of British Columbia, Canada) and Fleur Tynan for illuminating discussions, and/or help and comments on the manuscript; Patty Lovelace and Garry Nolan for invaluable help with the pFlow assays; Yael Rosenberg‐Hasson and Iris Herschmann for their great help in cytokine level determination; and the Stanford‐LPCH Vaccine Program staff: Sally Mackey MS, Sue Swope RN, Cynthia Walsh RN, Kyrsten Spann, Thu Quan and Michele Ugur who enthusiastically enrolled subjects into the study and obtained samples from the participants. This work has been supported by grants from the Ellison Medical Foundation (AG‐SS‐1788), the Howard Hughes Medical Institute and the National Institutes of Health (U19 AI057229 and U19 AI090019) to MMD, and grants for the Stanford CTRU (NIH contract M01 RR00070). PJU is the recipient of a Donald E. and Delia B. Baxter Foundation Career Development Award, and is supported by the National Institutes of Health HHSN288201000034C, U19‐AI050864 and U19 AI090019; and FP Grant number 261, Euro Commission, European Union Seventh Framework Program (FP7/2007‐2013) under grant agreement number 261382. DF was supported by the Stanford Center on Longevity.
Author contributions: DF performed pFlow experiments; DF, TT and PL performed antibody titer studies; JP and JJ performed peptide microarray assays; DF and HH performed vaccination experiments in mice; HTM coordinated gene expression, luminex and phenotyping experiments; DF, VJ, BK and SS analyzed data; PJU and BK provided support and guidance in peptide microarray experiments and data analysis, respectively; CLD coordinated and organized the clinical study and contributed to study design and manuscript review; DF and MMD wrote the manuscript.
Conflict of Interest
The authors declare that they have no conflict of interest.
Supplementary Information [msb201315-sup-0001.pdf]
Data set 1
Metadata, exposure history and antibody titer pre‐ and post‐vaccine [msb201315-sup-0002.csv]
Supplementary Table 1 [msb201315-sup-0003.xlsx]
Supplementary Table 2 [msb201315-sup-0004.xlsx]
Supplementary Table 3 [msb201315-sup-0005.xlsx]
Supplementary Table 4 [msb201315-sup-0006.xlsx]
This is an open‐access article distributed under the terms of the Creative Commons Attribution License, which permits distribution, and reproduction in any medium, provided the original author and source are credited. This license does not permit commercial exploitation without specific permission.
- Copyright © 2013 EMBO and Macmillan Publishers Limited