The Discovery of Novel BCR-ABL Tyrosine Kinase Inhibitors Using a Pharmacophore Modeling and Virtual Screening Approach

Chronic myelogenous leukemia (CML) typically results from a reciprocal translocation between chromosomes 9 and 22 to produce the bcr-abl oncogene that when translated, yields the p210 BCR-ABL protein in more than 90% of all CML patients. This protein has constitutive tyrosine kinase activity that activates numerous downstream pathways that ultimately produces uncontrolled myeloid proliferation. Although the use of the BCR-ABL tyrosine kinase inhibitors (TKIs), such as imatinib, nilotinib, dasatinib, bosutinib, and ponatinib have increased the overall survival of CML patients, their use is limited by drug resistance and severe adverse effects. Therefore, there is the need to develop novel compounds that can overcome these problems that limit the use of these drugs. Therefore, in this study, we sought to find novel compounds using Hypogen and Hiphip pharmacophore models based on the structures of clinically approved BCR-ABL TKIs. We also used optimal pharmacophore models such as three-dimensional queries to screen the ZINC database to search for potential BCR-ABL inhibitors. The hit compounds were further screened using Lipinski’s rule of five, ADMET and molecular docking, and the efficacy of the hit compounds was evaluated. Our in vitro results indicated that compound ZINC21710815 significantly inhibited the proliferation of K562, BaF3/WT, and BaF3/T315I leukemia cells by inducing cell cycle arrest. The compound ZINC21710815 decreased the expression of p-BCR-ABL, STAT5, and Crkl and produced apoptosis and autophagy. Our results suggest that ZINC21710815 may be a potential BCR-ABL inhibitor that should undergo in vivo evaluation.


