Strong Selection Significantly Increases Epistatic Interactions in the Long-Term Evolution of a Protein

Abstract

Epistatic interactions between residues determine a protein's adaptability and shape its evolutionary trajectory. When a protein experiences a changed environment, it is under strong selection to find a peak in the new fitness landscape. It has been shown that strong selection increases epistatic interactions as well as the ruggedness of the fitness landscape, but little is known about how the epistatic interactions change under selection in the long-term evolution of a protein. Here we analyze the evolution of epistasis in the protease of the human immunodeficiency virus type 1 (HIV-1) using protease sequences collected for almost a decade from both treated and untreated patients, to understand how epistasis changes and how those changes impact the long-term evolvability of a protein. We use an information-theoretic proxy for epistasis that quantifies the co-variation between sites, and show that positive information is a necessary (but not sufficient) condition that detects epistasis in most cases. We analyze the "fossils" of the evolutionary trajectories of the protein contained in the sequence data, and show that epistasis continues to enrich under strong selection, but not for proteins whose environment is unchanged. The increase in epistasis compensates for the information loss due to sequence variability brought about by treatment, and facilitates adaptation in the increasingly rugged fitness landscape of treatment. While epistasis is thought to enhance evolvability via valley-crossing early-on in adaptation, it can hinder adaptation later when the landscape has turned rugged. However, we find no evidence that the HIV-1 protease has reached its potential for evolution after 9 years of adapting to a drug environment that itself is constantly changing. We suggest that the mechanism of encoding new information into pairwise interactions is central to protein evolution not just in HIV-1 protease, but for any protein adapting to a changing environment.

Fig 1. Average per-site entropies at every position of the HIV-1 protease.

Untreated (top panel) and treated (bottom panel) datasets at the earliest (year 1998, red) and latest (year 2006, blue) time point of our analysis. 300 sequences are resampled from data for each year and average entropy for each position is calculated from the entropies in 10 resampled datasets. Site-specific variation generally increased across the protein following treatment. Entropy (variability) also increased from 1998 to 2006 for several positions. Error bars denote ±1 SE.

Fig 2. Changes in protease entropy and physico-chemical properties.

Changes in per-site entropies (top panel), residue isoelectric point (middle panel), and residue weights (bottom panel) due to treatment. The property difference at each site is obtained by subtracting property (entropy/pI/residue-weight) value of the untreated data from that of the treated data. Average values are obtained by sampling sequence data from all years (1998–2006, 10 subsamples/year of 300 sequences each). Error bars represent ±1 SE. Red dots represent positions known to be primary drug resistance loci, while black dots mark positions of compensatory or accessory mutations [52]. Resistance loci are shaded in red and accessory loci are shaded in black.

Fig 3. Estimates of the information content of the HIV-1 protease.

Filled black circles represent data from untreated subjects and blue triangles represent data from treated individuals. I₁ [see Eq (12)] is consistently low in treated sequence data over the years, indicating high sequence variability in the drug environment (top panel). The middle panel shows that the sum of pairwise mutual information significantly increases upon treatment (p ≤ 0.001). On adding the sum of pairwise mutual information to I₁, we obtain a comprehensive measure of information that considers pairwise interactions between residues [I₂, Eq (13)]. I₂ for both the treated and untreated data is comparable and unchanging over the years. We use data only for positions 15–90, as residues 1–14 as well as 91–99 have missing sequence data leading to error-prone estimates of entropy, as evidenced in S2 Fig. Error bars represent ±1 SD.

Fig 4. Increase in epistasis in HIV-1 protease over time.

Pairwise interactions in the HIV-1 protease are shown for years 1998, 2002, and 2006 in the drug-free (top row) and drug environment (bottom row). Each heatmap shows the mutual information for each pair of residues. Pairwise information (and thus epistatic effects) are fairly constant in the drug-free environment, but gradually increase in the treated group.

Fig 5. Epistatic interactions mapped onto the protease structure.

Epistatic interactions in the protease sequences in treated data from the year 1998 (red on the left chain) and 2006 (blue, right chain). The interacting residues are numbered. Only those interactions are shown where information is greater than 0.1 bits, indicating strong epistasis. Figure generated using Chimera [60].

Fig 6. Rates of mutation between four different genotypes.

The types are denoted as 0 = AA, 1 = Aa, 2 = aA, and 3 = aa.

Fig 7. Allele frequencies as fitness for type 3 (aa) is varied.

In this valley-crossing landscape, w₀ is always 1 and w₁ = w₂ = 0. Plot shows allele frequencies p_i at mutation rate μ = 0.1 as a function of w₃. The intermediate types aA and Aa occur only at the rate of mutation as they have zero fitness.

Fig 8. Two-loci two-allele model.

The left panel shows the fitness landscapes and epistasis given by Eq (9) in the first and second half of the simulation (updates 0–499: w₀ = 1 and w₁ = w₂ = w₃ = 10⁻⁵ ≈ 0; updates 500–1000: w₀ = w₃ = 1 and w₁ = w₂ = 10⁻⁵ ≈ 0). The xy-plane shows the four genotypes while the z-axis shows genotype fitness. The middle panel shows the genotype probabilities while the right panel shows the mutual information during the course of the simulation. Note that the increase in epistasis at the 500th update is reflected in the increase in mutual information. The mutation rate was 0.1 and starting population frequencies were p₀ = 1 and p₁ = p₂ = p₃ = 0.

Fig 9. Correlation between epistasis and information.

Each point corresponds to information and absolute value of epistasis calculated for one of the 10,000 combinations of w₀, w₁, w₂, and w₃. Here, w₀ is always 1, and other fitness values are uniformly randomly assigned between 0 and 1. The inset shows the percentage of points with negligible information (<0.001 bits) as a function of epistasis.