Analysis of HIV-1 envelope evolution suggests antibody-mediated selection of common epitopes among Chinese former plasma donors from a narrow-source outbreak

The HIV-1 envelope mutates rapidly to evade recognition and killing, and is a major target of humoral immune responses and vaccine development. Identification of common epitopes for vaccine development have been complicated by genetic variation on both virus and host levels. We studied HIV-1 envelope gp120 evolution in 12 Chinese former plasma donors infected with a purportedly single founder virus, with the aim of identifying common antibody epitopes under immune selection. We found five amino acid sites under significant positive selection in ≥50% of the study participants, and 22 sites consistent with antibody-mediated selection. Despite strong selection pressure, some sites housed a limited repertoire of amino acids. Structural modelling revealed that most of the variable amino acid sites were located on the exposed distal edge of the Gp120 trimer, whilst invariant sites clustered within the centre of the protein complex. Two sites, flanking the V3 hypervariable loop, represent novel antibody sites. Analysis of HIV-1 evolution in hosts infected with a narrow-source virus may provide insight and novel understanding of common epitopes under antibody-mediated selection. If verified in functional studies, such epitopes could be suitable as targets in vaccine development.

Five positions in Gp120 were under significant positive selection pressure in at least half of the study participants irrespective of HLA profile. Ten of the 12 study participants yielded sequences from two or more time-points and were subjected to evolutionary analyses (Fig. S1A). All time-points with available sequences for the ten study participants were included in the analyses. The median evolutionary rate ratio of positions 1 + 2 to position 3 among the participants was 0.861 (IQR: 0.779-0.892), indicating a general purifying selection over the HIV-1 env gp120 C2-V5 region in the majority of the study participants (the overall ratio was <1 in all study participants except SM021 [1.306]). Next, we evaluated the intrapatient dN/dS ratio of each codon by renaissance counting 41 . The majority of codons in Gp120 were under either significant positive or negative selection in seven of the ten study participants with longitudinal samples (Fig. 1), and neutral evolution was Scientific REPORTs | (2018) 8:5743 | DOI: 10.1038/s41598-018-23913-2 comparatively rare. Consistent with the 1 + 2:3 codon rate ratios, only one study participant (SM021) appeared to have more residues under significant positive than negative selection pressure.
Within the variable loops, substantial negative selection could be seen in V3, but not in V4 and V5 ( Fig. 2A). The negative selection in V3 corresponded with a marked increase in the density of neutralising antibody epitopes. Five codons, corresponding to positions T297, A337, S348, D415 and S468 in the HXB2 Gp120 (accession number K03455), were under significant positive selection in 50% or more of the study participants irrespective of HLA profile ( Fig. 1 and Table 1). These sites were further mapped to a homology-modelled structure of the SM cohort Gp120 consensus, and clustered either within or immediately flanking the variable loops on the distal exposed edge of the protein complex (Fig. 2B).
Significant selection pressure was exerted by the humoral response. We next considered how the virus evolves in response to significant positive selection pressure. For each position in the partial Gp120 sequence, deviations from the inferred founder were detected and 288 unique amino acid variants were recorded across 128 variable sites (Fig. 3). The remaining 51 positions were completely invariant (29%). Of the variants identified, many were present in a single or small collection of sequences and are likely of limited biological relevance. Major variants were therefore resolved, and to prevent overrepresentation by particular study participants, a single time-point was selected for each study participant (the selection aimed to get an equal distribution of samples from each year of cohort follow-up). Thirty-one major variants were detected in total, across 29 sites (Fig. 3).
In study participants with longitudinally-sampled sequences, the presence or absence of each major variant was recorded at each time-point. Whether these sites were under significant positive selection pressure in that study participant was also determined (Fig. 4A). Twenty-four of the 29 sites housing major variants were under significant positive selection pressure in at least one study participant, and of these sites, a higher proportion exhibited fluctuating patterns of variant emergence than a consistent pattern wherein the variant was present at all time-points sampled (p < 0.01, two-tailed Fisher's Exact Test).
Mapping these 24 sites to the homology-modelled structure of the SM cohort Gp120 consensus revealed that all but two were visible on the surface of the protein and likely accessible to antibodies (Fig. 4B). The exceptions were positions 345 and 424. Position 345 houses the major variant I345V, and is contained within the location of a known HLA-A11-restricted CTL epitope 42 . This position was found only to be under significant positive selection pressure in study participant SM176, who also expresses HLA-A11. Position 424 houses the CD4 binding site. The wholly invariant sites were similarly mapped, and the overwhelming majority clustered on the inner face of the quaternary structure. Despite significant antibody-mediated selection pressure, some sites housed a limited repertoire of amino acids. We next considered the biophysical diversity of amino acids in each of the 22 sites under significant positive selection pressure and visible on the surface of the protein (Fig. 5). In all positions, the biophysical properties of the amino acid in the inferred founder were preserved in approximately 50% or more of the sequences. Little divergence was seen in sites where the inferred founder residue was hydrophobic, with other properties being somewhat more variable. Tabulating the amino acids present in each site also reveals that seven positions flick back and forwards between just two or three amino acids.
Four novel sites were identified and consistent with antibody activity. The 22 sites under significant positive selection pressure and visible on the surface of the protein were also cross-checked for existing antibody epitopes in the LANL Immunology Database and Genome Browser (  previously been associated with neutralising antibody activity, whereof 18 were part of previously described antibody epitopes from different sources (human, mice, or both). In detail, five positions (T236, Q344, I345, V424, and D474) were contained within known human antibody epitopes, whilst 16 positions were contained within epitopes reported in mice. Two sites flanking the V3 hypervariable loop (E290 and L333) -which itself houses the majority of antibody epitopes -were identified against which no antibodies have yet been reported. The L333 site primarily switched between I/L.

Discussion
In line with previous observations 43 , mapping the ratio of non-synonymous to synonymous substitutions showed that the majority of sites within the C2-V5 region of Env Gp120 were under negative selection in all but one study participant out of ten. Whilst V4 and V5 loops exhibited a dearth of negative selection, the V3 hypervariable loop contained substantial negative selection. Of the five variable loops in Gp120, V3 is the most conserved with amino acid variation restricted to approximately 20% of the loop's residues 44 . It is also likely that V3 is subject to stronger functional constraints due to its important role in co-receptor binding [45][46][47][48] . Moreover, it has been shown that deletion of V3 abrogates viral infectivity 49 . Several sites within each study participant showed evidence of significant positive selection, and five of these were common to at least half of the study participants sampled. Structural modelling demonstrated that all but two of the 24 positively selected sites were found on exposed regions of the outer face of the protein complex. CTL epitopes in the HIV-1 Nef protein have been reported to cluster in hydrophobic regions 50 , whilst more recent evidence suggests that their distribution may be random across the genome 51 . Such strong clustering on the surface of the protein is therefore more consistent with antibody-mediated than CTL-mediated selection pressure 52 . The exceptions in terms of surface exposure were positions 345 and 424, which were buried within the protein. Notably, position 345 was found to be under significant positive selection pressure in only one study participant, SM176. This position is contained within a known HLA-A11-restricted CTL epitope, which is one of the HLA alleles expressed by study participant SM176. It is therefore feasible that this variant has emerged in response to CTL-mediated selection pressure in this study participant. Conversely, position 424 is important in Figure 1. Patient-specific selection pressure within the HIV-1 envelope Gp120 protein. Ten of the 12 study participants yielded sequences from two or more time-points and were subjected to evolutionary analyses (Fig. S1A). All time-points with available sequences for the ten study participants were included in the analyses. Mean dN/dS ratios for each codon within each study participant. A dN/dS estimate greater than 1 indicates positive selection whilst an estimate of less than 1 indicates negative selection. Sites under significant selection shown in blue whilst sites that have not reached statistical significance are shown in red and are assumed neutral. Differences in the number of codons across study participants are the result of length variation in the V4 and V5 hypervariable loops.
We were able to assign most positions to known antibody epitopes in humans and mice. We also identified two novel sites, which were not contained within any known antibody epitopes reported in the literature. Position 290 has, however, been associated with the CTL-restricted epitope AKTIIVQLTEPVE in the HIV-1 CRF02_AG lineage (https://www.hiv.lanl.gov/content/sequence/genome_browser/browser.html). Moreover, these sites flank the V3 hypervariable loop, which is the most epitope dense region of Gp120 (LANL Immunology Database; https://www.hiv.lanl.gov/content/immunology/), although this may be due to a bias in reporting stemming from the extensive study of V3 in vaccine design rather than a genuine increase in immune activity. Further characterisation by neutralisation experiments could help to confirm these novel surface-exposed epitopes as targets by humoral or cellular immune responses.
Whilst numerous antibodies against V3 have been described, the cross-neutralisation potential of these antibodies is generally low, reviewed by Hartley et al. 55 . Glycosylation, sequence variation, masking by V1-V2, and the specific amino acid make-up of the loop may contribute to this, reviewed by Pantophlet et al. 56 . In addition, more complex causes of mutations, such as blockade of the accessibility of antibodies to the real epitopes or compensatory mutations to primary changes induced by antibodies, cannot be excluded. However, some monoclonal and polyclonal antibodies specific to epitopes within V3 have been demonstrated to neutralise diverse HIV-1 strains in vitro [57][58][59][60] . Two of the sites exhibit particularly limited amino acid diversity and as such may warrant further investigation as potential components of vaccines targeting shared epitopes of very low diversity within V3.
Indeed, despite evidence for significant antibody-mediated selection pressure, some sites were relatively conserved in terms of their composition, containing only a limited number of biophysically similar amino acids. This could be due to functional or structural constraints on the protein, and may reduce the ability of the virus to successfully escape antibodies targeting these regions. We also identified sites containing biophysically diverse amino acids that may be contained within epitopes eliciting effective antibody responses which cycle continuously between a limited number of biophysically distinct forms throughout chronic infection, as predicted by a previous modelling study 21 . Consistent with this model, we found evidence within individuals that some sites contain major variants that appear and disappear over time, such as proline to leucine in position 369 and arginine to serine in position 444. Sites under significant positive selection in at least 50% of study participants are shown in purple. Variable loops V3, V4 and V5 are contained within the three boxes, respectively. The dotted lines denote 50% of study participants. Antibody epitope clustering is shown in grey, whereby intensity denotes number of epitopes spanning that residue as reported in the LANL Immunology Database (http://www.hiv.lanl.gov/content/immunology). Sequences have been aligned to the HXB2 Gp120 reference sequence (accession number K03455), and position is relative to this alignment. (B) Homologymodelled structure of the SM cohort consensus Gp120 sequence in surface representation. Variable loops V3, V4 and V5 are shown in grey and sites under significant positive selection in 50% or more study participants are shown in purple. Structure has been modelled on a glycosylated HIV-1 Gp120 trimer (RCSB PDB 3J5M) 74  In summary, our data provide insight into how and where the surface of Gp120 is mutating over the course of clinically latent infection. Our detailed analysis of HIV-1 evolution and selection allowed us to identify amino acid sites under positive selection that was likely attributed to host factors (since the study participants were infected with a purportedly identical founder virus). In line with this, it is also possible that site-specific selection may be virus-specific. Importantly, this comprehensive mapping resulted in the identification of both previously described and novel constrained antibody and T-cell epitopes. These sites may be crucial to viral envelope function, and if verified in functional studies, provide suitable targets for future drugs and vaccines.

Methods
Cohort characteristics and sampling. The SM cohort comprises HIV-1 infected study participants from a small rural community in Henan province, China, as described previously 38 . Between 1993 and 1995, the participants were infected with a narrow-source HIV-1 subtype B virus during a blood plasma donation scheme. Few individuals knew that they were infected until HIV screening programmes were implemented in China in 2004. Study participants were then recruited to the cohort.
Study participants with samples collected longitudinally over the course of HIV-1 infection permit detailed evolutionary analysis and identification of sites under selection. Twelve study participants from the SM cohort  (Table 1). All study participants gave informed consent for their samples to be used for research purposes.
HIV-1 env gp120 sequencing and sequence assembly. Viral RNA was isolated from cryopreserved plasma samples and purified using the QIAamp Viral RNA Extraction Kit (Qiagen) followed by reverse transcription using the SuperScript III Reverse Transcriptase System (Invitrogen). The Env gp120 C2-V5 (approximately 799 bp, HXB2 [accession number: K03455] positions 6883-7681) was amplified by nested touchdown PCR (primers listed in Table S1). Amplified DNA was purified using the MinElute Gel Extraction Kit (Qiagen), and ligated into a pCR4-TOPO sequencing vector using the TOPO TA Cloning Kit for Sequencing (Invitrogen). Chemically competent One Shot MAX Efficiency DH5α E. coli (Invitrogen) were transformed with the prepared plasmids, and cultured overnight at 37 °C. Eighteen colonies were selected for colony PCR (M13F and M13R primers, Table S1), and the resulting products were purified using ExoSAP-IT (Affymetrix) and sequenced to generate forward and reverse reads (Source BioScience). Contigs were assembled and controlled manually using Geneious v9.0.5 61 (HXB2 gp120 positions 661-1455, accession number K03455). Sequences were multiple aligned using

Inference of infecting founder strain.
To infer the founder HIV-1 strain of the infected study participants, a consensus sequence was generated for each study participant for each time-point, with an ambiguity threshold of 10%. One sequence was selected per study participant to generate a dataset with sequences evenly distributed across the sampling period. The sequences were aligned and codon-stripped to a final alignment length of 759 nucleotides. Study participant SM007 was excluded from this analysis because it was not possible to conclusively rule out dual-or superinfection as preliminary data exploration demonstrated that sequences from this study participant did not resolve monophyletically.
Bayesian Markov Chain Monte Carlo phylogenetic inference -implemented through BEAST v1.8.2 (http:// beast.bio.ed.ac.uk/beast/) 64 -was used to estimate the time to the most recent common ancestor (tMRCA). Divergence time was estimated using the SRD06 model 65 , and an uncorrelated lognormal relaxed molecular clock 66 with a rate prior of 0.001 substitutions per site per year. The Markov chain was run for 100,000,000 generations to allow for adequate mixing, with posterior samples extracted every 10,000 generations. A burn-in period of 10% was applied, and convergence of posterior probabilities was assumed once the effective sample size (ESS) of each parameter exceeded 200, as determined in Tracer v.1.5 (http://tree.bio.ed.ac.uk/software/tracer/). Three runs were combined in LogCombiner v1.8.2 (http://beast.bio.ed.ac.uk/LogCombiner/) and the mean root height of the tree was calculated. Phylogenetic trees were annotated in FigTree v.1.4.2 (http://beast.bio.ed.ac.uk/figtree).
As the SM cohort study participants were infected with a narrow source of virus, the sequence of the MRCA was used as a surrogate for the infecting founder. The sequence of the reconstructed ancestor at the tree root from each run following burn-in was extracted, and a consensus sequence was generated from an alignment of these sequences (26,000 sequences) with an ambiguity threshold of 10%. Viral subtyping. An unambiguous consensus sequence was generated from the sequences of each time-point for each study participant. The sequences were then aligned to the LANL 2005 gp120 reference dataset (http:// www.hiv.lanl.gov/), and viral subtyping was performed by Bayesian phylogenetic inference (BEAST v1.8.2). The generalised time-reversible (GTR) nucleotide substitution model plus invariant sites and gamma-distributed rate heterogeneity (GTR + I + G) was used, with a constant size coalescent tree prior, estimated base frequencies, and a strict molecular clock. Following exclusion of burn-in, TreeAnnotator v1.8.2 was used to determine the maximum clade credibility tree (MCC).  This site has been described as a supersite of vulnerability for antibody neutralisation 83 . The site has also been associated with the loss of transmitted-founder sequence, suggested to represent antibody-driven selection (Hraber et al. 82

Variant characterisation.
A variant was defined as any amino acid in any position in the alignment that differed from that present in the inferred founder. Owing to the extensive degree of variation, the hypervariable loops were conservatively stripped from the alignment prior to analysis (final length 179 amino acids). Major variants were defined as variants found in greater than 15% of the sequences. Whilst major variants are canonically defined as those present at a frequency greater than 5%, this value was conservatively tripled as the amplicon was approximately three times more variable than the full-length HIV-1 genome 69 .
Structural modelling. Homology modelling was implemented through SWISS-MODEL [70][71][72][73] to map the translated SM cohort Gp120 consensus sequence to a cryo-electron microscopy (cryo-EM) crystal structure of a glycosylated HIV-1 envelope trimer (RCSB PDB 3J5M) 74 . The sites of interest were annotated on the modelled structure in PyMOL v1.7.4 (Schrödinger, LLC. Availability of Data and Materials. The datasets generated and analysed during the current study are available in the Genbank repository (accession numbers: MF078678-MF079252). Custom R scripts used in the analysis of these data are available from the authors on request.  Table 2. Antibody epitope sequences corresponding to the sites identified on the surface of the HIV-1 Envelope. Data made available by the Los Alamos National Laboratory (LANL) Immunology Database (http:// www.hiv.lanl.gov/content/immunology). Position is relative to HXB2 Gp120 (accession number K03455). Additional information about the neutralising antibody positions and associations can be found with references in the LANL HIV Genome Browser (https://www.hiv.lanl.gov/content/sequence/genome_browser/browser. html). 1 Name of monoclonal antibody or "polyclonal" (if a general response is being studied) as listed in the LANL Immunology database. 2 Y = Yes; N = No. 3 NA = Not available. 4 Novel site identified in the current study.