The emergence of HIV/AIDS in the Americas and beyond

Under the “Haiti-first” model, non-Haitian subtype B strains are expected to be phylogenetically nested within an older and hence more extensive range of Haitian genetic variation, with Haitian lineages branching off closest to the B subtype ancestor. To test whether it is this or an alternative pattern that characterizes this HIV-1 subtype [supporting information (SI) Fig. 3], we conducted a detailed Bayesian Markov chain Monte Carlo (MCMC) phylogenetic analysis using an alignment of the Haitian archival sequences plus 117 previously published subtype B env sequences from a total of 19 countries. Five African strains of subtype D (the closest relative to subtype B) served as the outgroup.

On the env gene phylogeny, the archival sequences occupy basal positions within subtype B (Fig. 1 and SI Fig. 4); there is extremely strong statistical support for a Haitian origin of subtype B, in that the probability of a U.S. or other non-Haitian origin [i.e., the posterior probability of any non-Haitian sequence(s) occupying the most basal position in the B clade] is <0.001. These results indicate, with high statistical significance, that subtype B arrived and began spreading in Haiti before it spread elsewhere and was not originally introduced to Haiti from the United States. In addition, the unequivocal support for the monophyly of subtype B as a whole (P = 1.0) supports the previous suggestion of a single (epidemically successful) introduction of this subtype from central Africa to Haiti (10).

Drummond et al. (15) found that relaxed-clock models are more accurate and precise at estimating phylogenetic relationships than unrooted methods. In other words, when data have evolved in a somewhat clock-like fashion, incorporating knowledge of the tempo of evolution can make topology estimation more reliable compared with methods that ignore this information. Therefore, we were also interested in what the relaxed-clock results (described below) revealed with respect to where subtype B originated (i.e., the topological relationships among subtype B sequences, regardless of timing/branch lengths). Under the relaxed-clock model, the posterior probability of a U.S. origin was 0.00003, and the probability of a Haitian origin was 0.9979.

The inference that subtype B reached Haiti before spreading to other countries does not depend on a dating analysis. One of the advantages of a Bayesian statistical framework is that it yields direct estimates of the probability of phylogenetic hypotheses, and in this case there is strong evidence to reject a U.S. or other non-Haitian origin of subtype B. This means that even if there is some uncertainty regarding precisely when HIV-1 entered Haiti or the United States (see below), there is little doubt about the sequence of events; the clear-cut topological information implies that the entry to Haiti occurred first. Moreover, our sampling bias in favor of non-Haitian subtype B makes the “Haiti-first” inference conservative; the Haitian strains occupy the basal positions within subtype B even though there were many more opportunities for recovering non-Haitian basal strains, if they existed.

We next investigated how many times the older Haitian HIV-1 epidemic has seeded detectable secondary outbreaks elsewhere. We found evidence for only three such events, despite the fact that our data set comprised 109 non-Haitian subtype B strains from the Caribbean, North and South America, Europe, Africa, Asia, and Australia. One instance is confined to Trinidad and Tobago (Fig. 1). All of the HIV-1 sequences from this country form a distinct monophyletic grouping with the strongest possible support (posterior probability of 1.0). Previous studies have asserted homosexual/bisexual contact with North American foreigners in the late 1970s or early 1980s as the alleged route of introduction (16, 17). Our results indicate this was not the case; the Trinidad and Tobago clade is unequivocally nested among Haitian strains, not North American ones, clear evidence that the predominantly heterosexual epidemic in this country can be explained by a single introduction linked to Haiti.

The next, and most important, secondary epidemic accounts for all but three of the remaining non-Haitian strains, encompassing 96 sequences representing every other country in our data set (Fig. 1). This “pandemic” clade forms another unequivocally monophyletic cluster (P = 1.0) nested within the basal Haitian strains. As for the Trinidad and Tobago clade, the most parsimonious explanation for this pattern is that all these subtype B infections from across the world emanated from a single founder event linked to Haiti. This most likely occurred when the ancestral pandemic clade virus crossed from the Haitian community in the United States to the non-Haitian population there.

