Molecular Docking and Molecular Dynamics Simulation Studies of Quinoline-3-Carboxamide Derivatives with DDR Kinases–Selectivity Studies towards ATM Kinase

: Quinoline-3-carboxamides are an essential class of drug-like small molecules that are known to inhibit the phosphatidylinositol 3-kinase-related kinases (PIKK) family kinases. The quinoline nitrogen is shown to bind to the hinge region of the kinases, making them competitive inhibitors of adenosine triphosphate (ATP). We have previously designed and synthesized quinoline-3-carboxamides as potential ataxia telangiectasia mutated (ATM) kinase inhibitors to function as an adjuvant treatment with DNA damaging agents. This article discusses the molecular docking studies performed with these derivatives with the DNA damage and response (DDR) kinases-ATM, ataxia telangiectasia and rad3 related (ATR), and DNA dependent protein kinase catalytic subunit (DNA-PKcs) and highlights their selectivity towards ATM kinase. Docking studies were also performed with mTOR and PI3K γ , which are close homologs of the DDR kinases. Molecular dynamics simulations were performed for one of the inhibitors against all the enzymes to establish the stability of the interactions involved. Finally, the absorption, distribution, metabolism, and excretion (ADME) properties of the inhibitors were predicted using the QikProp manual in Maestro. In conclusion, the molecules synthesized showed high selectivity towards the ATM kinase in comparison with the other kinases, though the sequence similarity between them was relatively high.


Introduction
Synthesizing selective inhibitors for kinases belonging to the same family can be challenging, owing to the similarity in their sequences. One such family is the phosphatidylinositol 3-kinase-related kinases (PIKK) family of kinases, which comprises the DNA damage and response (DDR) pathway of kinases [1,2]. These kinases play a crucial role in mediating the detection and repair of damages imposed on DNA, thereby helping in the faithful transfer of genetic information from the parent cell to the daughter cell [3,4]. On the other hand, cancer cells exploit these kinases to repair their DNA damages, thereby eluding cell death. It is known that radiotherapy and most of the chemotherapeutic agents destroy cancer cells by causing DNA damage. Thus, the overexpression of the DDR kinases in cancer cells creates resistance to DNA damaging treatments [5][6][7]. Therefore, it is essential that these kinases are inhibited to avoid resistance towards the existing modes of treatment.
The three major kinases involved in the DDR pathway are the ataxia telangiectasia mutated (ATM), the ataxia telangiectasia and rad3 related (ATR), and the DNA dependent protein kinase catalytic subunit (DNA-PKcs) [1]. These three kinases belong to the PIKK  The kinase domain of ATM, ATR, DNA-PKcs, mTOR, and PI3Kγ share 25-40% of similarity in their sequence. Therefore, we decided to carry out molecular docking analysis for the synthesized molecules against all the mentioned kinases. In this article, we report the molecular docking and molecular dynamics simulation (MD simulation) studies of the synthesized molecules against the kinases of the DDR pathway (ATM, ATR, and DNA-PKcs), mTOR and PI3Kγ. The multiple sequence alignment (MSA) analysis for the proteins was also performed. Compound 6f showed the highest docking score against ATM, and thus, MD simulation was performed for the same with all the mentioned proteins. QikProp analysis was also done using Maestro 12.6 to quantify the absorption, distribution, metabolism, and excretion (ADME) properties of the molecules. The properties were all in the agreeable region. Finally, a brief analysis of the selectivity of the molecules  The kinase domain of ATM, ATR, DNA-PKcs, mTOR, and PI3Kγ share 25-40% of similarity in their sequence. Therefore, we decided to carry out molecular docking analysis for the synthesized molecules against all the mentioned kinases. In this article, we report the molecular docking and molecular dynamics simulation (MD simulation) studies of the synthesized molecules against the kinases of the DDR pathway (ATM, ATR, and DNA-PKcs), mTOR and PI3Kγ. The multiple sequence alignment (MSA) analysis for the proteins was also performed. Compound 6f showed the highest docking score against ATM, and thus, MD simulation was performed for the same with all the mentioned proteins. QikProp analysis was also done using Maestro 12.6 to quantify the absorption, distribution, metabolism, and excretion (ADME) properties of the molecules. The properties were all in the agreeable region. Finally, a brief analysis of the selectivity of the molecules towards the ATM kinase was also accomplished using the values obtained from the XP Visualizer module of Maestro 12.6.

