In Silico Investigation of Potential Src Kinase Ligands from Traditional Chinese Medicine

Src kinase is an attractive target for drug development based on its established relationship with cancer and possible link to hypertension. The suitability of traditional Chinese medicine (TCM) compounds as potential drug ligands for further biological evaluation was investigated using structure-based, ligand-based, and molecular dynamics (MD) analysis. Isopraeroside IV, 9alpha-hydroxyfraxinellone-9-O-beta-D-glucoside (9HFG) and aurantiamide were the top three TCM candidates identified from docking. Hydrogen bonds and hydrophobic interactions were the primary forces governing docking stability. Their stability with Src kinase under a dynamic state was further validated through MD and torsion angle analysis. Complexes formed by TCM candidates have lower total energy estimates than the control Sacaratinib. Four quantitative-structural activity relationship (QSAR) in silico verifications consistently suggested that the TCM candidates have bioactive properties. Docking conformations of 9HFG and aurantiamide in the Src kinase ATP binding site suggest potential inhibitor-like characteristics, including competitive binding at the ATP binding site (Lys295) and stabilization of the catalytic cleft integrity. The TCM candidates have significantly lower ligand internal energies and are estimated to form more stable complexes with Src kinase than Saracatinib. Structure-based and ligand-based analysis support the drug-like potential of 9HFG and aurantiamide and binding mechanisms reveal the tendency of these two candidates to compete for the ATP binding site.


