Population Diversity among Bordetella pertussis Isolates, United States, 1935–2009

Resurgence of pertussis was not directly correlated with changes in vaccine composition or schedule.

Since the 1980s, pertussis notifi cations in the United States have been increasing. To determine the types of Bordetella pertussis responsible for these increases, we divided 661 B. pertussis isolates collected in the United States during 1935-2009 into 8 periods related to the introduction of novel vaccines or changes in vaccination schedule. B. pertussis diversity was highest from 1970-1990 (94%) but declined to ≈70% after 1991 and has remained constant. During 2006During -2009.6% of the strains encoded multilocus sequence type prn2-ptxP3-ptxS1A-fi m3B, and 64% were multilocus variable number tandem repeat analysis type 27. US trends were consistent with those seen internationally; emergence and predominance of the fi m3B allele was the only molecular characteristic associated with the increase in pertussis notifi cations. Changes in the vaccine composition and schedule were not the direct selection pressures that resulted in the allele changes present in the current B. pertussis population. P ertussis, or whooping cough, is caused by the bacterium Bordetella pertussis and is the most frequently reported bacterial vaccine-preventable disease in the United States (1). Vaccination against pertussis began in the 1940s in the United States, using a whole-cell formulation (wP) that resulted in a dramatic decrease in infections and deaths (2). Acellular pertussis vaccines (aP) were licensed for the fourth and fi fth doses of the childhood booster series in 1991 and were recommended for all 5 doses of the childhood series by 1997; in 2005, a single-dose adolescent and adult booster (tetanus-diphtheria-aP, or Tdap) was recommended ( Figure 1). Despite a successful US childhood vaccination program with high coverage, the number of re-ported pertussis cases has increased since the early 1980s, with 27,550 cases reported in 2010 (3).
Before the current study, US B. pertussis isolates from 1935-1999 were characterized by pulsed-fi eld gel electrophoresis (4), and a subset of isolates was analyzed for 2 genes, prn and ptxS1 (5). Genetically, the B. pertussis population was largely homogeneous during this period, and only a few strain types caused most disease in the United States (4). Recently, the molecular typing methods multilocus variable number tandem repeat analysis (MLVA) and multilocus sequence typing (MLST) have been used to assess B. pertussis population trends in other countries (6)(7)(8)(9); used together, these methods offer discriminatory power similar to that of pulsed-fi eld gel electrophoresis (9,10). We used MLVA and MLST to type a large selection of B. pertussis isolates from the United States and examined a selection of molecular changes that occurred over time and how these changes related to increases in pertussis notifi cations or changes in vaccine policy.

Strain Selection
We selected 661 B. pertussis isolates of US origin from the Centers for Disease Control and Prevention (CDC) collection by using random sampling stratifi ed by geography (US states and territories) and period. The strains were divided in advance as follows: period 1 (prevaccine era), 1935-1945, n = 3; period 2 (early wP era), 1946-1969, n = 16; period 3 (late wP era), 1970-1990, n = 76; period 4 (aP transition for 4th and 5th dose of childhood series), 1991-1996, n = 86; period 5 (early aP), 1997-1999, n = 159; period 6 (middle aP), 2000-2002, n = 98; period 7 (late aP), 2003-2005, n = 98; and period 8 (early Tdap booster), 2006-2009, n = 125 ( Figure 1). Stratifi cation was used to ensure that all states and territories with isolates in the strain bank were represented in the random sample. The geographic distribution of strains and information regarding location and year of isolation are provided in the online Technical Appendix (wwwnc.cdc.gov/EID/eid-static/ spreadsheets/12-0082-Techapp.xlsx). We could not correct for the lack of representativeness of the isolates in the collection because CDC does not receive an isolate for every report of illness in the United States. After 3 days of incubation at 37°C, DNA was extracted by heat-lysis preparation from each isolate and stored at -20°C until ready for use in PCR.

