Abstract

   Despite robust antibody responses, immunity induced by acellular
   pertussis vaccine (DTaP) wanes over time and risk of pertussis seems to
   be lower in children who receive whole-cell vaccine (DTP) as their
   first dose. To interrogate the early immunologic response to pertussis
   vaccine, we enrolled 56 healthy infants who received either DTP or DTaP
   at 2-, 4-, 6-, and 18-months of age. RNA-sequencing and ribosome
   profiling of PBMC were performed prior to vaccination (Day 1) and on
   either Day 2 or Day 8. Pathway enrichment analysis on Days 2 and 8
   showed enrichment of TLR-signaling and FcϒR-mediated phagocytosis among
   DTP recipients. DTP also led to increases in IRAK-4 and IL-1ß. After
   booster vaccination, a higher frequency of PT-specific B-cells was
   observed in DTP- vs. DTaP recipients. These data provide insights into
   the early immunologic responses to pertussis vaccine and may guide
   next-generation pertussis vaccine development.

   Subject terms: Adaptive immunity, Infectious diseases, Paediatric
   research, Bacterial infection, Translational research

Introduction

   Bordetella pertussis is a Gram-negative bacterium that colonizes the
   upper respiratory tract and is the causative agent of whooping cough,
   also known as pertussis. Current global data estimate that
   approximately 24 million cases and 160,000 deaths occur annually in
   children less than 5 years of age, making pertussis the tenth leading
   cause of death in the world^[54]1.

   Two types of pertussis vaccines are available to reduce the incidence
   and severity of pertussis disease, each of which are given with
   diphtheria and tetanus toxoids and often other antigens in combination
   vaccines. Whole-cell pertussis vaccine has been the mainstay of
   pertussis prevention since the late 1940s and consists of inactivated,
   whole B. pertussis organisms combined with diphtheria and tetanus
   toxoids (DTP). Acellular pertussis vaccines (aPV) consist of one or
   more purified pertussis antigens combined with diphtheria and tetanus
   toxoids (DTaP). Both vaccines require a three dose priming series
   during infancy with boosters delivered in older children and adults.
   Acellular vaccines have been associated with fewer adverse events than
   commonly seen after administration of DTP, including reduced rates of
   fever, febrile seizures, hypo-responsive episodes, local swelling, and
   injection site redness^[55]2. Large randomized clinical trials
   demonstrated that the efficacy of DTaP was comparable to DTP, and DTaP
   has largely replaced DTP in high income countries because of its
   improved reaction profile^[56]2. However, over time, waning of immunity
   conferred by DTaP has been observed, particularly in high income
   countries^[57]3. The World Health Organization (WHO) continues to
   recommend the preferential use of DTP over DTaP in low- and
   middle-income countries^[58]3–[59]6. More recently, data from nonhuman
   primates challenged with pertussis organisms suggested that DTaP
   reduces the incidence and severity of pertussis disease shortly after
   vaccination but does not prevent nasopharyngeal colonization or
   transmission to cage mates^[60]7. It has been postulated that the
   reduced protection afforded by DTaP may be due, in part, to its
   inability to generate a robust Th1/Th17 response, with the
   lipo-oligosaccharide and multiple additional antigens contained in the
   whole-cell vaccine stimulating a more protective, longer lasting immune
   response^[61]8–[62]10. Data from adults primed with either whole-cell
   or acellular pertussis vaccine and boosted with acellular vaccine show
   important differences in the expression of pro-inflammatory genes (such
   as IL6, IL1B, and NFKBIA) and IgG subsets that may reflect the initial
   priming vaccine^[63]11.

   To define the differences in immune responses to whole cell and
   acellular pertussis containing vaccines more precisely in infants, we
   used a systems-biology approach to compare changes in gene expression
   following primary vaccination in infants with either whole cell or
   acellular pertussis containing vaccines. RNA-Sequencing (RNA-Seq) and
   ribosome profiling (RP) were used to identify early transcriptional and
   translational signatures in peripheral blood mononuclear cells (PBMC).
   These signatures were correlated with pertussis-specific antibody
   responses, pertussis-specific B-cell responses, serum cytokines, and
   local and systemic adverse events seen after vaccination.

Results

   Fifty-six infants were enrolled in the study; 12 were enrolled in the
   pilot study and 44 in the main study. The median age at enrollment was
   9 weeks. Infants received pertussis vaccines as standard of care and
   none of the infants were born to mothers who had received maternal
   pertussis immunization. In the pilot study DTP group, 5/5 participants
   had samples available for RNA-Seq and RP; 3/5 were included in the
   RNA-Seq and RP analysis based on sample integrity; 3/5 received all 4
   doses of DTP-containing vaccine; and 2 were lost to follow-up. In the
   main study DTP group, 21 received first dose DTP and one did not have
   sufficient Day 2 blood for analysis. Therefore, 20 remaining
   participants were available for and included in RNA-Seq and RP
   analysis. One participant received only 3 doses of DTP and 19
   participants (95%) received all 4 doses of DTP-containing vaccine as
   planned. In the pilot study DTaP group, 7 participants were enrolled
   but 2 were excluded due to lack of blood sample on Day 2. Therefore 5
   participants were available for RNA-Seq and RP analysis, and 4
   participants were included based on sample integrity. Of those, all but
   1 participant received Hexaxim (3 doses) plus Infanrix booster (both
   DTaP-containing vaccines); one participant received Hexaxim for all
   four doses. In the main study DTaP group, 23 participants were
   enrolled; three were excluded to lack of sample on Day 2. Therefore, 20
   participants were available for RNA-Seq and RP analysis, but one was
   excluded due to sample integrity. All but one participant received
   Hexaxim (3 doses) with an Infanrix booster; one participant received 3
   doses of Hexaxim with a DTP booster and was excluded from subsequent
   analysis. Deviations from the protocol that occurred after collection
   of samples for RNA-Seq and RP analysis (which were collected only after
   the first dose) did not prevent participants from being included in the
   primary transcriptomic analysis; however, since antibody concentrations
   could be affected using a different product, these participants were
   excluded from the final integrated analysis dataset for timepoints
   after deviations from the typical vaccination schedule (Fig. [64]1).

