Quantum Chemical Study of Molecular Recognition between Etravirine/Rilpivirine and the HIV-1 Reverse Transcriptase

Second generation non-nucleoside reverse transcriptase inhibitors (NNRTIs) etravirine and rilpivirine are essential components in the highly active antiretroviral therapy for the treatment of patients infected with human immunodeficiency virus type 1 (HIV-1). They are highly potent drugs against wild-type viruses and have exhibited excellent antiviral activities against some NNRTIs-resistant HIV-1 variants. In order to understand the underlying mechanism behind their robust resistance profile in comparison with the first generation NNRTIs nevirapine and efavirenz, it is necessary to quantitatively analyze their binding pockets in the wild-type HIV-1 reverse transcriptase (RT) and various HIV-1 RT mutants at the molecular level. Therefore, a high-level ab initio quantum chemical analysis was performed to decipher the molecular determinants for recognition of etravirine and rilpivirine by the wild-type RT and some RT mutants (K103N, K103N/Y181C, and K103N/L100I) of clinically important virus strains. Pair wise intermolecular interaction analysis determined the contribution of individual intermolecular interactions to the binding affinities between the second generation NNRTIs (etravirine or rilpivirine) and several variants of RTs, including the wild-type RT, and clinically relevant K103N, K103N/Y181C, and K103N/L100I mutant RTs. This quantitative analysis led to the identification of drug-protein interactions that persist despite mutations as well as to the evaluation of stabilization energy losses upon mutations. The results of this study enhanced our understanding of the molecular level mechanisms by which the second generation NNRTI drugs maintain their strong binding to mutant RTs. It is hoped that findings of this work would have a direct impact on designing new NNRTIs that are even more resilient to mutations in future.


Introduction
Human immunodeficiency virus type 1 reverse transcriptase (HIV-1 RT) represents an attractive target for the development of anti-AIDS therapies. This enzyme plays a pivotal role in the HIV replication cycle by converting single-stranded HIV RNA genome into double-stranded HIV-1 DNA that is subsequently integrated into the host cell's chromosome [1,2]. The current inhibitors of HIV-1 RT can be categorized into two main classes: nucleoside reverse transcriptase inhibitors (NRTIs) and non-nucleoside reverse transcriptase inhibitors (NNRTIs). The NRTIs are analogous of normal nucleotides but they lack the 3' hydroxyl group. They are competitive substrate inhibitors that incorporate Like the vast majority of drugs, NNRTIs bind to HIV-1 RT via non covalent interactions (e.g. hydrogen bonding, π-π stacking interactions, CH-π interactions and cation-π interactions). The binding pose of an inhibitor into its receptor, is determined by the collection of all the attractive and repulsive forces between the inhibitor and its surrounding interaction partners. Therefore, to gain an in-depth insight of the molecular recognition between a NNRTI and its receptor, it is of great importance to quantify the specific contributions of each amino acid in the binding pocket to the binding affinity between the NNRTI and its receptor. Deciphering the energetic contribution of the individual residues in the ligand binding pocket toward protein-NNRTI binding would have a great influence on designing new inhibitors. Binding affinity data determined experimentally represents the collective contributions of all protein residues that form the ligand binding pocket. Quantum chemical calculations, on the other hand, offer a way to rigorously quantify the magnitude and characterize the nature of individual intermolecular interactions between ligands and their surrounding protein residues.
The aim of the present study is to determine the underlying reason behind the origin of high binding affinities of second generation NNRTIs etravirine and rilpivirine to the wild-type HIV-1 RT and some clinically relevant mutant variants of the virus at the molecular level. For this purpose, data mining and high level quantum chemical calculations were performed to investigate the molecular recognition between etravirine/ rilpivirine and wild-type RT and mutated RTs. The energetic contributions of the individual intermolecular interaction between etravirine or rilpivirine and their surrounding amino acids to the binding affinities between these inhibitors and the various RTs were systematically analyzed in a pair wise manner. The results of this study will not only shed light on the molecular determinants for the molecular recognition of NNRTIs in HIV-1 RT, but also enhance our understanding of molecular recognition of drugs in protein in general. Therefore, this work is expected to have a direct impact on designing the next generation of NNRTIs with a higher potency and even more resilience to mutations