Materials and Methods
A Swiss modeling tool was used to model the structure of the kinase domain of DNA-PKcs. The template used to model DNA-PKcs was 7k0y.1, which corresponds to the Cryo-EM structure of the activated form of DNA-PK (complex VI) [16]. The model was validated using the Ramachandran plot. The kinase domains of ATM and ATR were used from previously available data [14,17]. The protein data bank (PDB) files for mTOR and PI3Kγ are 5WBY and 5G55, respectively [14,17]. Molecular docking and the QikProp studies for all the synthesized inhibitors ( Figure 1) were performed using Maestro 12.6 (Schrödinger, LLC, New York, NY, 2020). Molecular dynamics simulations were performed using NAMD 2.13 [18]. The key features of the proteins are highlighted in Table S1 in the Supplementary Materials.
MSA analysis and image generation-Clustal Omega and Jalview.

Molecular Docking
Protein structures and crystal structures, both modeled (details of the proteins included in Table 2), were subjected to the protein preparation module to fix the structural defects in the basic structures and make them compatible for further applications. The structures of the ligands were sketched using 2D Sketcher, and they were then subjected to the LigPrep module to obtain their energy minimized forms. Following this, SiteMap analysis was performed on every protein to identify the best druggable site. The best druggable site came out to be either the hinge region (where ATP binds) and/or a site in proximity to the hinge region. Thus, the best druggable site was chosen for every protein of interest for docking analysis. Further on, the 3D space where the docking was to be performed was confined to the site obtained. This was carried out in the receptor grid generation module. The ligands were docked in the grid generated using the outputs received from the LigPrep. The docking scores, penalties and rewards were analyzed in the XP Visualizer platform. Additionally, the LigPrep output file was subjected to ADME analysis by using the QikProp module.

Molecular Dynamics Simulation
Molecular dynamics simulations were performed on ATM, ATR, PI3Kγ, DNA-PKcs, and mTOR. The docked structures of the protein with an inhibitor were taken as a starting point for simulations. Simulations were conducted in a periodic water box for 100 ns using the CHARMM36 force field the and NAMD package version 2.13. The force field for ligands was generated from the CHARMM-GUI server. The water box (including 150 mM NaCl) was created by adding water for 20 Å in the positive and negative x, y, and z directions around the protein, yielding a cuboidal box. LJ cutoff was defined at a distance of 12 Å, with a switching distance of 10 Å. Long-range electrostatic interactions were handled using the particle mesh Ewald (PME) method. Before the production run, the systems were minimized for 5000 steps using a conjugate gradient algorithm. The simulations were performed in an NPT ensemble; the temperature and the pressure of the system were fixed at 300 K and 1 bar, respectively, using a Langevin thermostat and barostat. Post simulation analyses were performed using VMD (root-mean square deviation (RMSD) and secondary structure timeline analyses) [19].

Multiple Sequence Alignment Analysis
Multiple sequence alignment (MSA) of the kinase domains of the ATM(Q13315), ATR(Q13535), DNA-PKcs(P78527), mTOR(P42345), and PI3Kγ(P48736) proteins were performed using Clustal Omega [20]. The sequence of the kinase domain of the proteins was taken for analysis from the UniProt database. The results of Clustal Omega were further analyzed using Jalview software (version 2.11.1.4) [21].

