QSAR Study, Molecular Docking and Molecular Dynamic Simulation of Aurora Kinase Inhibitors Derived from Imidazo[4,5-b]pyridine Derivatives

Cancer is a serious threat to human life and social development and the use of scientific methods for cancer prevention and control is necessary. In this study, HQSAR, CoMFA, CoMSIA and TopomerCoMFA methods are used to establish models of 65 imidazo[4,5-b]pyridine derivatives to explore the quantitative structure-activity relationship between their anticancer activities and molecular conformations. The results show that the cross-validation coefficients q2 of HQSAR, CoMFA, CoMSIA and TopomerCoMFA are 0.892, 0.866, 0.877 and 0.905, respectively. The non-cross-validation coefficients r2 are 0.948, 0.983, 0.995 and 0.971, respectively. The externally validated complex correlation coefficients r2pred of external validation are 0.814, 0.829, 0.758 and 0.855, respectively. The PLS analysis verifies that the QSAR models have the highest prediction ability and stability. Based on these statistics, virtual screening based on R group is performed using the ZINC database by the Topomer search technology. Finally, 10 new compounds with higher activity are designed with the screened new fragments. In order to explore the binding modes and targets between ligands and protein receptors, these newly designed compounds are conjugated with macromolecular protein (PDB ID: 1MQ4) by molecular docking technology. Furthermore, to study the nature of the newly designed compound in dynamic states and the stability of the protein-ligand complex, molecular dynamics simulation is carried out for N3, N4, N5 and N7 docked with 1MQ4 protease structure for 50 ns. A free energy landscape is computed to search for the most stable conformation. These results prove the efficient and stability of the newly designed compounds. Finally, ADMET is used to predict the pharmacology and toxicity of the 10 designed drug molecules.


