Computational Biology and Chemistry

Title: CoMFA, CoMSIA, Topomer CoMFA, HQSAR, Molecular Docking and Molecular Dynamics Simulations Study of Triazine Morpholino Derivatives as mTOR Inhibitors for the Treatment of Breast Cancer

Authors: Dhara M. Chhatbar, Udit J. Chaube, Vivek K. Vyas, Hardik G. Bhatt

CoMFA, CoMSIA, Topomer CoMFA, HQSAR, Molecular Docking and Molecular Dynamics Simulations Study of Triazine Morpholino Derivatives as mTOR Inhibitors for the Treatment of Breast Cancer
Dhara M. Chhatbar#, Udit J. Chaube, Vivek K. Vyas, Hardik G. Bhatt*

Department of Pharmaceutical Chemistry, Institute of Pharmacy, Nirma University, Ahmedabad 382 481. India.
# Current affiliation: Torrent Research Center, Bhat, Gandhinagar. India.

*Address for correspondence: Institute of Pharmacy, Nirma University, S.G. Highway, Chharodi, Ahmedabad 382 481. India.
Phone no.+91 79 30642727; Fax no. +91 2717 241916; Email Id: [email protected]

Graphical abstract

Highlights
mTOR has become a promising target for the treatment of Breast Cancer
Triazine Morpholino derivatives as mTOR inhibitors were selected for 3D-QSAR study
3D-QSAR contour maps were validated by MD/MS assisted molecular docking study
Contour maps and molecular docking were used for the designing of mTOR inhibitors

Abstract
mTOR has become a promising target for many types of cancer like breast, lung and renal cell carcinoma. CoMFA, CoMSIA, Topomer CoMFA and HQSAR were performed on the series of 39 triazine morpholino derivatives. CoMFA analysis showed q2 value of 0.735, r2cv value of 0.722 and r2pred value of 0.769. CoMSIA analysis (SEHD) showed q2 value of 0.761, r2cv value of 0.775 and r2pred value of 0.651. Topomer CoMFA analysis showed q2 value of 0.693, r2 (conventional correlation coefficient) value of 0.940 and r2pred value of 0.720. HQSAR analysis showed q2, r2 and r2pred values of 0.694, 0.920 and 0.750, respectively. HQSAR analysis with the

combination of atomic number (A), bond type (B) and atomic connections showed q2 and r2 values of

0.655 and 0.891, respectively. Contour maps from all studies provided significant insights. Molecular docking studies with molecular dynamics simulations were carried out on the highly potent compound 36. Furthermore, four acridine derivatives were designed and docking results of these designed compounds showed the same interactions as that of the standard PI-103 which proved the efficiency of 3D-QSAR and MD/MS study. In future, this study might be useful prior to synthesis for the designing of novel mTOR inhibitors.

Keywords

mTOR; CoMFA; CoMSIA; HQSAR; Topomer CoMFA; Molecular Dynamics Simulation

1. Introduction

Cancer is one of the dangerous diseases, characterized by uncontrolled and undesirable cell division. Worldwide deaths due to cancer are increasing year by year. Breast cancer is the most common malignant tumor and the most common cause of cancer deaths in the female population.[1,2] Among the thousands of targets available for the treatment of breast cancer, mTOR (Mammalian Target of Rapamycin) is one of the promising target as the PI3K/AKT/mTOR signaling pathway get deregulated predominantly in breast cancer.[3-6] mTOR is present downstream of the PI3K/AKT/mTOR pathway and plays a decisive role in cancer therapy.[7-9] mTOR integrates intracellular and extracellular signals through this pathway and regulates cellular growth, proliferation, survival and protein synthesis. mTOR consists of two different multi-molecular complexes viz. mTORC1 and mTORC2. The difference between these two complexes is only of raptor and rictor, respectively.[10,11] Most of the first generation mTOR inhibitors like Everolimus, Sirolimus, and Ridaforolimus acts on TORC1 as shown in Fig.1(A), while recently developed second generation selective

mTOR inhibitors and dual PI3K/mTOR inhibitors act on both mTORC1 and mTORC2 as shown in Fig.1(B).[12-14] These second generation mTOR inhibitors have the added advantage of lower molecular weight as compared to the first generation and showed a binding affinity for ATP binding site of mTORC1 and mTORC2. Hence, in the present research work, the second-generation ATP competitive mTOR inhibitors triazine morpholino heterocyclic derivatives were selected for 3D-QSAR, HQSAR, Topomer CoMFA (Comparative Molecular Field Analysis). MD/MS (Molecular Dynamics Simulations) and molecular docking studies were successfully employed for the assistance of the designing of new small molecule drug candidate.[15,16]
The aim of the present research work is to design novel mTOR inhibitors for the treatment of breast cancer. CoMFA model was generated based upon the steric and electrostatic field while, on the other hand, 32 different combinations of steric, electrostatic, donor, acceptor, hydrophobic were optimized for the generation of the best CoMSIA (Comparative Molecular Similarity Index Analysis) model. Different contour maps were obtained from the generated CoMFA and CoMSIA models, which were further utilized for the identification of favored and disfavored regions of the most active compound of the series. HQSAR analysis suggested thousands of fragments and these fragments were segregated in the form of positive and negative contributions which gave insight about the important fragments responsible for the activity of the molecules. Further the results of the contour maps were validated by molecular docking study. Additionally, for the accurate molecular docking study, molecular dynamics simulations of mTOR protein along with the docked ligand were carried out which showed stable complex and accurate docking results.