Molecular Docking Studies
Prior to performing molecular docking, we carried out the in silico analysis of the pharmacokinetic properties of the synthesized small molecules (Figure 1). This highlighted some of the critical parameters (like bioavailability and the aqueous solubility of the molecules). Table 3 lists some of the vital ADME properties of the molecules. As evident from the table, the values for the selected properties were well within range, and the molecules showed an excellent percentage of human oral absorption. QPlogS-predicted aqueous solubility, S in mol/L (acceptable range; −6.5-0.5); QPPCaCo-predicted apparent Caco-2 cell permeability in nm/sec. Caco-2 cells are a model for the gut-blood barrier. QikProp predictions are for nonactive transport (<25 poor, >500 great); Percentage of human oral absorption (<25% is poor and >80% is high); QPlogKhsa-prediction of binding to human serum albumin (acceptable range; −1.0-1.5); QLogBB-prediction of brain/blood partition coefficient (acceptable range; −3.0-1.2); molecular weight (<500 Da); HBD-Hydrogen bond donor (<5); HBA-Hydrogen bond acceptor (<10); QPlogP(o/w)-predicted octanol/water partition co-efficient logP (acceptable range; −2.0-6.5).
After performing the ADME analysis, we carried out the docking studies for the molecules. Table 4 displays the docking score of the molecules with all the proteins of interest. From Table 4, it can be concluded that all of the molecules have a significantly lower docking score with ATM compared to those of other proteins (the more negative the docking score, the better the binding). Since molecule 6f had the highest docking score against ATM, the docking poses of 6f with ATM, ATR, DNA-PKcs, mTOR, and PI3Kγ are illustrated in Figure 2. Similar docking poses of the ligands for all the proteins have been attached in the Supplementary Materials ( Figures S1-S5). The π-π stacking interactions, halogen bonding, hydrogen bonding, and aromatic hydrogen bonding are the typical interactions being observed. Further analysis in order to understand the differences in the docking scores was carried out using MD simulations and reward/penalty analysis. Table 4. Docking scores for the ligands with the proteins.

Ligands
ATM ATR DNA-PKcs mTOR PI3Kγ ing score with ATM compared to those of other proteins (the more negative the docking score, the better the binding). Since molecule 6f had the highest docking score against ATM, the docking poses of 6f with ATM, ATR, DNA-PKcs, mTOR, and PI3Kγ are illustrated in Figure 2. Similar docking poses of the ligands for all the proteins have been attached in the Supplementary Section (Figures S1-S5). The π-π stacking interactions, halogen bonding, hydrogen bonding, and aromatic hydrogen bonding are the typical interactions being observed. Further analysis in order to understand the differences in the docking scores was carried out using MD simulations and reward/penalty analysis.   The interactions of molecule 6f with ATM, ATR, mTOR, DNA−PKcs, and PI3Kγ are abridged in Table 5. Similar tables for all the ligands with the remaining proteins have been added in the Supplementary Materials (Tables S1-S5). As seen from Table 5, either tyrosine or tryptophan is shown to have π−π stacking interaction in all four proteins, except for DNA−PKcs, which exhibited only aromatic hydrogen bonding with the ligand. MD simulation analyses to establish the stability of the protein−6f complexes have been performed, and the stabilities of the most common interactions between the proteins and compound 6f were also analyzed. Halogen bonding, π−π stacking interaction, aromatic hydrogen bonding, hydrogen bonding