Data Mining
The availability of crystal structures of etravirine and rilpivirine bound to the wild-type HIV-1 RT and some clinically important single and double point mutants from Protein Data Bank (PDB) forms the basis of our quantum chemical analysis present here. Table 1 lists the six PDB files of NNRTIs (etravirine and rilpivirine) bound to RTs that were selected and analyzed in this study. For all the calculations considered in this study, the atomic coordinates of non-hydrogen atoms in NNRTIs (etravirine and rilpivirine) and their interacting protein residues were extracted from the X-ray crystal structures of the complexes. The missing hydrogen atoms were added and their positions were determined by ab initio geometry optimization at the HF/6-31+G* level. During the geometry optimization processes, the positions of heavy (non-hydrogen) atoms were fixed. All geometry optimization calculations were carried out using the Gaussian09 software [23].

Intermolecular Interaction Energies in the Gas Phase
For a given level of theoretical method, the quality of quantum mechanical calculation depends on the level of theory used and the choice of the basis set. For an accurate description of the intermolecular interactions between ligands and their surrounding protein residues, the electron correlation must be taken into consideration. Intermolecular interactions between neutral molecules are usually dominated by dispersion interactions [24,25] which is caused by the mutual correlation of electrons that belong to the interacting partners. The magnitude of correlation energy can be as high as the magnitude of the intermolecular interaction itself.
Madridge J Nov Drug Res. ISSN: 2641-5232 Therefore, the inclusion of electron correlation is of a great importance for a proper description of non-bonded interactions in biological systems. The Møller-Plesset perturbation theory at the second order (MP2) represents a feasible and popular method that recovers a significant portion of electron correlations [26]. In our study, pair wise intermolecular interaction energies, in the gas phase, between etravirine or rilpivirine and their surrounding protein residues were calculated by ab initio electronic structure calculations at the MP2 level combined with the triple-zeta cc-pVTZ Dunning's correlation basis set for an adequate treatment of electron correlation.
The interaction energies between inhibitors and their interacting protein partners were estimated at the MP2/cc-pVTZ level using the super molecular approach in which the interaction energy between the two molecules (e.g. A and B) is computed as the potential energy difference between the energy of the interacting dimer E AB and the energies of the monomers E A and E B : The basis set superposition error (BSSE) was corrected by the Boys and Bernardi Counter Poise Method [27]. All calculations were carried out using the Gaussian09 software [23].

Free Energies of Solvation
Since all biological processes take place in an aqueous environment in living organisms, the solvation effect has to be taken into account in order to obtain a reliable estimate of binding strengths of the intermolecular interactions in biological systems. A number of different solvation models have been developed to deal with this issue. Depending on the way solvent molecules are treated, these models can be divided into two main classes: implicit solvation models and explicit solvent models [28,29]. In explicit solvation models, all of the interactions in the system are fully described at the atomic level, which means that all the solvent molecules are explicitly incorporated in quantum mechanics or molecular mechanics calculations. Unfortunately, explicit treatment of solvent molecules when simulating biological systems at the quantum chemical level is not feasible due to the prohibitively high computational cost. Instead, a common way to calculate the free energy of solvation when dealing with biological systems is using the implicit solvation models. In these models, the solvent is described as a continuum, a polarizable medium with a fixed dielectric constant equals to the bulk value for pure solvent. Here, we employed the implicit solvation model SM5.42R developed by Truhlar and his coworkers [30,31] at the HF/6-31+G* level as implemented in the GAMESS-version 11 Apr. 2008 (R1) software package [32] to quantify the free energy of solvation ΔG sol of the different species considered in this study.
The standard-state free energy of solvation ΔG sol in the SM5.42R model is composed of two components: the electrostatic term and the first-solvation-shell (nonelectrostatic) term. When the solute is solvated, the charge distribution of the solute polarizes the surrounding homogenous dielectric medium of solvent and subsequently, the polarized field of the solvent acts back on the solute electric charge distribution to be self-consistent with the solvent electric polarization. The electrostatic term includes solute-solvent favorable interactions, solvent rearrangement cost as well as the energy cost resulting from the distortion of solute charge distribution. The first-solvation-shell term includes contributions arising from cavity creation, dispersion and change in the solvent structure. Before the protein-NNRTI complex is formed in solution, both the protein and NNRTI are solvated and have to lose part of their solvation shells upon binding which incurs an energy cost commonly referred to as dehydration energy. For a general reaction of A combines with B to form the AB complex, the dehydration energy for the complex formation is given by the following equation: where, ΔG i sol = AB, A, B represents the free energies of solvation for the complex AB, and the monomers A, B, respectively.