2. Materials and Methods

2.1. Data Set

The triazine morpholino derivatives investigated in this report were synthesized by Venkatesan et al.[17-19] Among the 62 reported compounds, 23 were omitted due to their

low and intermediate potency or selectivity. The remaining 39 molecules were divided into the training set (31 compounds) and test set (8 compounds) for determination of CoMFA, CoMSIA, Topomer CoMFA, and HQSAR models. Training and test sets were made in such way that the test set was truly representative of the training set and both datasets cover whole activity range from 0.3 nM to 22.5 nM.[20-22] The structures of all compounds were constructed from the template molecule using SKETCH function in Sybyl X.[23] Partial atomic charges were calculated by Gasteiger-Huckel method. Tripose force field was
employed for energy minimization of these molecules.[24, 25] Structures of all compounds with their biological activity viz. IC50 and pIC50 are listed in Table 1.
2.2. Molecular Alignment

Molecular alignment of the dataset molecules is an important prerequisite of good 3D- QSAR models. Importantly statistical results of the generated CoMFA and CoMSIA models depends upon molecular alignment [26-28]. Usually, the molecular alignment brings the proper orientation and conformation of all the molecules in the space which consequently gives the good 3D-QSAR model. Hence in the present study, two different alignments methods viz. distill- based and docking-based were carried out to find out the best method for the generation of the 3D-QSAR model.
2.2.1. Alignment-1 (Distill – Rigid Body Alignment)

This alignment was carried out by the distill module implemented in Sybyl X in which common core scaffold was identified among all the 39 molecules which is shown in Fig. 2(A). Showing the highest activity, compound 36 of the series with the least conformational flexibility was chosen as the template molecule and the remaining molecules of the dataset were aligned on this highest active compound based upon the common substructure similarity shown in Fig. 2(B).
2.2.2. Alignment-2 (Docking Based Alignment)

This alignment was carried out by the Surflex-Dock module of Sybyl X, in which protein of
mTOR (PDB ID: 4JT6) was selected from RCSB protein data bank. Binding pocket was generated by the protomol generation technique, and all the molecules were docked in the active binding pocket as shown in Fig. 2(C).

2.3. CoMFA and CoMSIA Field Descriptors

In CoMFA, steric and electrostatic fields were calculated by Lennard-Jones (eq. 1) and Columbic potentials (eq. 2), respectively. 3D-lattice was generated around the series of the aligned molecules and steric and electrostatic fields were calculated using the following equations.[29]
where Evdw is van der Waals energy, rij is distance between atom i of the molecule and the grid point j where probe atom is located, Aij and Cij are constants, EC is the coulomb energy, qi is the partial charge of the atom of the molecule, qj is charge of the probe atom, D is the dielectric constant.
The reliability of the CoMFA model depends upon the spatial orientations of the aligned molecules in the grid box. Therefore, all orientation search (AOS) procedures optimize the field by rotating the molecular aggregate and select orientation that produces the best CoMFA model. The interaction energy of steric and electrostatic fields were calculated using a sp3 carbon atom and for probes +1 charge applied for steric and electrostatic fields, respectively. In CoMSIA analysis, descriptors like steric, electrostatic, hydrophobic, hydrogen bond donor, acceptor fields were obtained. This analysis includes a similar type of grid box which was generated earlier in CoMFA analysis. This grid box having different lattice points with grid spacing of 2 Å and
each point was considered as an sp3 hybridized carbon.

2.4. Topomer CoMFA

Topomer CoMFA is a combination of steric, electrostatic and topological properties of the molecules, hence its statistical results were independent of an alignment of the molecules. In this technique, topomer alignment represented 2D topology of the molecules which further utilized for the 3D-molecular subunit generation. In this methodology, molecules were broken into two fragments viz. R1 indicated with blue colored as while R2 indicated with red colored as shown in Fig. 3. Absolute fragmentation of any model in topomer CoMFA was obtained through the overlapping of the structural fragments on the common feature which may be in the form of open valence or attached bond. For each fragment, one standard topomer 3D model was constructed.
Steric and electrostatic fields were calculated for each set of topomer at a regular space grid of 2 Å, and were fixed into 1000 point cube. For steric field contribution, sp3 hybridized carbon atom was used, and for electrostatic field, oxygen atom was used. [30-32]
2.5. Hologram QSAR (HQSAR)

In HQSAR analysis, all possible fragments were generated with the help of training set molecules and these fragmented structures encoded through the hashing algorithm in the form of bins of the hologram. The hologram generated bins then were correlated with biological activity to generate the predictive HQSAR models. Usually, for the generation of the hologram; atomic number (A), bond type (B), atomic connection (C), hydrogen bond donor/ acceptor (DA), hydrogen (H) and chirality (C) were considered. In the present research work, three fragments viz. atomic number, bond type, and bond connections were selected for the HQSAR model generation. For the best model generation, all important structural features in the dataset were
being analyzed. Input molecules were broken into fragments. The value of minimum (M) and

maximum (N) number of fragments were defined appropriately. Length of the hologram was kept in the range of 50-500. To reduce fragment collisions, the value of L was selected as prime numbers. The default values of L were 53, 59, 61, 71, 83, 97, 151, 199, 257, 307, 353 and 401.
PLS analysis is a tool useful to correlate the hologram descriptors and biological activity. The obtained results could be employed for the designing of a new drug candidate with improved biological activity.[33, 34]