Molecular Dynamics Simulation
Molecular dynamics simulation was performed to provide insight into the protein−ligand stability and protein structural flexibility of the docked complexes. The simulations for the 6f−bound forms of ATM, ATR, DNA−PKcs, mTOR, and PI3Kγ were carried out using the docked structure as a starting geometry. Figure 3a shows the RMSD value of the protein−ligand, protein, and ligand for ATM−6f complex. The visual analysis of the protein−ligands and protein trajectory confirms that the protein's secondary structure is stable throughout the simulation (100 ns). Figure 3b presents the timeline analysis secondary structure of ATM protein during the simulation, and it was observed that the secondary structure does not show a significant variation protein throughout the simulations. However, some protein regions manifest structural fluctuation and interchange from one form to another in the simulation process, thus adding to the RMSD value observed (Figure 3b). In the ATM−6f complex, the critical binding site residue (Tyr (73) and Asp (207)   The RMSD and secondary structure timeline analysis for ATR−6f, DNA−PKcs−6f, and mTOR−6f complexes are shown in Figure S6−S8 (Supplementary Materials). In all three cases, the trajectory analysis suggests that the protein was stable during simulations. Likewise, the secondary structure analysis (frame by frame) showed that the secondary structure remains the same for most of the proteins' regions (all three). Simultaneously, some flexible regions manifest secondary structure changes and thus contribute to the global RMSD of the protein. Interestingly, the RMSD of the ligand (6f) in all these cases shows similar behavior to that of those observed in the ATM−6f complex (i.e., a periodic fluctuation in the RMSD value resulted because of the conformational switching due to the ro-  The RMSD and secondary structure timeline analysis for ATR−6f, DNA−PKcs−6f, and mTOR−6f complexes are shown in Figures S6-S8 (Supplementary Materials). In all three cases, the trajectory analysis suggests that the protein was stable during simulations. Likewise, the secondary structure analysis (frame by frame) showed that the secondary structure remains the same for most of the proteins' regions (all three). Simultaneously, some flexible regions manifest secondary structure changes and thus contribute to the global RMSD of the protein. Interestingly, the RMSD of the ligand (6f) in all these cases shows similar behavior to that of those observed in the ATM−6f complex (i.e., a periodic fluctuation in the RMSD value resulted because of the conformational switching due to the rotation around amide N−C bond (Figure 3d)).
For the PI3Kγ−6f complex, the trajectory analysis suggests that the protein is stable throughout the simulation (Figure 4a). Moreover, the RMSD value was observed to be in the range of 2.5 Å, half of that observed in other protein−6f complexes. The secondary structure timeline analysis suggests that the secondary structure of the protein remains the same during the simulation, with some minor fluctuations. The trajectory visualization of the ligand (6f) shows a nearly constant RMSD throughout the simulation process, unlike the other protein−6f complexes, wherein periodic fluctuations of the ligand RMSD were observed, as discussed above. A constant ligand RMSD in the case of the PI3Kγ−6f complex suggests no conformational switching of the ligand during simulations. The observation is opposite to that seen in other protein−ligand complexes, where conformational interconversion was seen due to free rotation around the amide N−C bond (Figure 3d). The restriction in the free rotation could possibly be due to the rigid packing of the inhibitor−binding site and, hence, an enhanced energy barrier of free rotation around the amide N−C bond of the ligand.

Multiple Sequence Alignment (MSA) Analysis
MSA analysis was conducted to understand the differences in docking score of all the synthesized inhibitors with the proteins under consideration. These proteins share an amino acid sequence similarity of 25−40% in their kinase domain. MSA analysis helped us to understand the conservation of residue across all the proteins. Furthermore, it also gave us an idea of whether the inhibitors are causing an interaction with the conserved or the nonconserved residue in the protein. Though all the proteins share a good sequence similarity, they have differences in preference for the druggable site in their kinase do-

Multiple Sequence Alignment (MSA) Analysis
MSA analysis was conducted to understand the differences in docking score of all the synthesized inhibitors with the proteins under consideration. These proteins share an amino acid sequence similarity of 25−40% in their kinase domain. MSA analysis helped us to understand the conservation of residue across all the proteins. Furthermore, it also gave us an idea of whether the inhibitors are causing an interaction with the conserved or the nonconserved residue in the protein. Though all the proteins share a good sequence similarity, they have differences in preference for the druggable site in their kinase domain. The result showed us that most of the interactions made by the inhibitors do not target the conserved residues of all the mentioned proteins ( Figure 5). There is a difference in the inhibitor binding sites and the type of interactions for all the proteins with the ligands. This variance might be one of the key reasons contributing to the differences in their binding affinity. This was shown in the docking and MD simulation studies for all the proteins with the inhibitor 6f.

Multiple Sequence Alignment (MSA) Analysis
MSA analysis was conducted to understand the differences in docking score of all the synthesized inhibitors with the proteins under consideration. These proteins share an amino acid sequence similarity of 25−40% in their kinase domain. MSA analysis helped us to understand the conservation of residue across all the proteins. Furthermore, it also gave us an idea of whether the inhibitors are causing an interaction with the conserved or the nonconserved residue in the protein. Though all the proteins share a good sequence similarity, they have differences in preference for the druggable site in their kinase domain. The result showed us that most of the interactions made by the inhibitors do not target the conserved residues of all the mentioned proteins ( Figure 5). There is a difference in the inhibitor binding sites and the type of interactions for all the proteins with the ligands. This variance might be one of the key reasons contributing to the differences in their binding affinity. This was shown in the docking and MD simulation studies for all the proteins with the inhibitor 6f.

Discussion
From the docking poses of 6f with ATM and PI3Kγ in Figure 2, it is evident that PI3Kγ makes better contact with the ligand. However, the docking score for all the ligands is higher than those for ATM. Analyzing the values of the rewards and penalties obtained from the XP visualizer panel of Maestro yielded a possible cause for the molecules to have

Discussion
From the docking poses of 6f with ATM and PI3Kγ in Figure 2, it is evident that PI3Kγ makes better contact with the ligand. However, the docking score for all the ligands is higher than those for ATM. Analyzing the values of the rewards and penalties obtained from the XP visualizer panel of Maestro yielded a possible cause for the molecules to have a better binding affinity towards ATM. Table 6 gives a comparison of both the rewards and penalties for ATM and PI3Kγ. Similarly, the penalties and rewards for the remaining three proteins with 6f have been included in Table S6 in the Supplementary Materials. Even though the rewards for the protein−ligand Van der Waal energy and hydrogen bond interactions that are hydrophobically packed are almost similar for these two proteins, the penalties for buried polar atoms, intraligand contacts, and exposed hydrophobic ligand groups are much higher for PI3Kγ. Therefore, the ligands are more stabilized when interacting with ATM, leading to a significantly higher docking score. LipophilicEvdW-ChemScore lipophilic pair term and fraction of the total protein−ligand vdW energy; PhobEnHB-Reward for hydrophobically packed hydrogen bond; HBond-ChemScore H-bond pair term; Electro-Electrostatic reward; Penalties-Polar atom burials and desolvation penalties and penalty for intra−ligand contacts; ExposPenal-Penalty for exposed hydrophobic ligand groups.
As discussed above, due to the large penalties resulting from intraligand contact, PI3Kγ shows a reduced binding affinity for the ligand, which results in a lower docking score. We analyzed the RMSD fluctuations of key binding site residues that interact with ligands to understand them further. In the PI3Kγ−6f complex, our docking study suggests that Trp(812), Tyr (867), and Asp(964) are the three key residues interacting with the ligand. The RMSD plot of these residues, along with ligand (6f), is shown in Figure 6a. Interestingly, in all three cases, the RMSD remains constant at around 1.5 Å, implying a stronger intermolecular interaction between these residues and ligands. This observation is consistent with our docking results, where the rewards (stabilization energy) for PI3Kγ are nearly the same as those seen in ATM (Table 6). However, such a shorter distance between active site residues and ligands may likely result in steric clashes between the atoms present, resulting in large penalties. Moreover, constant RMSD between key residues and ligands, with a magnitude of around 1.5 Å, exemplifies the rigidity of the PI3Kγ binding site. It is probably because of this rigidness that ligand molecules cannot undergo conformational interconversion, which was generally observed in the other four proteins. LipophilicEvdW-ChemScore lipophilic pair term and fraction of the total protein−ligand vdW energy; PhobEnHB-Reward for hydrophobically packed hydrogen bond; HBond-ChemScore H-bond pair term; Electro-Electrostatic reward; Penalties-Polar atom burials and desolvation penalties and penalty for intra−ligand contacts; ExposPenal-Penalty for exposed hydrophobic ligand groups.
As discussed above, due to the large penalties resulting from intraligand contact, PI3Kγ shows a reduced binding affinity for the ligand, which results in a lower docking score. We analyzed the RMSD fluctuations of key binding site residues that interact with ligands to understand them further. In the PI3Kγ−6f complex, our docking study suggests that Trp(812), Tyr (867), and Asp(964) are the three key residues interacting with the ligand. The RMSD plot of these residues, along with ligand (6f), is shown in Figure 6a. Interestingly, in all three cases, the RMSD remains constant at around 1.5 Å, implying a stronger intermolecular interaction between these residues and ligands. This observation is consistent with our docking results, where the rewards (stabilization energy) for PI3Kγ are nearly the same as those seen in ATM (Table 6). However, such a shorter distance between active site residues and ligands may likely result in steric clashes between the atoms present, resulting in large penalties. Moreover, constant RMSD between key residues and ligands, with a magnitude of around 1.5 Å, exemplifies the rigidity of the PI3Kγ binding site. It is probably because of this rigidness that ligand molecules cannot undergo conformational interconversion, which was generally observed in the other four proteins. PI3Kγ−6f ATM−6f In the ATM−6f complex, two interacting residues are Tyr (73) and Asp (207), forming π−π stacking and halogen bonding with the ligand, respectively. RMSD fluctuations for Tyr(73)−6f (ligand) and Asp(207)−6f(ligand) are shown in Figure 6b. For Tyr(73)−6f (ligand), the RMSD value shows a variation in the range of 1.5−2 Å, indicating a strong π−π stacking which remains stable during the simulation. Meanwhile, in the case of Asp(207)−6f(ligand), the RMSD shows larger variation, implying the dynamic formation and breakage of weaker halogen bonds, probably due to conformational interconversion of the ligand. Overall, the RMSD fluctuations between these residues and the ligand suggest a flexibility of binding sites, large enough to allow for the free conformational interconversion of the ligand. It is due to the flexibility of the binding site; the penalties due to intraligand contacts are negligible in the ATM−6f complex.
Similarly, for the ATR−6f complex and the mTOR−6f complex, the critical binding site residue shows similar fluctuations in RMSD ( Figure S9), implying binding site flexibility. This observation is also consistent with the ligands' conformational switching observed in all four proteins (except PI3Kγ, as discussed earlier). Furthermore, MSA analysis demonstrated that though the binding sites of all the proteins are confined to the kinase domain and its nearby region, the interacting residues were not a part of the conserved region. This factor can contribute to the differential selectivity of the inhibitors even though the nature of some of the residues were similar.
In summary, the docking affinity of 6f (ligand) with PIKK family kinases is determined by a subtle interplay of the natures of active site residues, stabilizing interactions, and penalties resulting from the rigidity/flexibility of the binding site.

Conclusions
One of the major problems in kinase drug discovery is the off−target kinase inhibition. Therefore, predicting the off−target binding of the inhibitors would help in fast−tracking the process of drug development. In this study, we analyzed the binding affinity of compound 6f against five different kinases sharing a reasonable sequence similarity. 6f was shown to bind selectively to ATM in comparison with the other proteins of interest. MD simulation showed that all of the protein−ligand complexes were stable and that secondary structure variations during the simulation were minimal. Interestingly, the trajectory analysis of the ligand (6f) in all cases suggested the ligand's conformational interconversion (except for in the case of the PI3Kγ−6f complex). The conformational changes observed in the ligand were attributed to the free rotation of the C-N bond and the binding site flexibility of the proteins. In the case of the PI3Kγ−6f complex, the binding site was found to be very rigidly packed, thus prohibiting the free rotation of the C-N bond. Moreover, the binding site rigidity in the PI3Kγ−6f complex resulted in steric clashes between active site residues and ligands, thus yielding a higher penalty and a lower docking score. Finally, MSA analysis also aided us in understanding the residue conservation across the proteins and rationalizing the selectivity of the molecules towards ATM.
To conclude, the molecules synthesized have higher binding affinity and selectivity towards ATM compared to the other kinases of interest. The proteins share a reasonable sequence similarity in the kinase domain where the designed inhibitors are speculated to bind. The interacting residues do not fall in the conserved region, thereby resulting in selectivity towards ATM. This result can be considered as a preliminary investigation, which warrants further validation using in vitro assays.