Intermolecular Interaction Energies in Solution
The binding energy for complex formation in solution was indirectly obtained as a sum of intermolecular interaction energies in the gas phase (Equation 1) and the dehydration energies for the complex formation in solution (Equation. 2) as follows: It is worth noting that a similar scheme was used to calculate solution phase binding affinities for ligand-protein complexes previously [33,34,35].

Results and Discussion
Since the main objective of this study is to identify molecular determinants for molecular recognition of second generation NNRTIs etravirine and rilpivirine in the wild-type RT and its mutants, we carefully examined crystal structures of etravirine or rilpivirine bound to the wild-type HIV-1 RT and mutant strains that are listed in Table 1. Despite the extensive experimental efforts on studying the mechanism behind the persistent potency of the secondary generation NNRTIs etravirine and rilpivirine against wild-type RT and mutated viruses, quantum chemical analysis is a very useful tool to quantify the contributions of the individual amino acids in the non-nucleoside inhibitor binding pocket (NNIBP) to the overall binding affinities between these drugs and RTs. The binding environments of both drugs in the wild-type RT and mutant strains were systematically examined using the Visual Molecular Dynamic (VMD) program [36] that is capable of three-dimensional stereographic display. Based on the 3D Madridge J Nov Drug Res. ISSN: 2641-5232 analysis of the crystal structures of RT···NNRTI complexes, amino acids residues of RTs that form the binding pockets of etravirine and rilpivirine were systematically identified. Subsequently, gas phase intermolecular interaction energies between these inhibitors and their interacting amino acid residues were calculated by means of the super molecular approach at the MP2/cc-pVTZ level, followed by solvation energy correction based on the implicit solvation model SM5.42R at the HF/6-31+G* level.