2.6. Partial Least Square (PLS) Statistical Analysis

CoMFA, CoMSIA and HQSAR models were derived using PLS regression analysis. These calculated CoMFA, CoMSIA and HQSAR descriptors were used as independent variables and mTOR inhibitory activity (pIC50) as the dependent variable in the PLS analysis. PLS analysis was performed and different statistical parameters viz. q2, r2, r2 , and r2 were determined.
The q2 value (Leave one out cross-validation correlation coefficient) was determined for the internal validation of the CoMFA and CoMSIA model. The value of q2 gives an idea about the predictive power of the CoMFA and CoMSIA model. This can be done with the aid of the reduced dataset of the molecules which are undertaken for the 3D-QSAR study. During this analysis, each time, one molecule from the data set was excluded and the new model would be generated. This newly generated model was used for prediction of the activity of the excluded molecule. This procedure continued for n number of times until all compounds have been excluded. The outcome of this procedure is termed as q2 value, which provides information about the robustness and predictive ability of the CoMFA and CoMSIA model.
The r2 value (person’s correlation coefficient) gives an idea about the relationship between the predicted activity and actual activity. The r2 value can be obtained by a plotting a graph of actual activity vs predicted activity.
The r2pre (r2predicted) was used for the CoMFA and CoMSIA model validation, as this value indicate8d

how well the model predicted the activity of external test set molecules. Generally during this process, training set molecules were used for the model building and this generated model was used for the prediction of the activity of the external test set molecules. This r2pre is very much similar to the q2 value. A 100-cylce bootstrap analysis (r2boot) was performed to assess the quality of generated model. The quality of 3D-QSAR model mainly depends upon the LOO cross-validated coefficient (q2) and predicted correlation coefficient (r2pred) which must be greater than 0.5 and can be obtained by equation 3 and equation 4, respectively.[35]
where Ya is an actual value of the activity, Yp is a predictive value of the activity, Yc and Ym are the average values of the observed activities, PRESS is a predictive sum of squares, and SD is a sum of squared deviations.

2.7. Molecular Docking and Molecular Dynamics Simulations

To validate the contour maps of CoMFA and CoMSIA, the molecular docking study between mTOR protein (PDB ID: 4JT6) and most active compound 36 of the series were carried out. The SurFlex-Dock module of Sybyl X was used for molecular docking studies. The X-ray crystallographic structure of mTOR (PDB ID: 4JT6) was retrieved from the RCSB Protein Data
Bank and modified for docking calculations.[36] For initializing molecular docking, water

molecules were removed from crystal structures then after polar hydrogen atoms were added, side chains were fixed during protein preparation and, at the end, Gasteiger-Huckel charge were assigned to protein atoms. Usually, molecular docking study can be done by applying three modes viz. automatic mode, ligand mode and residual mode. During the Automatic mode, docking usually acquired in the largest binding cavity of the protein. In case of ligand mode, docking usually occurred in the space vacated by the active co-crystallized ligand of the protein while on the other hand in case of residual mode, docking could be done on the specifically selected residue of the protein. For the present work, the ligand-based mode was selected where active co-crystallized ligand of mTOR protein (4JT6) was extracted, and protomol was generated which was considered as a binding pocket. This protomol is usually considered as vacant space
generated due to the extraction of the standard co-crystallized ligand of the protein where putative ligand could be placed. Protomol is the representation of an idealized co- crystallized ligand already present in the crystal structure of the protein. The docked ligand usually shows the docking score based upon the hydrogen bond interactions and orientation of the ligand in the binding pocket of the protein. Molecular dynamics simulations of the protein structure, in complex
with compound 36, was carried out to verify the stability of the docked pose of the compound 36 with the help of Sybyl X.

During this whole process of simulations Gaster Huckel charges with Assisted Model Building with Energy Refinement (AMBER7FF9) force field were applied to the protein-ligand complex. The AMBER7FF9 force field usually considers the bonding and non-bonding parameters for the simulations viz. charge, angle, torsion, and some van-der Waals parameters. Overall, molecular dynamics simulations of the protein-ligand complex was carried out which is based upon the principle of Boltzmann initial velocity. Protein was simulated up to 10,000 fs at a

time interval of 5 fs with respect to potential energy, kinetic energy and temperature.

3. Results and Discussion
3.1. Statistical Results of CoMFA

CoMFA model, generated by distill, was found to be better as it was predicted with good statistical results compared to docking based alignment in term of the q2, r 2pre, r 2cv, r 2ncv, F and SEE. The values of the q2, r 2pre, r 2cv, r 2ncv, F and SEE for distill-based alignment were found to be 0.735, 0.769, 0.722, 0.974, 152, and 0.088, respectively (Table 2). Among all the statistical parameters, r 2pre is one of the important parameter as it gave the predictive ability of the CoMFA model. During this process, 8 molecules among 39 molecules removed from the data set which acted as an external test molecule and their activity were predicted with the help of the generated CoMFA model. The value of r 2pre = 0.769 which is more than 0.5 validates the predictive capability of the CoMFA model. While on the other hand for docking-based alignment, the values of the q2, r2pred, r 2cv, r 2 , F and SEE were found to be 0.512, 0.500, 0.504, 0.802, 70 and 0.121, respectively.
Results obtained through the docking based alignment were not up to the mark which clearly indicated the efficiency of distil based alignment.

