Structural investigations on mechanism of lapatinib resistance caused by HER-2 mutants

HER-2 belongs to the human epidermal growth factor receptor (HER) family. Via different signal transduction pathways, HER-2 regulates normal cell proliferation, survival, and differentiation. Recently, it was reported that MCF10A, BT474, and MDA-MB-231 cells bearing the HER2 K753E mutation were resistant to lapatinib. Present study revealed that HER-2 mutant K753E showed some contrasting behaviour as compared to wild, L768S and V773L HER-2 in complex with lapatinib while similar to previously known lapatinib resistant L755S HER-2 mutant. Lapatinib showed stable but reverse orientation in binding site of K753E and the highest binding energy among studied HER2-lapatinib complexes but slightly lesser than L755S mutant. Results indicate that K753E has similar profile as L755S mutant for lapatinib. The interacting residues were also found different from other three studied forms as revealed by free energy decomposition and ligplot analysis.


Introduction
The epidermal growth factor receptor (EGFR), HER-2, HER-3, and HER-4 form the human epidermal growth factor receptor (HER) family. Via different signal transduction pathways. HER family regulates normal cell survival, proliferation and differentiation [1]. The intracellular domain of HER-2 is composed of three parts which contains~500 residues, namely a cytoplasmic juxta membrane linker, a tyrosine kinase (TK) domain and a carboxyl-terminal tail [2,3].
In a recent study on HER-2-negative tumors, novel mutations L768S and V773L were detected, whereas in HER2-positive another novel mutation K753E was found. The cells bearing the mutations L768S and V773L exhibited more rapid growth while cells bearing the K753E mutation were resistant to lapatinib [21]. Extensive research has been conducted regarding the molecular mechanisms leading to trastuzumab and lapatinib resistance [18], including the overexpression or hyper expression of other HER family receptors and their ligands [21]. In recent years, molecular dynamics (MD) simulations based studies have demonstrated great potential as a substitute for investigations in understanding the dynamical aspects of protein-ligand interaction, protein folding, and also protein-protein interaction [22][23][24][25][26][27][28][29][30]. In the present study, to gain insights into the effects of mutations (K753E, L768S and V773L) on the binding strength of lapatinib, MD simulations of the kinase domain of HER-2 (wild and mutants) in complex with the lapatinib were performed. We quantitatively evaluated the structural differences in HER-2 mutants with respect to initial structures using the root mean square deviation (RMSD) of the backbone of a protein. RMSD of lapatinib was also calculated in binding with all HER-2 mutants. To study change in molecular size over MD simulation time, the gyration radius (Rg) of the studied HER-2 mutants in lapatinib binding state was calculated. The Rg value represents the mean squared distance of the atoms from the centre of mass of the molecule [31]. To determine which parts of HER-2 mutants in lapatinib binding state are influenced by lapatinib, we compared the fluctuation of residues of studied HER-2 mutants in the binding state. We also calculated the free energy between HER-2 forms and lapatinib which provides valuable insights into the recognition and binding of the drug and elucidate the effect of these mutations on drug recognition and binding. To further explore the conformational space of HER-2 mutants stabilized by the lapatinib, the probability-based free energy landscapes (FEL) as a function of the first two principal components (PC1 and PC2) were relatively investigated.

Molecular dynamics simulations
Classical MD simulations were performed for HER-2 wild and mutant forms with lapatinib using the GROMACS package [32,33]. The initial structure of HER-2 was taken from PDB ID:3PP0. Discovery studio was used to generate HER-2 mutants (Discovery studio client). To generate the topologies for proteins and for the simulations, the Gromos43a1 force field [34] was used. Energy minimization, position restrain dynamics and final production MD simulations for 50 ns were carried out as described in our previous research papers [25,35]. Chimera was used to visualize the complexes [36].

Binding free energy calculations
Binding free energy calculations were performed from the snapshots of MD trajectory using the molecular mechanics Poisson Boltzmann surface area (MM/PBSA) method [37]. The binding free energy of the HER-2-lapatinib complexes were analyzed by taking snapshot at an interval of 1.5 ps from 45 to 50 ns MD simulation during equilibrium phase, using g_mmpbsa tool of GROMACS [38].

Results and discussion
The RMSD profiles were found to be below 0.50 nm for wild and mutant HER-2 backbone during the entire simulation period indicating the suitability of all systems for further analyses ( Fig 1A). The RMSD of lapatinib bound to wild and mutant HER-2 was also calculated to analyse the fluctuation from initial position. The lapatinib bound to K753E showed the highest RMSD. In case of wild and L768S, RMSD of lapatinib was found very close to each other. V773L bound lapatinib showed the lowest RMSD ( Fig 1B). This revealed the highest stability of lapatinib in V773L binding site in initial orientation. RMSF profile revealed that the residue fluctuation in K753E was more as compared to other forms (Fig 2A) and similarly, K753E has higher radius of gyration ( Fig 2B).