Binding of Etravirine/Rilpivirine to the Wild-Type HIV-1 RT
Due to the similar structural frames of the second generation NNRTIs etravirine and rilpivirine, both inhibitors adopt, to a large extent, very similar binding modes to RTs. Their interacting protein residues are located in the p66 subunit with the exception of Glu138 which belongs to the p51 subunit. The following amino acids Gly190, Val106, Val179, Leu100, Leu234, Pro236, Phe227, Tyr181, Tyr188, Tyr318, Trp229, His235, Lys101, Lys103, and Glu138 are located within 4.8 Å of both inhibitors and have the potential, based on 3-D analysis of structures, to be involved in favorable interactions with the inhibitors. The only two noticeable differences are that, etravirine seems to have two additional favorable interactions with Val189 and a water molecule (residue number 470 in the p51 subunit). The latter appears to bridge etravirine to RT enzyme via a mediated hydrogen bond. In addition, rilpivirine appears to be interacting with Pro95. Figure 1 shows stereographic displays of the binding pockets of etravirine and rilpivirine in the wild-type HIV-1 RT. Hydrogen bonds are shown as red dotted lines. Based on the above identification of intermolecular contacts, a high level quantum chemical analysis was employed to quantify the magnitude of intermolecular interaction energies between etravirine and rilpivirine and their interacting protein residues. Listed in Table 2 and Table 3 are the pair wise intermolecular interaction energies between etravirine/rilpivirine and the protein residues calculated at the MP2/cc-pVTZ level in the gas phase (ΔE MP2 ) and after solvation correction (ΔE sol ). Intermolecular interaction energy calculations after taking solute-solvent interaction into account yields the complete picture of interactions in a cellular environment. Therefore, the interaction energies after solvation correction presented in this study represents the complete picture of the binding strength of intermolecular interactions between these drugs and their interacting residues of RTs. The specific type of non covalent intermolecular interactions between etravirine/ rilpivirine and amino acids of RT are identified and shown in parentheses. structure of the complex (3MEE, [21]). This plot is generated with the program PyMOL [37].
As stated in the method section, prior to complex formation, both the HIV-1 RT and NNRTI are solvated and as they approach each other both molecules lose parts of their solvation shell upon binding to each other, which incurs an energy cost. Therefore, the interactions energies after solvation correction is usually lower than gas phase interaction energies especially for hydrogen bonding interactions, electrostatic interactions (salt-bridge) and cation-π interactions. The strengths of hydrophobic types of interactions such as π-π stacking interactions and CH-π interactions are less affected by solvation correction. As shown in Table 2 and Table 3, the magnitudes of the intermolecular interaction energies between etravirine/rilpivirine and RT after solvation correction are substantially large. This, to a large extent, can be attributed to the substantial contributions from π-π interactions and CH-π interactions between etravirine or rilpivirine and the aromatic and aliphatic residues of RT. Lys101 forms strong hydrogen bonding interactions with both inhibitors. Lys103 interacts strongly with etravirine and rilpivirine in the wild-type RT via CH-π and appears to be also forming cation-π interaction with rilpivirine.

Binding of Etravirine/Rilpivirine to the K103 Mutant
The K103N mutation is the most frequent mutation observed in viruses isolated from patients treated with NNRTIs [38,39]. It has been found that the K103N single point mutation is sufficient to confer high level of resistance to the first generation NNRTIs evirapine and efavirenz, reducing their susceptibility by 46 and 19-fold, respectively [40,41]. The second-generation NNRTIs etravirine and rilpivirine, on the other hand, show good profiles of activity against the single point K103N mutant, which allow them to retain potency against this mutant virus [42,43]. To understand the molecular basis of high binding affinities between etravirine/rilpivirine and K103N RT, we have quantitatively analyzed binding pockets of both drugs in the K103N RT mutant.
Following the same procedure as before, the intermolecular interactions between etravirine/rilpivirine and residues of K103N RT were energetically analyzed. Figure 2 Table 2 and Table 3 list the resulting intermolecular interaction energies in the gas phase (ΔE MP2 ) and after solvation correction (ΔE sol ) between etravirine and rilpivirine and residues of K103N mutant. Table 2 also presents a pair wise comparison of the intermolecular interaction energies between etravirine and residues of K103N mutant RT with those between etravirine and the wild-type RT residues. It can be seen that most residues in the K103N RT mutant maintain strong intermolecular interaction energies with etravirine in magnitudes that are roughly comparable to their binding strength with etravirine in the with wild-type RT. However, mutation of Lys103 to an asparagine results in a reduction of the solution phase interaction energies between etravirine and residue #103 from -3.32 kcal/molto -1.53 kcal/mol. Although the hydrogen bond interaction between etravirine and Lys101 is lost upon the K103N mutation, etravirine interacts favorably with Asn103 via NH-π and CH-π interactions. This single point mutation does not have any other noticeable effect on the intermolecular interactions between etravirine and the remaining amino acids in the binding pocket of etravirine. The favorable intermolecular interactions between etravirine and its binding pocket residues in the K103N mutant add up to a total of -35.12 kcal/mol stabilization energy. This strong binding affinity between etravirine and the K103N mutant is comparable to the binding affinity between etravirine and the wild-type RT.
Intermolecular interaction energies between rilpivirine and its binding pocket residues in the K103N mutant of HIV-1 RT are shown in Table 3, in comparison with those between rilpivirine and the wild-type RT residues. Mutation of lysine residue at position 103 to an asparagine resulted in a reduction of interaction energies from -5.90 kcal/mol to -0.39 kcal/mol, which represents a stabilization energy loss of -5.51 kcal/mol. Additionally, there are two significant stabilization energy losses at residue Lys101 (-5.59 kcal/mol stabilization energy loss) and at residue Val179 (-3.01 kcal/mol stabilization energy loss). Based on a careful examination of the optimized structures, it appears that Lys101 in wild-type RT interacts with rilpivirine via two hydrogen bonds and strong cation-π interactions meanwhile there is only one hydrogen bond formed between the carbonyl group of Lys101 and its adjacent NH group of rilpivirine in the K103N mutant of RT. In addition, there is change in the position of amino group of the Lys101 residue upon K103N mutation, which might have weakened cation-π interactions between Lys101 and rilpivirine. It was also observed that, the side chain of Val179 change its position with respect to the rilpivirin ephenyl ring nearby. One of the CH 3 groups of Val179 was pointing to the centroid of the rilpivirine aromatic ring in the wild-type RT, which results in strong CH-π interaction. Upon K103N mutation, the side chain of Val179 rotates so that CH 3 is no longer pointing toward the centroid of the aromatic ring of rilpivirine nearby. However, despite these stabilization energy losses, the remaining residues of the rilpivirine binding pocket in K103N RT maintain strong intermolecular interactions with rilpivirine. The total stabilization energy between rilpivirine and K103N RT residues is -25.46 kcal/mol. Despite K103N mutation, both etravirine and rilpivirine maintain strong interactions with a large number of RT residues, which enables these inhibitors to overcome drug resistance and retain potency against the K103N mutant.

