{'fraction_non_alphanumeric': 0.030040983606557378, 'fraction_numerical': 0.037090163934426226, 'mean_word_length': 5.219984705582462, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Constrained patterns of covariation and clustering of HIV-1 non-nucleoside reverse transcriptase inhibitor resistance mutations

Objectives We characterized pairwise and higher order patterns of non-nucleoside reverse transcriptase inhibitor (NNRTI)-selected mutations because multiple mutations are usually required for clinically significant resistance to second-generation NNRTIs. Patients and methods We analysed viruses from 13 039 individuals with sequences containing at least one of 52 published NNRTI-selected mutations, including 1133 viruses from individuals who received efavirenz but no other NNRTI and 1510 viruses from individuals who received nevirapine but no other NNRTI. Of the 17 reported etravirine resistance-associated mutations (RAMs), Y181C/I/V, L100I, K101P and M230L were considered major based on published in vitro susceptibility data. Results Efavirenz preferentially selected for 16 mutations, including L100I (14% versus 0.1%, P < 0.001), K101P (3.3% versus 0.4%, P < 0.001) and M230L (2.8% versus 1.3%, P = 0.004), whereas nevirapine preferentially selected for 12 mutations, including Y181C/I/V (48% versus 6.9%, P < 0.001). Twenty-nine pairs of NNRTI-selected mutations covaried significantly, including Y181C with seven other mutations (A98G, K101E/H, V108I, G190A/S and H221Y), L100I with K103N, and K101P with K103S. Two pairs (Y181C + V179F and Y181C + G190S) were predicted to confer >10-fold decreased etravirine susceptibility. Seventeen percent of sequences had three or more NNRTI-selected mutations, mostly in clusters of covarying mutations. Many clusters had Y181C plus a non-major etravirine RAM; few had more than one major etravirine RAM. Conclusions Although major etravirine RAMs rarely occur in combination, 2 of 29 pairs of covarying mutations were associated with >10-fold decreased etravirine susceptibility. Viruses with three or more NNRTI-selected mutations often contained Y181C in combination with one or more minor etravirine RAMs; however, phenotypic and clinical correlates for most of these higher order combinations have not been published.


Introduction
The non-nucleoside reverse transcriptase inhibitors (NNRTIs) bind to a hydrophobic pocket in the reverse transcriptase (RT) enzyme close to, but not contiguous with, the polymerase active site. These compounds inhibit HIV-1 replication allosterically by blocking translocation of the template-primer following nucleotide incorporation, or by displacing the catalytic aspartate residues relative to the polymerase-binding site. Mutations responsible for NNRTI resistance occur in this hydrophobic inhibitor-binding pocket. 1 A single mutation in this pocket can lead to high-level resistance to one or more first-generation NNRTIs, including efavirenz and nevirapine. However, two or more mutations are required to cause high-level resistance to etravirine, rilpivirine and other second-generation NNRTIs. 2 -5 Therefore, we decided to explore the landscape of NNRTI resistance mutations in combinations of two or more, as viruses with these combinations are likely to be most troublesome for NNRTIs with higher genetic barriers to resistance.
We used publicly available HIV-1 RT sequence data from the HIV Drug Resistance Database (HIVDB) 6 to identify patterns of mutations associated with efavirenz and nevirapine, the two most commonly used first-generation NNRTIs, and examined pairwise and higher order correlations among NNRTI-selected mutations. Our results provide new insights into the genetic mechanisms of NNRTI resistance and have important implications for the use of second-generation NNRTIs.

Virus RT sequences
HIV-1 Group M RT sequences were compiled from published studies in the HIVDB. These included a set of 47350 sequences from 39687 individuals and a set of 9349 sequences from a large reference laboratory. 7 To use no more than a single sequence per individual, we selected the sequence with the greatest number of NNRTI-selected mutations from each of the individuals in the first set and excluded all those (n¼2548) from the second set that differed from each other at ,5.0% of their nucleotide positions.
We studied only sequences obtained by direct PCR sequencing of plasma HIV-1 RNA because these are the sequences obtained in clinical settings. Sequences of molecular clones were excluded because of their higher probability of containing PCR errors. 8 Sequences obtained from other sources, particularly proviral peripheral blood mononuclear cell DNA, were excluded because they are at higher risk of containing defective genomes such as viruses with G to A hypermutation. 9,10 To assess the patterns in which NNRTI-selected mutations occur, we analysed only those sequences containing one or more NNRTI-selected mutations. This enabled us to study sequences obtained from individuals who are likely to be NNRTI experienced, regardless of the availability of a complete treatment history or the possibility of non-adherence.

Covariation analysis
We used the Jaccard similarity coefficient (J) to identify pairwise correlations among the 52 NNRTI mutations in 13039 sequences from unique individuals. Covariation analysis was run on the full set of sequences and on the subsets corresponding to sequences from individuals with a history of receiving either efavirenz or nevirapine but no other NNRTI.
We created 2×2 contingency tables to indicate the number of times a sequence contained both mutations, only one of the two mutations or neither of the two mutations. Overall, there were 1288 contingency tables representing all possible pairs of mutations except for permutations of mutation pairs at the same position: 1326 (52 choose 2) pairs minus 38 pairs of mutations at the same position. For each comparison, sequences containing an electrophoretic mixture at either of the two positions were excluded because of the inability to determine whether the two NNRTI-selected mutations were present in the same genome.
For a pair of mutations X and Y, the Jaccard similarity coefficient is calculated as J ¼N XY /(N XY +N X0 +N 0Y ) where N XY represents the number of sequences containing mutation X and mutation Y, N X0 represents the number of sequences containing X but not Y, and N 0Y represents the sequences containing Y but not X. This statistic therefore examines only the sequences containing both X and Y in the context of sequences containing at least one of the mutations, thereby eliminating the heavy weighting given to the double-negative category in such statistical analyses as the x 2 test. This is particularly important so as not to exaggerate the significance of pairs of rare mutations, where the double-negative category is very large. 15 To test whether the observed Jaccard similarity coefficients were statistically significant, we calculated the expected value of the Jaccard similarity coefficient (J RAND ) for each pair of mutations as if the two mutations (X and Y) occurred independently. J RAND was calculated as the mean Jaccard similarity coefficient after 2000 random rearrangements of the X or Y vector (containing 0 or 1 for the presence or absence of a mutation, respectively). J SE , the standard error of J, was calculated using a jackknifed procedure, which removed one sequence at a time, repeatedly for each sequence. 15 The standardized score Z, where Z¼(J2J RAND )/J SE , indicates a significant positive association when Z. 1.96 or a significant negative association when Z is less than 21.96, at an unadjusted P, 0.05. The case of N XY ¼0 results in a Z-score of negative infinity; Fisher's exact one-tailed test, using the R command fisher.test with a left-tailed alternative hypothesis, was used to approximate the P value for those comparisons with N XY ¼0.
Holm's correction was used to control the familywise error rate for multiple pairwise comparisons. 16 The P values of the n pairwise comparisons are ranked in ascending order, where p(1) is the smallest (ranked first) and p(i) is the ith ranked. Starting from p(1), one compares each p(i) with p/(n2i+1) where p is the false discovery rate (in this case, 0.05) and n is the total number of pairwise comparisons, until p(j). p/(n2j+1). All correlations between pairs for p(1),. . ., p(j 21) are considered statistically significant.
We modelled the distributions of mutations in each sequence using Poissonness plots, following the method of Hoaglin, 17 and conducted a parametric bootstrap using a multivariate Poisson simulation to assess the distribution of positive and negative correlations obtained from the Jaccard covariation analysis. This step was important in evaluating the possible biological significance of the negative correlations.
We also examined the extent to which the 52 NNRTI-selected mutations covaried with mutations at 12 major nucleoside reverse transcriptase inhibitor (NRTI) resistance positions, including M41L, K65R, D67N, T69S_SS, K70E, K70R, L74V, L74I, V75M/T, Y115F, Q151M, M184V/I, L210W, T215F and T215Y. 18 Sequences that contained at least one NNRTI-selected mutation and at least one of the above NRTIselected mutations were analysed using the same methods described for assessing covariation within the NNRTI mutation dataset.

Multidimensional scaling
We performed multidimensional scaling on the pairwise association data using a matrix D of dissimilarity coefficients (J D ¼12 J, where J is the Jaccard similarity coefficient) for the 22   . For a list of mutations (X 1 , X 2 , . . ., X n ), multidimensional scaling constructs points in two-dimensional (2-D) space such that the Euclidean distances between these points approximate the entries in the dissimilarity matrix. 19 For a given matrix, it computes points X 1 , is the Euclidean distance between X i and X j , and d ij is the dissimilarity between X i and X j in the matrix D. This was performed using the R function cmdscale (classical multidimensional scaling).

Cluster analysis
We also identified all sequences containing three or more NNRTI-selected mutations and characterized the patterns of mutations by their frequency in the dataset and by whether they were composed of mutations that covaried with one another. Sets of three or more mutations in which each mutation covaried with each of the other members in the set were referred to as clusters or cliques (a term derived from graph theory). A lower threshold for covariation was used for this analysis (uncorrected P,0.05) than in the pairwise comparison, to facilitate the identification of clusters. As for all previous covariation analyses, we used no more than one sequence per patient.

Covariation of NNRTI-selected mutations
To examine the extent of covariation among the NNRTI-selected mutations, we performed separate Jaccard analyses of the complete set of sequences from 13039 individuals with viruses containing at least one NNRTI-selected mutation, the 1133 sequences from individuals who had received only efavirenz and the 1510 sequences from individuals who had received only nevirapine. Among the 13039 sequences that comprised the full dataset, 66 of the 1288 possible mutation pairs had positive covariances, and 29 of these covariances were statistically significant after the Holm's adjustment for multiple comparisons (P,0.05). Table 1 shows the results of this analysis on the 29 pairs of significantly covarying mutations. G190A (n ¼9), Y181C (n¼ 8), V108I (n¼ 6), K101E (n¼ 4), H221Y (n ¼4) and F227L (n¼ 4) were the most common mutations among these 29 pairs. Tables S1 to S4 [available as Supplementary data at JAC Online (http://jac.oxfordjournals.org/)] show the results of the analysis on all 1288 mutation pairs for the complete dataset, the efavirenz-alone dataset and the nevirapine-alone dataset.
Of the 29 pairs of mutations that were significantly correlated in the full dataset, 8 pairs covaried (uncorrected P,0.05) in the efavirenz-alone subset, 5 pairs covaried in the nevirapine-alone subset and 4 pairs covaried in both subsets (Table 1). In 13 of the 29 significantly covarying pairs, both mutations were preferentially selected by the same NNRTI.
The major etravirine RAM Y181C was significantly associated with G190A, H221Y, V108I, A98G, G190S and V179F. L100I was significantly associated with K103N, and K101P was significantly associated with K103S. M230L was most strongly associated with K103N (P ¼ 0.003 uncorrected, not significant following correction for multiple comparisons).
We also found 225 negatively correlated mutation pairs, including 116 that had statistically significant negative correlations after adjusting for multiple comparisons. However, a parametric bootstrap using the multivariate Poisson analysis of a random dataset containing the same number of mutations per sample and the same prevalence of each mutation indicates that a large number of negative correlations is to be expected Constrained NNRTI mutational patterns 1479 JAC because our analysis included only isolates containing at least one NNRTI resistance mutation.

Multidimensional scaling
We used the 22 mutations that contribute to ≥1 of the 29 significantly covarying pairs in the full dataset for multidimensional scaling. Figure 2a, b and c shows the mutations on axes representing the first versus second, first versus third and second versus third principal coordinates. The first principal coordinate accounts for 36% of the total inertia and separates the predominantly nevirapine-selected mutations A98G, K101E, V108I, Y181C, G190A and H221Y from the other NNRTI resistance mutations. The second principal coordinate accounts for 25% of the total inertia and separates out V106A and F227L. The third principal component accounts for 15% of the total inertia and separates K103N and its associated mutations (L100I, P225H and K238T) from G190S and Y188L. The remaining principal coordinates did  not significantly improve the goodness of fit of the multidimensional scaling and were therefore not plotted. Although V179F covaried with Y181C, it is not shown in Figure 2 because it occurred too uncommonly to be accurately placed by multidimensional scaling.

NNRTI-NRTI correlations
To assess the covariation of NNRTI-and NRTI-selected mutations, we performed a Jaccard analysis on sequences from 9564 individuals with viruses containing at least one mutation from each category. After adjusting for multiple comparisons at P,0.05, we identified 39 positive and 40 negative statistically significant correlations between 52 NNRTI-selected mutations and the 17 major NRTI-selected mutations. The strongest positive correlations were between the NRTI mutations L74V/I and the NNRTI mutations Y181C and L100I. Another notable, though weaker, positive correlation exists between K65R and Y181C. The NNRTI-associated mutation A98G was also strongly positively correlated (corrected P,0.01) with the thymidine analogue mutations M41L, D67N, L210W and T215F/Y. The full results of this analysis are available in Tables S1 to S4 [available as Supplementary data at JAC Online (http:// jac.oxfordjournals.org/)].

Higher order clustering of NNRTI mutations
Approximately 17% (n¼ 2233) of RT sequences with one or more NNRTI-selected mutation contained three or more NNRTIselected mutations. Table 2 shows the 20 most commonly Indicates whether the pair of mutations also covaried significantly (uncorrected P,0.05) within the subsets of sequences from individuals experienced in efavirenz (EFV) or nevirapine (NVP). b Occurrences in 13039 sequences comprising the full dataset; only one isolate per patient was included in the analysis and sequences with mixtures at a given position were excluded from comparisons at that position. XY, isolates containing both mutations; X only, isolates containing mutation X but not mutation Y; Y only, isolates containing mutation Y but not mutation X. c P values are uncorrected values from the Jaccard covariation analysis. All mutation pairs shown above have correlations that are significant, with a Holm's correction for multiple pairwise comparisons at P,0.05.
Constrained NNRTI mutational patterns occurring mutational triplets, of which 14 meet the definition of a mutational cluster, or clique (each mutation covaries with each other mutation in the triplet). The remaining six triplets contain K103N in combination with a pair of covarying mutations. The 20 triplets listed in Table 2 occurred in 48% (n¼1069) of sequences with three or more mutations. In addition to the 14 mutational clusters included in the 20 most commonly occurring triplets, we identified 43 less   Table 1) in sequences containing at least one NNRTI-selected mutation: (a) compares the first and second principal coordinates; (b) compares the first and third principal coordinates; (c) compares the second and third principal coordinates. The distance between any two mutations is measured by their Jaccard dissimilarity coefficient, J D , where J D is equal to 1 minus the Jaccard similarity coefficient. The R command cmdscale was used to compare the first three principal coordinates from a table of J D s for each pairwise comparison. V179F is not included in this graphic; despite a significant positive correlation with Y181C, the comparison produces a high Jaccard dissimilarity coefficient, causing a misleading placement in multidimensional scaling.

Discussion
Prior to the licensing of etravirine in 2008, only a few studies examined the patterns of NNRTI resistance mutations arising following NNRTI treatment failure because failure with one NNRTI generally presaged failure with a second NNRTI. 20 -24 Since the licensing of etravirine, however, several studies have examined the likelihood of virological response to etravirine among individuals in whom previous NNRTI therapies were unsuccessful. 4,25 -29 Moreover, the introduction of etravirine has led to the identification of novel NNRTI-selected mutations that often decrease susceptibility to efavirenz and nevirapine as well as to etravirine or other second-generation NNRTIs. 2,30,31 In this study, we sought to identify and quantify associations among NNRTI-selected mutations. Our analysis differs from previous covariation analyses 15,22 in that we examined an expanded set of mutations-including recently identified non-polymorphic NNRTI-selected mutations 10 and the complete list of etravirine RAMs-in a very large number of sequences. In addition, we sought to control for associations that might reflect similar treatment exposure rather than a potential biophysical mechanism. Toward this end, we investigated only sequences containing one or more NNRTI-selected mutations, to minimize the likelihood that a history of NNRTI exposure (NNRTI naive versus NNRTI experienced) was responsible for the association. We then systematically examined the preferential selection of NNRTI resistance mutations by efavirenz or nevirapine and the covariation of mutations within subsets of efavirenz-or nevirapine-experienced individuals, to assess whether pairwise covariations could reflect the selection pressure of the NNRTI received.
Although the viruses that emerge during treatment with efavirenz or nevirapine are often resistant to both drugs, our analysis shows that efavirenz and nevirapine select for NNRTI resistance mutations in very different proportions. For example, nevirapine selects for V106A, Y181C/I/V and G190A significantly more frequently than efavirenz, whereas efavirenz preferentially selects for L100I, K101P, V106M, Y188L, G190S and P225H. These data suggest that patients developing virological failure while receiving nevirapine may be at increased risk of developing virological failure with etravirine because Y181C/I/V-which occurs in 45% of individuals receiving nevirapine-forms the foundation for high-level etravirine resistance. 4,12,13 Conversely, although the two efavirenz-selected mutations L100I and K101P are major etravirine RAMs, combined they occur in ,20% of patients receiving efavirenz. Nevirapine-associated etravirine RAMs, including Y181C, G190A and E138A, were also among the most commonly observed mutations in sequences from 782 women treated with a single dose of nevirapine.
With the exception of the uncommon mutations Y181I/V, multiple mutations are required for moderate to high-level etravirine resistance (.10-fold reduction in susceptibility). 4, 26 We identified 29 pairs of significantly covarying NNRTI-selected mutations in viruses from 13 000 individuals whose sequences contained at least one NNRTI-selected mutation. Of these 29 pairs, 22 contained at least one etravirine RAM. Nine of these 22 pairs contained major etravirine RAMs: Y181C (n¼ 7), L100I (n¼ 1) or K101P (n¼ 1), but only 2 pairs are predicted to confer moderate to high-level etravirine resistance: Y181C+ V179F and Y181C+ G190S. Although V179F and G190S do not reduce etravirine susceptibility when they occur individually, they are highly synergistic with Y181C; Y181C + V179F reduces etravirine susceptibility .100-fold and Y181C +G190S reduces etravirine susceptibility 20-fold. 2,4,12,13 Approximately 17% (n¼ 2233) of individuals with one or more NNRTI-selected mutations had viruses with at least three NNRTIselected mutations. Nearly half of these viruses contained at least 1 of 20 common patterns of mutation triplets, suggesting that these patterns of mutations are particularly compatible with enzymatic function and drug resistance. Most of these patterns consisted solely of mutations that covaried with one another (mutational clusters), and the remainder consisted of K103N in combination with one or more pairs of covarying mutations. Y181C was the only major etravirine RAM among the 20 most common mutational triplets. However, eight of these triplets contained 2 of the 17 etravirine RAMs, and six triplets had 3 etravirine RAMs. According to published phenotypic data, only 1 of the 20 triplets is predicted to have high-level resistance to etravirine-K101E, Y181C and G190S-but intermediate resistance to etravirine has previously been associated In a previous study, Rhee et al. 15 found 11 significant pairwise correlations between NNRTI-selected mutations (P,0.05 corrected), of which two pairs, K103N/Y181C and K103N/V108I, were not positively correlated in our analysis. This is almost certainly because we included only sequences containing one or more NNRTI-selected mutation, reducing the likelihood that the co-occurrence of mutations was confounded by a history of NNRTI exposure. The following nine pairs were still significantly correlated in our analysis: Y181C/G190A, K101E/G190A, Y181C/ H221Y, V108I/H221Y, L100I/K103N, V108I/Y181C, V106A/ F227L, K101E/Y181C and K103N/P225H. Likewise, Ceccherini-Silberstein et al. 22 reported five significant associations among the NNRTI-selected mutations that we studied: Y181C/H221Y, V108I/H221Y, L100I/K103N, K103N/P225H and K103N/K238T. Although we found that four of the pairs described by these studies consist of mutations selected preferentially by nevirapine (K101E, Y181C, G190A and H221Y) and that K103N and its associated mutations are all selected preferentially by efavirenz, we also found that the covariance of every pair was independent of the particular NNRTI received, suggesting a potential biophysical interaction rather than treatment pressure.
In conclusion, viruses from patients with virological failure following nevirapine treatment are more likely to have major etravirine RAMs than viruses from patients who failed treatment with efavirenz. This is because Y181C/I/V occur more commonly with nevirapine than with efavirenz, and because mutations at this position provide the foundation for high-level etravirine resistance. Few patterns among covarying pairs of mutations and the most common mutational triplets are predicted to be associated with moderate to high-level etravirine resistance. However, most patterns contained one or more minor etravirine RAMs in combination with a major etravirine RAM, particularly Y181C. These relatively constrained higher order mutational patterns constitute a potential target for future studies of the phenotypic and clinical significance of etravirine RAMs.

Funding
This work was supported by the National Institutes of Health (AI068581 and AI46148 to E. C. R, S. Y. R. and R. W. S., and R01GM086884 to S. H.).

Transparency declarations
None to declare.