Free energy calculations
The MM/PBSA calculations were performed for all the HER2-lapatinib complexes and binding free energy components ( Table 1). All mutants were found to be at highly favourable state than the wild type for lapatinib. The final binding energy for wild, K753E, L768S and V773L Molecular mechanism of lapatinib resistance caused by HER-2 mutants were -554.48, -811.89, -600.80 and -668.68 kJ/mol respectively. The electrostatic energy contribution in MM energy was found highest in case of K753E among the four forms. All mutant forms were found to have more vdW interactions as compared to wild (Table 1). These results indicate that studied mutations in HER2 increased the MM energy for lapatinib which is the most significant term in calculation of binding energy. To determine the contribution of each HER-2 residue which is involved in any type of interaction with inhibitor, we performed free energy decomposition and calculated per residues binding energy (Fig 3). The K753E form of HER-2 was found to have a very different profile of contributing residues as compared to rest of other three forms. The major contributing residues were Glu766, Asp769, Glu770, Asp871,  and Asp880. These residues were not involved in interaction with lapatinib in other HER-2 forms. In case of wild, L768S and V773L residues Asp808 and Glu812 were found to be the most significant contributor and were found absent in K753E. These results revealed that lapatinib binding orientation in K753E is different from rest of three forms.
To gain insight into the average binding pose of lapatinib, we extracted the average structure of HER-2-lapatinib complexes. We found that the orientation of lapatinib in K753E is reverse in the binding site while compared with other complexes (Figs 4-7). This suggests that in K753E, lapatinib binding is energetically more favourable in reverse orientation which leads to more binding energy. Ligplot analysis of average structure revealed that the solubilizing

Fig 5. Average structure of K753E HER-2-lapatinib complex (left top), surface representation of K753E HER2 binding site with lapatinib (left bottom), Interaction plot of K753E HER-2-lapatinib using average structure by LIGPLOT (right) (Dotted green and dotted red lines represent hydrophobic contacts and hydrogen bonds respectively).
https://doi.org/10.1371/journal.pone.0190942.g005  Molecular mechanism of lapatinib resistance caused by HER-2 mutants group interacts with residues Asp808, His809 and Glu812 were common residues in wild, L768S and V773L forms while in case of K753E Glu766, Asp880 and Arg868 were interacting residues (Figs 4-7). These results indicate the reverse orientation of lapatinib in K753E stabilized the complex in more favourable manner.

Free energy landscape (FEL)
To elucidate the sub conformational patterns of HER-2 forms and lapatinib complexes, we studied the FEL against first two principal components PC1 and PC2 which revealed ΔG value 0 to 18.70, 16.60 (Fig 8), 18.0 and 17.5 kJ/mol (Fig 9) for wild, K753E, L768S and V773L respectively. The size and shape of the minimal energy area (in blue) reveals the stability of a complex. Smaller and more concentrated blue areas suggest more stability of the corresponding complex [45]. All the HER-2 forms showed more centralized and concentrated minima except K753E. K753E HER-2 showed several minima which indicate that protein-lapatinib complex switch to various conformations to achieve minima. This is further support the change in orientation of lapatinib in binding site.
The cells bearing K753E mutation were found to be resistant to lapatinib. However the reason of resistance may be the overexpression of mutant HER-2. Same study have shown that the phosphorylation of downstream signalling protein PCLγ decreased as the concentration of lapatinib increased in MCF10A cell lines [21]. In previous studies it was reported that L755S mutation bearing cells showed resistance towards lapatinib might be due to overexpression of protein [46][47][48]. To compare with the L755S, we performed molecular dynamics simulation of HER-2 L755S-lapatinib complex and MM/PBSA calculations. Binding energy of HER-2 L755S-lapatinib complex was found -926.74 kJ/mol (Fig 10A) which is the highest in all studied complexes and close to K753E. Further free energy decomposition revealed that the contributing residues profile is very similar to the K753E. Glu766, Asp769, Glu770 and Asp871were the foremost contributing residues which were interacted with lapatinib ( Fig  10B). Average structure revealed the orientation of lapatinib is similar to K753E (Fig 10C). The FEL analysis of L755S complex showed several minima as found in case of K753E ( Fig  10D). Endogenous expression of L755S mutant HER-2 in models that do not have overexpression not sufficient to produce resistance to lapatinib due to lower levels of the mutant protein and show sensitivity towards lapatinib [49]. Zuo et al. [21] confirmed that the HER-2 K753E mutation may induce lapatinib resistance due to its proximity to the lapatinibresistant mutation L755S and cell lines with the HER-2 K753E mutation overexpression were certainly resistant to lapatinib. Our findings show consistency with the previous cell line studies and indicates that K753E lapatinib drug resistant mutant have more affinity towards lapatinib as compared to wild and other mutants so the loss of interaction might not be the reason behind the lapatinib resistance. Overexpression of mutant K753E HER-2 might be the cause of lapatinib resistance as found in case of L755S.

Conclusion
Present study revealed that mutant K753E was found to have some contrasting behaviour as compared to wild, L768S and V773L. Lapatinib showed stable reverse orientation in binding site of K753E which may be reason behind the highest binding energy among studied HER-2-lapatinib complexes but slightly lesser than L755S mutant. The interacting residues were also found different from other three studied forms as revealed by free energy decomposition and ligplot analysis. Results indicate that K753E has similar profile as L755S mutant for lapatinib. Overexpression of mutant K753E HER-2 might be the cause of lapatinib resistance as found in case of L755S.