MLVA
Analysis was performed by using a 6-target multiplex similar to that described (8) with some modifi cations. Using the HotStarTaq kit (QIAGEN, Valencia, CA, USA) yielded a fi nal reaction volume of 20 μL. Master mix 1 consisted of fl uorescently labeled oligonucleotides for variable number tandem repeats (VNTRs) 1 (0.13 μmol/L each primer), 5 (0.09 μmol/L each primer), and 6 (0.09 μmol/L each primer) and was supplemented with 1 mol/L betaine (Sigma-Aldrich, St. Louis, MO, USA) to facilitate primer-template interaction. Master mix 2 consisted of fl uorescently labeled oligonucleotides for VNTRs 2 (0.08 μmol/L each primer), 3 (0.23 μmol/L each primer), and 4 (0.08 μmol/L each primer) and was supplemented with 10% dimethyl sulfoxide (Sigma-Aldrich). PCR was performed in single-target reactions (6) for some targets that were not effi ciently amplifi ed by using the multiplex assay format. Amplifi ed products were diluted 1:50 and 1:100 and mixed with 0.5 μL MapMarker X-Rhodamine labeled 400-bp ladder (BioVentures, Murfreesboro, TN, USA). Sizes were determined by using the ABI Prism 3130xl Genetic Analyzer (Applied Biosystems, Foster City, CA, USA); VNTR sizes were determined by using GeneMapper version 4.0 software (Applied Biosystems). Sizing data for all strains were compared with those found for the B. pertussis prototype strain, Tohama I, to determine the repeat count for each locus. The assignment of an MLVA type was based on the combination of repeat counts for VNTRs 1, 3a, 3b, 4, 5, and 6 and was consistent with international nomenclature. Novel MLVA combinations were submitted to the laboratory of Frits Mooi (National Institute for Public Health and the Environment, Bilthoven, the Netherlands) for MLVA type designation.

Population Analysis
Typing data among strains were compiled by using BioNumerics software version 5.01 (Applied Maths, Sint-Martins-Latem, Belgium). Minimum spanning trees (MSTs) were generated by using default settings and the Manhattan coeffi cient. The Simpson index of diversity (DI) and 95% CIs were calculated as described by Hunter and Gaston (12) and Grundmann et al. (13), respectively. DI was calculated by using a combination of MLVA + MLST to defi ne types. For example, MLVA 27-prn2-ptxP3-ptxS1A-fi m3A was considered a unique type from MLVA 27-prn2-ptxP3-ptxS1A-fi m3B. DI is represented as 1 -D × 100 so that the level of diversity is proportional to the percentage. The Pearson correlation coeffi cient (r) was used to detect linear dependence between pertussis notifi cations and predominant molecular changes.