Introduction
Cell division is a biological marker, which is strictly regulated by a large number of proteins.The occurrence of cancer is due to the uncontrolled division, growth and spread of cancer cells.These cells go beyond their usual boundaries, then attack the adjacent parts of the body and spread to other organs.This can occur in any part of the body, leading to the death [1].Using scientific methods to prevent and control cancer has become one of the most important public health targets in the world [2].In the regulatory protein network of cell division, Aurora kinase plays an important role in cell division by controlling chromatin separation, which can lead to genetic instability, often accompanied by cancer [3].
Aurora kinase, a member of serine/threonine family, is responsible for cell cycle regulation.Aurora kinase is involved in many mitotic processes, such as centrosome maturation and separation, spindle formation and stability, chromosome separation and cytokinesis, so as to coordinate the orderly mitosis.Abnormal expression or function of Aurora kinase Aurora kinase, a member of serine/threonine family, is responsible for cell cycle regulation.Aurora kinase is involved in many mitotic processes, such as centrosome maturation and separation, spindle formation and stability, chromosome separation and cytokinesis, so as to coordinate the orderly mitosis.Abnormal expression or function of Aurora kinase can lead to a disorder of genetic material distribution, which can affect the stability of the genome [3].In recent years, relevant studies have confirmed that Aurora kinase is highly expressed in a variety of tumor cells and participates in the formation of tumors [4,5].Aurora A is a subtype of Aurora kinase.The Aurora A gene is located in the chromosomal region with active amplification.At present, many Aurora A compounds with different inhibitory activities have been found and show effective anticancer activities [6,7].They involves a wide range of cancers, such as colorectal cancer, prostate cancer, ovarian cancer, breast cancer and glioma [8][9][10][11].Therefore, Aurora A is considered to be a good marker of tumor progression and prognosis, and the inhibition of its kinase activity helps to reduce tumor invasiveness.(Figure 1 shows several kinds of Aurora kinase inhibitors entering clinical trials, including the first Aurora kinase inhibitor ZM-447439 (Figure 1(1)), the first Aurora B inhibitor, hesperidin (Figure 1(2)), and a new type of Aurora inhibitor, MLN-8237 (Figure 1(3)).Despite the great importance of Aurora kinases in antitumor research, the specific processes and mechanisms of the involvement of Aurora kinases in tumor formation are still unclear.Due to safety concerns, there is still no marketed class of Aurora kinase inhibitors.Therefore, research on the more potent and safer class of Aurora kinase A inhibitor anticancer drugs is still challenging.A series of imidazo [4,5-b]pyridine derivatives have been reported to possess excellent potency as orally bioavailable Aurora kinase A inhibitors [12,13].
Quantitative structure-activity relationship (QSAR) is a common method of drug design.It is usually used as a common research method to explore the relationship between molecular structure and biological activity.Since its inception, QSAR has been widely used in drug design.The mechanism is to develop the relationship between chemical structure, physicochemical properties and biological activity of compounds in order to obtain a reliable statistical model for predicting the activity of new chemical entities [14].The quantitative relationship between the chemical molecular structure and molecular biological activity is analyzed using a mathematical model.The quantitative relationship expression between activity and structure is established with molecular structure parameters as independent variables and molecular biological activity as dependent variables [15,16].Based on the established model, through the use of virtual screening, molecular docking, ADMET, molecular dynamic simulation and other means, reasonable drug design can be carried out.Despite the great importance of Aurora kinases in antitumor research, the specific processes and mechanisms of the involvement of Aurora kinases in tumor formation are still unclear.Due to safety concerns, there is still no marketed class of Aurora kinase inhibitors.Therefore, research on the more potent and safer class of Aurora kinase A inhibitor anticancer drugs is still challenging.A series of imidazo [4,5-b]pyridine derivatives have been reported to possess excellent potency as orally bioavailable Aurora kinase A inhibitors [12,13].
Quantitative structure-activity relationship (QSAR) is a common method of drug design.It is usually used as a common research method to explore the relationship between molecular structure and biological activity.Since its inception, QSAR has been widely used in drug design.The mechanism is to develop the relationship between chemical structure, physicochemical properties and biological activity of compounds in order to obtain a reliable statistical model for predicting the activity of new chemical entities [14].The quantitative relationship between the chemical molecular structure and molecular biological activity is analyzed using a mathematical model.The quantitative relationship expression between activity and structure is established with molecular structure parameters as independent variables and molecular biological activity as dependent variables [15,16].Based on the established model, through the use of virtual screening, molecular docking, ADMET, molecular dynamic simulation and other means, reasonable drug design can be carried out.

Dataset
The 3D conformations of all compounds in molecular construction are constructed using the SkechTool module in SYBYL2.0.The standard Tripos molecular force field and Powell energy gradient algorithm are used for energy optimization.All molecules are loaded with Gasteger-Huckel charges, the maximum number of iterations is 1000, the energy convergence limit is set to 0.005 kcal mol −1 , and the other parameters adopt the system default values [17,18].
The dataset consisted of 65 imidazo [4,5-b]pyridine derivatives, all of which were synthesized and evaluated for biological activity by Vassilios Bavetsias and colleagues [19][20][21].Compounds with similar structural framework and inhibitory ability were selected.The IC 50 was the inhibitory activity of imidazo [4,5-b]pyridine derivatives against Aurora kinase A. For the convenience of calculation, these IC 50 values were converted into corresponding pIC 50 values to characterize their biological activities and used as the dependent variable of QSAR research.The conversion formula is: pIC 50 = −lg IC 50 .In order to ensure the randomness of the data, we extracted 16 molecules from 65 molecules according to the ordinal number (one out of four) as the test set for external verification, and the other 49 molecules as the training set for model construction.The ratio of the test set to training set is 1:4.The main function of the training set is to build a QSAR model, and the test set is mainly used to externally verify the results of the QSAR model and evaluate the prediction ability of the model.In addition, when selecting the training set and test set, we should comprehensively consider the activity distribution of the whole molecular set, so as to ensure the effectiveness of the data and avoid its specificity (Figure 2).The structures and activities of all the compounds in the dataset are presented in Table 1.

HQSAR
Hologram Quantitative Structure-Activity Relationship (HQSAR) is a technique between 2D-QSAR and 3D-QSAR, which determines the relationship between biological activity and structural fragments and does not require molecular superimposition, 3D structure or choice of active conformation.It relies on 2D chemical database storage and the linear notation that defines the chemical structure, a process which involves the generation of fragments hash into an array called molecular holograms [22,23].The ability of molecular alignment and conformational specification is achieved by converting the representation of the chemical structure of a molecule into the corresponding molecular hologram [24,25].According to the principle of the holograms phase diagram, the effect of any group or atom in the drug molecule on the drug activity can be precisely determined, and the operation process is relatively fast [26].
The HQSAR method usually consists of several main steps.First, the molecule under study is cut into fragments of the appropriate size.The fragment structure is mainly determined by two parameters, Fragment Size and Fragment Distinction.Then, the encoded molecular fragment is converted into molecular holograms.Finally, after the molecular holograms are obtained, a quantitative model of the relationship between molecular holograms and compound properties is established by partial least squares regression.

Dataset
The 3D conformations of all compounds in molecular construction are constructed using the SkechTool module in SYBYL2.0.The standard Tripos molecular force field and Powell energy gradient algorithm are used for energy optimization.All molecules are loaded with Gasteger-Huckel charges, the maximum number of iterations is 1000, the energy convergence limit is set to 0.005 kcal mol −1 , and the other parameters adopt the system default values [17,18].
The dataset consisted of 65 imidazo [4,5-b]pyridine derivatives, all of which were synthesized and evaluated for biological activity by Vassilios Bavetsias and colleagues.[19][20][21].Compounds with similar structural framework and inhibitory ability were selected.The IC50 was the inhibitory activity of imidazo [4,5-b]pyridine derivatives against Aurora kinase A. For the convenience of calculation, these IC50 values were converted into corresponding pIC50 values to characterize their biological activities and used as the dependent variable of QSAR research.The conversion formula is: pIC50 = -−lg IC50.In order to ensure the randomness of the data, we extracted 16 molecules from 65 molecules according to the ordinal number (one out of four) as the test set for external verification, and the other 49 molecules as the training set for model construction.The ratio of the test set to training set is 1:4.The main function of the training set is to build a QSAR model, and the test set is mainly used to externally verify the results of the QSAR model and evaluate the prediction ability of the model.In addition, when selecting the training set and test set, we should comprehensively consider the activity distribution of the whole molecular set, so as to ensure the effectiveness of the data and avoid its specificity (Figure 2).The structures and activities of all the compounds in the dataset are presented in Table 1.In our study, the atoms of these fragments are encoded in holograms in different colors to reflect their contribution to biological activity, with negative, intermediate, and positive contributions represented in red, white and green, respectively.Each fragment is defined by its unique characteristic parameters, which are: atom (A), chemical bond type (B), atom connection (C), hydrogen atom (H), chirality (CH), and as donor or acceptor (DA).All the feature fragments of each molecule are mapped into a certain length integer string, the process is called molecular holography, and the length of integer string is called the molecular holographic length parameter.In HQSAR, the module provides 12 prime numbers (53,59,61,71,83,97,151,199,257,307,353 and 401) as a holographic length.Multiple HQSAR models with different predictive abilities can be obtained by combining different parameters.Therefore, we try to adopt as many combinations as possible to seek the optimal HQSAR model.
The method translates the molecular structure into characteristic molecular fingerprints, which consist of different types of molecular fragments and are labelled with 12 prime numbers (53, 59, 61, 71, 83, 97, 151, 199, 257, 307, 353 and 401).These labelled numbers were used as QSAR descriptors for modelling the biological activity prediction of the compounds.Meanwhile the HQSAR method combines different fragment descriptors such as atom (A), bond (B), linkage (C), chirality (Ch), hydrogen atom (H), donor and acceptor (DA), etc., and by combining different parameters, multiple HQSAR models with different predictive capabilities can be obtained.Therefore, in our study we consider both dimensionality and holographic length adjustments and try to use as many combinations as possible in order to build HQSAR models with good predictive ability.

CoMFA, CoMSIA and TopomerCoMFA
Comparative Molecular Field Analysis (CoMFA) is one of the most used methods to study the 3D-QSAR between drug and receptor.CoMFA can fully consider the 3D structural information of molecules, characterize the 3D field and static electric field through the molecular structure, conduct regression analysis of drug activity with these parameters as variables, and obtain the relationship between these 3D characteristic information and compound activities by partial least-squares method [27].By comparing with the series of molecular space near the point of spatiality, by selecting a public skeleton to artificial molecule composite, electrostatic potential, and the physical and chemical parameters of molecules with biological activities can be calculated according to the molecular force field around the probe atoms.In this way, structure-activity relationship analysis of drugs can be achieved that can guide the design of the new compounds [28].
The Comparative Molecular Similarity Index Analysis (CoMSIA) method is better than the CoMFA method.In CoMSIA, the molecular similarity index calculated from the modified SEAL similarity field is used as a descriptor to consider simultaneously the spatial, electrostatic, hydrophobic, and hydrogen bond donor and hydrogen bond acceptor fields [29].These indices are estimated indirectly by comparing the similarity of each molecule in the dataset to the common probe atoms located at the surrounding grid/lattice intersections [14].The Gaussian-type function is used to describe the distance between probe atoms and atoms aligned with molecules in the dataset to calculate the similarity of all grid points [15].
Comparative molecular field analysis (CoMFA) and comparative molecular similarity index analysis (CoMSIA) are commonly used 3D-QSAR methods.Firstly, the compounds are compared using the distillation module implemented in SYBYL 2.0 to identify a common core scaffold in all molecules, and then based on the SYBYL molecular comparison, comparative molecular field analysis (CoMFA) and comparative molecular similarity index analysis (CoMSIA) studies are performed to explore the contribution of different interactions.The CoMFA and CoMSIA fields of each molecule are calculated for compounds in each network point computational dataset using a carbon probe hybridized with sp 3 of charge +1.0, and then the relationship between the information of these threedimensional features and the activity of the compounds is obtained by the partial least squares method [27].In this way, the conformational relationship analysis of drugs can be achieved to guide the design of new compounds [28].Specifically, the CoMFA method calculates spatial (S) and electrostatic (E) fields.In contrast to CoMFA, CoMSIA uses Gaussian functions and Gaussian-distributed similarity indices to avoid mutations in lattice-based probe-atom interactions [14,29], taking into account spatial, electrostatic, hydrophobic, and hydrogen-bond donor and hydrogen-bond acceptor fields [15].
Topomer Comparative Molecular Field Analysis (Topomer CoMFA) was proposed by Cramer et al. [16].It is a 3D-QSAR tool which combines CoMFA and Topomer technology and can quickly predict the biological activities or properties of compounds.Topomer CoMFA method does not require manual superposition of molecules, can quickly and accurately predict the biological activity of compounds, and can build a relatively reliable model in a relatively short time [30].Topomer CoMFA uses topological molecular technology to cut ligand molecules into two or more small fragments, while preserving the common skeleton as much as possible.It automatically builds a standard 3D model of topologies for each fragment and generates a set of spatial and static electric fields for each set of topologies.Finally, the parameters obtained are taken as independent variables and the biological activity value as dependent variables.The relationship between molecular structure and biological activity is described.Finally, the partial least square regression analysis method is used to establish the QSAR model [16,31,32].
In this study, in the CoMFA method, the aligned molecules with the best orientation are in the 3D cube lattice, with the grid spacing of 2.0 Å in the x, y and z directions, and an extension radius of 4.0 Å around the aligned molecules in the Cartesian direction.The Van der Waals radius and charge of sp 3 hybrid carbon probe atom are 1.52 Å and +1.0, which are used to calculate the CoMFA space field and static electric field of each molecule.The Coulomb and Lennard-Jones potential functions are used to estimate the electrostatic and spatial interactions, respectively.The energy cut-off values for both the spatial field and the static field were set at 30 kcal mol −1 .The CoMSIA method calculates the similarity index descriptor using the same lattice used in CoMFA.Five physical and chemical properties of spatial field, electrostatic field, hydrophobic field, hydrogen bond donor field and acceptor field were evaluated by Å probe atom with charge of +1.0, radius of 1 Å, hydrophobicity of +1.0, hydrogen bond donor field and hydrogen bond acceptor field, attenuation coefficient of 0.3 and grid spacing of 2.0 Å.

External Validation of the QSAR Model
The QSAR model usually adopts the Left-One-Out (LOO) method as the internal validation.The cross-validation coefficient q 2 and non-cross-validation coefficient r 2 are calculated through the formula as the basis for judging the quality of QSAR model [33]; see Formulas (1) and (2).
Among the parameters of the model, the larger the cross-validation coefficient q 2 and non-cross-validation coefficient r 2 mean, the better the correlation of the model and the stronger the prediction ability.Generally, the cross-validation coefficient q 2 greater than 0.5, the non-cross-validation coefficient r 2 greater than 0.6, and the difference between r 2 and q 2 less than 0.3 can prove that the model has high prediction ability [34,35].Furthermore, ŷi and y i is the predicted value and the experimental value of the test set, respectively, and y and ŷ are the average activity values of the experimental value and the predicted value of the training set, respectively.The "n" is the number of molecules in the test set, and the "i" is the ordinal number of molecules in the test set. (1) However, only the interior of the QSAR model validation parameters obtained through the LOO method does not directly determine the stand or fall of model.Good internal validation can only show that the values of cross-validation coefficient q 2 and non-crossvalidation coefficient r 2 of compound training set are high, but it does not necessarily indicate that the prediction ability of the established model is good.Therefore, it is necessary to adopt external validation based on the test set [36].
External test set validation is an effective method to evaluate correction model prediction ability.Usually, the dataset is randomly divided into a training set and test set.The training set establishes a correction model, and the test set is used as an independent subset to test the prediction ability of the model.The external test set compares the predicted value with the experimental value to make an accurate evaluation of the prediction ability of the model.
In this study, the external validation of QSAR model can be verified and evaluated by the following parameters, and the prediction ability of QSAR model can be further verified by calculating the biological activity of compounds in the test set [37,38].
Prediction correlation coefficient q 2 pred based on test set.∑ n i=1 (y i − ŷi ) 2 is the sum of squares of the actual molecular activity of the test set and the mean molecular activity of the training set (SD).∑ n i=1 ( ŷi − y i ) 2 is the sum of squares of the deviation between the predicted value of the test set and the actual activity value (PRESS).Formula (3).
In order to obtain the best prediction model for the test set, additional validation is performed on the model, where r and r 0 are the regression correlation coefficients between the actual activity value and the predicted activity value, respectively, and k and k ′ represent the slopes of the model.Other validation statistic parameters r 2 m and ∆r 2 m are used to further evaluate the model.The parameters are as follows in Formulas ( 7)-( 13) [40].

Virtual Screen
Virtual screening is a better method for the discovery of lead compounds in drug research and development.It plays a great role in the identification and optimization of lead compounds in the early stage of drug research and development.It can effectively reduce the cost of drug research and development and improve the speed of drug research and development [41].In this study, we mainly used Topomer search technology in SYBYL2.0 for molecular virtual screening.Topomer search is a virtual screening method based on ligand drug design that can be used for structural optimization and skeleton transition of lead compounds [42].Moreover, Topomer search technology can be used in combination with receptor-based molecular docking, and can also be used as a preliminary screening tool when there is no receptor structure [43].
Usually, the Topomer CoMFA model will give the contribution value of each R group fragment in each molecule.The contribution value of the R group fully combines the 3D field and electrostatic field of the Topomer CoMFA model to intuitively show the good or bad degree of each fragment in the form of numerical size.In most cases, the higher the contribution value of the fragment, the higher the activity of the whole molecule.The basic principle of Topomer search is based on the contribution value of small fragments given by the Topomer CoMFA model.Using the R group in the model as the question formula, through the comparison of Topomer similarity, search the molecular fragments with high similarity in the compound database [44], and then select the fragments with high activity contribution value.The reasonable combination of the searched fragments and the basic skeleton is used to optimize and transform the structure of the original compounds, so as to obtain compounds with high biological activity [45], so as to achieve the purpose of rational drug design.

Molecular Docking
Molecular docking studies the interactions between the design drugs based on the characteristics and their receptors.This is a theoretical simulation method to predict the binding pattern and affinity through the geometric, energy matching and recognition between ligand and proteins [46], which is an important technology in the field of computeraided drug research.In this study, the Surflex-Dock module in the SYBYL2.0-Xsoftware package is selected to perform the simulated docking [47].Surflex-Dock is a semi-flexible docking with rigid protein structure, which uses a Protomol-based method to dock the ligand into the binding site of the receptor.It can avoid the tedious process of searching for the active site, and can compile specific docking methods, extract the ligands bound at the original site to generate a docking protocol, and dock the target ligand into the previous specific site.In addition, the conformations of the generated ligand molecules are compared with the original ligands, and are used as a reference to determine whether the docking method is applicable [48].The docking results of original ligands can also be compared with molecules of potential ligands to further screen drugs that are more suitable as inhibitors.
The protease used for this docking was derived from the PDB database (Protein Data Bank) [49], with the ID 1MQ4.The 1MQ4 macromolecular protein was pretreated before performing molecular docking to remove the required small molecule ligand from the macromolecular complex, the unwanted small molecule ligand and all water molecules, to hydrogenate the protein, and to add Gasteiger-Huckel charges [50,51].The docking regions were determined by analyzing the interactions of ligands and active residues, and the active sites for docking were determined according to the small molecule ligands.
The docking results were judged by the values of the compounds' scoring functions total score, crash and polar.The score function value of total score indicates the affinity between small molecule ligands extracted from macromolecular proteins and receptors, and the higher value means higher affinity.The crash represents the incompatibility between ligands extracted from macromolecular proteins and receptors, and the smaller crash value indicates the more compatible of molecular docking.Polar is the score of polarity function.A higher value is expected to be obtained when the binding site is on the protein surface, and a smaller value is expected when the binding site is inside the protein [52,53].

Molecular Dynamics Simulation
To explore the structural and energetic status of the newly designed compound, molecular dynamics (MD) simulation of 50 ns was performed using GROMACS 2020.4 software.The protein topology generation was used CHARMM36-February 2021 force-field, the ligand topology and the force field parameter files were generated from the CHARMM General Force Field server (CGenFF).These topology files were merged in GROMACS.Subsequently, the protein-ligand system was placed in a 12-sided box composed with a Tip3p water model.To neutralize the system charge, six Cl − were added to the system.Then the entire system was minimized by the steepest drop minimization algorithm.Following this, a 50 ns simulation was run under NPT (or NVT) ensemble at constant temperature (300 K) and pressure (1 atm) with the integration of trajectories at every two femtoseconds (fs), the atomic coordinates of simulated structures were recorded at every 10 ps.Root mean square deviation (RMSD), root mean square fluctuation (RMSF), radius of gyration (Rg) and number of hydrogen bonds were analyzed over the full trajectory.Furthermore, the free energy landscapes (FEL) were described to verify the stability of newly designed compound in complex with 1MQ4 protein.
The protein-ligand binding free energies were calculated by Molecular Mechanics-Poisson Boltzmann Surface Area (MM/PBSA).The last 10 ns of all MD simulations trajectories were subjected to calculate the binding free energy through MM/PBSA analysis.For calculations, 100 snapshots were extracted at every 100 ps from the last 10 ns of MD trajectories.Herein, the calculation of binding free energy ∆G bind was performed using the following equation: E vdw means the Van der Waals energy, E elec means electrostatic energy, E polar means electrostatic of salvation, and E SASA means the non-electrostatic energy of salvation.In addition, T is the temperature, S is the entropy, and TS is the entropic contribution in a vacuum.Because the variation in TS terms does not improve the anticipated results, the TS term is ignored.

ADMET Prediction
Currently, the research and development of drugs is a high-risk investment, which often faces unexpected and even catastrophic failures at different stages of drug discovery [54].One of the main reasons for the failures is the lack of efficacy and safety, which is related to the absorption, distribution, metabolism and excretion (ADME) properties and various toxicities (T) of the human body.Rapid ADMET evaluation to reduce failures in drug discovery is essential [55].ADMET (Absorption, Distribution, Metabolism, Excretion, and Toxicity) are prerequisites, and the nature of the molecule plays a key role in the preclinical phase.It is therefore necessary to pretest the ADMET properties of the designed compounds to ensure drug suitability for the human body [56].The ADMET properties acquired in this study are derived from the web flat ADMET-lab [57].

Result Analysis of HQSAR Model
The HQSAR model can be optimized by changing the parameters, including holographic length, Fragment Size and fragment characteristic parameters.In this article, the hologram length was the default (97, 151, 199, 257, 307, 353, and 401) and the Fragment Size set as the default (4-7).Different HQSAR models were established for 49 compounds in the training set using different fragment types.According to the combination of different fragmentary features, the best 32 combination methods selected are presented in Table 2.
As shown in Table 2, the best model could be obtained when the Fragment Distinction set to "A/CH" with the condition of all Fragment Size of "4-7" (Model 1-07), whose main parameters q 2 = 0.872, r 2 = 0.961, N = 6, HL = 151.Based on the results of Table 2, the Fragment Dispersion is set as "A/CH" again, with different Fragment Size chosen, 10 HQSAR models are built and shown in Table 3.The results show that the best model can be obtained when the Fragment Size is set as "2-4".The main parameters are shown as q 2 = 0.892, r 2 = 0.948, N = 5, HL = 257.Therefore, the best HQSAR model (Model 2-03) is built when the Fragment Dispersion is "A/CH" and the Fragment Size is "2-4".
The HQSAR results could be used to discover the potential effects of fragments and atoms on the activity of compounds, which is helpful to reveal the molecular mechanism affecting the activity of compounds.From the color code diagram, fragments and atoms that might be the key contributors for the activity of Aurora kinase inhibitors can be identified.Among them, the color code from green to red represents the activity contribution from high to low.

Result Analysis of CoMFA and CoMSIA Models
Generally, in the process of research, the molecule with the highest pIC50 value in the whole molecule set has the best molecular structure and is suitable to be used as the template molecule.Through the analysis of template molecules, the relationship between molecular structure and biological activity can be found quickly and accurately.Structure alignment is one of the most important input variables in 3D-QSAR analysis, and the reliability of the contour map based on molecular structure alignment greatly influences the predictive ability of the model.In this study, compound 52 with the highest activity is suitable for analysis as a template, and the common skeleton is selected as the overlap site.
As shown in Table 4, for the CoMFA model, the optimal composition number was 6, q 2 was 0.866, r 2 was 0.983, SEE was 0.156, and F was 403.587.This proves that this established model could strongly and stably predict the activity value of test set compounds.The contribution of spatial field and static electric field were 48.6% and 51.4%, respectively, which indicates the contribution of electric field to the CoMFA model is relatively large.
For theCoMSIA model, the results of the CoMSIA model are listed in Table 4 from 3-01 to 3-30.Model 3-30 is chosen as the research object for further discussion and analysis.The optimum composition was 9, q 2 was 0.877, r 2 was 0.995, SEE was 0.091, and F was 802.161.The contribution values of space field, electrostatic field, hydrophobic field, hydrogen bond donor field and hydrogen bond acceptor field were 18.3%, 31.6%,22.1%, 8.9% and 19.2%, respectively.The results show that the contribution value of electrostatic field to CoMSIA model is larger.

Result Analysis of CoMFA and CoMSIA Models
Generally, in the process of research, the molecule with the highest pIC 50 value in the whole molecule set has the best molecular structure and is suitable to be used as the template molecule.Through the analysis of template molecules, the relationship between molecular structure and biological activity can be found quickly and accurately.Structure alignment is one of the most important input variables in 3D-QSAR analysis, and the reliability of the contour map based on molecular structure alignment greatly influences the predictive ability of the model.In this study, compound 52 with the highest activity is suitable for analysis as a template, and the common skeleton is selected as the overlap site.
As shown in Table 4, for the CoMFA model, the optimal composition number was 6, q 2 was 0.866, r 2 was 0.983, SEE was 0.156, and F was 403.587.This proves that this established model could strongly and stably predict the activity value of test set compounds.The contribution of spatial field and static electric field were 48.6% and 51.4%, respectively, which indicates the contribution of electric field to the CoMFA model is relatively large.For theCoMSIA model, the results of the CoMSIA model are listed in Table 4 from 3-01 to 3-30.Model 3-30 is chosen as the research object for further discussion and analysis.The optimum composition was 9, q 2 was 0.877, r 2 was 0.995, SEE was 0.091, and F was 802.161.The contribution values of space field, electrostatic field, hydrophobic field, hydrogen bond donor field and hydrogen bond acceptor field were 18.3%, 31.6%,22.1%, 8.9% and 19.2%, respectively.The results show that the contribution value of electrostatic field to CoMSIA model is larger.
Figure 4 shows the 3D contour map of the spatial field and static electric field of the CoMFA model of template molecule compound 52.The green part indicates that increasing the volume of substituents is beneficial to the improvement the compound activity, while the yellow part indicates that reducing the volume of substituents is beneficial to the improvement of compound activity.The red area indicates that increasing the electronegativity of the group would increase the activity of the compound, and the blue area indicates that decreasing the electronegativity of the group would increase the activity of the compound.In Figure 4, the proportion of the green group and blue group are obviously larger.In order to improve the molecular activity, it is necessary to increase the volume of substituents in the corresponding parts or reduce the electronegativity of substituents in the blue region.Figure 5 shows the 3D contour map of CoMSIA model of compound 52.In Figure 5a, the green part indicates that increasing the volume of substituents is beneficial to improve the activity of the compound, and the yellow part indicates that reducing the volume of substituents is beneficial to improve the activity of the compound.In Figure 5b, the red region indicates that increasing the electronegativity of the group is beneficial to improve the activity of the compound, and the blue region indicates that reducing the electronegativity of the group is beneficial to improve the activity.In Figure 5c, the white part indicates that the increase of hydrophilic group is beneficial, and the yellow part indicates that the increase of hydrophobic group is beneficial.In Figure 5d, the cyan part indicates that increasing hydrogen bond donor is beneficial, and the purple part represents increasing hydrogen bond donor is not favorable.In Figure 5e, the magenta part indicates that increasing hydrogen bond receptor is beneficial, while the red part indicates that increasing hydrogen bond receptor is not conducive to activity.According to the 3D contour map of CoMSIA model, we can increase or decrease specific groups in corresponding parts to improve the activity of compounds.

The Result and Analysis of the Topomer CoMFA Model
In the Topomer CoMFA study, the most active compound (52) was analyzed as a template, and all the compounds were divided into four fragments according to the segmentation method shown in Figure 6.After cutting, the compounds were automatically Figure 5 shows the 3D contour map of CoMSIA model of compound 52.In Figure 5a, the green part indicates that increasing the volume of substituents is beneficial to improve the activity of the compound, and the yellow part indicates that reducing the volume of substituents is beneficial to improve the activity of the compound.In Figure 5b, the red region indicates that increasing the electronegativity of the group is beneficial to improve the activity of the compound, and the blue region indicates that reducing the electronegativity of the group is beneficial to improve the activity.In Figure 5c, the white part indicates that the increase of hydrophilic group is beneficial, and the yellow part indicates that the increase of hydrophobic group is beneficial.In Figure 5d, the cyan part indicates that increasing hydrogen bond donor is beneficial, and the purple part represents increasing hydrogen bond donor is not favorable.In Figure 5e, the magenta part indicates that increasing hydrogen bond receptor is beneficial, while the red part indicates that increasing hydrogen bond receptor is not conducive to activity.According to the 3D contour map of CoMSIA model, we can increase or decrease specific groups in corresponding parts to improve the activity of compounds.Figure 5 shows the 3D contour map of CoMSIA model of compound 52.In Figure 5a, the green part indicates that increasing the volume of substituents is beneficial to improve the activity of the compound, and the yellow part indicates that reducing the volume of substituents is beneficial to improve the activity of the compound.In Figure 5b, the red region indicates that increasing the electronegativity of the group is beneficial to improve the activity of the compound, and the blue region indicates that reducing the electronegativity of the group is beneficial to improve the activity.In Figure 5c, the white part indicates that the increase of hydrophilic group is beneficial, and the yellow part indicates that the increase of hydrophobic group is beneficial.In Figure 5d, the cyan part indicates that increasing hydrogen bond donor is beneficial, and the purple part represents increasing hydrogen bond donor is not favorable.In Figure 5e, the magenta part indicates that increasing hydrogen bond receptor is beneficial, while the red part indicates that increasing hydrogen bond receptor is not conducive to activity.According to the 3D contour map of CoMSIA model, we can increase or decrease specific groups in corresponding parts to improve the activity of compounds.

The Result and Analysis of the Topomer CoMFA Model
In the Topomer CoMFA study, the most active compound (52) was analyzed as a template, and all the compounds were divided into four fragments according to the segmentation method shown in Figure 6.After cutting, the compounds were automatically

The Result and Analysis of the Topomer CoMFA Model
In the Topomer CoMFA study, the most active compound (52) was analyzed as a template, and all the compounds were divided into four fragments according to the segmentation method shown in Figure 6.After cutting, the compounds were automatically divided into R a group (red), R b group (blue), R c group (green) and common skeleton (black).
divided into Ra group (red), Rb group (blue), Rc group (green) and common skeleton (black).Taking compound 52 as the template molecule, the relationship between the structure and activity of the compounds in training set was studied using the Topomer CoMFA method.The results are shown as follows: In Table 5, q 2 of the best model obtained is greater than 0.5, r 2 greater than 0.6, and residual between r 2 and q 2 less than 0.3.This indicates that the model established by this method is an ideal Topomer CoMFA model, and its statistical results have high predictive ability.At the same time, the predicted values of SEE, q 2 stderr, r 2 stderr and F are 0.199, 0.36, 0.20 and 369.402, respectively.These results indicate that the error of the model is quite small, which further demonstrates that the established model is reliable.Figure 7 shows the 3D contour map of the Topomer CoMFA model of the template molecule compound 52. Figure 7a1,a2 are the spatial contour map and electrostatic contour map of Ra group, respectively.Figure 7b1,b2 are the spatial contour map and electrostatic contour map of Rb group, respectively.Figure 7c1,c2 are the spatial contour map and electrostatic contour map of Rc group, respectively.In Figure 7, the green part indicates that increasing the volume of substituent is beneficial to improve the activity of the compound, and the yellow part indicates that reducing the volume of substituent is beneficial to improve the activity of the compound.The red region indicates that increasing the electronegativity of the group is beneficial, and the blue region indicates that reducing the electronegativity of the group is beneficial.Taking compound 52 as the template molecule, the relationship between the structure and activity of the compounds in training set was studied using the Topomer CoMFA method.The results are shown as follows: In Table 5, q 2 of the best model obtained is greater than 0.5, r 2 greater than 0.6, and residual between r 2 and q 2 less than 0.3.This indicates that the model established by this method is an ideal Topomer CoMFA model, and its statistical results have high predictive ability.At the same time, the predicted values of SEE, q 2 stderr, r 2 stderr and F are 0.199, 0.36, 0.20 and 369.402, respectively.These results indicate that the error of the model is quite small, which further demonstrates that the established model is reliable.Figure 7 shows the 3D contour map of the Topomer CoMFA model of the template molecule compound 52. Figure 7a 1 ,a 2 are the spatial contour map and electrostatic contour map of R a group, respectively.Figure 7b 1 ,b 2 are the spatial contour map and electrostatic contour map of R b group, respectively.Figure 7c 1 ,c 2 are the spatial contour map and electrostatic contour map of R c group, respectively.In Figure 7, the green part indicates that increasing the volume of substituent is beneficial to improve the activity of the compound, and the yellow part indicates that reducing the volume of substituent is beneficial to improve the activity of the compound.The red region indicates that increasing the electronegativity of the group is beneficial, and the blue region indicates that reducing the electronegativity of the group is beneficial.

Comparison of HQSAR, CoMFA, CoMSIA and Topomer CoMFA Models
The comparisons of the predicted activity values obtained by HQSAR, CoMFA, CoM-SIA and Topomer CoMFA models and the residual values obtained by subtracting the experimental values from the predicted values are listed in Table 6.The activities of these compounds predicted by the four methods are basically the same as the experimental values, which indicates that the results of the model established by the four methods are accurate and have a high predictive ability.

Comparison of HQSAR, CoMFA, CoMSIA and Topomer CoMFA Models
The comparisons of the predicted activity values obtained by HQSAR, CoMFA, CoMSIA and Topomer CoMFA models and the residual values obtained by subtracting the experimental values from the predicted values are listed in Table 6.The activities of these compounds predicted by the four methods are basically the same as the experimental values, which indicates that the results of the model established by the four methods are accurate and have a high predictive ability.Based on the results of training set, the pIC 50 values of the test set are predicted by correspond models and the residual values are also calculated (Table 6).According to the contents of Table 6, the linear regression analysis diagram (a 1 ) and residual analysis diagram (a 2 ) of the HQSAR model, the linear regression analysis diagram (b 1 ) and residual analysis diagram (b 2 ) of the CoMFA model, the linear regression analysis diagram (c 1 ) and residual analysis diagram (c 2 ) of the CoMSIA model and the linear regression analysis diagram (d 1 ) and residual analysis diagram (d 2 ) of the Topomer CoMFA model are drawn using Origin Pro (edu) software.In these figures, all samples in the linear regression diagram are evenly distributed around the 45 • line, and the data in the residual analysis diagram are mainly concentrated around the 0-tick line.This demonstrates that all four models built in this study have favorable predictive ability.
By comparing and analyzing the prediction results of HQSAR, CoMFA, CoMSIA and Topomer CoMFA models (Table 6) with the linear regression and residual analysis diagram (Figure 8), the predicted values obtained by the four established models are quite similar to the experimental values.This indicates that the Fragmentation Interval and Size in the HQSAR model, the overlapping position of CoMFA and CoMSIA model, and the choice of cutting method in the Topomer CoMFA model are all the optimum.Therefore, the model established in this study could be used to predict the biological activities of these types of new compounds and their derivatives more accurately.

External Validation Results of the QSAR Model
Through the analysis and comparison of the models built by QSAR, good results achieved.However, these results were theoretically tested based on the training set, wh might be accidental.In order to further verify the predictive ability of QSAR model, in

External Validation Results of the QSAR Model
Through the analysis and comparison of the models built by QSAR, good results are achieved.However, these results were theoretically tested based on the training set, which might be accidental.In order to further verify the predictive ability of QSAR model, independent test sets not used for model generation are used for further analysis (Table 7).In Table 7, combined with the relevant parameters and corresponding standards introduced above, the results of external verification of the four models are basically within a reasonable range, which further shows that the model established in this study has high reliability.

The Result and Analysis of Molecular Design
Based on the analysis of the Topomer CoMFA model results, suitable R a , R b and R c groups are designed based on the template compound 52.Taking the R a , R b and R c groups as the query structures, the Topomer search technology is used to search the conformation in the ZINC database.Among the structural fragments, the R a group (contribution value is 0.80), R b group (contribution value is 1.54) and R c group (contribution value is 1.77) with topological distance close to 185 and activity contribution value higher than the template molecule are selected.Finally, 10 new imidazo [4,5-b]pyridine derivatives are designed and screened as Aurora kinase A inhibitors.According to the previous method, the structure of the new compound is constructed in SYBYL2.0-X,and the molecular optimization and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8.
As shown in Table 8, the predicted activity values (pIC 50 ) of all the newly designed compounds are higher than the template molecules (pIC 50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For R a group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a 1 ,a 2 ).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b 1 ,b 2 ).For the R b group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8.  and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8. and naming are carried out by the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table 8.As shown in Table 8, the predicted activity values (pIC50) of all the newly designed compounds are higher than the template molecules (pIC50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For Ra group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a1,a2).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b1,b2).For the Rb group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH4, and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively ().The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.
Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be As shown in Table 8, the predicted activity values (pIC50) of all the newly designed compounds are higher than the template molecules (pIC50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For Ra group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a1,a2).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b1,b2).For the Rb group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH4, and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively ().The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.
Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be As shown in Table 8, the predicted activity values (pIC50) of all the newly designed compounds are higher than the template molecules (pIC50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For Ra group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a1,a2).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b1,b2).For the Rb group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH4, and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively ().The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.
Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be As shown in Table 8, the predicted activity values (pIC50) of all the newly designed compounds are higher than the template molecules (pIC50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For Ra group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a1,a2).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b1,b2).For the Rb group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH4, and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively ().The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.
Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be As shown in Table 8, the predicted activity values (pIC50) of all the newly designed compounds are higher than the template molecules (pIC50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For Ra group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a1,a2).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b1,b2).For the Rb group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH4, and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively ().The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.
Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be As shown in Table 8, the predicted activity values (pIC50) of all the newly designed compounds are higher than the template molecules (pIC50 = 8.70).The structural analysis of the newly designed compound molecules are carried out.For Ra group, when the -CN group is replaced by -Br group, the volume of the substituent group is increased and the electronegativity is enhanced.The corresponding criteria are that enlarging the volume around the green contour is favorable to the activity of inhibitors (Figure 8a1,a2).Increasing the electronegativity around the red group is conducive to the increase in activity (Figure 8b1,b2).For the Rb group, the newly designed compounds increased the volume near C-33 of the Rb group, and C-36 position is deleted, corresponding to the principle that increasing the volume around the green group and reducing the volume around the yellow group are beneficial to increase activity.These results are consistent with the analysis results of the QSAR models.

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH4, and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively ().The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.
Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be 0.80 1.54 1.91 8.96

