Prediction of High Level of Multiple Drug Resistance Mutations in HIV-1 Subtype C Reverse Transcriptase Gene among First Line Antiretroviral-Experienced Virological Failure Patients from North India Using Genotypic and Docking Analysis

ISSN:2155-6113 JAR, an open access journal J AIDS Clinic Res Drug Resistance: HIV Prediction of High Level of Multiple Drug Resistance Mutations in HIV-1 Subtype C Reverse Transcriptase Gene among First Line AntiretroviralExperienced Virological Failure Patients from North India Using Genotypic and Docking Analysis Jaideep S Toor1, Rajender Kumar2, Prabha Garg2, Aman Sharma3 and Sunil K Arora1* 1Department of Immunopathology, Postgraduate Insitute of Medical Education & Research, Chandigarh, India 2Department of Pharmacoinformatics and Computer center, National Institute of Pharmaceutical Sciences, Mohali (Punjab), India 3Internal Medicine, Postgraduate Insitute of Medical Education & Research, Chandigarh, India


Introduction
Acquired Immunodeficiency Syndrome (AIDS), caused by the human immunodeficiency virus (HIV), is one of the leading causes of death with major medical and economic impact on the society. Due to high mutation rates associated with RNA replication and retrotranscription, there is spontaneous emergence of a pool of mutant viruses, some of which may be associated with drug resistance (DR) [1][2][3]. Keeping in mind the increasing demand of improving the drug efficacy, there is an urge to improve the methodology of identification of these mutations and know how they affect the susceptibility to a particular drug.
Subtype C of the HIV-1 is accountable for over 50% of the HIV-1 infections in Southeast Asia and African Countries while subtype B predominates in Western Europe and North America [4][5][6][7]. An increasing body of experimental evidence suggests that different HIV-1 subtypes exhibit disparate biological behaviors, and might respond differently to diagnostic, immunologic and therapeutic interventions [8,9]. Recent studies have identified subtype specific differences in HIV susceptibility to specific anti-retroviral drugs [10,11] and signature mutations selected by treatment [12][13][14]. An increased amount of resistance in subtype C HIV has been documented as compared to other subtypes in response to single-dose NVP therapy [15,16]. This indicates that subtype of HIV might have influence on the level of resistance development against a particular anti-retroviral drug under that drug pressure and deserves further attention to prove this fact. For genotypic analysis of isolates from infected individuals the sequences need to be analyzed using one of the few drug-resistance databases available which and binding energy calculations were performed. All drugs were built using SYBYL7.1 (SYBYL7.1. Tripos Inc., St. louis, MO 63144, USA.) and minimized using Tripos force field, GasteigerHuckel, partial atomic charges [20] and Powell's conjugate gradient method with energy gradient convergence criteria of 0.05 kcal/mol [21]. AutoDock Tool was used for preparation of all mutated modeled protein structures for docking. All water molecules were removed from protein structures. Also, polar hydrogen's were added and non polar hydrogen's were merged, finally atom type parameter and Kallman united atom charge was added. For docking study Lamarckian genetic search algorithm was employed and docking run was set to 30. All other parameters were set to default value. Energy calculations were performed by placing the receptor protein inside a user defined 3D grid and allowing various probes to systematically visit every grid point. Each probe consisted of an atom type found in the ligand whose dispersion/repulsion energies and, if appropriate, hydrogen-bonding energies were calculated for all receptor atoms within 8 Å of each grid point. Desolvation energy values were obtained by calculating burial of polar atoms, whilst conformation and torsion energies were calculated according to AMBER force field terms.

Baseline characteristics of the population
All subjects were from north India with age ranging from 7 to 60 years (Median age 32 years) and median CD4 count of 226 cells/µl (Interquartile range: 209.5 to 244.7 cells/µl). More than 92% individuals revealed a history of heterosexual transmission.
The patients who failed the first-line ART had received a combination of 2NRTI+NNRTI treatment for 3-4 years (median time of therapy 42 months; interquartile range 34.3 to 51.0months) and showed clinical signs of treatment failure with baseline median CD4 count of 84.5 cells/µl (interquartilerange 55.5 to 139.4 cells/µl) and viral loads ranging from 1.74E4 to 1.83E6 copies/ml (Table 2).