Fig. 1. Study design and schema.

   [65]Fig. 1
   [66]Open in a new tab

   Panel A displays the overall study schema and study procedures; Panel B
   displays the CONSORT-style disposition of participants in the study and
   vaccines received. Participants for whom RNA-Seq or RP data did not
   pass quality control or were not available were excluded from the
   RNA-Seq/RP analysis. Created in BioRender. Creech, B. (2025)
   [67]https://BioRender.com/n33d820.

DTP generated higher antibody responses to pertussis toxin (PT) while DTaP
led to higher antibody responses to filamentous hemagglutinin (FHA)

   We began by assessing serologic responses to PT, FHA, FIM, and
   pertactin using an antigen-specific liquid-bead-based assay calibrated
   to WHO standards and expressed as International Units/mL (IU/mL)
   (Supplementary Fig. [68]1). Formal comparisons of immunogenicity are
   limited to PT and FHA as only these antigens are included in both
   vaccines. One month after completion of the three-dose primary
   vaccination series (Study Day 150), IgG antibody responses directed
   against pertussis toxin (PT) increased from pre-vaccination for both
   vaccine groups, though antibody concentrations were significantly
   higher in the DTP group compared to DTaP group (Geometric Mean
   International Units [GMIU] = 153.0 vs. GMIU = 66.0, p = 0.023). In
   contrast, DTaP recipients had significantly higher antibody
   concentrations to FHA at Day 150 (GMIU = 138.5 vs. GMIU = 17.2,
   p < 0.001. One month following the booster dose of vaccine (Day 510),
   the response to PT was higher in the DTP group than the DTaP group
   (GMIU 389.9 vs. GMIU 86.7, p < 0.001). As before, responses to FHA were
   significantly higher in the DTaP group (GMIU 311.9 vs. GMIU 67.3,
   p < 0.001). In addition, DTP recipients also mounted IgG responses to
   FIM 2/3 and pertactin, antigens that are not present in the bivalent
   DTaP vaccine formulation used in Peru (Table [69]1).

Table 1.

   Components of Hexaxim (DTaP) and DTP-containing vaccine administered to
   study population
   Component Content
   Hexaxim^a Bordetella pertussis antigens
    Pertussis toxoid 25 mcg
    Filamentous hemagglutinin 25 mcg
   Diphtheria toxoid ≥20 IU^b
   Tetanus toxoid ≥40 IU
   Poliovirus (inactivated)
    Type 1 (Mahoney) 40 D antigen units
    Type 2 (MEF-1) 8 D antigen units
    Type 3 (Saukett) 32 D antigen units
   Hepatitis B surface antigen 10 mcg
   Haemophilus influenzae type b polysaccharide conjugated to tetanus
   22–36 mcg
   Infanrix (GSK)^c Bordetella pertussis antigens
    Pertussis toxin 25 mcg
    Filamentous hemagglutinin 25 mcg
    Pertactin 8 mcg
   Diphtheria toxoid 25 Lf^d
   Tetanus toxoid 10 Lf^d
   Diphtheria, Tetanus, Pertussis, Hepatitis B, and Haemophilus influenzae
   type B conjugate vaccine^e
   Bordetella pertussis (whole-cell) ≥4 IU
   Diphtheria toxoid ≥30 IU
   Tetanus toxoid ≥40 IU
   Hepatitis B surface antigen ≥10 mcg
   Haemophilus influenzae type b polysaccharide conjugated to tetanus

   Polysaccharide: 10 mcg

   Tetanus toxoid protein: 19–33 mcg
   [70]Open in a new tab

   ^aAdsorbed onto aluminum hydroxide.

   ^bIU international units.

   ^cAdsorbed onto aluminum hydroxide.

   ^dLimits of flocculation units.

   ^eAdsorbed onto aluminum phosphate with thiomersal 0.005% preservative.

Among measured cytokines, IL-6 cytokine concentration was most increased
following both DTP and DTaP

   Next, we assessed serum cytokine concentrations for 21 cytokines on
   Days 1, 2 and 8 using a multiplex assay. Radar plots that contrast
   median fold changes in cytokine concentrations from baseline between
   groups are presented in Supplementary Fig. [71]2. Following DTP
   administration, there were no increases in serum cytokine
   concentrations on Day 2; however, on Day 8, DTP recipients had ≥2 fold
   increased concentrations of GM-CSF, IL-1β, IL-6 (>5-fold), IL-13,
   IL-21, and IL-23, compared to baseline. Following DTaP administration,
   we observed ≥2-fold increases in IL-2, IL-4, IL-5, and IL-6 on Day 2;
   on Day 8, DTaP recipients continued to have increased levels, compared
   to baseline, of IL-6. Given the modest sample size, these differences
   in post-vaccination median fold changes between DTP and DTaP recipients
   were not statistically significant.

PT-specific B-cell frequencies were higher in the DTP group

   The frequency of antigen-specific B-cells was determined one month
   after completion of the primary series (Day 150) and one month after
   the booster dose of vaccine (Day 510), Fig. [72]2. On Day 150,
   PT-specific B-cell responses were significantly higher in the DTP group
   compared to the DTaP group (0.112% vs. 0.045%, p = 0.011) and the
   differences increased by Day 510 (0.65% vs. 0.10%, p < 0.001). B-cell
   responses to FHA were similar after both vaccines at Day 150 (0.117%
   vs. 0.061%, p = 0.477), with both vaccine groups showing similar
   increases in FHA antigen-specific reactive B-cell frequencies by Day
   510 (0.34% vs. 0.23%, p = 0.374).

Fig. 2. Reactive B-cell frequencies for FHA, PT, FIM 2/3, and Pn.

   [73]Fig. 2
   [74]Open in a new tab

   Percent reactive B-cells one-month following the primary series (Day
   150) and booster immunization (Day 510). In orange: DTP; in purple:
   DTaP. Extreme outliers are removed from display, inclusive of 2 values
   for DTaP at Day 510 for FHA (3.5 and 5.6 percent reactive B cells), and
   1 value for DTP at Day 510 for FHA (2.3 percent reactive B cells).
   Extreme outliers were identified as having outlying profiles in one or
   all of the principal component analysis plots, multidimensional scaling
   plots, and Euclidean distance based clustering dendrograms. Boxplots
   are presented as the median and the 25–75th percentile (box) and
   whiskers (5–95th percentile).

PBMCs showed more pronounced transcriptomic and translational changes at
early time points following immunization

   To better understand the mechanisms underlying immunologic responses at
   the transcriptional and translational level, we carried out RNA-Seq and
   RP assessments on PBMCs collected pre-vaccination (Day 1) and Days 2
   (Supplementary Fig. [75]3) and 8 post-vaccination (Supplementary Fig.
   [76]4). Most differentially transcribed (DTSC) and differentially
   translated (DTRL) genes from baseline were identified on Day 2 for both
   vaccine groups, as demonstrated by highest values in the upper barplot
   portion of the UpSet plot where leftmost 5 bars show highest numbers of
   DTSC and DTRL genes exclusively for Day 2 timepoints (2171 unique DTSC
   genes overall, 644 unique for DTP alone and 551 unique for DTaP alone;
   1021 unique DTRL genes, 527 unique for DTP alone and 147 unique for
   DTaP alone), Supplementary Fig. [77]5. On Day 2, most DTSC genes were
   shared between DTP and DTaP, but an overall stronger fold change
   response was seen for DTP. This was more pronounced on the translation
   level with 3.6 times more unique DTRL genes for the DTP group compared
   to the DTaP group (527 vs. 147 unique DTRL genes, respectively and
   larger fold changes (1.05 and 1.00 LFC for DTP and DTaP for RP,
   respectively, compared to 0.89 and 0.81 LFC for DTP and DTaP,
   respectively).

Key innate immune signaling pathways were uniquely activated based on RP
analysis following DTP vaccination

   Next, we characterized differential translation on the pathway level
   using KEGG as the basis for pathway enrichment^[78]12,[79]13. When
   comparing KEGG pathway enrichment by RP between the two vaccine groups
   on Day 2, 25 pathways were shared between both, 16 were unique to DTP,
   and 3 were unique to DTaP [Supplementary Tables [80]1–[81]3]. Both the
   DTP and DTaP arms were enriched for NOD-like receptor signaling,
   phagosome, complement and coagulation cascade, neutrophil extracellular
   trap formation, B-cell receptor signaling, TNF signaling, and
   osteoclast differentiation-related pathways, demonstrating increased
   activation of these response pathways (Supplementary Fig. [82]6). In
   addition, several infectious disease pathways were activated for both
   groups, including the pertussis disease pathway. Pathways that were
   only enriched following DTP vaccination included TLR-signaling pathway,
   FcγR-mediated phagocytosis, and NK-cell mediated cytotoxicity. The top
   20 KEGG pathways identified across study groups and study days are
   shown in Fig. [83]3. In general, differences between DTP and DTaP were
   most pronounced on Day 2 (RNA-Seq).

Fig. 3. Top 20 KEGG pathways identified across all assays, arms, and days
shaded by enrichment score and color−coded by fold change direction.

   [84]Fig. 3
   [85]Open in a new tab

   Enrichment score is defined as −1 ×log[10] (FDR-adjusted p-value). Top
   pathways are defined based on the greatest change in enrichment score
   across assays, post-vaccination days, and treatment groups. Circle size
   is determined based on the mean of mean log[2] fold change relative to
   pre-vaccination for genes included in the pathway. RP Ribosomal
   Profiling. In red: increased compared to pre-vaccination, in blue:
   decreased compared to pre-vaccination. For each of the 9 gene sets and
   differentially transcribed or translated gene sets, p-values were
   adjusted using the Benjamini–Hochberg procedure to control the FDR.
   Gene sets with an FDR-adjusted p-value < 0.1 were considered to be
   significantly enriched.

   Next, we inspected differential responses between groups within
   enriched pathways. The pertussis KEGG pathway maps indicated that
   complement receptor 3 (CR3), the cellular receptor for FHA and the
   non-vaccine antigen adenylate cyclase toxin (ACT)^[86]14, was
   significantly increased in DTP recipients, but not DTaP recipients, on
   Day 2. The C1 complex, specifically C1Q, was differentially upregulated
   in DTP recipients vs. DTaP, as were genes encoding IL-1, IL-1β, and
   IRF1. Analysis of TLR-signaling pathway maps showed statistically
   significant activation of translation of TLR5 and TLR6 proteins that
   recognize flagellin and lipoprotein following DTP but not following
   DTaP. TLR4-mediated signaling, which is involved in the recognition of
   lipopolysaccharide (LPS), showed differential responses via the
   activation of MD2 (which binds to TLR4) in the DTP but not the DTaP
   group. In addition, the transcription of IRAK-4, a gene involved in
   mediating TLR4, TLR5, and TLR6 signaling within cells, was
   significantly increased from pre-vaccination in the DTP but not the
   DTaP group. Within both pathways, the translation of the
   pro-inflammatory cytokine IL-1ß was significantly increased while the
   RANTES chemokine, which is involved in migration and honing of effector
   and memory T-cells, was significantly reduced in the DTP group vs. the
   DTaP group.

DTP vaccination led to activation of gene expression signatures present in
mucosal-associated invariant (MAIT) T-cells

   To understand the distributions of cell populations within PBMCs based
   on cell-specific gene transcription signatures, we carried out
   exploratory analyses using cell deconvolution based on the
   transcriptomics data (Fig. [87]4). On Day 2, we observed an increase in
   monocytes and neutrophils following vaccination. Among DTP recipients,
   we observed increases, from baseline, in mucosal-associated invariant T
   (MAIT) cells on Day 8. At the same time, there was a decrease in the
   number of CD4 Memory T cells. This change in MAIT cells was not
   observed in DTaP recipients, though formal comparisons were not made.

Fig. 4. Cell deconvolution data.

   [88]Fig. 4
   [89]Open in a new tab

   Radar plot of the median of the fold change relative to Day 1,
   estimating the proportion of cell types in the PBMC fraction for DTaP
   and DTP on Days 2 and 8. In orange: DTP; in purple: DTaP; solid lines:
   Day 2; dashed lines: Day 8. MAIT Mucosal-associated invariant T cells,
   mDCs Myeloid dendritic cells, C Monocytes classical monocytes, NC
   Monocytes non-classical monocytes, LD Neutrophils low-density
   neutrophils, NK natural killer, pDCs plasmacytoid dendritic cells.

Correlation between gene expression, translation, and antibody responses to
pertussis toxin

   On Day 8, translation for IGLV10-54, an immunoglobulin-encoding gene,
   was increased 3.7-fold in the DTP group relative to the DTaP group. The
   overall translational efficiency for IGLV10-54 was also 3.3-fold higher
   in the DTP group compared to the DTP group (DTP vs. DTaP, RP vs.
   RNA-Seq), suggesting that the gene was more highly translated vs.
   transcribed in the DTP vs. the DTaP group. The HP gene, which encodes
   haptoglobin, showed the strongest differential signal between groups
   with 9.6-fold higher translation in the DTP vs. DTaP group on Day 2.
   The overall relative differential translational efficiency for this
   protein was 2.6-fold (DTP vs. DTaP, RP vs. RNA-Seq), suggesting that
   the gene was more highly translated vs. transcribed in the DTP group
   when compared to the DTaP group. Boxplots of participant-specific
   values for gene translation and translational efficiency for IGLV10-54
   and HP are shown in Supplementary Fig. [90]7.

   To assess which genes on Day 2 or 8 were correlated with Day 150
   antibody responses to PT and FHA, we performed regularized linear
   regression analysis in combination with bootstrapping to identify
   robust gene predictions. The predictor gene variable set was based on
   log2 fold change in LCPM (log counts per million) and included genes
   with an average absolute baseline fold change of ≥1.5 in either
   treatment group at either of the two post-vaccination visits (Day 2 or
   8). The analysis was performed for each combination of assay (RP and
   RNA-Seq), treatment group, day, and antibody response measurement (PT
   and FHA). Overall, 15 genes were identified robustly, as being genes
   identified for most bootstrap iterations for at least one assay, day,
   and treatment group combination. A listing of overlap in all genes
   whose responses at treatment arms and post-vaccination days are
   correlated with FHA and/or PT values at Day 150 and have at least 50%
   selection among bootstrap fitted model results is presented in
   Supplementary Table [91]4.

   Supplementary Fig. [92]8 displays genes of interest measured via
   RNA-Seq that correlated with PT antibody responses. Changes in CXCL10
   gene expression on Day 2 correlated with PT antibody response following
   both DTaP and DTP vaccination. CXCL10 is activated following interferon
   signaling and is a member of the KEGG TLR-signaling pathway. Its gene
   product, IP-10, is a potent immune cell chemoattractant. In both
   vaccine groups, CXCL10 gene expression was increased from
   pre-vaccination on Day 2 (1.49-fold for DTP, 1.98-fold for DTaP)
   (Supplementary Table [93]5, Supplementary Table [94]6). While for DTP
   this increase in CXCL10 expression on Day 2 was associated with an
   increased antibody responses against pertussis toxin (r = 0.40), the
   inverse was observed for DTaP, in which increased CXCL10 expression was
   negatively correlated (r = −0.52). In addition, Day 2 and Day 8
   expression of two long non-coding RNAs were predictive of Day 150
   antibody levels to PT following DTP vaccination. These genes included
   NRIR, a long non-coding RNA that acts as a negative regulator of
   interferon response; expression of NRIR was increased on Day 2
   post-vaccination (1.61-fold) and Day 8 (1.15-fold) (Supplementary Table
   [95]5, Supplementary Table [96]7) and this increase was associated with
   an increase in antibody responses against PT (r = 0.77 and r = 0.40,
   respectively). Secondly, ENSG00000259926 (lnc-CHD9-1:1), a long
   non-coding RNA novel transcript, was increased on Day 2 (4.56-fold) and
   Day 8 (2.61-fold) (Supplementary Table [97]5, Supplementary Table
   [98]7) and this increase was negatively associated with later antibody
   responses against PT (r = −0.73 and r = 0.16, respectively). These
   findings imply different regulatory functions of these long non-coding
   RNA genes.

Reactogenicity correlation analysis

   Supplementary Table [99]8 displays the frequency of participants
   experiencing any solicited adverse events (reactogenicity) following
   dose 1 of DTaP or DTP-containing vaccine. Though not statistically
   significant, the frequency of injection site erythema and swelling was
   slightly more common in the DTP group. Fever occurred in 5/26 (19.2%)
   DTP recipients and 5/30 (16.7%) DTaP recipients; fussiness occurred in
   26/26 (100%) DTP recipients and 17/30 (56.7%) DTaP recipients. Gene
   transcription fold changes for 12 genes on Day 2 were associated with
   the development of injection site redness. The strongest positive
   correlations were seen with TRAV14DV4 (encoding a T-cell receptor
   protein), NBEA (neurobeachin, involved with membrane protein
   trafficking), and HP (haptoglobin, an acute phase protein) gene
   expression. The strongest negative correlations were observed with
   transcription of EGR1 (involved in cytokine signaling), OTOF
   (otoferlin) and LGR6, each associated with reduced likelihood of
   redness. Boxplots summarizing gene expression for genes whose responses
   best predict injection site erythema at Day 2 as measured via RNA-Seq
   are displayed in Supplementary Fig. [100]9.

Discussion

   In this study, we used RNA-Seq and RP, in combination with serologic
   and cellular assays, to characterize the immune response to DTP and
   DTaP vaccination more completely in infants. Genes responsible for
   TLR-signaling, innate immune responses, and energy metabolism were most
   highly upregulated in DTP recipients, as was antibody mediated
   phagocytosis. These innate immune responses were most apparent on Day 2
   but persisted through Day 8. DTP was also associated with more robust
   PT antibody responses and PT-specific B-cell responses after the
   primary series and booster dose of vaccine.

   Previous work has demonstrated that the type of vaccine administered
   during infancy is clinically and immunologically important, likely due
   to antigenic imprinting^[101]8,[102]15–[103]17. Though
   pertussis-specific antibody concentrations generated after DTaP booster
   are generally considered to be similar to those seen after DTP^[104]15,
   previous data indicate that individuals primed with DTaP demonstrate
   decreased IFNγ and IL-17 production, decreased expansion of memory
   cells, and decreased T-cell proliferation in vitro in response to DTaP
   boost compared to DTP primed individuals. The fact that we found that
   DTP priming resulted in greater translation of genes involved in innate
   immune signaling pathways in PBMCs soon after vaccination suggests that
   DTP leads to improved antigen presentation and processing and
   immunologic crosstalk between effector cells. Da Silva Antunes and
   others have shown that adults primed with DTP vs. DTaP differ
   substantially in T-cell phenotype and polarization, with changes
   persisting for decades after priming. Upon boosting, acellular
   pertussis vaccine in adults results in increased cell activation;
   upregulation of TLR4; induction of CCL3, CXCL8, and C-type lectins; and
   increased frequency of antigen secreting cells^[105]11,[106]18. In
   contrast, individuals primed with DTP retain the capacity for a
   Th1/Th17 response, IFNγ production, and higher T-cell proliferation
   after adult acellular pertussis vaccine
   boosting^[107]11,[108]15,[109]18. A direct comparison between the
   infant data presented in this study and data from adults is challenging
   given that our ‘omics dataset characterizes only the initial response
   to vaccination rather than characterizing the response to later doses
   or to adult boosting in the context of DTP or DTaP priming.

   It is important to note that pertussis vaccine efficacy likely does not
   result merely from serum antibody (indeed there is no antibody-based
   correlate of protection) but from a combination of innate, cellular,
   and mucosal immunity^[110]19,[111]20. Kapil and Merkel have reviewed
   this topic^[112]8 and suggest that priming with acellular pertussis
   vaccine leads to upregulation of IL-1, IL-4, IL-5, IL-9 and TGF-ß
   (inhibitor of T-cell proliferation) within the CD4 T cell population
   while whole-cell pertussis vaccine leads to generation of IL-1ß, IFNγ
   and IL-17^[113]15,[114]21. Analysis of circulating serum cytokines in
   our study did not reveal statistically significant differences, but we
   did observe a 5-fold increase in IL-6 among DTP recipients on Day 8,
   compared to 3.5-fold increase among DTaP recipients, a finding
   consistent with infection models of pertussis in nonhuman
   primates^[115]22. IL-1ß was increased 2.6-fold among DTP recipients
   compared to a decrease among DTaP recipients on Day 8 and IFNγ was
   increased nearly 3-fold in DTP recipients on Day 8. In contrast, IL-4,
   a Th2 cytokine, was 2.5-fold higher than baseline on Day 2 among DTaP
   recipients, compared to a reduction from baseline among DTP recipients.

   IL-1ß, also known as lymphocyte activation factor, is a
   pro-inflammatory cytokine that induces expression of IL-17 by γδ-T
   cells^[116]23. We observed increased translation of IL-1ß, as well as
   increased abundance of the cytokine in the serum following DTP
   vaccination. Place et al. have shown that IL-1ß is required for
   clearance of Bordetella pertussis in a murine model of
   pertussis^[117]24; our study extends these findings to suggest that
   increased IL-1ß was likely mediated through translation of several
   genes within the TLR and NOD-like receptor signaling pathways. Our
   results support the hypothesis that DTP-specific activation of
   TLR-signaling through TLR4 (via MD2), TLR2, and TLR5 receptors, as
   quickly as 24 h after vaccination, contributes to both increased
   reactogenicity (through the elaboration of pro-inflammatory cytokines)
   and more durable antibody responses in an IL-1ß-mediated manner. To
   test this hypothesis, TLR4, TLR5, TLR6 receptor agonists (such as
   monophosphoryl lipid A (MPLA), flagellin, and lipoprotein,
   respectively) could be evaluated as potential vaccine adjuvants for
   DTaP to assess whether vaccine performance of the acellular vaccine
   could be improved^[118]25–[119]27. In addition to the TLR-signaling
   pathway, antibody mediated phagocytosis (FcγR-mediated) was also
   enriched in DTP recipients, which may also lead to improved antibody
   and cellular responses. Our deconvolution data analysis also suggests
   that cell populations may be altered post-immunization distinctly after
   DTP vs. DTaP, though additional prospective evaluations of cell
   population dynamics should be pursued given the nature of this
   analysis.

   Haptoglobin translation was approximately 10-fold higher in the DTP
   group at Day 2; moreover, the translational efficiency for this protein
   was 2.6-fold higher in DTP vs. DTaP. Thus, not only was the HP gene
   transcribed in more abundance, but the transcript was translated to
   protein more rapidly in the DTP group based on ribosome profiling
   analysis. In plasma, haptoglobin binds to free hemoglobin, which is
   released into the plasma by damaged red blood cells. Haptoglobin
   binding to hemoglobin allows degradative enzymes to access hemoglobin
   while preventing iron loss. More relevant to infectious diseases,
   haptoglobin sequesters iron within hemoglobin, preventing
   iron-utilizing bacteria from benefitting from hemolysis; as such,
   haptoglobin is considered an acute phase reactant, likely protecting
   the host from bacterial growth early in the course of disease. In
   addition, the gene cluster CRSDTaPD2D8-068, which contains 6 genes, was
   the most differentially transcribed gene cluster on Day 2 based on
   RNA-Seq; this cluster includes the inflammatory mediators
   S100A8/A9^[120]28. Though haptoglobin, S100A8, and S100A9 are
   differentially transcribed and translated in DTP recipients, they may
   also represent non-specific, general inflammatory responses, as many
   inflammatory processes may increase levels of these proteins. Of note,
   haptoglobin (Hp) has high affinity for PT, which, along with FHA,
   serves as a hemagglutinin^[121]29. This may be relevant to vaccine
   development as a current, live attenuated pertussis vaccine candidate,
   BPZE1^[122]30,[123]31, contains mutations in the enzymatically active
   moiety, S1, that abrogate haptoglobin binding in vitro^[124]32.

   Our findings should be considered in the context of certain
   limitations. Samples for RNA-Seq and RP were obtained one day or one
   week after vaccination; as a result, there are likely biologically
   important changes occurring in between these timepoints that were not
   assessed. Future studies should consider additional timepoints,
   particularly Day 4 (~72 h post-immunization). Variability in DTP and
   DTaP vaccine manufacturing may lead to variability in immunogenicity
   and reactogenicity; therefore, differences seen between DTP and DTaP
   may not be generalizable to all pertussis-containing vaccines. This is
   particularly important given that the DTaP preparation (Hexaxim)
   contains inactivated polio vaccine (IPV) while the DTP preparation
   (pentavalent vaccine) does not; thus, IPV was given as a separate
   injection at 2 and 4 months of age. While the sample size is robust for
   transcriptomic and ribosome profiling analyses, integration with
   antibody concentrations and cytokine measurements was more challenging
   given the limited number of participants and the small volume blood
   samples that were obtained. Moreover, the use of serum cytokine
   measurements (rather than plasma) may have modified the cytokine
   profiles that we identified, though this was necessary due to
   optimization of RNA-Seq and RP workflows. Due to the use of different
   DTaP formulations for booster doses during the second year of life, we
   could only compare ELISA antibody responses against PT and FHA and not
   Fim 2/3 and Pn antigens between the two vaccine groups. Moreover, the
   bivalent acellular vaccine used in Peru precludes a direct comparison
   to the response that may occur after quadrivalent vaccines that are
   used in the US. Nevertheless, we identified differential gene
   translation signatures that correlated with differential antibody
   levels against PT, corroborating the importance of molecular processes
   that were activated to shape later humoral immune responses. Future
   work will focus on the applicability of modulating specific innate
   immune pathways such as TLR signaling via novel adjuvants or alternate
   routes of vaccine administration to improve the immunogenicity,
   durability, and protection of pertussis vaccines.

   With pertussis remaining a leading cause of childhood mortality
   globally, updated strategies are needed what will either improve the
   immunogenicity of DTaP or reduce the reactogenicity of DTP.
   Understanding the immunophysiology of these vaccines, as well as
   vaccines in development, is critical to the control of pertussis.

Methods

Study design

   The study was conducted in Lima, Peru, where both DTP and DTaP vaccines
   are available and there is no preferential recommendation for their
   use. The study was conducted in accordance with the study protocol,
   International Council on Harmonization Good Clinical Practices
   standards, the Declaration of Helsinki, and applicable regulatory
   requirements. Both the Vanderbilt University Medical Center (VUMC)
   Human Research Protection Program and the Ethics Review Board of the
   Instituto de Investigación Nutricional (IIN) approved the study
   protocol and amendments.

Participants

   Following written informed consent by a parent or legal guardian,
   healthy male and female infants ≥2 months of age were enrolled. A pilot
   phase of the study (n = 12 participants) was conducted first to
   optimize workflows and establish sufficiency of blood volume for
   subsequent analysis. Pilot study participants were included only in
   analyses from Day 150 (pertussis-protein reactive B-cell frequency and
   pertussis serology), Day 480 (pertussis serology), and Day 510
   (reactive-B cells and pertussis serology) and did not contribute to
   correlation analyses involving RNA-Seq, RP, or cytokine data. The pilot
   study was followed by the main study phase (n = 44 participants).
   Infants from whom blood could not be obtained on Days 2 or 8 were
   excluded and replaced. Participants included infants receiving care in
   either the Peruvian Ministry of Health system who received DTP by
   standard of care or infants receiving care at private pediatric clinics
   who received DTaP as routine vaccination. Infants were eligible to
   participate if they were ≥37 weeks’ gestation at birth, 50–89 days of
   age at the time of enrollment, and considered to be in good general
   health by medical history and physical examination. Participants were
   excluded if they were of low birth weight (<2.5 kg); anemic (hemoglobin
   <7 mg/dL); had evidence of a primary immunodeficiency or HIV; were
   receiving systemic corticosteroids at enrollment or had received
   corticosteroids for more than 14 days since birth; or had a history of
   receiving any blood products since birth. Children with acute illnesses
   remained eligible, though enrollment was delayed until all symptoms
   were resolved for >48 h. Participants who received a vaccine that
   differed from the expected sequence, or participants for whom data were
   not available due to loss to follow-up/parent withdrawal, were excluded
   from analyses after the relevant event. Maternal pertussis immunization
   was not available during this time for pregnant individuals in Peru.

Study procedures

   The overall study scheme and consort diagram are illustrated in Fig.
   [125]1. Briefly, 12 evaluable participants were enrolled during the
   pilot phase of the study. Following successful venipuncture, 7 infants
   received a hexavalent vaccine containing DTaP (Hexaxim, Sanofi Pasteur,
   containing PT/FHA, diphtheria, tetanus, hepatitis B, HiB, and
   inactivated poliovirus) and 5 infants received a pentavalent vaccine
   containing DTP (Serum Institute of India, containing whole-cell
   pertussis, diphtheria, tetanus, hepatitis B, and Hib), Table [126]1.
   Infants returned to the clinic on Day 2, approximately 24 h after
   vaccination, for physical examination and repeat venipuncture. Blood
   samples from Days 1 and 2 in the pilot study were used for RNA-Seq, RP,
   and cytokine measurements. Infants subsequently received doses 2 and 3
   of DTP or DTaP per standard of care at 4 and 6 months of age and
   returned on Day 150 (approximately 1 month following completion of the
   primary series) for measurement of pertussis antibodies and enumeration
   of pertussis-specific B-cells. Where possible, the same vaccine
   formulation was administered for each vaccine dose. At the time of the
   pertussis vaccine booster (Day 480), pre-immunization blood was
   obtained for serum pertussis antibodies and on Day 510, approximately 1
   month following the booster dose, blood was obtained for measurement of
   serum pertussis antibodies and enumeration of pertussis-specific
   reactive B-cells^[127]33.

   In the main study phase, infants were randomly allocated at enrollment
   to either a Day 2 or Day 8 blood sample to ascertain early and late
   transcriptomic changes in response to vaccine. As with the pilot study,
   infants in whom paired blood samples pre- and post-vaccination could
   not be obtained were withdrawn from the study and replaced. The
   remainder of the study schedule was unchanged during the main study
   phase from the pilot study.

   All participants were provided a thermometer and memory aid to measure
   solicited adverse events (reactogenicity) for the 7 days following the
   first dose of vaccine. Unsolicited adverse events, deemed by the
   investigator to be due to vaccine, were also recorded throughout the
   study.

Sample preparation, RNA-Seq, and RP assays

   Detailed descriptions of laboratory procedures are described
   elsewhere^[128]34 accompanying manuscript, by Leguia et al.). Briefly,
   whole blood was collected on Day 1 (pre-vaccination) and on either Day
   2 or Day 8 and immediately processed to maintain fidelity of
   transcripts for sequencing. In the pilot study, plasma was separated
   first, and peripheral blood mononuclear cells (PBMC) were recovered
   using SepMate and Lymphoprep (STEMCELL Technologies) according to
   manufacturer’s instructions. Due to insufficient recovery of ribosome
   footprints for RP in the pilot phase, we modified the approach in the
   main study and added 100 ug/mL cycloheximide (CHX), which freezes
   translating ribosomes on their in vivo positions on mRNA, to blood
   collection tubes prior to obtaining blood samples. As a result, only
   sera were available for cytokine analyses in the main study. PBMCs were
   frozen after separation. At the time of analysis, cells were removed
   from liquid nitrogen, thawed, and resuspended in lysis buffer; 80% of
   the cell lysate of each sample was used to prepare ribosome profile
   libraries, as described elsewhere^[129]34 accompanying manuscript, by
   Leguia et al., and 20% of the lysate was used for RNA-Seq, using
   standard methods^[130]35.

Serology

   The humoral immune response to vaccine was assessed by measuring serum
   IgG antibody responses by liquid bead-based assay to pertussis toxin
   (PT), pertactin (PRN), fimbriae 2/3 (FIM, List Biologicals, Campbell,
   CA), and filamentous hemagglutinin (FHA, Enzo Biochem, Framingdale,
   NY)^[131]36. Antigens were individually and covalently conjugated to
   carboxylate modified Luminex beads, each having a distinct spectral
   signature. The reportable value of the assay is the mean of two tests
   and expressed as the serum concentration of antigen-specific IgG in
   IU/mL; the assay is calibrated to the WHO International Standard
   Pertussis Antiserum WHO 06/140.

   Pre- and post-vaccination differences in antibody concentrations for
   each pertussis antigen and each time point were evaluated for each
   group using a two-sided Welch’s t test on the log scale adjusting for
   unequal variance.

Cytokine assay

   A commercially available 21-plex inflammatory cytokine panel
   (Milliplex, Millipore) was used to measure pro-inflammatory cytokines
   from plasma (pilot study) and sera (main study). Briefly, 25 µl of neat
   sera was added to 25 µl of assay buffer in a prepared flat bottom
   96-well plate with standard and controls. All wells had 25 µl of the
   bead mixture added and incubated overnight at 4 °C. After washing,
   50 µl of detection reagent was added for 1 h followed by 50 µl of
   streptavidin-phycoerythrin for 30 min. The plate was then acquired on a
   Luminex FlexMap 3D. Inflammatory cytokine concentrations, as well as
   their post-vaccination fold changes, were summarized using the median,
   95% CI of the median, first and third quartile, minimum and maximum.
   The 95% CI of the median concentration and median fold change from
   pre-vaccination were determined using the bootstrap method (1000
   replications) and were visualized using time trend plots. To identify
   inflammatory cytokines that showed a differential response from
   pre-vaccination (Day 1), a Wilcoxon signed-rank test was carried out
   for each cytokine, vaccine group, and post-vaccination time point (Days
   2 and 8).

Reactive B-cell frequency assay

   Methods for the enumeration of antigen-specific B cells have been
   published previously^[132]37 and includes EBV-transformation of cells,
   formation of lymphoblastoid cell lines (LCL), screening for PT
   antibody, and calculation of the number of LCL’s in each well. Briefly,
   previously frozen PBMC were thawed rapidly in a 37 °C water bath and
   washed with a warm medium. Each well was seeded with 50 μL of PBMC
   resuspended in a medium with B95.8 cell supernatant at a concentration
   of 0.5 ×10^6 cells/mL and CpG oligonucleotides to stimulate
   transformation. Plates were incubated at 37 °C and checked daily for
   the formation of LCLs. After 6–12 days, cell supernatants were
   harvested and screened for antibodies to PT in a plate-based ELISA,
   resulting in a positive, indeterminate, or negative result.
   Indeterminate results occur due to a lack of cell viability or failure
   to form LCLs. The RBF was calculated as the number of positive ELISA
   wells divided by the number of LCLs. Significant differences between
   the pertussis-specific B-cell responses were evaluated using a Wilcoxon
   rank-sum test.

RNA-Seq and ribosome profiling analysis

   Raw sequence reads determined using RNA-Seq and RP were pre-processed
   by removing adapters, low-quality reads, and reads mapping to known
   rRNA and tRNA genes. For each sample, sequence reads were mapped to the
   latest human reference transcriptome/genome using HISAT2 (Version
   2.2.0). Gene expression quantification was carried out using Subread
   (Version 2.0.1) against the latest Ensembl gene model reference
   collection (Version 100; April 2020). Systematic differences in
   sequence coverage between samples were normalized using the TMM method
   as implemented in the edgeR R package (Version 3.18.1). Post
   normalization, genes located on the X or Y chromosomes were excluded to
   avoid sex-specific effects. Moderated log counts per million (LCPM) for
   the remaining genes were determined using edgeR. LCPM was used to
   identify technical bias (e.g., batch effects or GC-content bias) and
   outlying samples (e.g., sample labeling error). More specifically,
   principal component analysis, non-metric multidimensional scaling, and
   correlation heatmaps based on LCPM were used to identify potential
   global outliers and systematic effects. Samples with outlying profiles
   based on both principal component analysis and non-metric
   multidimensional scaling (4 samples for RNA-Seq, 1 sample for RP) were
   considered global outliers. Reverse cumulative distribution function
   plots of average LCPM across samples were used to identify a suitable
   minimum average LCPM filter cut off to exclude lowly expressed genes.

   Negative binomial models as implemented in edgeR were used to determine
   DE genes for each assay and vaccine group after removal of lowly
   expressed genes and outliers under the assumption that counts are
   negative binomial distributed. Models included fixed effects for
   subject, study visit day (pre-vaccination or post-vaccination day). The
   statistical significance of the study day effect (post -
   pre-vaccination) was evaluated for each gene using a likelihood ratio
   test. To compensate for testing multiple genes, p-values were adjusted
   by calculating false-discovery rates (FDR) based on the
   Benjamini–Hochberg procedure. Genes with an FDR-value < 0.05 and a fold
   change of ≥1.5-fold (up or down compared to pre-vaccination) were
   deemed to be significant DE genes. The following number of participants
   were included in each DE gene evaluation: n = 10 for RNA-Seq, Day 2,
   DTP, n = 9 for RNA-Seq, Day 2, DTaP, n = 9 for RNA-Seq, Day 8, DTP,
   n = 7 for RNA-Seq, Day 8, DTaP, n = 9 for RP, Day 2, DTP, n = 9 for RP,
   Day 2, DTaP, n = 7 for RP, Day 8, DTP, n = 7 for RP, Day 8, DTaP.
   Results for each DE gene including FDR-adjusted p-values are detailed
   in Supplementary Tables [133]9–[134]15.

   Pathway enrichment analysis was carried out separately for each
   differentially transcribed or translated gene set identified using
   published gene set information obtained from KEGG (KEGG Pathways
   Version 98.0; April 2021 and KEGG Modules Version 98.0; April 2021).
   The analysis was conducted using the goseq R package (Version 1.28)
   which adjusts for gene length bias when modeling the null distribution.
   To compensate for testing multiple gene sets, p-values were adjusted
   using the Benjamini–Hochberg procedure to control the FDR. Gene sets
   with an FDR-adjusted p-value < 0.1 were considered to be significantly
   enriched.

   To estimate changes in cell type proportions that were present in bulk
   PBMC RNA-seq data, in silico deconvolution was performed using the R
   package granulator (Version 1.6.0). Deconvolution algorithms used for
   this analysis were based on the following methods: ordinary least
   squares, non-negative least squares, quadratic programming with
   non-negativity and sum-to-one constraint, quadratic programming without
   constraints, re-weighted least squares, support vector regression, and
   linear mixing model, as implemented via the granulator R
   package^[135]38. Cell proportions were then estimated by combining the
   results from each deconvolution method using the median.

   Regularized linear regression analysis in combination with
   bootstrapping was performed to gene expression fold change responses
   that predicted antibody response measurements using the glmnet R
   package (Version 2.0-13). Separate models were fit for each day and
   assay. To avoid overfitting (n«p and collinearity among gene responses)
   and to facilitate variable selection, an elastic net regularization
   step (combination of L1 Lasso and L2 ridge penalization, a = 0.5) was
   included as part of the fitting procedure. Cross validation was used to
   determine the optimum regularization parameter. The predictor gene
   variable set was based on log[2] fold change in LCPM and included genes
   with an average absolute baseline fold change of ≥1.5 in either
   treatment group at either of the two post-vaccination visits (Day 2 or
   8). The analysis was performed for each combination of assay (RP and
   RNA-Seq), treatment group, day, and antibody response measurement (PT
   and FHA).

   Regularized logistic regression models were fit to determine gene
   expression fold change responses that predicted reactogenicity group
   membership. Analysis was conducted as described for regularized linear
   regression analysis aside from (1) the absence of bootstrapping and (2)
   the predictor gene variable set included genes with an average absolute
   baseline fold change of ≥1.5 in either reactogenicity group at either
   of the two post-vaccination visits (Day 2 or 8).

Supplementary information

   [136]Supplementary Information^ (12MB, pdf)

Acknowledgements