Analysis of Molecular Docking Results
In order to understand the inhibitory protein mechanism of the designed compound and the binding between the newly designed compound and the receptor protein, it was necessary to perform the molecular docking.According to the three probes NH, CH 4 , and CO, which form binding pockets with small molecule ligands, the three probes represent the hydrogen bond donor, hydrogen bond acceptor, and hydrophobic sites, respectively.The cocrystal ligand in protein macromolecules is then re-accessed into the crystal structure by docking technique to observe its docking in the interface pocket to verify the reliability.The ligand molecules before docking are contrasted with those after docking (Figure 9).The red sticks in the Figure represent the ligands after docking, and the green sticks are the original ligands.The conformations of crystal ligands almost completely overlap with those after ligand docking (similarity of 0.903), and their rotation tendency is basically similar.According to Figure 9, the docking conformation of ligand molecules basically matched the original conformation, and the docking method in this study is reasonable and reliable.seen from Figure 9, the docked conformations of the ligand molecules basically match the original conformations, and the docking method in this study is reasonable and reliable.The template molecule and all the newly designed compounds are put into the docking pocket for molecular docking, and the template molecule is used as a control, with the value of the scoring functions, total score, crash and polar, and the number of hydrogen bonds formed as the criteria.A higher value of the total score indicates that the docking results meet the requirements of the analysis results (Table 9).Molecular docking was necessary to understand the mechanism of protein inhibition by the designed compounds and the binding of the newly designed compounds to the receptor proteins.The eutectic ligand (ADP) in the protein macromolecule was rewired into the crystal structure by docking technique, and the ligand molecules before docking were compared with the ligand molecules after docking (Figure 9).The red bars in the figure represent the docked ligand and the green bars are the original ligand.The conformation of the crystalline ligand overlaps almost completely with that of the ligand after docking (similarity of 0.903), and their rotational trends are essentially similar.As can be seen from Figure 9, the docked conformations of the ligand molecules basically match the original conformations, and the docking method in this study is reasonable and reliable.
The template molecule and all the newly designed compounds are put into the docking pocket for molecular docking, and the template molecule is used as a control, with the value of the scoring functions, total score, crash and polar, and the number of hydrogen bonds formed as the criteria.A higher value of the total score indicates that the docking results meet the requirements of the analysis results (Table 9).In the 2D analysis diagram of molecular docking (Figure 10), the ball-and-sticks are ligands, and the round spheres are the amino acid residues forming interactive forces.Hydrogen bonding is shown as the green dotted line, and hydrophobic interaction is shown as the pink dotted line.Among them, the hydrogen bond is the main force to maintain the ligand molecules of protein and compound, which made the binding between them more stable.The hydrophobic bonds are mainly used to enhance the affinity between ligand and amino acids around active site.
Figure 10a shows a docking analysis diagram of ligand extracted from 1MQ4 protein crystals.In this Figure, the ligand mainly forms four hydrogen bonds with A/ASN261, A/GLU260, A/ASP274, A/GLU211 and A/ALA213 residues in protein crystals, and hydrophobic interactions formed with amino acid residues such as A/VAL147, A/ALA160, A/LEU263 and A/LEU139.Total score, crash and polar are 9.3093, −0.5178 and 1.0743, respectively.The selection of ligands and protein crystals used in this docking is appropriate, and the docking method is reasonable and reliable.
Figure 10b shows the 2D analysis diagram of the docking 1MQ4 protein crystal with compound 52 in the training set as the template molecule.Compound 52 mainly forms a hydrogen bond with A/LYS162 residue in the crystal structure, and forms a hydrophobic interaction with amino acid residues such as A/LYS143, A/PHE144, A/LEU210 and A/VAL147.The total score, crash and polar are 6.3356, −0.5557, and 1.1914, respectively.Figure 10c shows that the newly designed compound N1 forms four hydrogen bonds with A/GLN177, A/ALA213, A/VAL174 and A/GLU181 residues, and forms a hydrophobic interaction with amino acid residues such as A/LEU139, A/VAL147 and A/LYS143, with total score, crash and polar as 6.3210, −1.7537 and 1.2747, respectively.Figure 10d shows that the newly designed compound N2 forms two hydrogen bonds with A/ALA213, A/ARG137 and A/THR260 residues, and forms a hydrophobic interaction with amino acid residues such as A/ALA160, A/LEU263, A/VAL147 and A/LEU139.The total score, crash and polar are 7.8842, −1.4843 and 2.2739, respectively.Figure 10a shows a docking analysis diagram of ligand extracted from 1MQ4 protein crystals.In this Figure, the ligand mainly forms four hydrogen bonds with A/ASN261,A/GLU260,A/ASP274,A/GLU211 and A/ALA213 residues in protein crystals, and hydrophobic interactions formed with amino acid residues such as A/VAL147, A/ALA160, A/LEU263 and A/LEU139.Total score, crash and polar are 9.3093, −0.5178 and 1.0743, respectively.The selection of ligands and protein crystals used in this docking is appropriate, and the docking method is reasonable and reliable.
Figure 10b shows the 2D analysis diagram of the docking 1MQ4 protein crystal with compound 52 in the training set as the template molecule.Compound 52 mainly forms a hydrogen bond with A/LYS162 residue in the crystal structure, and forms a hydrophobic interaction with amino acid residues such as A/LYS143, A/PHE144, A/LEU210 and A/VAL147.The total score, crash and polar are 6.3356, −0.5557, and 1.1914, respectively.Figure 10c shows that the newly designed compound N1 forms four hydrogen bonds with A/GLN177, A/ALA213, A/VAL174 and A/GLU181 residues, and forms a hydrophobic interaction with amino acid residues such as A/LEU139, A/VAL147 and A/LYS143, with total score, crash and polar as 6.3210, −1.7537 and 1.2747, respectively.Figure 10d shows that the newly designed compound N2 forms two hydrogen bonds with A/ALA213, A/ARG137 and A/THR260 residues, and forms a hydrophobic interaction with amino acid The results show that the newly designed compounds formed strong hydrogen bonds with amino acid residues such as A/ALA213, A/VAL174 and A/LEU263 and forms hydrophobic interactions with amino acid residues such as A/VAL147, A/LYS143 and A/LEU263.These interactions enhance the binding strength of ligands and receptors, and their scoring functions are all ideal.Thus, the designed compound docking results are reliable.