3.2. Statistical Results of CoMSIA

Similar to CoMFA, this model was also generated using two different alignment techniques viz. rigid alignment (distill) and docking-based alignment; among which, rigid alignment (distill) was found to be the best as it gave good statistical results compared to docking-based alignment.
Rigid alignment (Distill) attempts to align molecules in a database to a template molecule on a common backbone or core scaffold. In addition to this, different combinations of Steric, Electrostatic, Hydrophobic, Hydrogen Bond Donor and Hydrogen Bond Acceptor (SEHDA) were
carried out for the determination of efficient CoMSIA model generation as shown in Table 3. Statistical results revealed that the combinations of Steric, Electrostatic, Hydrophobic, and Hydrogen Bond Donor (SEHD), were found to be efficient among the 31 different combinations shown in Table 3. For CoMSIA (SEHD) analysis the values of q2, r2pre, r2 r2 and r2 and SEE were found to be 0.761, 0.946, 0.775, 0.976 and 0.085 respectively which is shown in Table 2.
Similar to the CoMFA model, CoMSIA (SEHD) gave r2pre = 0.946 which is more than 0.5 proving the predictive capability of the CoMSIA (SEHD) model.

3.3. Statistical Results of Topomer CoMFA

In Topomer CoMFA analysis, the cross-validated coefficient (q2) and noncross-validated coefficient (r2) were found as 0.693 and 0.940 respectively; which suggested that the model has good predictive ability. The generated model also showed good predictive power as it showed r 2pre value of 0.720.

3.4. Statistical Results of HQSAR

For statistically better HQSAR model generation, first combinations of the parameters like atomic number (A), bond type (B), atomic connection (C), hydrogen bond donor/ acceptor (DA), hydrogen (H) and chirality (Ch) with minimum 4 and maximum 6 atom counts were performed which are listed in Table 4. The results suggested that the parameters such as atomic number (A), bond type (B) and atomic connection (C) were important for a model generation. Among 14 generated models, model 4 was selected for further improvement of q2, and for this purpose, the fragment size of 1 to 10 was explored. A significant difference was observed in the statistical parameters as shown in Table 5 where q2 and r2 were found to be 0.655 and 0.891, respectively.
This was done to check the influences of the fragment size on the q2 and r2 values, which are the key statistical parameters responsible for model robustness. The experimental and predicted pIC1502 with residual values of CoMFA, CoMSIA, Topomer CoMFA and HQSAR are shown in Table 6, and graphically depicted in Fig.4

3.5. Contribution Maps of 3D-QSAR Analysis

The most active compound 36 (IC50: 0.5 nM) was chosen as a representative molecule for visualization of contour maps. The models resulted from the 3D-QSAR analysis (CoMFA, CoMSIA, and Topomer CoMFA) were graphically represented by color contour maps, where each color represented favorable and disfavorable regions while for HQSAR analysis, models were graphically represented by color coded structures in which the color of each atom of a molecule reflected the contribution of the atom towards the activity of the molecules. Contour maps of CoMFA (steric and electrostatic fields) with favorable and disfavorable regions are shown in Fig. 5(A) and 5(B). Topomer CoMFA with R1 and R2 fragments with steric and electrostatic fields are shown in Fig. 5(C), 5(D), 5(E), and 5(F). CoMSIA hydrophobic, and hydrogen bond donor contour maps with their favorable and unfavorable regions are shown in Fig. 5(G) and 5(H). HQSAR contribution with their color coded structure for positive and negative contributions are shown in Fig. 6(A) and 6(B), respectively.

3.5.1. CoMFA Contour Maps

The green color of the contour map as shown in Fig. 5(A), represented sterically favorable groups with 80 % contribution while yellow color represented unfavorable for sterically bulky groups with 20 % contribution. The green colored portion found near R1 substitution indicated that sterically bulky group at this site increased the biological activity. The above statement validated by comparing the activity of compound 30 (IC50 = 22.5 nM) and compound 32 (IC50 = 1.3 nM). Compound 30 possess methyl group at the R1 position, while compound 32 possess an

isopropyl side chain at the R1, which clearly indicated that there was an increase in the activity as there was an increment in the bulkiness at R1 position. Importance of sterically bulky group at R1 could be proven with the several compounds of the series but most prominent change was occurred in the compounds 1 (IC50 = 3.9 nM), 2 (IC50 = 3.9 nM) and 3 (IC50 = 6.3 nM) as compared to compounds 18 (IC50 = 2 nM) , 19 (IC50 = 1.9 nM) and 20 (IC50 = 1.3 nM) where morpholine at R1 position was replaced by the 3-oxa-8-azabicyclo [3.2.1]-octane which was more bulkier as compared to morpholine. On the other hand yellow region at R2 indicated that this position was favorable for less bulky groups and this is proved by the compound 10 (IC50 =
2.5 nM) and compound 11 (IC50 = 1.5 nM). The difference between these two compounds, 10 and 11, was only of methyl benzoate and benzoic acid at R2 position, respectively, which made compound 11 less bulky at R2 position compared to compound 10.
The CoMFA electrostatic contour map is shown in Fig. (5B). Blue colored region was favorable