Identifi cation of Strains using MLVA + MLST
The prevaccine era (period 1,[1935][1936][1937][1938][1939][1940][1941][1942][1943][1944][1945] is depicted in Figure 2, panel A, left side; all 3 strains encoded the same MLST profi le, prn1-ptxP1-ptxS1D-fi m3A. The strain identifi ed in Figure 2, panel A, as MLVA 167 is the 10536 strain used in the manufacture of the Sanofi -Pasteur (Swiftwater, PA, USA) aP in the United States (7). Neither the MLVA types found (167 or 205) nor the MLST profi le for Many MLVA and MLST types were found among the 76 strains in period 3 (1970-1990), the late wP era ( Figure  2, panel B). MLVA 10 was still present from period 2, but other types also dominated, including MLVA 29. MLVA 27, the dominant type among isolates from period 8, emerged in 2 strains from Ohio and 2 from Missouri isolated in 1989. Many of the strains characterized in period 3 differed from period 2 in ptxS1 by encoding the A allele, which was fi rst observed in a 1970 isolate from Colorado. In addition, the fi rst prn2 allele was found in a 1983 isolate from Washington, DC, whereas the ptxP3 allele was fi rst characterized in an Ohio isolate from 1989 ( Figure 1). The 1989 isolate from Ohio was the fi rst in our random selection of US isolates that encoded the combined typing data of MLVA 27 with prn2-ptxP3-ptxS1A-fi m3A ( Figure 2, panel B). This MLST pattern was dominant for the subsequent 2 periods and represents a single-locus intermediary to the prn2-ptxP3-ptxS1A-fi m3B MLST pattern (dominant in period 8).
The MST of 86 strains from period 4 (1991-1996) is shown in Figure 2, panel C. The aP vaccine was recommended for the fourth and fi fth doses of the childhood series in 1991. During this time, MLVA 27, ptxP3, and prn2 were dominant. The fi m3B allele was fi rst noted in an Idaho isolate from 1994. The prn1-ptxP1-ptxS1B-fi m3A MLST type that was widely distributed throughout periods 2 and 3 was restricted to 8% of the selected strains during period 4 and corresponded with MLVA 10 (also observed in periods 2 and 3).
Period 5 (1997-1999) is depicted in Figure 3, panel A. The aP was recommended for all 5 doses of the childhood immunization series in 1997. MLVA 27 increased to 73.6% of the strains compared with 62.1% for period 4 ( Figure 4). Diversity constricted to a total of 7 MLST types (10 were seen previously). The fi m3A* allele was identifi ed in 1999 in a New Hampshire isolate and is shown within the MLVA 28 circle in Figure 3, panel A. The proportion of the population that encoded prn2-ptxP3-ptxS1A-fi m3B (Figures 2, 3) increased from 6.9% in period 4 to 30.8% in period 5. Figure 3, panel B, shows the MST for the mid-aP era (period 6). MLVA 27 and the prn2-ptxP3-ptxS1A-fi m3B profi le continued to dominate. However, MLVA 27 represented 76.5% of selected strains, thus reaching a plateau in frequency. Meanwhile, the prn2-ptxP3-ptxS1A-fi m3B MLST profi le increased to represent 58% of strains. The MLST type prn1-ptxP3-ptxS1B-fi m3A that was previously found in periods 2-4 reappeared with a new MLVA type, 238.  The late-aP era (period 7) is shown in Figure 3, panel C. During this time, a novel pertactin allele (prn14; GenBank accession no. HQ165753) was identifi ed in an isolate from New York in 2004, shown in Figure 3, panel C. MLVA 27 decreased to 67.3% of the strains while the MLST profi le (prn2-ptxP3-ptxS1A-fi m3B) increased to 77.6%.
The MST for period 8 is shown in Figure 3, panel D. The Tdap booster for adolescents and adults was recommended for use in 2005. MLVA 27 frequency remained approximately the same as for period 7, 64%. The MLST profi le prn2-ptxP3-ptxS1A-fi m3B increased to 81.6% of the strains. Meanwhile, the prn1 (n = 6), ptxP1 (n = 6), and ptxS1B (n = 2) alleles reemerged; these alleles are also encoded by vaccine strain Tohama I. The last time these alleles were seen in multiple strains was in period 4, with 17 strains encoding prn1, 28 strains encoding ptxP1, and 8 strains encoding ptxS1B.

Trends in Typing Data during 74 Years of US History
DI values and 95% CIs are provided in the Table. During periods 1 and 2, the DI was in the mid to upper 60% range, but the CIs were large due to a low sample size. Period 3 (1970-1990) had a DI of 94.0% with a small CI that was distinct from previous time periods. To determine if the length of the time interval (20 years) was biasing the results, period 3 was subdivided into 5-and 10-year intervals; all DI values remained ≈90% with small 95% CIs. DI decreased to 75.7% in period 4 and remained relatively constant following the introduction of aP.
The frequencies of individual MLST alleles and MLVA types over time are shown in Figure 5. Changes in ptxS1 occurred fi rst with the transition of ptxS1B to ptxS1A beginning in the 1970s (fi rst observed in a 1970 isolate from Colorado). Later, changes within prn, ptxP, and MLVA 27 occurred at approximately the same time, with transitions  to MLVA 27, prn2, and fi nally ptxP3. In 1995, ptxS1A was encoded by 100% of tested strains, and in 1996, >90% of strains tested encoded prn2 or ptxP3, but more recently, frequency of prn2 and ptxP3 has declined (period 8). The transition within fi m3 observed in the early 2000s was a more gradual increase that did not approach 100% as with the other MLST alleles.

Comparison of Typing Data and Increase in US Pertussis Notifi cations
We aligned annual US pertussis notifi cations with vaccine coverage data, the trend lines for each of the dominant MLST alleles, MLVA 27, and DIs ( Figure 6). An inverse relationship was observed between DI and notifi cations as well as vaccine coverage. Increases among ptxS1A, prn2, ptxP3, and MLVA 27 were not signifi cantly correlated with the increase in notifi cations, whereas the proportion of fi m3B-encoding strains was signifi cantly correlated with pertussis notifi cations (r = 0.8608; p = 0.0277).

Discussion
Changes in the US B. pertussis population have followed trends that are largely consistent with other nations and yet were unique in the correlation of annual case counts with fi m3B. Our fi ndings regarding DI trends were similar to those found for the United Kingdom (9), but unlike a previous study that examined B. pertussis in the Netherlands (14), we did not fi nd a correlation between the ptxP3 allele and annual reporting of pertussis cases in the United States.
MLVA + MLST analysis showed that the US B. pertussis population changed genetically during the period covered by our study (Figures 2, 3). We found a high degree of diversity in period 3 (1970s and 1980s), with a DI of 94%, that was consistent with fi ndings (84%) for a similar period in the United Kingdom (7). The results from the United Kingdom were attributed to a decline in vaccine coverage rates below 80% during 1975-1989, with a low of 31% in 1978 (7). Parent-supplied data from the US Immunization Survey for 1962-1985 showed US coverage among 2-yearold children for >3 doses of wP declined from a peak of 77.9% in 1967 to a low of 63.6% in 1985; vaccination rates then increased dramatically, to >94% by 1994, and have remained high ever since ( Figure 6) (15).
Given that B. pertussis has no nonhuman hosts or environmental niche, vaccine-mediated immunity is the most likely selective pressure against B. pertussis. Therefore, vaccination coverage may have contributed to the increase in diversity during period 3. This hypothesis also supports the correlation between the decline in diversity as the rate of vaccination coverage (>3 doses) increased to ≈95% in the mid-1990s. The emergence and subsequent dominance of MLVA 27 and the MLST alleles prn2 and ptxP3 (in that order, temporally) in the United States occurred during the period with the highest diversity, period 3.The timing and order of these transitions are consistent with global trends (16) in which the emergence of nonvaccine-type alleles for ptxS1 and prn appeared 15-30 years after the introduction of pertussis vaccines. Despite increasing pertussis incidence in the United States, diversity has remained stable for B. pertussis during the past 20 years; a similar trend was observed in the Netherlands (6).   Regional clustering of B. pertussis MLVA types may be best exemplifi ed by comparing the United States to Australia. Kurniawan et al. found that a single MLVA type, 70, emerged in Australia's wP era and dominated during the aP transition (reaching ≈40% of isolates), then declined once the transition to aP was complete (8). A similar rise and fall of MLVA 70 surrounded the transition to aP in the United Kingdom, where MLVA 70 disappeared after 2004 (7). However, to our knowledge, MLVA 70 has only been detected once in the United States, in a Missouri isolate from 1997. Further, MLVA 64, which represents ≈10% of the B. pertussis population in Australia, was not detected in the United States. Our fi ndings, combined with the fi ndings from the United Kingdom, do not support the conclusion that the introduction of aP was the driving factor responsible for the emergence and dominance of MLVA 70 or 64. On the contrary, the data suggest that aP helped to eliminate MLVA 70. Wholegenome sequencing has implicated regional bottlenecks as a likely contributor to the geographic restriction of particular MLVA types (17,18), which may explain the exclusive prevalence of MLVA 70 in Australia (8). Furthermore, the timing of the emergence and dominance of the MLVA 27 and MLST alleles prn1, ptxS1A, and ptxP3 in the United States predate the completion of the wP to aP transition (1997) by ≈10 years (Figures 1-3, 5, 6). The emergence and dominance of the fi m3B allele ( Figure 6) probably coincides with the increase in notifi cations.
Many of the allele changes we found have been identifi ed in clinical isolates of B. pertussis throughout the world. However, any conclusions involving individual gene loci in clinical isolates are vulnerable to phenotypic results that arise from mutations elsewhere in the genome. The clearest way to identify the effect of allele mutations is to examine them alone and in combination in a genetically controlled bacterial background (19). Therefore, it is premature to associate allele changes with a phenotype, disease severity, or events of epidemiologic importance until they are functionally analyzed individually and cumulatively in a model system. Data related to the functional or clinical effects of allele changes for the MLST targets used in this study are limited.
Recently, an increase in pertussis notifi cations and a 1.41-fold increase in hospitalizations were correlated with the increasing presence of the ptxP3 allele in circulating B. pertussis isolates (20). Unfortunately, the ptxP allele was not characterized for B. pertussis among hospitalized patients, so a direct correlation could not be made between ptxP3 and disease severity. In addition, ELISA was used to demonstrate a modest increase (1.6-fold) in pertussis toxin production among ptxP3 strains relative to ptxP1-encoding strains when grown in vitro. In vivo experimentation using genetically controlled B. pertussis mutants for ptxP3 is needed to determine whether a 1.6-fold increase is suffi cient to cause more severe disease or to overwhelm vaccine-mediated anti-pertussis toxin antibody response. Bart et al. hypothesized that ptxP3 may be a "hitchhiker" mutation that benefi ted from advantageous mutations selected elsewhere in the genome (18); our fi ndings lend support to this hypothesis.
More information is known about the effects of the prn and ptxS1 allele changes, but studies assessing the fi m3 locus are lacking. The divergence among the pertactin alleles is proximal to the encoded RGD motif that is involved in eukaryotic cell binding and antigen presentation to B cells (21). In theory, such mutations could biochemically affect protein folding, host cell binding, or recognition by B cells (14); this hypothesis is supported by the fi nding that the pertactin variants 1-3 induce type-specifi c antibodies (22). However, the effects of these insertions/deletions on Vaccine coverage data were collected for the United States Immunization Survey (USIS, 1962(USIS, -1985, National Health Interview Survey (NHIS, 1991(NHIS, -1993, and National Immunization Survey (NIS, 1994(NIS, -2009. No data are available for 1986-1990 because USIS was cancelled (15). The fi m3B trend line was temporally and signifi cantly associated with the rate of increase for pertussis notifi cations. A color version of this fi gure is available online (wwwnc.cdc.gov/EID/ article/18/8/12-0082-F6.htm).
the pertactin protein product have not been determined experimentally. The ptxS1A allele encodes 3 amino acid changes relative to the 10536 vaccine type, ptxS1D, and 1 amino acid change compared with ptxS1B for Tohama I-based vaccine (23). The pertussis toxin remains biologically functional despite these changes (24). In vivo, mousederived anti-pertussis toxin antibodies tolerate numerous amino acid substitutions in ptxS1 with equal neutralization between wild-type and mutant isolates (25). Therefore, ptxS1 allele changes may not be clinically or immunologically relevant.
Little is known about the functional or in vivo effects of the fi m3B mutation on protein function, bacterial survival, and adherence. The fi m3B allele results in an alanine-toglutamic acid mutation at aa 87 (11). Biochemically, this is a potentially important residue change, and the mutation is located in a surface epitope of fi m3 (aa 79-91) that has been shown to interact with human serum (26). Given the signifi cant correlation between the increases in fi m3B and US pertussis case notifi cations, the effect of the fi m3 mutation needs to be functionally and clinically determined. Alternatively, this could be another example of a regional bottleneck (17,18); additional data regarding the prevalence of fi m3B in other countries is needed to rule out this possibility. Moreover, the strain collection available for this study may not be fully representative of the B. pertussis population in the United States over time. Efforts were made to ensure that the strains selected for this study were diverse in year of isolation as well as geography within the United States. According to Mouillot (27), DI can be infl uenced by selection bias and sample size, but almost all studies evaluating a historical collection of strains encounter this limitation, including the recent study in the United Kingdom (7).
In summary, the US B. pertussis population has evolved in the time since vaccinations were introduced in the 1940s (Figures 2-5). Our fi ndings demonstrate that the resurgence of pertussis in the United States was not correlated with the ptxP3 allele but with the presence of the fi m3B allele among the B. pertussis population. The commonly circulating strains of B. pertussis in the United States encode different alleles compared with the strains used for manufacture of the pertussis vaccines, but the relevance of these allele changes remains to be fully elucidated. Because B. pertussis has no nonhuman host, the selective pressures it encounters are limited to the human immune system and the vaccine, but the infl uence of this selection pressure versus natural evolution on the modern US B. pertussis population is unclear. For example, minor types are beginning to emerge, including the reemergence of vaccine-type alleles. In addition, the vaccine policies of other nations may have contributed to the makeup of the US B. pertussis population in ways that could not be measured by using this study of US-based isol ates.
As vaccine coverage rates improve among adolescents and adults, changes in the B. pertussis population should be monitored through molecular typing. The in vivo effects of MLST allele changes in a genetically controlled model for pathogen and host should be characterized to determine what effects, if any, these allele changes have with respect to vaccine-mediated immunity to circulating B. pertussis. More specifi c studies, such as genomic sequencing of particular strains and genetic expression of the multiple alleles in animal models, should be performed to determine the virulence and pathogenesis of these variants.