The Analysis of Molecular Dynamics Simulation
MD simulation analysis was used to study the nature of newly designed compound in dynamic states and the stability of the protein-ligand complex.All molecular dynamics simulations were carried out for template 52, N3, N4, N5 and N7 docked with 1MQ4 protease structure for 50 ns.Meanwhile, the interaction between compounds N3, N4, N5, and N7 with the receptor protein is depicted in Figure S1.
The RMSD is one of the pivotal parameters which can explain the overall stability and describe its structural conformation changes in the protein-ligand complex from MD simulation trajectories.The high stability of the complex can be obtained by the low and stable RMSD value.Figure 11a represents the RMSD plot for the original ligand and newly designed compounds (N3, N4, N5, N7) docked with 1MQ4 protein.Compared to the template ligandprotein complex, the RMSD fluctuation of new compounds system became progressively lower.After 40 ns, the RMSD fluctuation reaches the minimum and equilibration is expected to occur.This suggests the system of newly designed compound combined with 1MQ4 protein is more stable.The free energy landscape is computed to characterize the conformational changes in the four protein-ligand complexes during the MD simulation.Stability is color-coded from dark blue to red.The free energy landscapes plotted between RMSD and Rg coordinates for N3, N4, N5 and N7 are displayed in Figure 12 with the Gibbs free energy varying from 0 to 9.4 kJ/mol, 11.2 kJ/mol, 9.2 kJ/mol and 9.4 kJ/mol, respectively.
As shown in the free energy landscape, the ensemble of conformations of the newly designed compounds (N3, N4, N5, N7) in complex with protein 1MQ4 are restricted into a single conformation cluster and present favorable basins.Multiple favorable free energy basins of N3 in complex with protein 1MQ4 are formed at RMSD values fluctuating between 0.18 and 0.31 nm, and Rg values varying between 1.89 and 1.92 nm.The free energy basins of N4 in complex with protein 1MQ4 are located at RMSD values between 0.32 and 0.45 nm, and Rg values varying between 1.91 and 1.94 nm.The free energy basins of N5 in complex with protein 1MQ4 are located at RMSD values between 0.18 and 0.25 nm, and Rg values varying between 1.87 and 1.90 nm.The free energy basins of N7 in complex with protein 1MQ4 are relatively scattered, and located at RMSD values fluctuating between 0.20 and 0.35 nm, and Rg values varying between 1.88 and 1.92 nm.These free energy profiles clearly correlate with the previous RMSD and Rg analysis.The surface formats and snapshots of lowest-energy conformation from free energy basins are extracted and represented in Figure 12, which depict the stable conformational states of the complexes.The mean and standard deviation values of RMSD and Rg of every complex are shown in Table S1.
The RMSF study confirms the flexible regions of the protein and determines the fluctuation of protein during the MD simulation.A higher RMSF value indicates that the atoms in the molecule have a greater change in position during the simulation, suggesting that these atoms are more flexible.On the contrary, a lower RMSF value indicates that the atoms have less positional changes, suggesting that these atoms are more stable during the simulation.Figure 11b shows the fluctuations of the residue of newly designed compounds (N3, N4, N5, N7) docked with 1MQ4 protein.In the region of 175-190 and 210-270, the fluctuations of RMSF values are smaller, which is because the key residues, such as A/GLU181, A/ALA213, A/LEU263 and A/THR260 at the active sites of protein 1MQ4 are more stable than other residues during the whole MD simulation process.These results are consistent with the molecular docking results.Furthermore, from the lower RMSF value, designed compounds N3, N5 and N7 showed better stability with protein 1MQ4.
The radius of gyration (Rg) is calculated to describe the compactness changes in the system and the stability of the system.It denotes the folding and unfolding of the proteins during MD simulations.The values of Rg for newly designed compounds are shown in Figure 11c.Five complexes show stable compactness and radius of gyration.These results indicate the compounds designed are superior in stability.
The hydrogen bonds between protein and ligand are important factors in maintaining the stability of complex conformations.Figure 11d shows the number of hydrogen bonds of five protein-ligand complexes during the 50 ns simulation.The number of hydrogen bonds varies for N3 (0-6), N4 (0-5), N5 (0-6) and N7 (0-6).The averages of the hydrogen bonds are calculated for 0.66 (template 52), 2.21 (N3), 1.84 (N4), 2.21 (N5) and 1.51 (N7).These results show that these newly designed compounds are more stable than template compound 52 during the simulation process.
The free energy landscape is computed to characterize the conformational changes in the four protein-ligand complexes during the MD simulation.Stability is color-coded from dark blue to red.The free energy landscapes plotted between RMSD and Rg coordinates for N3, N4, N5 and N7 are displayed in Figure 12 with the Gibbs free energy varying from 0 to 9.4 kJ/mol, 11.2 kJ/mol, 9.2 kJ/mol and 9.4 kJ/mol, respectively.The Van der Waals energy (Evdw), electrostatic energy (Eelec), polar solvation energy (Epolar), SASA energy (ESASA) and binding free energy (Gbind) are calculated and listed in Table 10.The results show that the binding free energies for N3, N4, N5 and N7 complex structure are −68.515± 2.218, −54.869 ± 1.705, −62.317 ± 2.418 and −52.857 ± 4.201 kcal/mol whereas that of template compound 52 complex are −26.939± 5.302 kcal/mol.The negative values represent the binding interaction, while the positive values represent the opposite.As shown in the free energy landscape, the ensemble of conformations of the newly designed compounds (N3, N4, N5, N7) in complex with protein 1MQ4 are restricted into a single conformation cluster and present favorable basins.Multiple favorable free energy basins of N3 in complex with protein 1MQ4 are formed at RMSD values fluctuating between 0.18 and 0.31 nm, and Rg values varying between 1.89 and 1.92 nm.The free energy basins of N4 in complex with protein 1MQ4 are located at RMSD values between 0.32 and 0.45 nm, and Rg values varying between 1.91 and 1.94 nm.The free energy basins of N5 in complex with protein 1MQ4 are located at RMSD values between 0.18 and 0.25 nm, and Rg values varying between 1.87 and 1.90 nm.The free energy basins of N7 in complex with protein 1MQ4 are relatively scattered, and located at RMSD values fluctuating between 0.20 and 0.35 nm, and Rg values varying between 1.88 and 1.92 nm.These free energy profiles clearly correlate with the previous RMSD and Rg analysis.The surface formats and snapshots of lowest-energy conformation from free energy basins are extracted and represented in Figure 12, which depict the stable conformational states of the complexes.
The Van der Waals energy (E vdw ), electrostatic energy (E elec ), polar solvation energy (E polar ), SASA energy (E SASA ) and binding free energy (G bind ) are calculated and listed in Table 10.The results show that the binding free energies for N3, N4, N5 and N7 complex structure are −68.515± 2.218, −54.869 ± 1.705, −62.317 ± 2.418 and −52.857 ± 4.201 kcal/mol whereas that of template compound 52 complex are −26.939± 5.302 kcal/mol.The negative values represent the binding interaction, while the positive values represent the opposite.The newly designed compounds show negative values which strongly suggest that these compounds bind with Aurora protein (1MQ4) efficiently and provide a direction for further research.