The only other Haiti-linked outbreak we detected is comparatively insignificant, here comprising a single Brazilian sequence (BR020) and one American one (US4). Additional rare chains of transmission may have emerged from Haiti but remained undetectable in this study. Likewise, pandemic clade viruses may have reentered Haiti but were undetectable. Regardless, it is evident from this large international sample that the subtype B epidemics in most afflicted countries and the bulk of subtype B infections worldwide are caused by viruses belonging only to the pandemic clade.

Three additional nominally non-Haitian strains are noteworthy for falling among the basal Haitian lineages rather than within the three non-Haitian clades (Fig. 1). Each one provides further support for the notion that strains linked to Haiti occupy the deepest branches within subtype B. RF was sampled in 1983 from a Haitian immigrant to the United States. Like our five newly sequenced strains, it represents a Haitian strain that entered the United States via an immigrant host (10). WMJ2 was sampled in 1985 from a perinatally infected infant born to an HIV-positive Haitian immigrant mother (18). KR5086 was recovered in 1995 from a South Korean sailor infected in the Dominican Republic (19), a country directly linked to Haiti in terms of both geography and HIV-1 epidemiology (20). These cases, plus the additional ones considered in our study, show that the virus moved out of Haiti on many separate occasions but did so mostly as dead-end infections that evidently failed to ignite successful epidemics. A key question for future research will be to determine why Haiti-to-U.S. epidemic outbreaks have apparently become established so rarely since the initial introduction in or around 1969, despite presumably frequent movement of the virus because of rising HIV-1 prevalence, continued migration, and the once-thriving sex tourism industry linking Haiti and the United States.

Several additional analyses corroborated these findings. First, we reanalyzed our env gene data by replacing the African subtype D outgroup sequences with a variety of alternatives from different subtypes and geographical regions. This had no impact on the basal position of the Haitian strains within subtype B (SI Fig. 5a); their ancestral position cannot be dismissed as an artifact of convergent evolution with African subtype D sequences. Second, the gag sequences also unambiguously place the Haitian sequences in the most ancestral positions within subtype B (P = 0.9930) (SI Fig. 5b). We also inspected the nucleotide substitutions that mapped specifically onto the branch leading to the pandemic clade. For env, all eight such changes were silent at the amino acid level. For gag, only one of the six changes on the relevant branch indicated a change in amino acid between the Haitian strains and the pandemic clade strains. The paucity of amino acid substitutions along this clade-defining branch suggests that the ancestor of the pandemic clade probably possessed no selective advantage over other Haitian strains; its remarkable epidemic success may simply reflect ecological factors rather than evolutionary ones (chance colonization of a new population, as opposed to competitively superior transmission fitness). However, analysis of complete genomes would be necessary to definitively rule out selection.

We used a second Bayesian MCMC method that simultaneously estimates phylogenetic relationships and times of most recent common ancestors (15) to perform a supplementary phylogenetic analysis on a reduced data set. This method uses a “relaxed molecular clock” model, so-called because it relaxes the need to assume a constant rate of molecular evolution across the tree to obtain date estimates from gene sequences. We estimated the time of the most recent common ancestor (TMRCA) of subtype B at 1966 (1962–1970) (Fig. 2), a date that suggests its arrival in Haiti may have occurred with the return of one of the many Haitian professionals who worked in the newly independent Congo in the 1960s (21).

The TMRCA of the U.S. epidemic is estimated to be 1969 (1966–1972) (Fig. 2), suggesting that HIV-1 was circulating cryptically in the United States for ≈12 years before the initial recognition of AIDS in 1981. The evidence of a single origin of the pandemic variant of subtype B allowed us to date the beginning of the actual U.S. HIV-1 epidemic, rather than the ancestor of multiple viruses introduced from Haiti; if multiple introductions from Haiti had occurred, the “U.S. MRCA” would actually correspond to a Haitian virus that predated the initial entry of HIV-1 to the United States (11). Serological evidence of an ≈5% prevalence of HIV-1 by 1978 in men who have sex with men (MSM) populations in both San Francisco (22) and New York City (23) suggests that several thousand individuals in the United States would already have been infected by then. Even assuming the fastest-documented growth rates for HIV-1 (24), this implies that the virus had been spreading in the MSM population for several years before this point, consistent with our findings here. The virus may well have been spreading slowly for an extensive period, perhaps in the heterosexual population, before entering the highest-risk MSM subpopulation, where it spread explosively enough to finally be noticed. We contend that the phylogenetic estimate, with appropriate confidence intervals, provides more reliable information on the date of the origin of the U.S. epidemic than the available epidemiological data, which cannot resolve this question.

Nevertheless, although our relaxed-clock methods should be reasonably robust to rate variation among lineages and uncertainty in phylogenetic inference, some caution is always warranted when such inferences are made. For example, the relatively sparse sampling of both Haitian and Trinidadian sequences means it is conceivable that more intensive future sampling could recover deeper-branching lineages and push back these TMRCA estimates slightly. This is unlikely to apply to the U.S. TMRCA, on the other hand, because of the already dense sampling of this HIV-1 population. If obtainable, additional archival sequences should help clarify the early spread of subtype B with greater precision.

The three-decade gap between the estimated timing of the HIV-1 M group ancestor (9) and the earliest evidence of HIV/AIDS in Africa (25, 26) seems unexceptional in comparison to the U.S. cryptic period, especially because a good deal of tuberculosis-caused mortality in Africa must have gone unrecognized as AIDS-related, then as now. Taken together with our dating analysis, including the subtype B/D ancestor dated to 1954 (1946–1961) (Fig. 2), the extensive cryptic period in the United States, therefore, provides compelling corroboration of an early 20th century M group ancestor (9).

Peripheral blood mononuclear cell (PBMC) samples were collected in 1982 and 1983 at Jackson Memorial Hospital in Miami, FL, during one of the first investigations establishing that Haitians in Haiti and elsewhere were at risk for AIDS (6). One of the six PBMC samples obtained for this study failed to yield any amplifiable HIV-1 PCR products. As described in Pitchenik et al. (6), all of the patients were Haitian immigrants who had entered the United States after 1975 and progressed to AIDS by 1981 and hence were presumably infected with HIV-1 before entering the United States. The position of these sequences on the subtype B phylogeny (distinct from and basal to the dominant U.S. variant of subtype B) is consistent with this sequence of events.

DNA was extracted from 10 μl of peripheral blood mononuclear cells by using QIAamp DNA micro kits (Qiagen, Valencia, CA), following the manufacturer's instructions for extractions from blood. After extraction, DNA was eluted into 100 μl of elution buffer AE and stored frozen at −20°C until required for DNA analyses. DNA was PCR-amplified from the extracts by using a nested PCR approach (28). First-round PCRs were undertaken in 25 μl of final volume reactions, using 0.1 μl of Platinum Taq HiFi enzyme (Invitrogen, Carlsbad, CA)/0.1 μl of 25 mM dNTP mix/2.5 mM (final concentration) MgSO₄/10× PCR buffer/1–5 μl of DNA extract. Second-round amplifications were performed on 1 μl of the first-round PCR product by using the same reagent concentrations. Enzyme activation, dissociation, and extension temperatures followed the manufacturer's guidelines, with an extension time of 3 min. Annealing temperatures varied by extract in response to initial amplification success rates.

After amplification, the PCR products were visualized on 0.8% agarose gels stained with ethidium bromide and then purified by using QIAquick spin columns (Qiagen). Purified products were sequenced by using several overlapping primer pairs (28) by the University of Arizona Genomic Analysis and Technology Core Facility with ABI Big Dye 3.1 chemistry (Applied Biosystems, Foster City, CA) on Applied Biosystem 3730xl DNA Analyzers. Each sample was extracted, PCR-amplified, and sequenced twice to ensure that the sequences generated were not modified through low template copy number. We recovered five full-length env sequences and five partial (0.7- to 1.2-kb) gag sequences.

We used the Los Alamos National Laboratory HIV sequence database (29) to download all full-length published env and gag gene sequences of subtypes B and D. We then subjected the resulting sequence set to strict quality control measures to remove (i) incomplete sequences; (ii) sequences not published in a peer-reviewed journal; (iii) multiple sequences from the same patient; (iv) sequences suspected a priori of possibly anomalous evolutionary patterns (from long-term nonprogressors with nef deletions, laboratory workers infected accidentally, sequences exhibiting evidence of hypermutation, etc.); (v) sequences with midpeptide stop codons, frame-shift mutations, or nonnucleotide characters; or (vi) sequences for which there was any uncertainty regarding which subtype they belonged to. To search for such sequences, which might include unidentified intersubtype recombinants, we screened all sequences using the REGA HIV-1 Subtyping Tool (30) and then removed any sequence with bootstrapping support <100% or bootscanning support <1.0 for clustering with subtype B or D.

To ensure that no very-early-diverging subtype B strains were removed by this procedure, we inferred additional phylogenies that included the “cleaned” sequences (data not shown). For env, the only such strains that were positioned basal to the pandemic clade were from Trinidad and Tobago, and these clustered with the other sequences from the Trinidad and Tobago clade. Similarly, for gag, the only nonpandemic clade sequences that were removed (US4 and RF) were ones that had already been identified as basal in the env analysis; all other sequences with bootstrap or bootscan support <100% or 1.0 were from the pandemic clade. Unlike intersubtype recombinants, the possible nonexclusion of intrasubtype recombinants is not expected to affect inferences regarding the geographical origin and emergence of subtype B because intrasubtype recombination cannot plausibly lead to strains from one locality systematically falling basal to all of the others. Moreover, although unidentified intrasubtype recombination might increase the variance of dating estimates, it is unlikely to systematically bias these dates in one direction or the other in an exponentially growing population (31).

The resulting data sets were codon-aligned and then adjusted by eye in Squint Ver. 1.0 (M. Goode, University of Auckland, Auckland, New Zealand), and regions of ambiguous alignment were removed. We constructed three additional env alignments replacing the D subtype outgroup with a subtype C sequence from India, a subtype A sequence from Kenya, or a CRF01 (A/E) sequence from Thailand. All previously published sequences (see the taxon labels in SI Figs. 4 and 5b) are available from the Los Alamos National Laboratory HIV sequence database. All alignments are available from the authors upon request.

We used MrBayes, Ver. 3.1 (32), to perform two independent runs of 20 million steps (for env). Examination of the MCMC samples with Tracer, Ver. 1.3 [A. Rambaut (University of Edinburgh) and A. J. Drummond (University of Auckland, Auckland, New Zealand); http://beast.bio.ed.ac.uk], indicated adequate mixing of the Markov chain. We discarded the first 2 million steps from each run as burn-in and combined the resulting MCMC samples (n = 36,002) for subsequent estimation of posteriors. Fewer steps (5 million) were required for convergence and adequate mixing in the gag analyses.

We used tree filtering in PAUP, Ver. 4.0b10 (33), to calculate the posterior probabilities of a Haitian, non-Haitian, or effectively simultaneous Haitian/non-Haitian origin of subtype B. Briefly, we removed from the posterior sample any tree with a Haitian sequence(s) in the most basal position. Only 4 of 36,002 env trees had a non-Haitian basal sequence (P_{non-Haitian-origin} = 0.00011), and 24 others placed the Haitian and non-Haitian sequences in reciprocally monophyletic clades (P_{simultaneous-origin} = 0.00066), which left 35,974 with a Haitian sequence or group of sequences, in the most-ancestral position(s) within subtype B (P_{Haitian-origin} = 0.9992). A similar approach was followed for the estimation of the probability of a Haitian or non-Haitian origin of subtype B for the gag and the relaxed-clock analyses.

For both the env and gag data sets, we used a parsimony approach as implemented in MacClade, Ver. 4.08 (34), to identify the nucleotide substitutions that mapped onto the branch leading to the pandemic clade of subtype B. We then determined whether these changes were synonymous or nonsynonymous.

To infer the timescale of HIV-1 group M subtype B evolution, we used a Bayesian molecular clock method (15), as implemented in BEAST (http://beast.bio.ed.ac.uk), under an uncorrelated log-normal relaxed molecular clock model with a Bayesian Skyline coalescent tree prior. For this analysis to run in a reasonable time, we considered only subtype B sequences from Haiti, Trinidad and Tobago, and the United States (plus one from Canada), and we removed some pandemic clade sequences from overrepresented years and localities. The maximum a posteriori tree from the MrBayes analysis (available upon request) revealed that the sequences were scattered across the entire pandemic clade, consistent with a U.S. entry of the founding virus. We ran 10 independent MCMC analyses, with each run consisting of 100 million steps, and then discarded the first 5 million steps from each run as burn-in and combined the resulting postburn-in MCMC samples for subsequent estimation of posteriors.