Binding of Rilpivirine to the K103N/Y181C double Mutant
The double drug resistance mutations K103N and Y181C are among the most commonly observed mutations in patients treated with NNRTIs. The double mutation K103N/Y181C severely impaired the effectiveness of the first generation NNRTIs nevirapine and efavirenz. It has been found that K103N/Y181C double-point mutation reduces nevirapine and efavirenz potencies by 625 and 1176 fold, respectively [22]. Second generation NNRTIs etravirine and rilpivirine potently inhibit K103N/Y181C double mutant with EC 50 values less than one nanomolar [44]. Only the X-ray crystal structure of the K103N/Y181C RT double mutant bound to rilpivirine was available from the Protein Data Bank (PDB). Therefore, we were able to conduct a quantitative analysis of the intermolecular interactions between rilpivirine and its binding pocket residues in the K103N/Y181C RT double mutant. Following the same computational procedure as described earlier, intermolecular interactions between rilpivirine and the residues of K103N/Y181C RT were systematically analyzed. Figure 3 shows the stereographic drawing of the rilpivirine binding pocket in the K103N/Y181C mutant. Potential hydrogen bonds are shown as red dotted lines. structure of the complex (3BGR, [27] ).This plot is generated with the program PyMOL [37].
Listed in Table 4 are pair wise intermolecular interaction energies in gas phase (ΔE MP2 ) and after solvation correction (ΔE sol ) between rilpivirine and its interacting residues in the K103N/Y181C double mutant of HIV-1 RT. For the purpose of comparison, the interaction energies between rilpivirine and amino acids of the wild-type are also listed in Table 4. Comparing the intermolecular interaction energies between rilpivirine and the K103N/Y181C mutant with those between rilpivirine and the wild-type RT, we can see that the doublepoint mutation K103N/Y181C severely impaired the intermolecular interactions between rilpivirine and enzyme residues at position 101 and 103, in a similar manner as in the case of the K103N mutant. The pair wise intermolecular interaction energy, in solution, between rilpivirine and Lys101 changes from a favorable -3.29 kcal/mol in wild-type RT to an unfavorable 1.83 kcal/mol in the K103N/Y181C mutant. This is a -5.12 kcal/mol stabilization energy loss that can be attributed to the loss of hydrogen bonding interaction and cation-π interaction at residue 101 as a result of the K103N/Y181C mutation. Another great loss of intermolecular interaction energy is found at residue 103. Upon mutation of Lys103 to an asparagine, the interaction energy, in solution, changes from a favorable -5.90 kcal/mol for rilpivirine-Lys103 (in wildtype) to an unfavorable 1.21 kcal/mol for rilpivirine-Asn103 in the double mutant RT. This is a -7.11 kcal/mol stabilization energy loss that can be attributed to the loss of cation-π, CH-π interactions at residue 103 as a result of the K103N/ Y181C mutations. Additionally, it was found that Asn103 forms a repulsive interaction with rilpivirine after solvation correction. By examining intermolecular interactions between rilpivirine and the other residues in the binding pocket, we found that rilpivirin retains strong interactions with residue 181 after the Y181C mutation. Some other points worth mentioning here are that, interaction energy after solvation correction for the rilpivirine-Phe227 pair changes from an unfavorable 1.28 kcal/molto a favorable -1.65 kcal/mol as a result of the K103N/Y181C mutation. In addition, there is an enhancement of interaction energies at residue W229. As can be seen in Table 4, the interaction energy between rilpivirine and W229 changes from -2.61 kcal/mol to -4.27 kcal/mol upon the double mutation. It has been also found that the binding of rilpivirine to the K103N/Y181C mutant is enhanced by forming two additional attractive interactions with residues Tyr183 and a water molecule (residue number 780 in the p51 subunit). The total stabilization energies gained from these two new interactions are -1.39 and -2.52 kcal/mol, respectively. The total stabilization energy for the formation of drugprotein complex between rilpivirine and the K103N/Y181C double mutant, calculated as the summation of favorable intermolecular interactions after solvation correction, is -34.07 kcal/mol. Tyr181 (π-π, CH-π) -6.06 -3.98 Cys181 (CH-π, S-π) -4.02 -3.38 -0.60 Tyr188 (π-π, CH-π) -7. 29 -5.43 Tyr188 (π-π, CH-π) -5.34 -3.33 -2.10 Tyr318 (π-π) -3.47 -1.95 Tyr318 (π-π) -3.24 -2.01 0.06 Trp229 (π-π, CH-π) -5.51 -2.61 Trp229 (π-π, CH-π) -6.

Binding of Rilpivirine to the K103N/L100I double Mutant
The K103N/L100I double mutation is another common combination of mutations that has a severe impact on the effectiveness of the first generation NNRTIs nevirapine and efavirenz. Although, the K103N/L100I double-point mutation has a larger effect on the potency of rilpivirine than K103N and K103N/Y181C mutations, it has been found that second generation NNRTI rilpivirine retains excellent activity against the K103N/L100I double mutant virus with an EC 50 value of ~ 8nM [22]. The availability of X-ray structure of rilpivirine bound to the K103N/L100I mutant of HIV-1 RT enabled us to quantify the magnitude of intermolecular interaction energies between rilpivirine and its binding pocket residues in the K103N/L100I double mutant. Following a similar manner for analyzing drug-protein interactions, the intermolecular interactions between rilpivirine and residues of K103N/L100I HIV-1 RT were systematically analyzed. Figure 4 shows a stereographic drawing of the rilpivirine binding pocket in the K103N/L100I mutant. Hydrogen bond between rilpivirine and Lys101 is shown as red dotted line. Table 5 lists the intermolecular interaction energies in gas phase (ΔE MP2 ) and after solvation correction (ΔE sol ) between rilpivirine and its interacting residues in the K103N/L100I mutant. For an easy comparison of intermolecular interaction energies between rilpivirine and the K103N/L100I double mutant residues with those between rilpivirine and residues of the wild-type RT, the intermolecular interaction energies between rilpivirine and its binding pockets residues in both RTs are listed side by side. It was found that the attractive interaction between rilpivirine and Lys101 via hydrogen and cation-π interaction is lost upon the K103N/L100I mutation. After solvation correction the intermolecular interaction energy between rilpivirine and Lys101 changes from a favorable -3.29 kcal/mol (attraction) to an unfavorable 2.38 kcal/mol (repulsion). structure of the complex (2ZE2, [22]). This plot is generated with the program PyMOL [37]. 85 Tyr181 (π-π, CH-π) -6.06 -3.98 Tyr181 (π-π, CH-π) -6.75 -4.73 2.77 Tyr188 (π-π, CH-π) -7.29 -5.43 Tyr188 (π-π, CH-π) -6.73 -4.71 1.3 Tyr318 (π-π) -3.47 -1.95 Tyr318 (π-π) -2.36 -0.80 0.41 Trp229 (π-π, CH-π) -5.51 -2.61 Trp229 (π-π) -5. Another significant stabilization energy loss occurred at residue 103. The interaction energy changes from a favorable -3.29 kcal/mol between rilpivirine and Lys103 to an unfavorable 7.51 kcal/mol between rilpivirine and Asn103. This can be attributed to the loss of cation-π and CH-π interactions at residue 103 as a consequence of the K103N/L100I double mutation. A new strong favorable interaction is formed between rilpivirine and Phe227 in the K103N/L100I mutant. Mutation at residue 100 in which the leucine residue is replaced by an isoleucine has a subtle effect on the magnitude of the interaction energy between rilpivirine and residue 100. Both Leu100 in the wild-type and Ile100 in the K103N/L100I mutant form strong CH-π interactions with rilpivirine. The magnitudes of these intermolecular interaction energies, in solution, are-3.50 kcal/ mol and -2.04 kcal/mol, respectively. Most of the remaining amino acids in the binding pocket of rilpivirine in the K103N/ L100I mutant form favorable interactions that are comparable to those between rilpivirine and wild-type RT residues. The total stabilization energy between rilpivirine and residues of the K103N/L100I double mutant RT add up to -25.07 kcal/mol. This strong stabilization energy explains the good activity of rilpivirine against the K103N/L100I mutant of RT.

Conclusion
In summary, a high level quantum chemical analysis was carried out to analyze molecular determinants for recognition between the second generation NNRTIs etravirine and rilpivirine in HIV-1 RTs. To understand the origin of the potency of these drugs against wild-type and mutant RTs, the contributions of the individual amino acid residues in the NNRTI-binding pocket to the overall binding forces between etravirine/rilpivirine and both the wild-type RT and some clinically important mutants RTs were quantum chemically evaluated. In addition to the large flexibility of these drugs that allows them to adopt different conformations in response to the different mutations, non-bonded interaction involving the aromatic moieties of these drugs are critical for their excellent activity against the wild-type RT and its mutant strains. The pair wise intermolecular interaction calculations conducted in our study have deciphered the molecular level drug-protein interactions that help these drugs to overcome drug resistance. It has been discovered that, unlike the first generation of NNRTIs nevirapine and efavirenz that rely heavily for their binding to RT on their interactions with Tyr181 and Tyr188, both etravirine and rilpivirine attain their tight binding to RT and its mutants by forming strong favorable interactions with a wide range of their binding pockets residues including the highly conserved residues Trp229, Leu234, Phe227 and Tyr318. Therefore, any stabilization energy losses resulting from single or double point mutations of the non-nucleoside inhibitor binding pocket (NNIBP) have less effects on etravirine and rilpivirine potencies. The results of our quantum chemical calculations presented in this study, enhanced our understanding of molecular recognition between the second generations NNRTI etravirine/rilpivirine and RTs, unravelling at the molecular level the mechanisms that used by these drugs to achieve their high binding affinities to the clinically relevant mutants of RTs. It is hoped that the insights gained thorough this rigorous comparative theoretical investigation of molecular determinants for binding of NNRTIs with their targeted wild type and mutant RTs will lay the foundation for designing the next generation of NNRTIs that are even more resilient to mutations.