The Analysis of ADMET Prediction Results
From the results of ADMET property prediction in Table 11, it is clear that all the new compounds showed better human intestinal absorption (HIA), most of them had high oral bioavailability (F) and blood-brain barrier (BBB) distributions, but the plasma protein-binding (PPB) indices were comparatively low, so that the dose of the drug needs to be taken into account when designing the drug.Inhibitory capacity prediction showed that all compounds could act as substrates for CYP3A4 and CYP2C9 inhibitors.Drug clearance prediction showed that most of the compounds exhibited moderate clearance (5-15 mL/min/kg).On this occasion, compound toxicity prediction assessment revealed that all compounds passed Ames toxicity and were not carcinogenic in rats, which is a positive feature for drug development.In conclusion, the ADMET prediction results suggest that the new compounds have potential medicinal value for drug development.However, they still need to pass rigorous testing to ensure drug safety and efficacy that can be guaranteed in subsequent studies and clinical trials.
Based on the results of the model built by four methods, HQSAR, CoMFA, CoMSIA and Topomer CoMFA, virtual screening is performed using Topomer search technology, and 10 novel compounds belonging to the Aurora kinase A inhibitor class are designed based on the screening results.However, the pharmacology and safety of the designed compounds are uncertain.It is necessary to pretest the ADMET properties of the designed compounds to ensure drug suitability in humans.
The results are shown in Table 11.Each value of the designed drugs basically meets the requirement without rat carcinogenicity, which indicates that the designed drug molecules have high safety and obvious pharmacological effects.