for electropositive groups (electron donating group) with 80% contributions while the red color portion was favorable for electronegative groups (electron withdrawing groups) with 20% contribution. A red contour found near R1 substitution could be suggested that electronegative groups were required for the biological activity. This can be easily recognized after viewing the structures of the compound 14 (IC50 = 1.4 nM) and compound 30 (IC50 = 22.5 nM). This drastic decrease in the activity of compound 30 might be due to the absence of the electronegative atoms like oxygen and nitrogen in the form of a morpholine ring at R1 position. Small red color contour was found near to the urea linker of the scaffold which proved the importance of carbonyl group for the activity. R2 position was found to be favorable for both combinations of electron donating and withdrawing groups. This could be validated by the compound 1 (IC50 = 3.9 nM) and 12 (IC50 = 0.7 nM). Compound 1 possess benzene ring at R2 position while compound 12
possess N-methyl benzamide at R2 position which was a balanced combination of –NH and –C=O.

3.5.2. CoMSIA Contour Maps

CoMSIA (SEHD) was found to be the best among all other combinations of the fields as it predicted good statistical results. Contour maps were generated based upon this model. CoMSIA steric and electrostatic contour maps are shown in Fig. 5(C) and 5(D) which showed a similar type of trend as that of the CoMFA steric and electrostatic. Hence, this outcome additionally validated the CoMFA contour maps and their interpretation.
CoMSIA hydrophobic contour map is shown in Fig. 5(G); where hydrophobic field represented by yellow color contour map with 80 % contribution, and white color contour map with 20 % contribution. R1 substitution was found to be favorable for the hydrophobic group as it was surrounded by yellow color portion. Similarly, R2 position was also found to be favorable for the hydrophobic groups. Therefore, substitution with hydrophobic groups at R1 and R2 may increase the mTOR inhibitory activity of the molecules. Hence in the whole series, most of the compounds consisted of bulky groups at both R1 and R2 positions. It can also be revealed by the compounds 28 (IC50 = 7.7 nM) and compound 29 (IC50 = 14.5 nM). This reduction in the activity of the compound 29 was might be due to the replacement of one methyl group of dimethylene diamine side chain by hydrogen at R2 position which made compound 29 less hydrophobic as compared to the compound 28. Similarly, the importance of the bulky group at R1 position could be explained by looking at compounds 6 and 22. More potency of the compound 22 (IC50 = 0.43 nM) compared to compound 6 (IC50 = 1.9 nM) might be due to the presence of 3-oxa-8- azabicyclo [3.2.1]-octane instead of the simple morpholine ring.

The hydrogen bond donor (HBD) contour map is shown in Fig. 5(H) where cyan color contour map with 80 % contribution was favorable for hydrogen bond donor groups while purple color contour map with 20 % contributions was unfavorable for hydrogen bond donor groups. The R2 position of the scaffold was found to favorable for hydrogen bond donor and it could be proved by the compounds 11 (IC50 = 1.5 nM), 12 (IC50 = 0.7 nM), and 13 (IC50 = 0.5 nM) which were more potent compared to the compound 1 (IC50 = 1.9 nM) as there was no hydrogen bond donor at R2 position of the compound 1.
3.5.3 Topomer CoMFA Contour Maps

The contour maps of Topomer CoMFA were generated around R1 as shown in Fig. 5(C) and 5(D) while for R2 is shown in Fig. 5(E) and 5(F). Contour maps generation includes fragment based analysis. The steric contour map of R1 substitution shown in Fig. 5(C) surrounded by green color portion indicated that this position is favorable for the sterically bulky group while steric contour map of the R2 substitution shown in Fig. 5(E) was found to be disfavored for the sterically bulky group as it is surrounded by the yellow color region which is in line with the outcome of the CoMFA steric contour maps. Similarly, electrostatic contour map of the R1 substitution was found to be favorable for the electronegative groups as this region surrounded by small red color portion while on the other hand electrostatic contour map of R2 substitution shown in Fig. 5(F) surrounded by the blend of blue and red color indicated that this position was found to be favorable for both combination of electron donating and withdrawing groups. This results of the electrostatic contour maps of R1 and R2 substitutions were found to be in line with the outcome of the CoMFA electrostatic contour map. Overall, Topomer CoMFA showed a similar type of contours as that of the CoMFA steric and electrostatic. Hence additionally topomer CoMFA contour maps validated the CoMFA steric and electrostatic contour maps.
3.5.4. HQSAR Contribution Maps

In HQSAR contribution maps, atoms or fragments with one end of the spectrum (yellow, yellow-

green, and green) represented favorable or positive contribution towards the inhibitory activity, while another end of the spectrum (red, red-orange, and orange) represented unfavorable or negative contribution towards the activity. White colored atom or fragments indicated intermediate contribution to the activity. For molecular optimization, HQSAR serves as an important tool because the fragments, which were responsible for the negative contribution of the molecules, could be modified to improve the biological activity. The fragments which were part of the most active compound and least active compounds occurred in the spectrum as per their contribution to the activity. Fourteen HQSAR models were generated among which model 4 was found to be the best, as listed in Table 3. Further, this model was optimized with different combinations which suggested model 3 as the best model as shown in Table 4, and the spectrum obtained through this model is depicted in Fig. 6(A) and 6(B). The examination of the contours associated with compound 36 showed large areas of the molecule covered by green and yellow contours. This indicated that a large percentage of the molecular fragments made a positive contribution towards the biological activity as shown in Fig. 6(A). On the other hand, compound 30 was found at the negative end of the spectrum indicating that the fragments of this compound diminished the activity of the compound 30. Green color of urea linker shown in Fig. 6(B) suggested that this linker might play an important role for the mTOR inhibitory activity.
Most of the green colored atoms were found in the benzene ring of R2 substitution proving the importance of R2 fragment. Fragments showed as F1 and F2 in Fig. 6(C) contributing positively to the mTOR inhibitory activity as it possesses the positive coefficient value of 0.023 and 0.020 respectively.

3.5.5. Molecular Docking and Molecular Dynamics Simulations (MD/MS)
To validate the contour maps, MD/MS study was performed with respect to kinetic energy,

potential energy, and temperature (data shown in the supplementary file as Table ST-1). Initially, the compound 36 was docked and followed by molecular dynamics simulations (MD) which stabilized the protein-ligand complex and suggested more precise and accurate results.
Compound 36 in complex with simulated mTOR protein (shown in Fig. 7) indicated that tetrahydrofuran substitution at R1 of the compound 36 surrounded by TRP2549 and LYS2171 which were hydrophobic and positively charged, respectively. These clearly proved that the oxygen atom of tetrahydrofuran and ethereal linkage at R1 were favorable for electronegative groups as LYS2171 was positively charged amino acid. Since TRP2549 was hydrophobic in nature, it proved that the R1 position was favorable for bulky and hydrophobic groups. The -NH of the urea linker of compound 36 forming a hydrogen bond with ASP2195 validated the contour map of hydrogen bond donor where this position of compound 36 was found to be favorable for hydrogen bond donor groups. As per the CoMSIA hydrophobic contour map, few portion of the R2 position was also found to be favorable for the hydrophilic group, which was further validated by the presence of charged amino acids viz. ASP2191, ASP2195, and ASP2360 surrounded to the R2 position.
3.5.6. Designing

To validate the effectiveness of 3D-QSAR and MD/MS study, acridine scaffold was selected as the core moiety and substitutions were added in this acridine scaffold followed by docking which is shown in the Fig. 8. All these substitutions were added based upon the contour maps which is summarized in Fig. 9. The R1 substitution of triazine morpholino scaffold was found to be favorable for bulky, hydrophobic and electronegative groups. On the other hand, R2 was found to be favorable for less bulky group, electronegative polar functional group and hydrogen bond donor. In lieu with this the R1 position of the acridine, scaffold was substituted
with the benzenesulfonamide, nitrobenzene, benzene and fluorobenzene while R2 position

substituted with the methoxy group (-OCH3) group which resulted in the new acridine derivatives shown in Table 7 along with their predicted pIC50. These pIC50 values of the newly designed acridine molecules predicted from the generated CoMSIA (SEHD) model revealed the better r2pred (0.946) compared to the generated CoMFA model during the external test validation. These newly designed acridine derivatives shown predicted pIC50 values very close to the highest active compound 36 of triazine morpholino derivatives which further validated by the molecular docking study. Designed molecules shown in Table 8 were docked with the mTOR protein (PDB ID: 4JT6) using GOLD 5.2 software where these molecules showed the same interactions as that of the standard PI-103 shown in Fig.8. The outcome of the in silico results of the designed molecules proved that 3D-QSAR and MD/MS are very useful methods for the designing of novel mTOR inhibitors.

4. Conclusion

mTOR is proven to be a promising target for the treatment of cancer. 3D-QSAR techniques (CoMFA, CoMSIA, Topomer CoMFA) and HQSAR were applied on the series of morpholino triazine derivatives as mTOR inhibitors. All the models CoMFA, CoMSIA, Topomer CoMFA, and HQSAR; were found satisfactory according to the statistical parameters. Additionally, MD/MS study was performed to validate the contour maps results. Based upon these studies, the structural requirement of mTOR inhibition is summarized in Fig. 9, where R1 position of triazine morpholino scaffold was found to be favorable for bulky, hydrophobic and electronegative groups. On the other hand, R2 was found to be favorable for less bulky group, electronegative and hydrogen bond donor. The obtained features were incorporated in the acridine scaffold and four acridine derivatives were designed and docked with mTOR. The results of the molecular
docking study of all of the designed molecules proved that the conclusions of the present research

work, in the form of important features required for mTOR inhibition, could be explored in the future for the discovery of novel potent inhibitors.
Acknowledgement

The authors are thankful to Nirma University, Ahmedabad, India for providing necessary facilities to carry out the research work. Udit Chaube is thankful to the Department of Science and Technology (DST), Govt. of India for providing INSPIRE fellowship [IF140932].

Supplementary Information

Supplementary Information file contains details about MD/MS study with respect to kinetic energy, potential energy, and temperature. Data are shown in the supplementary file as Table ST-1.

References

1. Ali I, Wani W, Saleem K. Cancer Scenario in India with Future Perspectives. Cancer Ther; 2011; 8, 56-70.
2. WHO website. (www.who.int.) assessed on 17 September 2015.

3. Ciruelos GEM. Targeting the PI3K/AKT/mTOR Pathway in Estrogen Receptor-Positive Breast Cancer. Cancer Treat Rev; 2014;40,862–71.
4. Polivka J, Janku F. Molecular Targets for Cancer Therapy in the PI3K/AKT/mTOR Pathway. Pharmacol Ther; 2014;142:164-75.
5. Guertin DA, Sabatini DM. An Expanding Role for mTOR in Cancer. Trends Mol Med; 2005;11:353–61.
6. Sun SY. mTOR Kinase Inhibitors as Potential Cancer Therapeutic Drugs. Cancer Lett; 2013;340:1–8.
7. De P, Miskimins K, Dey N, Leyland JB. Promise of Rapalogues versus mTOR Kinase Inhibitors in Subset Specific Breast Cancer: Old Targets New Hope. Cancer Treat Rev; 2013;39:403–12.
8. Husseinzadeh N, Husseinzadeh HD. mTOR Inhibitors and Their Clinical Application in Cervical, Endometrial and Ovarian Cancers: A Critical Review. Gynecol Oncol;2014; 133:375–81.
9. Efeyan A, Sabatini DM. mTOR and Cancer: Many Loops in One Pathway. Curr Opin Cell Biol; 2010;22:169–76.
10. Zhang YJ, Duan Y, Zheng XFS. Targeting the mTOR Kinase Domain: The Second Generation of mTOR Inhibitors. Drug Discov Today; 2011;16:325–31.
11. Zagouri F, Sergentanis TN, Chrysikos D, Filipits M, Bartsch R. mTOR Inhibitors in Breast

Cancer: A Systematic Review. Gynecol Oncol; 2012;127:662–72.

12. García EC. Allosteric and ATP-Competitive Kinase Inhibitors of mTOR for Cancer Treatment. Bioorg Med Chem Lett; 2010;20:4308–12.
13. Ran T, Lu T, Yuan H, Liu H, Wang J, Zhang W, Leng Y, Lin G, Zhuang S, Chen YA. Selectivity Study on mTOR/PI3Kα Inhibitors by Homology Modeling and 3D-QSAR. J Mol Model; 2012;18:171–86.
14. Yang H, Rudge DG, Koos JD, Vaidialingam B, Yang HJ, Pavletich, NP. mTOR Kinase Structure, Mechanism and Regulation. Nat; 2013;497:217–223.
15. Saqib U, Kumar B, Siddiqi MI. Structural investigations of anthranilimide derivatives by CoMFA and CoMSIA 3D-QSAR studies reveal novel insight into their structures toward glycogen phosphorylase inhibition. SAR QSAR Environ Res; 2011;22:411–49.
16. Gadhe CG, Thirumurthy M, Gugan K, Seung C. In Silico Quantitative Structure-Activity Relationship Studies on P-Gp Modulators of Tetrahydroisoquinoline-Ethyl-Phenylamine Series. BMC Struct Biol; 2011;11:1-15.
17. Dehnhardt CM, Venkatesan AM, Chen Z, Delos SE, Ayral KS, Brooijmans N, Yu K, Hollander I, Feldberg L, Lucas J. et al. Identification of 2-Oxatriazines as Highly Potent Pan-PI3K/mTOR Dual Inhibitors. Bioorg Med Chem Lett; 2011;21:4773–4778.
18. Venkatesan, AM, Dehnhardt CM, Delos SE, Chen Z, Dos SO, Ayral KS, Khafizova, G, Brooijmans N, Mallon R, Hollander I. et al. Bis(morpholino-1,3,5-Triazine) Derivatives: Potent Adenosine 5′-Triphosphate Competitive Phosphatidylinositol-3-kinase/Mammalian Target of Rapamycin Inhibitors: Discovery of Compound 26 (PKI-587), a Highly Efficacious Dual Inhibitor. J Med Chem; 2010;53:2636–45.
19. Venkatesan AM, Chen Z, Santos OD, Dehnhardt C, Santos ED, Ayral-KS, Mallon R, Hollander I, Feldberg L, Lucas J. et al. PKI-179: An Orally Efficacious Dual

Phosphatidylinositol-3-Kinase (PI3K)/mammalian Target of Rapamycin (mTOR) Inhibitor. Bioorg Med Chem Lett; 2010;20:5869–73.
20. Sybyl X Molecular Modeling Software, Tripos Associates, V. 1.2, St. Louis, USA, 2011; software available at http://www.tripos.com.
21. Gasteiger J, Marsili M. Iterative partial equalization of orbital electronegativity a rapid access to atomic charges. Tetrahedron; 1980; 36:3219–28.
22. Yu Z, Li X, Ge C, Si H, Cui L. et al. 3D-QSAR modeling and molecular docking study on Mer kinase inhibitors of pyridine-substituted pyrimidines. Molec Divers;2014;19: 135-47.
23. Vaidya A, Jain AK, Prashantha Kumar BR, Sastry GN, Kashaw SK, Agrawal RK. CoMFA, CoMSIA, kNN MFA and docking studies of 1,2,4-oxadiazole derivatives as potent caspase-3 activators. Arab J Chem; 2014:1-11
24. Chaube UJ, Bhatt H. 3D-QSAR, molecular dynamics simulations, and molecular docking studies on pyridoaminotropanes and tetrahydroquinazoline as mTOR inhibitors. Molec Divers; 2017;21:741-59
25. Vyas VK, Bhatt HG, Patel PK, Jalu J, Chintha C, Ghate M. CoMFA and CoMSIA studies on C-aryl glucoside SGLT2 inhibitors as potential anti-diabetic agents. SAR QSAR Environ Res; 2013;27:37–41.
26. Pandey G, Saxena AK. 3D QSAR studies on protein tyrosine phosphatase 1B inhibitors: Comparison of the quality and predictivity among 3D QSAR models obtained from different conformer- based alignments. J Chem Inf Model; 2006;46: 2579–90.
27. Yang W, Shu M, Wang Y, Wang R, Hu Y, Meng L, Lin Z. 3D-QSAR and Docking Studies of 3-Pyridine Heterocyclic Derivatives as Potent PI3K/mTOR Inhibitors. J Mol Struct; 2013:1054-55.