Pattern of mutations in HIV-1 subtype C reverse transcriptase gene
Therapy naïve: Mutations were defined as amino acid differences from the consensus subtype-B RT sequence (as per Stanford HIV RT sequence database) and polymorphisms were defined when they occurred in at least 1% of the sequences from untreated subjects. Among these, five mutations (121Y (56.3%), 173A (60.17%), 177E (83.18%), 207E (65.45%), 211K (66.36%)) occurred in more than 56% of the HIV-1 isolates and were considered as a local wild type for molecular modeling studies.
In the absence of a subtype-C based database, the amino acid sequence of each isolate in our study was compared to the subtype B consensus amino acid sequence (Stanford HIV RT sequence database) for mutations associated with resistance to RT inhibitors (http:// hivdb.stanford.edu). Of the 122 viral isolates from the therapy naïve interpretation of genotype variants in a subtype-C related infections with respect to their behavior in response to existing drug regimens, we have a proposed a computational model of subtype C specific wild type RT constructed after introducing the genetic variations. Mutated HIV-1 subtype C HIV-1 RT and inhibitor complexes were created after docking various nucleotide reverse transcriptase inhibitors (NRTIs) and binding energy was calculated for each. The data suggest that best drug combinations for a given set of mutations in subtype C infected individuals can be designed using the information generated from our proposed model, which was capable of calculating the combinatorial effect of all the mutations in that particular isolate, rather than effects in isolation for each mutation.

Patient selection and sample collection
For the analysis of Drug Resistance (DR) mutations and genetic diversity in the RT gene of HIV-1 subtype C strains in North India, 122 HIV-1 positive drug naïve individuals and 13 patients who had received first-line therapy (combination of 2 NRTIs and one NNRTI) for 4-5 years and shown clinical as well as virological signs of treatment failure (HIV-1 RNA >1000 copies/ml) were recruited. The study was approved by the Institutional Ethics Committee. After informed consent, 3 ml of venous blood was drawn in a K 3 EDTA vacutainer tube (Becton Dickinson, USA) for CD4 count, viral load and DR genotyping.

CD4 count and plasma viral load
The CD3/CD4+ cells were enumerated by a bead-based flow cytometry kit (Trucount TM , BD Biosciences Immunocytometry Systems, San Jose, CA, USA) and the plasma viral load was estimated using Cobas-Amplicor HIV-1 monitor test, version 1.5 (Roche, Branchburg, NJ) according to manufacturer's protocol.

Drug resistance genotyping
For DR genotyping, the viral RNA was isolated from plasma using QIAamp viral RNA isolation kit (QIAGEN, Valencia, CA) and reverse transcribed to make complementary DNA (cDNA) copies using a cDNA synthesis kit (MBI Fermentas, Vilnius, Lithuania). The RT region (68 to 231 codon) was amplified using RT specific nested primers as previously described [17]. The PCR products were sequenced in an automated DNA sequencer (ABI PRISM 3730 version 3.0, Applied Biosystems, Foster city, CA, USA). All sequences were subjected to quality assessment and determination of DR mutation profile using the Stanford DR Database (http://hivdb.stanford.edu). NCBI HIV subtyping tool using HIV-1 2009 reference set was used to identify the subtype of the virus.

Computational analysis of HIV-1 RT mutations
The computational analysis of the interactions of mutated and wild type RT gene with NRTIs was performed. In this study, the highly identical hits of Crystal structure of HIV-1 RTIs were recovered from Protein Data Bank (PDB). The PDB ID 1RTD for NRTI with tri-phosphate in the active site of RT with DNA template: primer complexwere taken [18]. Molecular models of the wild-type HIV-1 subtype-C RT with NRTI complexes were built from the template. The modeled proteins were optimized by using molecular dynamics (MD) approach together with simulated annealing (SA) methods in the Modeller program [19], and eventually validated on SAVES. The mutated models were created by including reported mutations and optimized. Modeled and mutated HIV-1 subtype C RT and inhibitor complexes were created by docking     HIV-1 individuals, eight showed mutations known to affect the drug susceptibility to RT inhibitors (Table 1). Among these, two isolates had G190A/V, a major NNRTI resistance mutation that is known to cause high level resistance to NVP, intermediate resistance to efavirenz (EFV), and low-level to etrinavir (ETR). Another important mutation, K103N/R detected in HIV isolate from another subject, is also known to cause high-level resistance to NVP (50-fold), DLV (50-fold), and EFV (25-fold). K103R occurred in about 1%-2% of isolates from untreated individuals and by itself has no effect on NNRTI susceptibility. However the combination of K103R with V179D, reduces the susceptibility to each of the NNRTIs about 15-fold. Other mutations detected were: for NNRTI: A98 (V/S) in 4 isolates, K101R in 1, V179D in 3, V179T/I in 6 and for NRTIs: K219Q in 1 isolate, V118I in 2, V75I in 1 and L210F in 1, which are known to be associated with very low-level resistance to certain RTIs.

Molecular modelling and docking studies of NRTI
Consensus sequence of RT gene obtained by sequencing 113 therapy naïve isolates was used for molecular modeling of wild type HIV-1 subtype C docked with NRTI using highly identical templates (PDB ID 1RTD). Two isolates from first-line therapy failure patients, SK-206 and SK-208, were selected for the modeling and docking study with difference of just one extra mutation, 75M in SK-206, to see the effect of this mutation in combination with other mutations (118I, 184V, 210W, 215Y) common in both (Figure 4) and compare it with interpretation from Stanford database. Both patients were on ART for 36 months and had high viral load and low CD4 count and revealed mutations affecting drug susceptibility to NRTI as per Stanford DR database ( Table 2). On verification and validation of all modeled proteins on SAVES (Structure Analysis and Verification Server) it was observed that all amino acids in the models were located in the same pattern as crystal structures in Ramachandran plots.  (Table 3). This decrease in the free energy is in concordance with Stanford DR database prediction. On the contrary, our results show the presence of V75M mutation along with other mutations (as found in isolate SK-208) in isolate from SK-206,seem to be partially neutralizing the effect of118I, 184V, 210W and 215Y mutations in terms of change in the binding energy for drugsD4T (∆G=-1.46 Kcal/mol), DDI (∆G=-1.52 Kcal/mol), and FTC (∆G=-0.57 Kcal/mol) where as in case of AZT there was further increase in binding energy (∆G=0.87) with the presence of same mutation (75M) which is likely to decrease the susceptibility and might cause high level of resistance (Table 3). On the other hand,in case of isolate SK-206 where Stanford DR database is predictinghigh to intermediate level resistance to D4T and DDI, our model predictsthis isolate to be susceptible to these inhibitors due to nullification effect of the drug resistance mutations by V75M, indicating the power of the model to take into consideration the combinatorial effect of all mutations instead of calculating the mutational effect in isolation for each mutation.

Discussion
Drug resistance is a major obstacle to achieve and maintain virus suppression. In this study, we have sequenced the amplified RT region of pol-gene of HIV isolates from the treatment-naïve as well as firstline ART failure patients in north India and identified drug resistance mutations as per Stanford DR database. There is very high prevalence of DR mutations in the RT gene in HIV-1 isolates from individuals who failed first-line therapy, re-emphasizing the importance of genotyping prior to the start of second line treatment for a rational selection of regimen to have better efficacy of the therapy.
Recent concerns have been raised regarding the first line of ART failures which comprised of RTI drugs. Due to the absence of PI in this regimen in India, a higher rate of DR mutations have been reported in the RT region as compared to PR region in therapy naïve individuals [22].
Among the 4 classes of antiretroviral drugs, the prevalence of resistance in antiretroviral drug-naive patients is often highest for NNRTIs [23,24]. In the United States, rates of primary NNRTI resistance have been observed to be as high as 16% in certain populations [23]. In our current study isolates from 8 out of 122 therapy naïve individuals had mutations that are known to affect susceptibility to RTIs. Amongst these individuals, 75% (n=6) had mutations affecting NNRTI susceptibility (Table 1). We found high occurrence (2.46%) of 179D mutations compared to the reported frequency of 1% in untreated individuals (http://hivdb.stanford.edu). Mutation V179D/E is known to cause low-level reduction in susceptibility to each of the NNRTIs. The combination of K103R and V179D is indicated to reduce the susceptibility of NVP, DLV, and EFV by 15 folds. The combined effect of this combination on ETR is not known. Other mutations found causing high level resistance to NNRTIs were 190V/A and K103N.
Consequently the DR genotyping results of isolates from all first line therapy failure patients in our study depicted high level of mutations which are known to cause resistance to RTIs.It was not surprising to find one or more DR mutations in RT region of all the isolates from 13 virological failure individuals ( Figure 1 and Figure 2) with G190A being most frequent (53.85%) followed by 103N/A, 98G, 181C (38.46% each) and 101E (30.77%) effecting NNRTI susceptibility while 184V being most frequent (92.31%) followed by 215Y(53.85%), 210W/S (46.15%) and 118I (30.77%) known to be affecting the susceptibility to NRTIs. All, except one isolates (92.3%) showed high levels of drug resistance to NRTI drugs 3TC and emtricitabine (FTC). Most isolates (46.15%) showed susceptibility to TDF (Figure 1), whereas all isolates (100%) showed mutations known likely to cause high level of drug resistance to NVP followed by mutations likely to affect susceptibility to EFV in 76.92%. Though most isolates had mutations associated with high level resistance to most of NNRTIs, etravirine (ETR) seems to be most acceptable drug among these according to our genotyping data ( Figure  2). With the rolling out of second line of therapy by Govt. of India, that comprises 2 or 3 NRTI drugs (Tenofovir (TDF), 3TC, AZT) along with one protease inhibitor (Lopinavir (LPV) or Ritonavir(RTV)), it becomes all the more important to carefully select the regimens, as it is observed that all the first line ART failures have not only RTI resistance mutations but also PI resistance mutations which can reduce the efficacy of the second line of treatment with faster emergence of PI drug resistance.
Subsequently we have proposed a new model based on the computational assessment of change in binding energy during docking of the drugs with mutated RT and finally interpreted the results in terms of prediction of best 2 nd line of therapy for these patients. In order to achieve this we generated a computational NRTI docks HIV-1 subtype C RT model from the consensus sequence of the wild type HIV-1 subtype C in this region, and further verified genotype findings on this model. The computation was verified/validated on SAVE and Ramachandran plot. The RT sequences from isolates from two 1 st line treatment failure patients were analyzed using the Stanford DR data base as well as our proposed model with constructed crystal RT structure, complexed with DNA template-RNA primer by adding mutations in the wild type subtype-C structure. On the basis of docking energy scores, the effects of mutations on the binding free energy have given very interesting results. In comparison of the binding energy scores obtained after docking of different RTIs with crystal structures of WT sequence or with mutated structure generated after the introduction of mutations detected in isolates SK-208 (V118I, M184V, L210W, T215Y) and SK-206 (V75M, V118I, M184V, L210W, T215Y), the analysis reveals that the presence ofV75M mutation partially neutralizes the resistance effect of mutations 118I, 184V, 210W and 215Y. This may be a combinatorial effect of all mutations occurring together in a particular sequence because of an increase in the number of hydrogen bonds with RTIcausing a change in conformation with presence of V75M along with other mutations, as seen in docking attempt of RTIs with structure generated from a sequence of isolated from SK-206 and SK-208, both of which have similar mutations except for the existence of V75M additionally in SK-206 isolate ( Table 3). Because of an increase in the number of hydrogen bonds with Lys65 due to presence of V75M mutation,it probably resulted in stabilizing the bound state of RT-NRTI, thus our model predicts patient SK-206 to be susceptible to d4T whereas the Stanford DR database predicts high level of resistance to