Conclusions
Cancer is a major disease that seriously threatens human life and social development.Utilizing scientific methods for cancer prevention and control has become one of the most important global public health issues.This study established reliable 3D-QSAR models based on 65 imidazo [4,5-b]pyridine derivatives, analyzing the relationship between inhibitor bioactivity and structure.It designed novel imidazo [4,5-b]pyridine-based Aurora kinase A inhibitor anticancer drugs, ultimately identifying 10 new molecules.Molecular docking was employed to explore the binding modes and targets between ligands and protein receptors, identifying interactions between amino acid residues forming hydrogen bond interactions and crystal structures.The newly designed compounds exhibited strong hydrogen bonding interactions with amino acid residues such as A/ALA213, A/VAL174, and A/LEU263, as well as hydrophobic interactions with residues like A/VAL147, A/LYS143, and A/LEU263, enhancing the binding strength between ligands and receptors.Molecular dynamics simulations were performed on four newly designed compounds, with results indicating that compound N3 bound to the 1MQ4 protein is the most stable.Snapshots were extracted from the lowest energy conformational state of the compounds in complex with protein 1MQ4.Finally, ADMET prediction was conducted for the 10 newly designed drug molecules, validating their high safety and pharmacological effects.This study provides insights for the design and development of novel Aurora kinase A inhibitor anticancer drugs, aiding in a better understanding of their inhibition mechanism, and offering a theoretical basis for synthesizing new anticancer drugs.