28. Lu P, Wei X, Zhang R. CoMFA and CoMSIA 3D-QSAR studies on quionolone caroxylic acid derivatives inhibitors of HIV-1 integrase. Eur J Med Chem; 2010;45: 3413–19.
29. Caballero J. 3D-QSAR (CoMFA and CoMSIA) and pharmacophore (GALAHAD) studies on the differential inhibition of aldose reductase by flavonoid compounds. J Mol Graph Model; 2010;29: 363–71.
30. Ke, Zhipeng, Tao L, Haichun L, Haoliang Y, Ting R, Yanmin Z, Sihui Y, et al. 3D-QSAR and Molecular Fragment Replacement Study on Diaminopyrimidine and Pyrrolotriazine ALK Inhibitors. J Mol Struct; 2014;1067:127–37.
31. Ajala, AO, Cosmas O. 3D-QSAR Topomer CoMFA Studies on 10 N-Substituted Acridone Derivatives. Open J Med Chem; 2012:43–49.
32. Cramer, Richard D. Topomer CoMFA: A Design Methodology for Rapid Lead Optimization. J Med Chem; 2003;46:374–88.
33. Ugarkar, Apoorva G, Premlata KA, Evans CC, Santosh N, Raghuvir RS, Pissurlenkar.

Extracting Structural Requirements for Activity of GPR119 Agonists: A Hologram Quantitative Structure Activity Relationship (HQSAR) Study. Can J Chem;2014;92, 670– 76.
34. Dixit A, Kashaw SK, Gaur S, Saxena AK. Development of CoMFA, advance CoMFA and CoMSIA models in pyrroloquinazolines as thrombin receptor antagonist. Bioorg Med Chem;2004;12:3591–98.
35. Roy KK, Dixit A, Saxena AK. An investigation of structurally diverse carbamates for acetylcholinesterase (AChE) inhibition using 3D-QSAR analysis. J Mol Graph Model; 2008;27:197–208.

36. Sridhar J, Foroozesh M, Stevens CLK. QSAR models of cytochrome P450 enzyme 1A2
inhibitors using CoMFA, CoMSIA and HQSAR. SAR QSAR Environ Res; 2011; 22:681– 97.

Figure Captions

Figure 1. Chemical structure of rapamycin along with novel second generation mTOR & mTOR/PI3K dual inhibitors

Figure 2. (A) Core structure of triazine morpholino derivatives with marked points used for alignments (B) Distill Rigid Body Alignment (C) Docking Based Alignment.

Figure 3. Fragmentation pattern of R1 and R2 substitution for all molecules of the dataset in Topomer CoMFA analysis. R1 fragment is represented by blue color, R2 fragment is represented by red color and maximum common moiety is represented by green color.

Figure 4. Correlation between the predicted versus experimental activities of CoMFA, CoMSIA, Topomer CoMFA, and HQSAR models.

Figure 5. Illustration of CoMFA, Topomer CoMFA, CoMSIA (std*coeff) contour maps. Template molecule 36 is shown inside the field, (A) CoMFA steric field shown in green contour (80% contribution) refer to sterically group favored regions; yellow contour (20% contribution) indicates disfavored areas. (B) CoMFA electrostatic field shown in blue contour (80% contribution) indicates regions where electropositive groups (electron donating groups) are favored and red contour (20% contribution) refer to the region where electronegative (electron withdrawing groups) are favored.

(C) Topomer steric contour maps for R1 fragments. (D) Topomer electrostatic contour maps for R1 fragments. (E) Topomer steric contour maps for R2 fragments. (F) Topomer electrostatic contour maps for R2 fragments. (G) CoMSIA Hydrophobic field shown in yellow contour (80 % contribution) indicates favored region and white contour (20 % contribution) indicates disfavored region. (H) CoMSIA Donor field shown in cyan contour (80 % contribution) indicates favored region and blue contour (20 % contribution) indicate disfavored region.

Figure 6. The contribution of HQSAR maps. (A) Contribution map for the most active compound 36 (IC50: 0.3 nM). Green color indicate most favored region, green-yellow color indicate favorable regions, and yellow color indicate the least favorable regions for the activity. (B) Contribution map for least active compound 30 (IC50: 22.5 nM). Red, red-orange, orange colors indicate negative regions for the activity. The white color code indicates a molecule’s atom with no contributions towards the activity. (C) The Fragments F1, F2 and F3 identified through this study in which fragments F1 and F2 are found to be important for the activity.

Figure 7. Molecular docking interactions of compound 36 with mTOR protein (PDB: 4JT6) after molecular dynamics simulations of compound 36 and mTOR protein complex.

Figure 8. A) Standard co-crystallized ligand PI-103 in a complex with mTOR protein (PDB: 4JT6). Molecular docking results of (B) Hydrogen bond interactions of acridine derivative-1 with amino acids ASP2195, TYR2225, and ASP2357 (C) Hydrogen bond interaction of acridine derivative-2 with amino acid TYR2225 (D) Hydrogen bond interaction of acridine PI-103 derivative-3 with amino acid TYR2225 (E) Hydrogen bond interaction of acridine derivative-4 with amino acid TYR2225