INTRODUCTION
Chronic myelogenous leukemia (CML) is a hematopoietic neoplastic disease that primarily results from the reciprocal translocation between chromosomes 9 and 22 (t9;22) (q34; q11) in a hematopoietic stem cell (HSC), resulting in the formation of the bcr-abl oncogene (Deininger et al., 2005;Ren, 2005;Kimura, 2006). Typically, the most common BCR-ABL proteins produced in CML patients, based on the loci of the chromosomal breaks, are p210, p190, and p230 (Winter et al., 1999;Arana-Trejo et al., 2002). The p210 protein occurs was observed in more than 90% of CML patients and has been hypothesized to be the initial step in the pathogenesis of CML (Arana-Trejo et al., 2002). The p210 protein is a constitutively active, non-receptor tyrosine kinase that phosphorylates its Tyr177 residue, which subsequently activates a number of downstream signaling pathways that causes HSC cells undergo significant abnormal proliferation and differentiation compared to normal cells, producing the pathogenic changes reported in CML patients (Deininger et al., 2000;Gora-Tybor and Robak, 2008;Cao et al., 2015;Wang et al., 2016).
Imatinib (Gleevec), which selectively targets the ATP binding site of the BCR-ABL protein, is a first generation BCR-ABL tyrosine kinase inhibitor (TKI) that was approved in 2001 by the FDA for treating CML (Gontarewicz et al., 2008;Wu et al., 2014). Initially, it was reported that imatinib significantly increased the survival time of CML patients (Machova Polakova et al., 2015;Deng et al., 2020). However, subsequent clinical trials indicated that patients with advancedstage CML experience relapse due to development of resistance caused by point mutations within the BCR-ABL domain (notably T315I), amplification of the BCR-ABL protein, overexpression of BCR-ABL protein or the ABC transporter, P-glycoprotein (Coutre et al., 2000;von Bubnoff et al., 2002;Apperley, 2007;Tanaka and Kimura, 2008).
The second generation TKIs, dasatinib, nilotinib, and bosutinib, have been shown to have efficacy in treating CML patients with certain BCR-ABL mutations, but not those with the T315I mutation (O'Hare et al., 2005;Abbas and Hsyu, 2016). The mutations of BCR-ABL are a recognized cause of resistance to second-generation, especially T315I (Baccarani et al., 2009;Khorashad et al., 2013). Ponatinib, a third-generation TKI, is also efficacious in a number of mutations in the BCR-ABL protein, including T315I (Okabe et al., 2014). Compound mutations of BCR-ABL (T315I/F359V or Y253H/T315I) are usually resistant against ponatinib therapy (Zabriskie et al., 2014). However, ponatinib has been reported to significant adverse cardiovascular effects, notably serious venous thromboembolic events (Bu et al., 2014). Therefore, there is still the need to develop highly efficacious and tolerable BCR-ABL inhibitors for the treatment of CML.
Hypogen pharmacophore model were built based on known inhibitors, the structures and activity data of known active compounds were analyzed to generate quantitative 3D pharmacophore models with common pharmacodynamic action models of the compounds, which can be used to predict the efficacy of unknown compounds (Rampogu et al., 2018;Musoev et al., 2019). The Hiphop pharmacophore models were qualitative pharmacophore models generated according to the common characteristics of known inhibitors without considering the activity data of compounds (Zheng et al., 2009). Pharmacophore models can provide guidances for drug design and virtual screening.
In this study, the quantitative Hypogen pharmacophore model and the qualitative Hiphop model were built, these models were then used to screen the ZINC database to identify potential BCR-ABL inhibitors. The identified hit compounds were further evaluated using Lipinski's rule of five, ADMET and CDOCKER docking. The molecules obtained from the Hypogen and Hiphop pharmacophore models mapping were selected to determine their anti-proliferative efficacy in specific leukemia cell lines. Subsequently, the skeleton of the most active compound was chosen to identify additional derivatives from the database to perform structural optimization. Finally, our results led us to evaluate the efficacy of the compound ZINC217108151 to: (1) inhibit the proliferation of certain leukemia cell lines; (2) inhibit the phosphorylation of the p210 BCR-ABL tyrosine; and (3) affect the activity of the downstream proteins, STAT5 and Crkl, which are activated by BCR-ABL.

Virtual Screening
The best Hypogen and Hiphop pharmacophore models (the method of building pharmacophore shown in Supplementary Material) were used as 3D structures to search the ZINC clean drug-like database (1,469,373 molecules) to obtain the compounds that matched the chemical features (Bhadauriya et al., 2013;Kumar et al., 2015). The workflow of each screening step and the number of compounds is shown in detail in Figure 1. Compounds that satisfied the fit value of Hypogen and Hiphop were retained. All hit compounds were evaluated using Lipinski's rules of five and absorption, distribution, metabolism, excretion, and toxicity (ADMET) screening process (Rampogu et al., 2018). Lipinski's rule of five was predicted by using calculate molecular properties in small molecules protocol in DS2.5. Lipinski's rule of five consists of molecular weight ≤500 Daltons, rotatable bonds number ≤10, hydrogen bond acceptors ≤10, and hydrogen bond donors ≤5 (Rampogu et al., 2017). ADMET was predicted by using ADMET descriptors in calculate molecular properties protocol in DS2.5. The cut-off value of the solubility, absorption, and the BBB were 3, 0, and 3, respectively (Rampogu et al., 2018). Finally, the database was further screened using CDOCKER. Docking analysis was used to determine whether small molecules can enter the active site of the protein (Rampogu et al., 2020). The three-dimensional structure of ABL tyrosine kinase receptor (PDB: 1IEP) was obtained from Protein Data Bank (Protein Data Bank, 1971;Machova Polakova et al., 2015). All ligands were minimized by Chemistry at Harvard Macromolecular Mechanics (CHARMm) force field and the protein was prepared by removing all of the water molecules, adding hydrogen atoms and filling the loop area (Al-Balas et al., 2013). The active site was defined as 12.9 Å around the cocrystal ligand. The prepared protein and ligands were imported to perform CDOCKER docking. Higher docking energy and interaction energy represent the favorable binding between the protein and the ligands (Rampogu et al., 2018).

Reagents and Cell Culture
After virtual screening, the most optimal compounds were purchased from TargetMol company (Boston, MA, United States) and dissolved in DMSO to make a stock solution (1 × 10 −1 mol/L) and stored at −20 • C for future use.
Antibodies C-Abl and p-C-Abl were purchased from Cell Signaling Technology (MA, United States), p-Crkl was purchased from Affinity Biosciences (OH, United States), p-STAT5A and β-actin were purchased from Bioss (MA, United States), LC3 was purchased from Sigma-Aldrich (United States), Beclin1 was purchased from Abcam (Cambridge, United Kingdom), Caspase3 was purchased from Cell Signaling Technology (MA, United States). HRP-conjugated Affinipure Goat Anti-Rabbit antibody was purchased from Proteintech (IL, United States).
3-(4,5-dimethylthiazol-2-yl)-2,5diphenyltetrazolium Assay K562, BaF3/WT, and BaF3/T315I leukemia cells were seeded in 96-well plates at a density of 5 × 10 3 cells/well and 50 µL of various drug concentrations were added. CCC-HEL-1 cells were grown in 96-well plates at the same density and exposed to the same drug concentrations. The plates were incubated at 37 • C for 72 h. Cell proliferation was evaluated using MTT [3-(4,5dimethylthiazol-2-yl)-2,5-diphenyltetrazolium] assay. Fifteen µL of the MTT solution was added to each well of 96-well and plates were continually incubated at 37 • C for 4 h and three liquid solutions. SDS, HCl, and isobutanol were added to each well and the plates were incubated at 37 • C overnight. The optical density (OD) of each sample was measured at a wavelength of 570 nm and the 50% inhibitory concentration (IC 50 ) was determined by from the concentrationresponse curve.

Cell Cycle and Cell Apoptosis Analysis
K562 cells were plated in 6-well plates with 1 × 10 6 cells/well, BaF3/WT and BaF3/T315I cells were plated in 6-well plates with 5 × 10 5 cells/well. After all of the cells were incubated with 0 (vehicle), 0.1, 1, or 10 µM concentration of the test compounds for 24 h, they were harvested and washed with PBS and 1 ml of 70% cold ethanol was added to the cells (Kim et al., 2015;Jiang et al., 2016). The ethanol-fixed cells were stored at −20 • C overnight. Subsequently, the cells were washed with PBS and centrifuged at 1,800 rpm for 8 min to obtain the pellet. The pellet was resuspended in 100 µL of RNase A, incubated at 37 • C for 30 min, and mixed with 400 µL of propidium iodide (PI) and incubated at 4 • C in the dark for 30 min. DNA content was determined using flow cytometry, at an emission wavelength of 488 nm and analyzed using ModFit LT software.
All cells were plated in 6-well plates with 5 × 10 5 cells/well and incubated with 0 (vehicle), 0.1, 1, or 10 µM of the test compounds for 24 h. The cells were harvested, washed twice with cold PBS, resuspended with binding buffer, stained using Annexin V-fluorescein isothiocyante (FITC) and PI and mixed with PBS to obtain a final volume of 500 µL. Cellular apoptosis was determined using flow cytometry method within 1 h of staining and dot plots were set up to detect FITC and PI. Cells stained with FITC and without PI (FITC+, PI-) were categorized as being in early apoptosis, and cells stained with both dyes (FITC+, PI+) were categorized as being in late apoptosis or necrotic (Peng et al., 2001).

Western Blot Analysis
Cells were plated in 6-well plates with 5 × 10 5 cells/well and incubated with 0 (vehicle), 0.1, or 1 µM of the test compounds for 48 h. The harvested cells were washed with PBS, and lysed in a RIPA lysis buffer and collected by centrifugation at 14,000 × g for 15 min at 4 • C. The cell extracts were separated by SDSpolyacrylamide gel electrophoresis (PAGE) and transferred to PVDF membranes. The membranes were blocked with 5% nonfat milk at room temperature for 75 min. The membranes were incubated with C-Abl, p-C-Abl, p-Crkl, p-STAT5A, β-actin, LC3, Beclin1 or caspase-3 antibodies overnight at 4 • C, then washed three times with TBST and incubated for 1 h with the HRP-conjugated Affinipure Goat Anti-Rabbit antibody at room temperature. The protein expressions were quantified by ImageJ software.
The chemicals used to construct the Hypogen pharmacophore models are shown in Supplementary Figure 1. The top 10 pharmacophore models containing HBA, AR and H chemical features are shown in Supplementary Table 1. Among all of these models, Hypogen1, composed of HBA and 4H chemical features, had the lowest total cost (91.20), highest cost difference (69.50), lowest RMS (0.50), best correlation coefficient (0.99), and highest fit value (11.89). Cost analysis was applied to determine the statistical significance of the Hypogen models (Ye et al., 2010). The cost difference of Hypogen1 (69.4973), the difference between the null cost and total cost indicated more than 90 instances of statistical significance of achieving an acceptable pharmacophore (Kavitha et al., 2015;Kumar et al., 2015). The low RMS and high correlation coefficient suggested that the model is stable and internally predictive (Kandakatla and Ramakrishnan, 2014). The fit value was indicative of the overall fitness of the training set compounds on a pharmacophore model during pharmacophore generation (Dube et al., 2012). Therefore, the Hypogen1 was selected based on the highest correlation coefficient, lowest total cost and RMS value (Aparoy et al., 2010), 3D spatial relationship and distance constraints of Hypogen1 are shown in Figure 2A. The Figures 2B,C show the most active and inactive compounds in the training set aligned with Hypogen1.
If the compounds were divided into three categories based on their efficacies as being highly active (IC 50 ≤ 0.1 µM, +++), moderately active (0.1 µM < IC 50 < 3 µM, ++), and least active (IC 50 ≥ 3 µM, +), this categorization could be used to rapidly estimate the predictive accuracy of the pharmacophore (Huang et al., 2012). Table 1 shows that the efficacies of the training set compounds could be accurately predicted.

Validation of Hypogen1
The test set, Fischer randomization and decoy set were used to verify the Hypogen pharmacophore model.

Test set prediction
Thirty eight compounds, with different scaffold and efficacy ranges, were used as a test set (Supplementary Figure 2) to assess whether the pharmacophore model can accurately predict the efficacy of the compound (Huang et al., 2012). The experimental and predicted efficacies of the test set are shown in Table 2. A significant correlation coefficient of 0.847 was obtained for the test set compounds using regression analysis. The error factors of all compounds were less than 5%. Figure 3 shows the scatter plot of the experimental and the predicted efficacies, indicating that the samples are clustered near the diagonal line, y = x. These results indicated that the Hypogen1 model has excellent predictive validity.

Fischer randomization test
The Fischer randomization test is a cross validation method that can be used to evaluate the statistical correlation between structures and activities in training set compounds (Pal et al., 2019). To achieve a confidence level of 95% for Fischer's randomization, 19 random hypotheses were generated (Sakkiah et al., 2012). Figure 4 clearly shows that the Hypogen1 pharmacophore model was not generated by chance as the total cost of these randomly generated hypotheses is much higher than the initial hypothesis.

Decoy set
A decoy set was used to validate the efficiency of Hypogen1 by computing various parameters such as false positive, false negative, enrichment factor (EF), goodness of fit score (GF), total number of hit molecules (Ht), and total number of active compounds (Ha) (Halgren et al., 2004). The decoy set database (D) contained 325 compounds, including 25 BCR-ABL inhibitors and 300 inactive compounds. Our results indicated that 26 compounds were identified as being active based on the Hypogen1 model and 21 compounds were known inhibitors of the BCR-ABL protein. The relevant parameters for the decoy set are shown in Table 3. An EF of 10.5 and goodness of GF of 0.8022 verifies the high efficiency of Hypogen1, indicating the high predictive ability of the model.

Construction of Hiphop Model
Eight compounds were selected as the training set to generate qualitative Hiphop models for common feature pharmacophore generation. The structures and bioactivities of the training set are shown in Supplementary Figure 3. The top 10 common feature hypotheses were produced and are shown in Supplementary Table 2, according to their ranking scores. Based on the ranking scores and feature similarities of hypotheses, ten hypotheses were clustered into two groups: the first group models contained 2RA,  Table 3, we can see that compared with Hiphop9, Hiphop1 has a more pronounced difference in the active and inactive compounds, and thus, Hiphop1 was chosen as the optimal Hiphop model ( Figure 5A). Figure 5B shows the most active compound in the training set aligned with Hiphip1.

Validation of Hiphop1
To validate the Hiphop1 model, 21 active compounds, with different scaffold and efficacies (0.0005-41 µM) (Supplementary Figure 4) and 220 compounds with unknown efficacies from Maybridge, were combined to create a new database (D 1 ). Hiphop1 was used to select the active compounds. As a result, 22 compounds were screened from the database by Hiphop1 and 21 of these compounds have been previously proven to be known effective compounds. The validation results ( Table 4) indicated an EF of 10.95 and a goodness of GF of 0.91, thus validating the high efficiency of Hiphop1 and demonstrating that Hiphop1 can discriminate active compounds from inactive compounds.

Virtual Screening
The best validated Hypogen1 and Hiphop1 pharmacophore models were used as 3D queries to screen the drug-like subset of the ZINC Database. ZINC is a free database of purchasable molecules, including multiple subsets such as drug-like and lead-like and it offers similarity searching on its website (Irwin and Shoichet, 2005;Awale et al., 2015). The work flow for the virtual screening is shown in Figure 1  ADMET to assess their drug-like and pharmacokinetic properties (Rampogu et al., 2018), and 3,865 and 162 compounds meet these standards, respectively. Subsequently, molecular docking was used to determine if these compounds interact with the BCR-ABL protein. In total, 3,856 and 161 compounds were screened using Lipinski's rule of five and ADMET, respectively, interacted with the active site of protein. Finally, identical compounds screened from the two pharmacophore models were selected for the in vitro efficacy assay. The docking scores for these compounds are shown in Supplementary Table 4.

Compounds Inhibit the Proliferation of K562 Cells
There were 13 identical compounds (Supplementary Table 4) and we purchased six of these for the proliferation experiments: ZINC36617838, ZINC30201139, ZINC65008391, ZINC45895251, ZINC36617852, and ZINC36617849. MTT assays were performed to determine the inhibitory efficacy of compounds on the proliferation of K562 leukemia cells. The cells were incubated with various concentrations of the six compounds for 72 h and the results are shown in Figure 6. The compound, ZINC36617838, had a greater efficacy in inhibiting the proliferation of K562 leukemia cells (IC 50 10.002 µM) compared to the other ZINC compounds. However, the efficacy of ZINC36617838 was approximately 200 times lower than that of imatinib (IC 50 0.047 µM). This finding led us to select ZINC36617838 as a potential hit compound for structural optimization. The structure of ZINC36617838 and the interaction between ZINC36617838 and 1IEP is shown in Supplementary Figure 5.

Structure Optimization
In order to find new compounds with higher inhibitory efficacy, the skeleton structure of ZINC36617838 (Figure 7), was used as a substructure to search for its derivatives in the ZINC database. The compounds that had 40% similarity with ZINC36617838 were selected (circled in red, green, and blue in Figure 7). As shown in Figure 7, the analysis identified 290 compounds. The validated Hypogen1 and Hiphop1 pharmacophores were used as 3D queries to screen the database and four and 153 compounds   matched all of the pharmacophore features of Hiphop1 and Hypogen1, respectively. Subsequently, four and 148 compounds meet the standards of the Lipinski's rules of five and ADMET filters, respectively, and were docked into active site of protein to remove the detected false positives (Rampogu et al., 2018). After screening, four compounds from the Hiphop screening were also identified from the Hypogen screening. However, the predicted IC 50 values or efficacy of 29 out of the 148 compounds were less than 0.1 µM in the Hypogen pharmacophore screening pathway, and the CDOCKER INTERACTION ENERGY of six compounds were greater than 60. Therefore, these ten compounds were selected for further experimental evaluation. We purchased ZINC21710815, ZINC36617889, and ZINC20617585. The binding modes of these three compounds are shown in Figures 8A-C.
We used the MTT assay to determine the efficacy of these compounds to inhibit the proliferation of K562, BaF3/WT, and BaF3/T315I leukemia cells (Figure 8). All cells were incubated with different concentrations of each compound for 72 h. As shown in Figure 8E, imatinib and ZINC21710815 inhibited the proliferation of K562 cells with IC 50 values of 0.047 and 0.531 µM, respectively. The IC 50 value for ZINC21710815 was approximately 10 times lower than that of imatinib. Imatinib and ZINC21710815 decreased the viability of BaF3/WT leukemia cells, with IC 50 values of 0.196 and 0.512 µM, respectively ( Figure 8F). In the BaF3/T315I leukemia cells, which have the T315I mutation that produces resistance to imatinib, ZINC21710815 significantly inhibited the proliferation of BaF3/T315I leukemia cells, (IC 50 = 0.88 µM), and its IC 50 value was significantly lower than that for imatinib (IC 50 = 10.11 µM), and its IC 50 for normal CCC-HEL-1 cells was 89.587 µM (Figure 8H). The selectivity indexes (SI) were 1, 168.71, 174.97, and 101.80, respectively, for the CCC-DEL-1, K562, BaF3/WT, and BaF3/T315I cells. These results suggest that ZINC21710815 was not cytotoxic in normal CCC-HEL-1 cells and it significantly inhibited the proliferation of wild-type BCR-ABL and T315I mutated BCR-ABL leukemia cells.

ZINC21710815 Inhibits the Growth of CML Cells
Based on the proliferation experiments, we next determined the effect of ZINC21710815 on the cell cycle of the leukemia cells. The K562, BaF3/WT, and BaF3/T315I leukemia cells were incubated with various concentrations of ZINC21710815 for 24 h. These cells were incubated and stained with PI, which stains DNA, allowing for the analysis of cellular DNA content using flow cytometry. ZINC21710815 significantly increased the accumulation of K562 leukemia cells in the G 2 phase, suggesting an arrest of the cells in the G 2 phase (Figures 9A,D). In contrast, ZINC217108155 significantly increased the accumulation of BaF3/WT and BaF3/T315I cells leukemia cells in the G 1 phase (Figures 9B,C,D). These findings suggested that ZINC21710815 may inhibit the proliferation of leukemia cells by interrupting the cell cycle in CML cells.

ZINC 21710815 Induces Apoptosis in CML Cells
We conducted experiments to determine if ZINC21710815 inhibits the proliferation of K562, BaF3/WT, and BaF3/T315I leukemia cells by inducing apoptosis. The cells were cultured with various concentrations of ZINC21710815 for 24 h and the cells were stained using Annexin V-FITC/PI and analyzed using flow cytometry. As shown in Figure 10A, ZINC21710815 produce a concentration-dependent increase in the apoptosis of K562, BaF3/WT, and BaF3/T315I leukemia cells. Apoptosis is regulated FIGURE 7 | The work flow chart of structure optimization and its screening process.
by specific caspases, and cleaved caspase-3 has been shown to be a reliable marker for detecting apoptotic cells (Crowley and Waterhouse, 2016). Therefore, we measured the expression levels of caspase-3 using western blotting ( Figure 10B). ZINC21710815 significantly increased the levels of cleaved caspase-3 in BaF3/WT and BaF3/T315I leukemia cells, suggesting that this compound induces apoptosis.

ZINC21710815 Induces Autophagy in BaF3/WT and BaF3/T315I Leukemia Cells
To determine whether ZINC21710815 induced autophagy, we determined the expression of the autophagy-related proteins, LC3-II, and Beclin1. LC3-II is a hallmark protein of autophagy and the level of LC3-II expression is significantly correlated with the number of autophagosomes (Mizushima and Yoshimori, 2007;Yoshii and Mizushima, 2017;Lu et al., 2018). In the early stages of autophagosome formation, Beclin1 is considered to be an essential component for the initiation of autophagy (Hamurcu et al., 2018). As shown in Figures 10C,D, the level of LC3-II/LC3-I and Beclin1 are significantly increased in BaF3/WT and BaF3/T315I leukemia cells following incubation with ZINC21710815 for 48 h. These data suggest that ZINC21710815 induces autophagy, in part, by increasing the expression of the levels LC3-II and Beclin1 in BaF3/WT and BaF3/T315I leukemia cells.

ZINC21710815 Inhibits Tyrosine
Phosphorylation of the BCR-ABL Protein and Its Downstream Protein Targets, Signal Transducer and Activator of Transcription 5 (STAT5) and CRK Like Proto-Oncogene, Adaptor Protein (Crkl) BCR-ABL activates multiple downstream signaling pathways, including STAT5 and p-Crkl (Gupta et al., 2016). The activation of the transcription factor STAT5 is crucial for the progression of CML (Warsch et al., 2012). STAT5 was necessary for cellular proliferation and overexpression of constitutively active STAT5 could stimulate cell proliferation (de Groot et al., 1999). Crkl is an adaptor protein and a major substrate of BCR-ABL in CML cells (Oda et al., 1996). Crkl is phosphorylated by BCR-ABL and it subsequently activates other pathways that play a role in leukemic cell transformation (Senechal et al., 1996). To confirm the connection between the compound and BCR-ABL and downstream proteins, the effect of ZINC21710815 on the tyrosine phosphorylation of BCR-ABL and its downstream protein, STAT5 and Crkl, was determined using western blotting. BaF3/WT leukemia cells were incubated with various concentrations of imatinib and ZINC21710815 for 48 h. ZINC21710815 significantly reduced the phosphorylation of BCR-ABL at a concentration of 0.1 µM, decreased the phosphorylation of STAT5 at a concentration of 0.1 µM and reduced the phosphorylation of Crkl at a concentration of 1 µM in BaF3/WT leukemia cells (Figure 11).

DISCUSSION
In this study, Hypogen pharmacophore models were constructed based on 21 BCR-ABL TKIs and the best quantitative pharmacophore model consisted of five features. The correlation coefficient of Hypogen1 with the training and test set were 0.973 and 0.847, respectively. Hypogen1 was further validated by Fischer's randomization method and by using a decoy set. Hiphop pharmacophore models were constructed based on five highly activity BCR-ABL inhibitors and three inhibitors with low efficacy. The best Hiphop pharmacophore model consists of eight features. The best pharmacophore models were selected to validate the test set and were highly efficient in distinguishing active from inactive compounds. The best Hypogen and Hiphop pharmacophore models were used as 3D queries to screen the ZINC database in order to find novel potential BCR-ABL inhibitors. The hit compounds were further screened by Lipinski's rule of five, ADMET and molecular docking. Finally, three compounds, ZINC21710815, ZINC36617889, and ZINC20617585, were subjected to further analysis, of which compound ZINC21710815 had the greatest efficacy in inhibiting the proliferation of K562, BaF3/WT, and BaF3/T315I leukemia cells.
ZINC21710815, ZINC36617889, and ZINC20617585, as well as imatinib, interact with the BCR-ABL protein via a carbonyl of amide, which forms an O-NH hypogen interaction with Asp381, and all of these compounds form π-π interactions with residue Lys271. Given that imatinib binds to the ABL kinase domain known as the DFG-out conformation (Asp381-Phe382-Gly383) (Pereira et al., 2012), we hypothesize that ZINC21710815, ZINC36617889, and ZINC20617585 can induce the DFG-out conformation. Our results indicated that ZINC21710815 and imatinib have a higher -CDOCKER ENERGY compared to ZINC36617889 and ZINC 20617585 (Figure 8) and this could potentially explain why ZINC21710815 and imatinib have similar efficacies in inhibiting the proliferation of the leukemia cell lines used in this study. Imatinib inhibited the proliferation of K562 (IC 50 = 0.047 µM), BaF3/WT (IC 50 = 0.196 µM), and BaF3/T315I (IC 50 = 10.11 µM) leukemia cells. ZINC21710815 was shown to be a novel compound that inhibited the proliferation of K562 (IC 50 = 0.531 µM), BaF3/WT (IC 50 = 0.512 µM), and BaF3/T315I leukemia cells (IC 50 = 0.88 µM) in a concentration-dependent manner. ZINC21710815 and imatinib decreased the viability of K562 and BaF3/WT leukemia cells and ZINC21710815 produced a greater inhibition of BaF3/T315I leukemia cell viability compared to imatinib ( Figure 8F). Although ZINC21710815 was the lead compound, it produced a lower magnitude of inhibition of BaF3/T315I leukemia cell proliferation compared to dasatinib in BaF3/WT and BaF3/T315T leukemia cells (IC 50 = 1 nM/300 nM) and ponatinib in BaF3/WT and BaF3/T315I leukemia cells (IC 50 = 1 nM/8 nM) (Cassuto et al., 2012), we will continue to optimize the structure of ZINC21710815 to improve its efficacy, particularly in the leukemia cells with the T315I mutation.
Our mechanistic experiments suggested that the ZINC21710815-induced decrease in the K562, BaF3/WT, and BaF3/T315I leukemia is due to it producing cell cycle arrest and apoptosis. ZINC21710815 increased the accumulation of K562 leukemia cells the G 2 phase in and increased accumulation of BaF3/WT and BaF3/T315I leukemia cells in the G 1 phase. As reported, imatinib induced an accumulation in the S phase, and dasatinib caused accumulation at the G 1 phase and were efficacious in inhibiting leukemia cell proliferation (Song et al., 2013;Wu et al., 2014). In contrast, ZINC21710815 significantly increased the accumulation of cells in the G 1 or G 2 phases. These results suggest that ZINC21710815, in part, decreases the proliferation of the leukemia cell lines used in this study by producing cell cycle arrest in either the of the G 1 or G 2 phases. The anti-proliferative efficacy of ZINC21710815 may also be due to it inducing the apoptosis of CML cells as this compound increased the expression of the apoptosis-inducing protein, cleaved capase-3 in BaF3/WT and BaF3/T315I. It has been reported that the activation of the enzyme caspase-3 plays a critical role in inducing cancer cell apoptosis, and imatinib, ponatinib, and dasatinib can increase the levels of the caspase-3 by increasing the proteolysis of its inactive zymogen, pro-caspase-3, thereby increasing apoptosis   Quintás-Cardama et al., 2008;Okabe et al., 2014;Wu et al., 2014). Similarly, ZINC21710815, at 0.1 uM, increased the level of cleaved caspase-3, producing leukemia cell apoptosis. Finally, our results indicated that ZINC21710815 induces autophagy in BaF3/WT and BaF3/T315I leukemia cells, based on the increase in the expression of the autophagy-related proteins LC3-II and Beclin1, in BaF3/WT, and BaF3/T315I leukemia cells. Numerous studies have reported that the induction of autophagy by various TKIs, including imatinib and dasatinib, may protect CML cells from death, thereby promoting their survival (Bellodi et al., 2009;Helgason et al., 2011). Indeed, the inhibition of the BCR-ABL kinase activity is a main trigger of autophagy that allows the cells to survive in a stressful environment (Bellodi et al., 2009;Calabretta and Salomoni, 2012). In vitro studies indicate that the inhibition of autophagy increases the efficacy of certain TKIs in CML cells (Bellodi et al., 2009;Salomoni and Calabretta, 2009;Crowley et al., 2011;Yu et al., 2012). In the current study, we do not know if ZINC21710815's induction of autophagy increases or decreases its anti-proliferative efficacy. Therefore, future experiments, where autophagy of the three leukemia cell lines is inhibited, must be done to ascertain if ZINC21710815's efficacy is increased or decreased. However, it is important to note that the induction of autophagy in leukemia cells by imatinib increases the sequestration of the BCR-ABL protein in autophagosomes, decreasing its levels, thereby attenuating its oncogenic efficacy (Elzinga et al., 2013). Furthermore, the inhibition of autophagy in wild type and T315I leukemia cells decreases the biodegradation of the BCR-ABL protein, thus increasing its stability (Shinohara et al., 2019).
In order to further determine the efficacy ZINC21710815 on the BCR-ABL kinase, we measured the effect of ZINC21710815 on levels of the proteins STAT5 and Crkl, which are substrates for the BCR-ABL kinase. BCR-ABL kinase phosphorylates and activates the protein transcription factor, STAT5, which is translocated to the nucleus, where it increases the transcription of genes that produce proteins that promote cell survival and proliferation (Huang et al., 2002). Our results indicated that ZINC21710815 and imatinib, at 0.1 µM, significantly decreased p-STAT5 level in BaF3/WT leukemia cells. Our imatinib results are congruent with previous studies in leukemia cells (Dong et al., 2018;Tu et al., 2018) and similarly, dasitinib also decreases p-STAT5 levels in BaF3/WT leukemia cells (Fiskus et al., 2006). We also determined the effect of ZINC21710815 on the levels of phosphorylated Crkl, a protein that is phosphorylated by the BCR-ABL (Bhat et al., 1997). Phosphorylated Crkl is one of main tyrosyl-phosphoproteins present in peripheral blood cells of CML patients, where it functions as a nuclear adaptor and transcriptional activator in BCR-ABL expressing cells (Nichols et al., 1994;Rhodes et al., 2000). Furthermore, the levels of phosphorylated Crkl can be used as an indirect marker of BCR-ABL function and its levels can be used to determine the status of BCR-ABL kinase activity due to the fact that Crkl is only phosphorylated by BCR-ABL (Nichols et al., 1994;Hamilton et al., 2006). Our in vitro results indicated ZINC21710815 significantly decreased the levels of phosphorylated Crkl in BaF3/WT leukemia cell lines, the magnitude of inhibition of phosphorylated Crkl levels by ZINC21710815 was similar to that of imatinib. Similarly, dasatinib has been reported to decrease p-Crkl (Fiskus et al., 2006) levels with an efficacy comparable to that of ZINC21710815. Based on these findings, we hypothesize that ZINC21710815 may decrease leukemia cell proliferation by decreasing the phosphorylation of BCR-ABL, STAT5, and Crkl.
Our experimental results showed that the compound ZINC21710815 inhibits the proliferation of leukemia cells, inhibits cell cycle, induces autophagy and apoptosis and inhibits the tyrosine phosphorylation of BCR-ABL target protein and downstream proteins. However, ZINC21710815 had no significant effect on the expression of BCR-ABL kinase in the T315I leukemia cells. Therefore, in the future, we will conduct experiments to optimize the structure of ZINC21710815 to increase its inhibitory efficacy in T315I leukemia cells.
The BCR-ABL fusion gene can be seen in CML patients, it is also found in other types of leukemia patients, such as acute myeloid leukemia (AML) patients (Bucur et al., 2013), B-cell acute lymphoblastic leukemia (B-ALL) patients (Zhou et al., 2018), and mixed phenotype acute leukemia (MPAL) patients (Alexander et al., 2018). ZINC21710815 could decrease the expression of BCR-ABL, which may indicate that the compound may also have an inhibitory effect in other cancers.
In conclusion, ZINC21710815 significantly decreased the proliferation of K562, BaF3/WT, and BaF3/T315I leukemia cells due to it (1) producing cell cycle arrest; (2) inducing apoptosis; and (3) inhibiting the phosphorylation of BCR-ABL kinase and the downstream targets, STAT5 and Crkl, in BaF3/WT leukemia cells. As stated above, it remains to be determined if inducing autophagy increases or decreases the anti-proliferative efficacy of ZINC21710815. Overall, our results suggest that ZINC21710815 is an inhibitor of BCR-ABL tyrosine kinase activity. Additional studies must be done to determine the in vivo efficacy and toxicity of ZINC21710815.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.

AUTHOR CONTRIBUTIONS
T-TH and XW designed research. T-TH, S-JQ, and Z-XW performed research. T-TH and Z-NZ analyzed data. T-TH and CA wrote the manuscript. CA and Z-SC revised the manuscript. J-ZL and Z-SC supervised the research. All authors contributed to the article and approved the submitted version.