Supplementary Materials:
The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/molecules29081772/s1, Figure S1: 2D interaction diagram of compounds N3, N4, N5, and N7 docked with receptor proteins; Table S1: The mean and standard deviation values of RMSD and Rg of every Complex.

Figure 1 .
Figure 1.Several Aurora kinase inhibitors that have entered clinical trials.

Figure 1 .
Figure 1.Several Aurora kinase inhibitors that have entered clinical trials.

Figure 2 .
Figure 2. Distribution of experimental inhibitory activities (pIC50) for the training and test sets compounds in the QSAR models.

Figure 2 .
Figure 2. Distribution of experimental inhibitory activities (pIC 50 ) for the training and test sets compounds in the QSAR models.

33 Figure 3 .
Figure 3. Atomic contribution maps of compound 52 (a), compound 32 (b), compound 21 (c) and compound 40 (d).(Cyan indicates common backbone, green or yellow indicates positive contribution, orange and red indicate negative contribution, and white indicates intermediate contribution to biological activity).

Figure 3 .
Figure 3. Atomic contribution maps of compound 52 (a), compound 32 (b), compound 21 (c) and compound 40 (d).(Cyan indicates common backbone, green or yellow indicates positive contribution, orange and red indicate negative contribution, and white indicates intermediate contribution to biological activity).

