We selected 661 B. pertussis isolates of US origin from the Centers for Disease Control and Prevention (CDC) collection by using random sampling stratified 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). Stratification 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 Technical Appendix. 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.

Analysis was performed by using a 6-target multiplex similar to that described (8) with some modifications. Using the HotStarTaq kit (QIAGEN, Valencia, CA, USA) yielded a final reaction volume of 20 µL. Master mix 1 consisted of fluorescently 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 fluorescently 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 efficiently amplified by using the multiplex assay format. Amplified 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.

Our algorithm consisted of 4 DNA targets: the pertactin (prn) gene, the first gene in the pertussis toxin operon and its respective promoter (ptxP-ptxS1), and the fimbrial protein-encoding gene (fim3). The prn and fim3 genes were amplified by using oligonucleotides and conditions as described (6,11). The ptxP-ptxS1 region was amplified by using oligonucleotides Ptox1Fpert (5′-CCCTCGATTCTTCCGTACATCC-3′) and Ptox2R (5′-CGCGATGCTTTCGTAGTACA-3′), resulting in an amplified product of 964 nt. Products were sequenced and analyzed as described (10).

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 coefficient. 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 define types. For example, MLVA 27-prn2-ptxP3-ptxS1A-fim3A was considered a unique type from MLVA 27-prn2-ptxP3-ptxS1A-fim3B. DI is represented as 1 – D × 100 so that the level of diversity is proportional to the percentage. The Pearson correlation coefficient (r) was used to detect linear dependence between pertussis notifications and predominant molecular changes.

The prevaccine era (period 1, 1935–1945) is depicted in Figure 2, panel A, left side; all 3 strains encoded the same MLST profile, prn1-ptxP1-ptxS1D-fim3A. The strain identified 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 profile for period 1 are seen again in later periods. The dotted line between MLVA circles 191 and 167 in Figure 2, panel A, indicates a distant genetic relationship between the respective clusters for periods 1 and 2. The predominant MLVA type during period 2 (1946–1969), the early wP era, was 10 (Figure 2, panel A, right side), and the MLST profile shifted to prn1-ptxP1-ptxS1B-fim3A, identical to Tohama I (identified as MLVA circle 83 with a star), representing a single-locus change to ptxS1B compared with period 1. The strain used for manufacture of the GlaxoSmithKline (Research Triangle Park, NC, USA) pertussis vaccine in the United States is Tohama I.

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 first observed in a 1970 isolate from Colorado. In addition, the first prn2 allele was found in a 1983 isolate from Washington, DC, whereas the ptxP3 allele was first characterized in an Ohio isolate from 1989 (Figure 1). The 1989 isolate from Ohio was the first in our random selection of US isolates that encoded the combined typing data of MLVA 27 with prn2-ptxP3-ptxS1A-fim3A (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-fim3B 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 fifth doses of the childhood series in 1991. During this time, MLVA 27, ptxP3, and prn2 were dominant. The fim3B allele was first noted in an Idaho isolate from 1994. The prn1-ptxP1-ptxS1B-fim3A 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 fim3A* allele was identified in 1999 in a New Hampshire isolate and is shown in green within the MLVA 28 circle in Figure 3, panel A. The proportion of the population that encoded prn2-ptxP3-ptxS1A-fim3B (Figure 2, Figure 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-fim3B profile continued to dominate. However, MLVA 27 represented 76.5% of selected strains, thus reaching a plateau in frequency. Meanwhile, the prn2-ptxP3-ptxS1A-fim3B MLST profile increased to represent 58% of strains. The MLST type prn1-ptxP3-ptxS1B-fim3A 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 identified in an isolate from New York in 2004, shown in light orange in Figure 3, panel C. MLVA 27 decreased to 67.3% of the strains while the yellow MLST profile (prn2-ptxP3-ptxS1A-fim3B) 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 profile prn2-ptxP3-ptxS1A-fim3B 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.

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 first with the transition of ptxS1B to ptxS1A beginning in the 1970s (first 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 finally 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 fim3 observed in the early 2000s was a more gradual increase that did not approach 100% as with the other MLST alleles.

We aligned annual US pertussis notifications 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 notifications as well as vaccine coverage. Increases among ptxS1A, prn2, ptxP3, and MLVA 27 were not significantly correlated with the increase in notifications, whereas the proportion of fim3B-encoding strains was significantly correlated with pertussis notifications (r = 0.8608; p = 0.0277).