Introduction
Src kinases are nonreceptor tyrosine kinases that are of physiological importance in cell survival, bone metabolism, angiogenesis, proliferation, migration, and invasion [1]. Overexpression of Src kinase has been linked to various cancers and is now a well-established proto-oncogene [2][3][4][5][6][7]. The physiological pathway involved in hypertension is also associated with Srcdependent signaling pathways, suggesting a potential link between hypertension and Src [8][9][10][11][12][13][14][15]. Figure 1 illustrates the components of Src kinase and its activation mechanism [16,17]. In general terms, the catalytic activity of Src is co-regulated by SH3 and SH2 domains. Src is locked in the closed conformation (inactive) when SH2 binds to the phosphorylated Tyr530, and SH3 binds with prolines on the linker domain ( Figure 1A). When Tyr530 is dephosphorylated, Src assumes an open conformation, achieving full activity when Tyr416 within the catalytic domain is autophosphorylated ( Figure 1B). This opening of the Src structure frees the SH2 and SH3 domains to interact with surface receptors such as focal adhesion FAK and initiate downstream signaling governing the aforementioned physiological pathways [1]. In this regard, inactivation of Src can be achieved through hindering disassembly of the regulatory SH2 and SH3 Src domains, or by inhibiting ATP binding to the Src catalytic site [16,18].
Many small molecular Src inhibitors have been identified due to the involvement of Src in cancer. Comprehensive reviews on such advancements are detailed elsewhere [16]. Most Src inhibitors discovered to date are Type I inhibitors that compete with ATP for binding at the ATP binding pocket [19,20]. Structures of the three most studied Type I inhibitors Bosutinib, Dasatinib, and Saracatinib and are shown in Figure 2 along with their respective status in clinical trials [21][22][23][24][25][26][27][28][29][30][31][32]. The varying efficacies of these commercial drugs highlight the need for novel compounds that can exhibit more consistent inhibition of Src.
The goal of this research is to investigate novel small compounds from traditional Chinese medicine (TCM) that may be potential Src kinase ligands. During the past decade, our laboratory has focused on constructing the most comprehensive TCM database (TCM Database@Taiwan, url:) (http://tcm.cmu.edu.tw/) [33]. In addition, we've also developed the first cloud-computing webserver based on TCM Database@Taiwan, url: (http://iscreen.cmu.edu. tw/) [34] and an integrative website combining TCM and systems biology (http://iSMART.cmu.edu.tw/) [35]. Utilizing these TCM computational resources, several novel lead compounds from TCM with application potential for different diseases have been successfully uncovered [36][37][38][39][40][41][42][43][44]. In the current research, we utilize the newly updated TCM Database@Taiwan to screen for novel, TCM-origin ligands with drug-like properties against Src kinase.

Results and Discussion
Docking Based on the DockScore, top ranking TCM candidates selected for further investigation were Isopraeroside IV, 9alpha-hydroxyfraxinellone-9-O-beta-D-glucoside (9HFG), and aurantiamide (Table 1). Isopraeroside IV is a coumarin isolated from the root of Angelica dahurica [45]. 9HFG originates from the root of Dictamnus dasycarpus, which has been used as a folk remedy for inflammation and skin ailments [46]. Aurantuamide is a composition of Curcuma wenyujin, a dried rhizome used for bleeding and menstrual disorders in traditional Chinese medicine. Structural comparisons of the TCM candidates with Saracatinib are illustrated in Figure 3.
DockScore is the negative sum of the ligand/receptor interaction and the ligand internal energy calculated by DS. For clarification purposes, the ligand internal energy consists of a van der Waals (vdW) term computed using a standard 9-6 (unsoftened) potential and an optional electrostatic term. Since TCM candidates and Saracatinib have similar predicted binding energies, differences in DockScore were primarily affected by the respective calculated ligand internal energies of each compound ( Table 1). The low DockScore observed for Saracatinib is mainly due to its high calculated ligand internal energy. The positive ligand internal energy value is indicative of a highly strained chemical structure, a phenomenon commonly observed among synthetic compounds. By contrast, low calculated ligand internal energies of the TCM compounds suggest energetically stable structures in which intermolecular atoms are not too close to exhibit replusion and result in high DockScore values. Based on DS2.5 calculations, the TCM candidates have higher DockScores and binding energies comparable to Saracatinib, which may lead to equal if not better binding to the Src kinase binding site compared to Saracatinib.
The ADMET properties of the TCM candidates are summarized in Figure 4 and Table 2. All candidates were estimated to be within the 99% adsorption ellipse ( Figure 4) and have good to moderate adsorption ( Table 2). The TCM candidates are more highly bound to plasma proteins, and less likely to inhibit cytochrome P4502D6 metabolism.
As our designated binding site overlapped with the Src kinase ATP binding pocket, it is referred to as the latter for clarification purposes. Docking poses of the candidates in the ATP binding pocket are illustrated in Figure 5. The docking model of Saracatinib indicates p-interactions with Leu273 and Lys295, and an H-bond with Met341 within the ATP binding pocket ( Figure 5A). Isopraeroside IV formed a total of five H-bonds with Met341, Ser345, Asp348 ( Figure 5B), possibly contributing to the high estimated binding energy in Table 1. The three hydroxyl groups on the cyclohexane moiety of 9HFG interacted with Lys295, Met341, and Asp404, forming a total of five H-bonds ( Figure 5C) in which three were formed with Lys295. Auraniamide interacted with Leu273, Lys295, Ser345, Asp348, forming two p interactions at the terminal benzenes with Leu273 and Lys295 and H-bonds with Ser345 and Asp 348 ( Figure 5D).
Ligplot diagrams illustrating hydrophobic and H-bond interactions are shown in Figure 6. In general, the higher ability of ligands to form hydrophobic interactions with hydrophobic amino acids of the binding site, the higher the binding affinity. As shown in Figure 6, TCM candidates are surrounded by hydrophobic amino acids, therefore are harder to remove from the ATP binding pocket under physiologically dynamic conditions. Less hydrophobic interactions were formed by Saracatinib. Minor variations in H-bond formation were observed in LigPlot diagrams compared to 3D-docking poses generated by DS 2.5 ( Figure 5). Additional Hbonds were formed with Tyr340 and Glu399 by Isopraeroside IV ( Figure 6B). 9HFG formed two additional H-bonds at Ala390 and Asn391 ( Figure 6C). In auraniamide ( Figure 6D), the H-bonds at Leu273 and Lys295 were replaced by hydrophobic interactions. These differences may be due to the DS cutoff distance of 2.5 Å or to the HBPLUS program used to calculate H-bonds in LigPlot. The number of H-bonds and hydrophobic interactions determined through LigPlot is summarized in Table 3. These simulation results suggest that TCM candidates have a higher number of stabilizing interactions than Saracatinib.
Key amino acids and important structural moieties can be summarized based on the docking simulation results. Leu273 is important for H-bonds and p-interactions, and ligands form interactions with its nonpolar, aliphatic alkyl group. Saracatinib, 9HFG, and auraniamide interact with Src kinase through Lys295. The importance of Lys295 is reinforced by the three H-bonds and two p-interactions formed by the aforementioned ligands. Other amino acids serving as anchor points in the Src kinase ATP binding pocket include Met341, Ser345, Asp348, and Asp390. The chemical structures of the candidates facilitate hydrophobic interactions with surrounding amino acids, and enable higher ligand-protein affinity and complex stability.

Bioactivity Prediction by SVM and MLR
The 75 descriptors relevant to the training set were determined and GFA further employed to determine the top ten descriptors that could most accurately describe the training set. The descriptors thus determined were: C_Count (the amount of carbon atoms within the ligand), ES_Sum_sCH 3 (sum of the electrotopological state (E-state) values for carbons with single bonds), E_ADJ_equ (molecule differentiation based on edge adjacency), CHI_2 (number of pairs of bonds within a molecule), CHI_3_P (number of phosphate molecules with three bonds), Kappa_2 (shape index of order 2 describing topological descriptors), Kappa_3 (shape index of order 3 describing topological descriptors), Jurs_FNSA_1 (sum of negative surface areas divided by total molecular solvent-accessible surface area), Jurs_PNSA_1 (sum of solvent-accessible surface areas of all negatively charged atoms), Shadow_nu (characterizes molecule shape using ratio of largest to smallest dimension).
Using these descriptors, the MLR model generated was: As shown in Figure 7, the square correlation coefficients for the MLR and SVM models were 0.837 and 0.7922, respectively, indicating reliable models. Predicted bioactivities (pIC 50 ) of the TCM candidates by the generated MLR and SVM models are listed in Table 2. Both models predicted higher pIC 50 values for Isopraeroside IV and 9HGF than aurantiamide. The predicted bioactivities for Saracatinib was 4.9410 using MLR and 7.3036 using SVM. These predicted results, when compared to the in vitro inhibition concentration of Saracatinib (pIC 50 = 8.57) [47], indicate that values predicted by the SVM model may be closer to actual inhibition concentrations observed in vitro. The higher accuracy of the SVM model is expected as its predictions are generated on a non-linear model. Bioactivity Prediction by 3D-QSAR Table 4 summarizes the PLS results of CoMFA and CoMSIA models. Steric fields are the dominant factor in CoMFA. Under an optimal number of components (ONC) of 6, cross-validated correlation coefficient (q 2 ) and non-cross-validated correlation coefficient (r 2 ) of 0.430 and 0.920 indicate a reliable model. Multiple CoMSIA models were generated using different factor combinations ( Table 4). Models SHD and SEHD were the best overall models with high q 2 (q 2 .0.4), high r 2 , low SEE, and large F values. However, addition of an electrostatic variable in model SEHD did not significantly increase model correlation or strength (F-value), thus model SHD was selected as the optimum CoMSIA model. Reliability of the selected CoMFA and CoMSIA models are validated by the q 2 and r 2 values ( Table 4).
The observed, predicted and calculated residual pIC 50 values of the training and test set compounds by CoMFA and CoMSIA are summarized in Table S1. Differences between observed and predicted pIC 50 ranged between 20.667 and 1.066. Residuals using the CoMSIA model ranged from 20.5 to 0.821, indicating high accuracy of the predictions from our generated model with actual biologically verified values. All ratio values calculated through CoMSIA were within 160.1 indicating high accuracy. The correlation coefficient between the observed and predicted pIC 50 using the CoMFA model was 0.8448 ( Figure 8A). With the CoMSIA model, the correlation coefficient was 0.9014 ( Figure 8B). Both models are highly reliable.
Superimposing CoMFA and CoMSIA contour maps on Saracatinib and the TCM candidates provide insights into the predicted bioactivity of the compounds. Saracatinib ( Figure 9A), Isopraeroside IV ( Figure 9B), and 9HFG ( Figure 9C) contoured well to the CoMFA model. No bulk structures were located in the steric disfavoring (yellow) region. In addition, bulky ring moieties in Saracatinib and 9HFG were located in the steric favoring (green) region. By comparison, a benzene ring of aurantiamide was located within the steric disfavoring region, suggesting a lower biological activity due to its deviation from the CoMFA contour map ( Figure 9D). These results agree with the MLR and SVM results in which aurantiamide had the lowest predicted bioactivity   among the tested compounds. With regard to CoMSIA maps, Saracatinib forms H-bonds with Met341 which is located near the region favoring hydrogen bond donors (purple) and hydrophobic interactions (cyan) ( Figure 10A). Isopraeroside IV was a smaller compound and its interactions with Src did not fall within the disfavoring regions of the CoMSIA map ( Figure 10B). A similar contour was observed in 9HFG ( Figure 10C). In aurantiamide, the benzene moieties were in proximity to the steric disfavoring regions, suggesting a lower bioactivity of the ligand ( Figure 10D). Results of the ligand-based studies are consistent, and suggest high bioactivity for Saracatinib, moderate to high bioactivity for Isopraeroside and 9HFG, and moderate to low bioactivity for aurantiamide.

Molecular Dynamics Simulation
More realistic interactions of the candidates with Src kinase in a biological system were simulated through MD. Interaction differences at key residues for each tested compound during docking and MD are summarized in Table 5. H-bond formation and occupancies during and the H-bond distance profiles are presented in Table 6 and Figure 11, respectively. The primary binding residue for Saracatinib was Lys295 where stable pinteractions and a high occupancy H-bond were formed ( Figure 11A). The H-bond observed at Met341 during docking was substituted by an H-bond with Thr338 throughout MD. The p interaction with Leu273 during docking was not observed during MD. As illustrated in Figure 12A, the head-on direction of the H atom and the near perpendicular spatial arrangements facilitate the formation of p-sigma interaction between Leu273 and Saracatinib. During MD, angle increase due to ligand and Leu273 fluctuations disrupted the formation of p-sigma interactions ( Figure 12A).
Isopraeroside IV formed stable H-bonds with Met341, Lys343, and Ser345 (Table 6). H-bonds formed with Tyr340 ( Figure 11B) and Asp348 ( Figure 11C) were not stable and could not contribute to stability of Isopraeroside when in complex with Src kinase. Sharp decreases in H-bond distances observed for Lys343 and Ser345 ( Figure 11B-C) are likely due to rotations of the residues which bring the ligand into closer proximity for H-bond formation ( Figure 13A).
The highest number of H-bonds was formed by 9HFG with Lys295 and Asp404 (Table 5-6) Rhythmic fluctuations of H-bond distances with the key binding residue Lys295 ( Figure 11D-F) are likely due to rotations of the amine group on Lys295 ( Figure 13B). 9HFG also formed two H-bonds with Asp404, but the H-bond with Ala390 during docking was lost ( Figure 11G).
Aurantiamide was primarily bound to the ATP binding pocket through the p interaction at Lys295 and H-bonds with Ser345 and Asp348 (Table 5). Increased distance fluctuations in Ser345 ( Figure 11H) could be due to the rotation of benzene moieties ( Figure 13C). Similar to Saracatinib, the loss of the p interaction at Leu273 involved torsion of the benzene moiety ( Figure 12B). Benzene ring torsion during MD reduced angles between the H       atom and the benzene normal plane to approximately 90 degrees, making the formation of p-interactions unlikely ( Figure 12B). In general, MD results were in agreement with docking (Table 5). Torsion angles also provide molecular insights into ligand stability ( Figure 14). Torsion angle fluctuations critical to the stability of Saracatinib were not observed ( Figure 14A). Torsion angles recorded for Isopraeroside IV were relatively stable with the exception of c and d ( Figure 14B). Notable torsion angle changes at c and d, may account for the decrease in H-bond distances observed for Tyr340 (6.34 ns) ( Figure 11B). Subsequent increase in bond distance and resulting loss of the Tyr340 H-bond could be due to sharp angle changes of d from 17.44 ns to the end of MD simulation. Torsion angles for 9HFG provide supporting evidence for the cause of H-bond fluctuations in Figure 11D-G. Since the stable torsion angles at c and e ( Figure 14C) imply a relatively stable ligand, fluctuations in bond distance ( Figure 11D-F) should be due to changes in Lys295. Further evidence is given by torsion angle changes at f. As f is associated with the binding of O9 with Lys295 and H39 with Asp404, should ligand stability be the cause of fluctuation, bond distances at both locations should be equally affected. However, H39 forms stable H-bond with Asp404, but O9     does not form stable bonding with Lys295, validating that the primary reason for the lack of stable binding with Lys295 is due to the instability of the amino acid itself. Among the torsion angles recorded for aurantiamide ( Figure 14D), only a, d, and l are relatively unstable with angle fluctuations greater than 50 degrees. Stable H-bonds are formed as a result of the stability of the ligand. Figure 15 summarizes the RMSD and total energy trajectories of the ligands during MD. All compounds, with the exception of 9HFG from 20-24 ns, had relatively stable ligand RMSDs. Whole molecule RMSDs and total energies reached equilibration after 16 ns. Among the four test compounds, Saracatinib had the highest complex RMSD and total energy. The TCM candidates exhibited similar complex RMSDs, but aurantuamide had higher total energy than Isopraeroside IV and 9HFG after stabilization. The higher energy and complex RMSD of Saracatinib may be associated with its high internal ligand energy (Table 1). Once bound, the internal energy is transferred to the complex, causing more fluctuations and higher energy levels. Isopraeroside IV and 9HFG are smaller structures compared to Saracatinib. The lower energy profiles and complex RMSD may be attributed to the compact structures and formation of multiple H-bonds. Flexible moieties in Isopraeroside IV are anchored to the ATP binding pocket through Hbonds, limiting the fluctuations of the complex and increasing stability. By contrast, H-bonds formed by 9HFG stabilizes the phenyl moiety, but the remaining structure is unbound and can freely rotate. This could be the reason for high ligand RMSDs observed for 9HFG ( Figure 15B). Aurantuamide is a larger compound that shares similarities with Saracatinib. As previously discussed, aurantuamide interacts with Src kinase ATP binding pocket through a p interaction with Lys295 and an H-bond at Ser345, leaving two benzene moieties free for rotation. The larger ligand RMSD and higher total energy is likely due to these rotations. In summary, all four test compounds formed stable complexes at the ATP binding pocket for the duration of MD. TCM candidates exhibited higher stability than Saracatinib, indicating advantages as substitutes for the commercial drug.

Assessment of TCM candidate mode of action
Despite promising drug-like characteristics, one must consider binding mechanisms of the ligand compounds to objectively speculate their effect on Src kinase. As described previously, most Src inhibitors, including Saracatinib, compete with ATP for binding to Lys295 (Type I inhibitors). Src kinase is activated when interactions between SH2 and SH3 are broken, and when Tyr416 becomes phosphorylated [18] (Figure 1). Among our TCM candidates, 9HFG and aurantiamide show Src Type I inhibitorlike docking characteristics. Aurantiamide interacts primarily with Lys295, and thus its mode of inhibition would also be through competitive binding against ATP for Lys295 ( Figure 16A). However, its lower total energy ( Figure 15) implies a more stable protein-ligand complex than Saracatinib. Intriguingly, 9HFG is a potential candidate that mimics the natural autoinhibitory mechanisms in Src kinases. Not only does 9HFG bind strongly to Lys295, but it simultaneously binds to Asp404 as well ( Figure 16A). This is an important feature as it mimics the bridging function of Mg 2+ between these two amino acids in native Src kinases and stabilizes the structural integrity of the catalytic cleft ( Figure 16B) [18]. The dual function of 9HFG as a binding site competitor and stabilizer of the catalytic cleft suggests good potential for Src kinase inhibition. On the contrary, Isopraeroside IV may not affect Src kinase activity despite its good binding affinity and predicted bioactivities. Isopraeroside IV binds to regions located deeper within the catalytic cleft than Lys295 (Table 3, Table 5) and does not interact directly with Lys295. Though Isopraeroside IV may limit the ability of ATP to bind to Lys295 through steric hindrance, there is no evidence from our structure-based analysis to support such a claim.

Conclusion
The potential of Isopraeroside IV, 9HFG, and aurantiamide as potential leads for Src-kinase was assessed through structure-and ligand-based in silico approaches. Docking and MD suggest that Isopraeroside IV, 9HFG, and aurantiamide have higher and more stable binding affinities with Src kinase than Saracatinib. Assessment of biological activity through ligand structures further implied that the TCM candidates were biologically active compounds. Considering the interaction between the TCM candidates and the Src ATP binding region, 9HFG, and aurantiamide may exert inhibition against Src kinase by competitive binding against ATP. Our results may be applicable to drug development in different ways. TCM compounds 9HFG and aurantiamide may be directly used as candidate lead compounds in biological studies based on their high stability and predicted bioactivities. In addition, structural insights such as key amino acids and important moieties which form stable interactions can be utilized for de novo synthesis of more stable compounds that interact with Src kinase.

Docking
The protein structure of human Src kinase (PDB: 2H8H) was obtained from Protein Data Bank [47]. The binding site used in this study was based on the space occupied by Saracatinib within 2H8H. The binding site was located within the cleft separating the carboxyl-terminal lobe (residues 345-523) and the amino-terminal lobe (residues 270-340), and surrounded by residues Lys295, Trp340, Met341, Lys343, Gly344, Ser345, Leu373, Leu393, and Asp404. Over 20,000 ligands from TCM Database@Taiwan were used for docking in the Src kinase binding site. All ligands were pre-treated with CHARMm [48] to attach missing H-atoms. Saracatinib was used as the control. TCM ligands were docked into the binding site using the LigandFit program [49] within Discovery Studio 2.5 and DockScore selected as the primary scoring function for assessing binding affinities. TCM ligands were further screened with Lipinski's Rule of Five [50,51] and absorption, distribution, metabolism, excretion and toxicity (ADMET) [52] in DS 2.5 to rule out potentially toxic derivatives.  Attraction forces (H-bond and hydrophobic interactions) between the ligand and protein were analyzed using LigPlot v.2.2.25 [53]. The 3D protein-ligand complexes generated during docking were flattened to 2D diagrams through the ligplot algorithm.
Bioactivity Prediction by Support Vector Machine (SVM) and Multiple Linear Regression (MLR) DS 2.5 was used to calculate individual molecular property descriptors of 53 Src inhibitors with known pIC 50 values [54]. Representative descriptors for the inhibitors were determined from a pool of 552 descriptors through genetic function approximation (GFA) [55,56]. A QSAR models were calculated by GFA and ranked by square correlation coefficient (R 2 ) corresponding to the fitness of each model. Descriptors from the model with the highest R 2 were used to construct MLR and SVM models for predicting the bioactivity of TCM candidates. The linear MLR model [57] was constructed with the representative descriptors using MATLAB [58].
In addition to the linear MLR model, a non-linear QSAR model utilizing support vector machine (SVM) [59,60] for the regression of continuous bioactivity data (pIC 50 ), also termed support vector regression (SVR), was also constructed. SVR uses the Kernel function to map the input data set into a highdimensional feature space in order to find the hyperplane that best predicts data distribution. The SVR formulation is min w,bj,jÃ s:t: The parameters (x 1 , y 1 ), …, (x l , y l ) represent data within a given training set and can also be expressed as f(x) = ,w,x.+b, The parameter w is calculated by performing e-  Table 5). Blue lines indicate the torsion angles recorded; the red and gray lines indicate the angle during docking and at time 0 of MD, respectively. doi:10.1371/journal.pone.0033728.g014 insensitive loss function. The performance of SVR is dependent on parameters C,e, the kernel function and kernel parameters. To further optimize this to nonlinear functions, Lagrange multipliers and kernel k(x,x9) = Sw(x),w(x 0 )T are used. The nonlinear functions and predictions are determined by: max: f (x)~X l i~1 (a i {a Ã i )k(x i ,x)zb s:t: 0ƒa Ã i ƒC,0ƒa i ƒC ð4Þ The SVM model was constructed by applying the LibSVM [61] within the gridregression.py program to determine key parameters C cost, epsilon and gamma.

Bioactivity Prediction by 3D Quantitative Structure-Activity Relationship (QSAR) Analysis
The 53 Src kinase inhibitors [54] used in this study were randomly divided into a training set of 43 compounds and a test set of 10 compounds. Alignment of the training set molecules was performed using atom-fit module of SYBYL-X 1.1 [62] Comparative force field analysis (CoMFA) and comparative similarity indices analysis (CoMSIA) models were constructed for 3D-QSAR. The steric and electrostatic field descriptors in CoMFA were calculated using Lennard-Jones potential and Coulombic potential, respectively. For CoMSIA models, the steric, electrostatic, hydrophobic and hydrogen bond donor and acceptors were calculated through Gaussian functions. Partial least squares (PLS) regression was utilized to analyze the 3D-QSAR models with descriptor used as an independent variable. Structural contour to the generated CoMFA and CoMSIA models were used to predict bioactivity of TCM candidates.

Molecular Dynamics (MD) Simulation
Src kinase complexes of the TCM candidates were subjected to molecular dynamics (MD) simulation under the force field of CHARMm [48] using DS 2.5. The energy of each complex was minimized using 500 steps of Steepest Descent and 500 steps of Conjugate Gradient method. The system was heated without constraint from 50 K to 310 K in 50 ps, and equilibrated for 200 ps. The final production was conducted for 40 ns in NVT ensemble with snapshots saved at 2.5 ps intervals. The time step was set to 2 fs. Hydrogen bond frequencies, energy trajectories, and hydrogen bond distances were calculated via DS 2.5 Analyze Trajectory module to analyze molecular interactions within the protein-ligand system.

Supporting Information
Table S1 Predicted and observed activity (pIC 50 ) of 53 known Src kinase inhibitors 1 using generated CoMFA and CoMSIA models. (DOC)

Acknowledgments
We are grateful to the Asia University cloud-computing facilities and Dr. Su-sen Chang for technical assistance in preparing the manuscript.