Figure 4 .
Figure 4. Three-dimensional (3D) (a)spatial field and (b)static electric field action diagrams for the CoMFA model with compound 52 as the template molecule.

Figure 4 .
Figure 4. Three-dimensional (3D) (a)spatial field and (b)static electric field action diagrams for the CoMFA model with compound 52 as the template molecule.

Figure 4 .
Figure 4. Three-dimensional (3D) (a)spatial field and (b)static electric field action diagrams for the CoMFA model with compound 52 as the template molecule.

Figure 6 .
Figure 6.Schematic representation of segmentation methods and fragments using compound 52 as a template molecule.

Figure 6 .
Figure 6.Schematic representation of segmentation methods and fragments using compound 52 as a template molecule.

Figure 7 .
Figure 7. Three-dimensional (3D) contour map of the Topomer CoMFA model based on template 52.(a1) the spatial contour map of Ra group; (a2)the electrostatic contour map of Ra group; (b1) the spatial contour map of Rb group; (b2)the electrostatic contour map of Rb group; (c1) the spatial contour map of Rc group; (c2)the electrostatic contour map of Rc group.

Figure 7 .
Figure 7. Three-dimensional (3D) contour map of the Topomer CoMFA model based on template 52.(a 1 ) the spatial contour map of R a group; (a 2 ) the electrostatic contour map of R a group; (b 1 ) the spatial contour map of R b group; (b 2 ) the electrostatic contour map of R b group; (c 1 ) the spatial contour map of R c group; (c 2 ) the electrostatic contour map of R c group.

Figure 9 .
Figure 9. Conformational superposition of docked ligand and extracted primary ligand (ADP).(The red bars in the figure represent the docked ligand and the green bars are the original ligand).

Figure 9 .
Figure 9. Conformational superposition of docked ligand and extracted primary ligand (ADP).(The red bars in the figure represent the docked ligand and the green bars are the original ligand).

Figure 11 .
Figure 11.RMSD values of ligands at the active site of 1MQ4 during 50 ns MD simulation (a).The RMSF values of chain A in five protein-ligand complexes (b).Radius of gyration values in five protein-ligand complexes (c).The hydrogen-bond numbers between 1MQ4 and the ligands (d).

Figure 11 .
Figure 11.RMSD values of ligands at the active site of 1MQ4 during 50 ns MD simulation (a).The RMSF values of chain A in five protein-ligand complexes (b).Radius of gyration values in five protein-ligand complexes (c).The hydrogen-bond numbers between 1MQ4 and the ligands (d).

Figure 12 .
Figure 12.Free energy landscape (FEL) plotted between RMSD and Rg coordinates, snapshot from lowest energy conformational state of the newly designed compounds ((a): N3, (b): N4, (c): N5, (d): N7) in complex with protein 1MQ4.The colored scale plot shows the free energy profile (kJ/mol), the red color region in plot represents the metastable state, and the dark blue region represents the lowest energy state.

Figure 12 .
Figure 12.Free energy landscape (FEL) plotted between RMSD and Rg coordinates, snapshot from lowest energy conformational state of the newly designed compounds ((a): N3, (b): N4, (c): N5, (d): N7) in complex with protein 1MQ4.The colored scale plot shows the free energy profile (kJ/mol), the red color region in plot represents the metastable state, and the dark blue region represents the lowest energy state.

Table 1 .
Compound structures and their activity values.

Table 1 .
Compound structures and their activity values.

Table 3 .
HQSAR model analysis with the same fragment type (A/CH) and different fragment length.

Table 4 .
Statistical Parameters of CoMFA and CoMSIA Models.
S, E, H, D and A represent space field, electrostatic field, hydrophobic field, hydrogen bond donor field and acceptor field, respectively.

Table 5 .
Results of the Topomer CoMFA model.
N: Best composition score; q 2 : Cross-validation coefficient; r 2 : Non-cross-validation coefficient; q 2 stderr: Standard error of cross-validation coefficient; r 2 stderr: Standard error of non-cross validation coefficient; F: Statistical value of F test; SEE: Estimate standard error.

Table 5 .
Results of the Topomer CoMFA model.
N: Best composition score; q 2 : Cross-validation coefficient; r 2 : Non-cross-validation coefficient; q 2 stderr: Standard error of cross-validation coefficient; r 2 stderr: Standard error of non-cross validation coefficient; F: Statistical value of F test; SEE: Estimate standard error.

Table 6 .
Comparison of predicted activity value and residual value of HQSAR, CoMFA, CoMSIA and Topomer CoMFA models.

Table 6 .
Comparison of predicted activity value and residual value of HQSAR, CoMFA, CoMSIA and Topomer CoMFA models.

Table 7 .
External verification of the QSAR Model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.carriedoutby the same method.The activities of newly designed compounds are predicted by Topomer CoMFA model as shown in Table8.

Table 8 .
Molecular structure of the new compounds and activity value prediction using the Topomer CoMFA model.

Table 9 .
Molecular docking results and scoring function.

Table 9 .
Molecular docking results and scoring function.Rank of the affinity of ligands that bind to the active site of the receptor and report the output of the total score.Crash: Crash score shows inappropriate penetration into the binding site.A crash score close to 0 is advantageous.A negative number means penetration.Polar: Polarity indicates the contribution of polar interactions to the total score. Total-score: