In Silico Approach Using Free Software to Optimize the Antiproliferative Activity and Predict the Potential Mechanism of Action of Pyrrolizine-Based Schiff Bases

In the current study, a simple in silico approach using free software was used with the experimental studies to optimize the antiproliferative activity and predict the potential mechanism of action of pyrrolizine-based Schiff bases. A compound library of 288 Schiff bases was designed based on compound 10, and a pharmacophore search was performed. Structural analysis of the top scoring hits and a docking study were used to select the best derivatives for the synthesis. Chemical synthesis and structural elucidation of compounds 16a–h were discussed. The antiproliferative activity of 16a–h was evaluated against three cancer (MCF7, A2780 and HT29, IC50 = 0.01–40.50 μM) and one normal MRC5 (IC50 = 1.27–24.06 μM) cell lines using the MTT assay. The results revealed the highest antiproliferative activity against MCF7 cells for 16g (IC50 = 0.01 μM) with an exceptionally high selectivity index of (SI = 578). Cell cycle analysis of MCF7 cells treated with compound 16g revealed a cell cycle arrest at the G2/M phase. In addition, compound 16g induced a dose-dependent increase in apoptotic events in MCF7 cells compared to the control. In silico target prediction of compound 16g showed six potential targets that could mediate these activities. Molecular docking analysis of compound 16g revealed high binding affinities toward COX-2, MAP P38α, EGFR, and CDK2. The results of the MD simulation revealed low RMSD values and high negative binding free energies for the two complexes formed between compound 16g with EGFR, and CDK2, while COX-2 was in the third order. These results highlighted a great potentiality for 16g to inhibit both CDK2 and EGFR. Taken together, the results mentioned above highlighted compound 16g as a potential anticancer agent.


Introduction
Cancer is still one of the leading causes of death in the world [1]. Among females, breast cancer is the most common type of cancers, whereas, in males, colon cancer was ranked third in the incidence in the USA in 2021. Due to the consistently high rates of ranked third in the incidence in the USA in 2021. Due to the consistently high rates cancer incidence and mortality, research in this field is always in demand and of contin ous interest.
Before thirty years, the antiproliferative effect of nonsteroidal anti-inflammato drugs (NSAIDs) was discovered [2]. Since then, there has been much greater interest fro many researchers around the world, hoping to explore the antiproliferative potential NSAIDs and investigate their mechanisms of action. The selective COX-2 inhibitors, co ibs 1-5 (Figure 1), attracted much attention in cancer research, although most of them we removed from the market [3]. Considering their chemical structures, we can observe th compounds 1-5 have similar pharmacophoric features which include three aromatic rin and hydrophobic atom/group attached to the five-membered ring. The aryl rings in the compounds are either unsubstituted or substituted at the para-position with hydropho groups (CH3/Br) or hydrogen bond-forming groups (SO2CH3/SO2NH2) which occupy t side pocket in COX-2. Among these drugs, celecoxib 2 was extensively studied for its antiproliferative, ther alone or in combination with anticancer drugs [4][5][6][7][8]. Celecoxib 2 showed antipro erative activity against breast and colon cancer cells [4][5][6]. It also showed a chemopreve tive effect and was approved by the FDA for familial adenomatous polyposis (FAP) [ Moreover, Celecoxib 2 exhibited anticancer and antimetastatic effects against seve types of ovarian cancer cells [7].
Substantial research was also conducted to investigate the mechanism of action NSAIDs [2,[8][9][10][11]. Based on these studies, the anticancer activity of NSAIDs can be tributed to the COX-dependent mechanisms [2]. However, COX-independent mech nisms were also reported to mediate the anticancer potential of celecoxib 2 [8][9][10][11]. T induction of apoptosis in HT29 cells by celecoxib 2 was attributed to the inhibition of phosphoinositide-dependent protein kinase-1 [9]. The inhibition of Akt kinase activati could also play an important role in the induction of apoptosis in prostate cancer ce [10]. Mechanistic studies also showed that celecoxib 2 blocks the activation of MAP p kinase and downregulates COX-2 [11].
To date, celecoxib is still the only NSAID that has been approved for FAP [6]. T may be attributed to the weak anticancer activities of NSAIDs or to their toxic side effe of these drugs [2,12]. According to this, extensive research in this field should be co ducted to develop new scaffolds with potent anticancer activities and good safety profil Among these drugs, celecoxib 2 was extensively studied for its antiproliferative, either alone or in combination with anticancer drugs [4][5][6][7][8]. Celecoxib 2 showed antiproliferative activity against breast and colon cancer cells [4][5][6]. It also showed a chemopreventive effect and was approved by the FDA for familial adenomatous polyposis (FAP) [6]. Moreover, Celecoxib 2 exhibited anticancer and antimetastatic effects against several types of ovarian cancer cells [7].
Substantial research was also conducted to investigate the mechanism of action of NSAIDs [2,[8][9][10][11]. Based on these studies, the anticancer activity of NSAIDs can be attributed to the COX-dependent mechanisms [2]. However, COX-independent mechanisms were also reported to mediate the anticancer potential of celecoxib 2 [8][9][10][11]. The induction of apoptosis in HT29 cells by celecoxib 2 was attributed to the inhibition of 3-phosphoinositide-dependent protein kinase-1 [9]. The inhibition of Akt kinase activation could also play an important role in the induction of apoptosis in prostate cancer cells [10].
Mechanistic studies also showed that celecoxib 2 blocks the activation of MAP p38 kinase and downregulates COX-2 [11].
To date, celecoxib is still the only NSAID that has been approved for FAP [6]. This may be attributed to the weak anticancer activities of NSAIDs or to their toxic side effects of these drugs [2,12]. According to this, extensive research in this field should be conducted to develop new scaffolds with potent anticancer activities and good safety profiles.
Investigation of the chemical structure of compounds 7-10 also revealed similar pharmacophoric features which include a central five-membered pyrrole ring attached to two (un)substituted phenyl rings. These features are also similar to those compounds 1-5 ( Figure 1). However, compound 10, which lacks the acidic group, exhibited higher antiproliferative activity compared to both celecoxib 2 and licofelone 7 [6,17].
The high antiproliferative activity of compound 10 may deserve further investigation and optimization. Accordingly, the current study aimed to optimize the antiproliferative activity of compound 10. This aim will be achieved by designing a series of new analogs and evaluate their antiproliferative activities. The design of the new analogs will use different techniques of computer-aided drug design [21][22][23], which proved successful in optimization of lead compounds, will be used to in the optimization of compound 10. Licofelone 7 (Figure 2) also displayed a strong anti-inflammatory activity mediated by inhibition of COX and 5-LOX enzymes [15,16]. Compounds 8-9 are also pyrrolizine derivatives which lack the carboxylic acid group of licofelone; however, they also exhibited inhibitory activities against COXs [15]. Biological evaluation of licofelone 7 also revealed cytotoxic activities against breast, colon, and prostate cancer cell lines [17][18][19]. In addition, compound 10 showed moderate anti-inflammatory activity compared to ibuprofen and cytotoxic activity against three cancer cell lines [20]. Mechanistic studies of compound 10 revealed weak inhibitory activity against COX-2.
Investigation of the chemical structure of compounds 7-10 also revealed similar pharmacophoric features which include a central five-membered pyrrole ring attached to two (un)substituted phenyl rings. These features are also similar to those compounds 1-5 ( Figure 1). However, compound 10, which lacks the acidic group, exhibited higher antiproliferative activity compared to both celecoxib 2 and licofelone 7 [6,17].
The high antiproliferative activity of compound 10 may deserve further investigation and optimization. Accordingly, the current study aimed to optimize the antiproliferative activity of compound 10. This aim will be achieved by designing a series of new analogs and evaluate their antiproliferative activities. The design of the new analogs will use different techniques of computer-aided drug design [21][22][23], which proved successful in optimization of lead compounds, will be used to in the optimization of compound 10.

Design of the Compound Library
To optimize the antiproliferative activity of compound 10, a small compound library of pyrrolizine-based Schiff bases was designed. The new analogs were obtained by modification of the chemical structure of compound 10. In this study, we will focus mainly on variation of the substituents on the two phenyl rings (A and B), Figure 3.

Design of the Compound Library
To optimize the antiproliferative activity of compound 10, a small compound libra of pyrrolizine-based Schiff bases was designed. The new analogs were obtained by mo ification of the chemical structure of compound 10. In this study, we will focus mainly variation of the substituents on the two phenyl rings (A and B), Figure 3. The compound library included 288 Schiff bases designed bearing electron-don ing/withdrawing substituents at the ortho/meta/para-positions of the two phenyl rings. T chemical structures of these compound library are provided in supplementary data (F ures S52-S61).

Pharmacophore Search
To identify the new analogs with the highest potential inhibitory activity agai COX-2, pharmacophore search of the compound library was performed using Pharm (http://pharmit.csb.pitt.edu) [24]. This software package is available online and provide free tool to perform pharmacophore/shape search of the large compound libraries (up 10 million compounds). In an attempt to build a valid pharmacophore to screen our co pound library, three of the crystal structures of COX-2 bound to the selective inhibit SC-558 1 (pdb: 1CX2) [25], celecoxib 2 (pdb: 3LN1) [26], and rofecoxib 3 (pdb: 5KIR) [ were analyzed to identify different types of binding interactions. Analysis of the bindi interactions was done using Discovery Studio Visualizer (DSV) [28]. Among the th compounds 1-3, celecoxib 2 displayed the highest number of hydrogen bonds (supp mentary data, Figure S1). However, rofecoxib 3 exhibited only one conventional hydrog bond with Arg513 and displayed higher inhibitory activity and selectivity against COX (IC50 = 0.53 μM, SR = 35.5) compared to celecoxib 2 [29]. In addition, SC-558 1 also show higher inhibitory activity and selectivity toward COX-2 compared to celecoxib [25], hough it formed a fewer number of hydrogen bonds. Moreover, analysis of the binding int actions of compounds 1-3 revealed multiple hydrophobic interactions with the amino ac in COX-2.
Based on the above results, a simple pharmacophore model A ( Figure 4) was gen ated based on the binding interactions of SC-558 with COX-2 (pdb: 1CX2). This mo consists of seven pharmacophoric features including three aromatic rings, hydropho groups and one hydrogen bond acceptor. The compound library included 288 Schiff bases designed bearing electron-donating/ withdrawing substituents at the ortho/meta/para-positions of the two phenyl rings. The chemical structures of these compound library are provided in Supplementary Data ( Figures S52-S61).

Pharmacophore Search
To identify the new analogs with the highest potential inhibitory activity against COX-2, pharmacophore search of the compound library was performed using Pharmit (http://pharmit.csb.pitt.edu) [24]. This software package is available online and provides a free tool to perform pharmacophore/shape search of the large compound libraries (up to 10 million compounds). In an attempt to build a valid pharmacophore to screen our compound library, three of the crystal structures of COX-2 bound to the selective inhibitors SC-558 1 (pdb: 1CX2) [25], celecoxib 2 (pdb: 3LN1) [26], and rofecoxib 3 (pdb: 5KIR) [27] were analyzed to identify different types of binding interactions. Analysis of the binding interactions was done using Discovery Studio Visualizer (DSV) [28]. Among the three compounds 1-3, celecoxib 2 displayed the highest number of hydrogen bonds (Supplementary Data, Figure S1). However, rofecoxib 3 exhibited only one conventional hydrogen bond with Arg513 and displayed higher inhibitory activity and selectivity against COX-2 (IC 50 = 0.53 µM, SR = 35.5) compared to celecoxib 2 [29]. In addition, SC-558 1 also showed higher inhibitory activity and selectivity toward COX-2 compared to celecoxib [25], although it formed a fewer number of hydrogen bonds. Moreover, analysis of the binding interactions of compounds 1-3 revealed multiple hydrophobic interactions with the amino acids in COX-2.
Based on the above results, a simple pharmacophore model A ( Figure 4) was generated based on the binding interactions of SC-558 with COX-2 (pdb: 1CX2). This model consists of seven pharmacophoric features including three aromatic rings, hydrophobic groups and one hydrogen bond acceptor. To perform the pharmacophore search, the compound library was first upload the Pharmit server. The pharmacophore search was performed using the pharmacop model A, Figure 4. The active site of COX-2 enzyme where each derivative of the pound library was located overlaid with SC-558 was illustrated in Figure 5. To perform the pharmacophore search, the compound library was first uploaded to the Pharmit server. The pharmacophore search was performed using the pharmacophore model A, Figure 4. The active site of COX-2 enzyme where each derivative of the compound library was located overlaid with SC-558 was illustrated in Figure 5. The Pharmacophore search was performed, and the results in the form o based on their root-mean-square deviation (RMSD) were presented in Table 1 ular weight (MWs) and the number of rotatable bonds (RBs) in each hit wer lated. 10ipyr 0.777 510 The Pharmacophore search was performed, and the results in the form of hits ranked based on their root-mean-square deviation (RMSD) were presented in Table 1. The molecular weight (MWs) and the number of rotatable bonds (RBs) in each hit were also calculated.
Analysis of the results of the pharmacophore search revealed that the top seven scoring hits 1-7 bear the methoxy group on the phenyl ring A, Table 1. Of these derivatives, hit 1 with the methoxy group on ring A and dimethylamino group at para-position of the ring (B) exhibited the lowest RMSD. The chemical structure of the top fifteen hits 1-15, overlaid with SC-558 were arranged based on their RMSD values as illustrated in Figure 6.
The results of the pharmacophore search were analyzed to identify the substituents on the phenyl rings of the top-scoring hits. The results revealed that the order of substituents on the phenyl ring (A) in the following sequence: OCH 3 (in the top 7 Hits), Cl, I, CH 3 , F, and Br. On the other hand, the substituents of phenyl ring (B) were in the following order: N(CH 3 ) 2 , OCH 3 , Cl, I, CH 3 , F, and Br, Figure 7.  The results of the pharmacophore search were analyzed to identify the substituents on the phenyl rings of the top-scoring hits. The results revealed that the order of substituents on the phenyl ring (A) in the following sequence: OCH3 (in the top 7 Hits), Cl, I, CH3, F, and Br. On the other hand, the substituents of phenyl ring (B) were in the following order: N(CH3)2, OCH3, Cl, I, CH3, F, and Br, Figure 7. Although hit 4, hit 9, and hit 13 have displayed low RMSD values (Table 1), they have high molecular weights (>500 daltons). Considering the limits of Lipinski's rule, the three hits will be excluded from the compounds selected for the chemical synthesis.  Although hit 4, hit 9, and hit 13 have displayed low RMSD values (Table 1), they have high molecular weights (>500 daltons). Considering the limits of Lipinski's rule, the three hits will be excluded from the compounds selected for the chemical synthesis.

Compounds Selected for the Synthesis
Based on the above results of the pharmacophore search, hit 1 was selected for the next step in this study. The three isomers of hit 1 (Figure 8) were evaluated for their binding affinities toward COX-2 in a preliminary docking study. The study was to evaluate the impact of the position of the methoxy group on the binding affinities toward COX-2. The results of this study are provided in the Supplementary Data ( Figure S2). The results revealed higher binding affinities for hits 1 m,p than hit 1 o which has a methoxy group in ortho-position. Investigation of the binding interactions of the three derivatives revealed that hit 1 o lacks any conventional hydrogen bonds with COX-2, while hits 1 m,p showed two conventional hydrogen bonds each.  In addition, logP (ilogP) and synthetic accessibility (SA) of the three hits 1o,m,p w calculated using SwissADME (http://www.swissadme.ch) [30], while their drug-liken scores (DLSs) were calculated by Molsoft (http://molsoft.com/mprop/) [31]. The resu are presented in Figure 8. The results revealed the lowest SA and highest DLS for hit These results are also matched with the substitution pattern of compounds 1-10, which clude phenyl groups substituted at the para-position (Figures 1 and 2). Accordingly, hit was used as the first target compound to be synthesized in this study.
To study the SAR of hit 1p, the electron-donating (methoxy) group was replaced w an electron-withdrawing flour-atom to investigate the impact of the electronic effect substituents on their antiproliferative activities, Figure 9. In addition, four other substi ents including two electron-donating (-N(CH3)2 and -CH3) and two electron-withdrawi (-F and -Cl) atoms. These substituents were selected based on the results of the pharm cophore search, Figure 7. In addition, logP (ilogP) and synthetic accessibility (SA) of the three hits 1 o,m,p were calculated using SwissADME (http://www.swissadme.ch) [30], while their drug-likeness scores (DLSs) were calculated by Molsoft (http://molsoft.com/mprop/) [31]. The results are presented in Figure 8. The results revealed the lowest SA and highest DLS for hit 1 p . These results are also matched with the substitution pattern of compounds 1-10, which include phenyl groups substituted at the para-position (Figures 1 and 2). Accordingly, hit 1 p was used as the first target compound to be synthesized in this study.
To study the SAR of hit 1 p , the electron-donating (methoxy) group was replaced with an electron-withdrawing flour-atom to investigate the impact of the electronic effect of substituents on their antiproliferative activities, Figure 9. In addition, four other substituents including two electron-donating (-N(CH 3 ) 2 and -CH 3 ) and two electron-withdrawing (-F and -Cl) atoms. These substituents were selected based on the results of the pharmacophore search, Figure 7.

Synthesis of the New Derivatives
In this section, the starting materials 12, 14a,b, and 15a,b (Scheme 1) were prepared according to the previously reported procedures [32][33][34]. In addition, the target Schiff bases 16a-h were prepared by refluxing compounds 15a,b with the appropriate aldehydes following the previous report [20].

Synthesis of the New Derivatives
In this section, the starting materials 12, 14a,b, and 15a,b (Scheme 1) were prepar according to the previously reported procedures [32][33][34]. In addition, the target Schiff b ses 16a-h were prepared by refluxing compounds 15a,b with the appropriate aldehyd following the previous report [20]. The IR spectra of 16a-h revealed an absorption band at the range of 2209-2214 cm −1 indicating the cyano groups. In addition, an absorption band was also observed at the range of 1656-1669 cm −1 indicating the carbonyl groups in 16a-h.
The proton magnetic spectrum of compound 16a,b revealed a singlet signal at 3.13 ppm indicating the six protons of the N(CH3)2 group. Four doublets at the range of δ 6.79-7.82 ppm indicating the para-substituted phenyl rings in compound 16a, while the aromatic protons of 16b was observed at the range of δ 6.86-7.95. Two singlet signals at the range of δ 8.77-11.00 ppm indicating the benzylidene (N=CH) and amide (CONH) protons of 16a,b.
The 13 C-NMR also revealed a signal at δ 40.27 and 40.67 ppm indicating the carbon atoms of N(CH3)2 group in 16a and 16b, respectively. DEPT C 135 spectrum showed a signal at δ 159.59 and 158.43 ppm indicating benzylidene carbon (N=CH) in 16a and 16b, respectively. The IR spectra of 16a-h revealed an absorption band at the range of 2209-2214 cm −1 indicating the cyano groups. In addition, an absorption band was also observed at the range of 1656-1669 cm −1 indicating the carbonyl groups in 16a-h.
The proton magnetic spectrum of compound 16a,b revealed a singlet signal at 3.13 ppm indicating the six protons of the N(CH 3 ) 2 group. Four doublets at the range of δ 6.79-7.82 ppm indicating the para-substituted phenyl rings in compound 16a, while the aromatic protons of 16b was observed at the range of δ 6.86-7.95. Two singlet signals at the range of δ 8.77-11.00 ppm indicating the benzylidene (N=CH) and amide (CONH) protons of 16a,b. Copies of the spectra data including IR, mass, 1 H-NMR, 13 C-NMR, and DEPT C 135 spectra of compounds 16a-h are provided in Supplementary Data ( Figures S3-S45).

Antiproliferative Activity Antiproliferative Activity Assay
The antiproliferative activities of the eight Schiff bases 16a-h against three cancer cell lines (MCF-7, A2780, and HT29) were evaluated using the MTT assay. The selection of these cell lines was done to compare the antiproliferative activities of the new compounds with those of compound 10 [20]. The assay was performed following the previous report [35]. These results of the MTT assay are presented in Table 2. Compounds 16a-h exhibited their antiproliferative activity at IC 50 values in the range of 0.01-40.50 µM against the three cancer cell lines compared to compound 10 (IC 50 = 0.33-0.44 µM). These results also indicated high antiproliferative activity for compounds 16a-h compared to lapatinib (IC 50 = 5.98-13.22 µM).
In addition, compounds 16b,c, 16e, and 16g showed higher activity against MCF7 cells than compound 10. Moreover, higher antiproliferative activities were observed for compounds 16a, 16d-g against HT29 cells compared to compound 10. However, only compound 16c was more active against A2780 cells than compound 10. Among the new derivatives, compound 16g was the most active against MCF7 cells, while 16c and 16d showed the highest antiproliferative activities against A2780 and HT29 cells, respectively, Table 2.

Evaluation of Antiproliferative Selectivity
Selectivity of the antiproliferative agents toward cancer cells plays a critical role in the development of these agents as potential anticancer drugs. Although the lead compound 10 exhibited high antiproliferative activity, its toxicity and selectivity were not evaluated [20]. In the current study, all of the new compounds 16a-h were evaluated for the antiproliferative activities against the normal MRC5 cells. The cells were treated compounds for 72 h, and the IC 50 values were calculated in Table 3. The aim of this study was to assess the toxicity of the new compounds against normal MRC5 cells. In addition, the IC 50 values were also used to calculate the selectivity index (SI) which measures the selective cytotoxicity of the new compounds. The results showed that compounds 16a-h exhibited their antiproliferative activities against MRC5 cells at IC 50 values in the range of 1.27-24.06 µM compared to lapatinib (IC 50 = 14.89 µM). In addition, the new compounds exhibited SIs in the range of 0.03-802 toward the three cancer cell lines compared to MRC5 cells. Compounds 16b,c,e,g showed SIs > 10 toward MCF7 cells, while 16c,d,g were selective to HT29 cells (SIs >10). On the other hand, only two compounds (16c,g) showed more than 10-fold higher selectivity toward A2780 cells compared to MRC5 cells, Table 2.
Among these new compounds, 16c showed the least toxicity toward MRC5 cells. Moreover, 16c showed the highest SI against MCF7 and A2780 cells, while 16g was the most selective to HT29 cell line. The results also revealed an SI in the range of 10.32-578.00 for compound 16g against the three cancer cell lines compared to MRC5 cells, Table 2.

Structure-Activity Relationship (SAR)
In an attempt to illustrate the impact of different substituents on the antiproliferative activity and selectivity of the new compounds 16a-h, the study of SAR of these compounds was illustrated in Figure 10. To evaluate the effect of different substituents on activity/selectivity, we started with compound 16a which exhibited IC 50 values in the range of 0.19-40.50 µM against the tested cancer cell lines, and SI in the range of 0.03-6.68. We found that the replacement of the methoxy group in 16a by fluoro atom caused a significant increase in the antiproliferative activities against MCF7 and A2790 cells, while activity against HT29 cells was decreased. An improvement in the selectivity of compound 16b toward MCF7 cells compared to 16a was also observed, Figure 10. Compound 16c also showed higher antiproliferative activity against MCF7 and A2790 cells compared to 16a. Among the new compounds, compound 16c was the most selective toward MCF7 and A2790 cells. However, a sharp decrease in antiproliferative activity and selectivity toward MCF7 and A2780 cells was observed upon replacement of the methoxy group in 16c by fluoro atoms, Figure 10.
Meanwhile, the replacement of the dimethylamino group in 16a by fluoro resulted in an increase in the antiproliferative activity and selectivity toward MCF7 and A2790 cells. However, the replacement of the methoxy group in 16e by F atom decreased the antiproliferative activity against MCF7, while the activities against A2780 HT29 cells were slightly improved. Similarly, the replacement of the dimethylamino group in 16a by chloro also increased the antiproliferative activity and selectivity toward MCF7 and A2790 cells. On the other hand, replacement of methoxy group in 16g by F also decreased antiproliferative activity and selectivity against the three cancer cell lines.
In conclusion, the antiproliferative activity and selectivity of compounds 16c,e,g toward MCF7 cells were decreased upon replacement of methoxy group by the flour atom. These results suggested that the methoxy analogs are more favored for high antiproliferative activity and selectivity.

Determination of Cell Cycle Perturbations
Cell cycle analysis of MCF7 cells was also performed to investigate the mechanism which mediates the antiproliferative activity of the new compounds. Compound 16g, the most active in MTT assay (Table 2), was selected for this assay. The assay was performed according to the previous report [36]. The results are represented in Figure 11. Compound 16c also showed higher antiproliferative activity against MCF7 and A2790 cells compared to 16a. Among the new compounds, compound 16c was the most selective toward MCF7 and A2790 cells. However, a sharp decrease in antiproliferative activity and selectivity toward MCF7 and A2780 cells was observed upon replacement of the methoxy group in 16c by fluoro atoms, Figure 10.
Meanwhile, the replacement of the dimethylamino group in 16a by fluoro resulted in an increase in the antiproliferative activity and selectivity toward MCF7 and A2790 cells. However, the replacement of the methoxy group in 16e by F atom decreased the antiproliferative activity against MCF7, while the activities against A2780 HT29 cells were slightly improved. Similarly, the replacement of the dimethylamino group in 16a by chloro also increased the antiproliferative activity and selectivity toward MCF7 and A2790 cells. On the other hand, replacement of methoxy group in 16g by F also decreased antiproliferative activity and selectivity against the three cancer cell lines.
In conclusion, the antiproliferative activity and selectivity of compounds 16c,e,g toward MCF7 cells were decreased upon replacement of methoxy group by the flour atom. These results suggested that the methoxy analogs are more favored for high antiproliferative activity and selectivity.

Determination of Cell Cycle Perturbations
Cell cycle analysis of MCF7 cells was also performed to investigate the mechanism which mediates the antiproliferative activity of the new compounds. Compound 16g, the most active in MTT assay (Table 2), was selected for this assay. The assay was performed according to the previous report [36]. The results are represented in Figure 11.

Annexin V FITC/PI Apoptosis Assay
To investigate the ability of compound 16g to induce apoptosis, MCF7 cells wer treated with three different concentrations (0.01, 0.05, and 0.10 μM) of the test compound following the previous report [37]. The results showed that compound 16g induced dose-dependent increase in the apoptotic events in MCF7 cells compared to the contro group (from 7.9-17.8%), Figure 13.

Annexin V FITC/PI Apoptosis Assay
To investigate the ability of compound 16g to induce apoptosis, MCF7 cells were treated with three different concentrations (0.01, 0.05, and 0.10 µM) of the test compound following the previous report [37]. The results showed that compound 16g induced a dose-dependent increase in the apoptotic events in MCF7 cells compared to the control group (from 7.9-17.8%), Figure 13. The increase in apoptotic events was associated with a dose-dependent decrea the number of the living cells, Figure 14. These results also indicated the ability of c pound 16g to induce apoptosis in MCF7 cells at much lower concentration than c pound 10 [20]. The increase in apoptotic events was associated with a dose-dependent decrease in the number of the living cells, Figure 14. These results also indicated the ability of compound 16g to induce apoptosis in MCF7 cells at much lower concentration than compound 10 [20].

Target Prediction
In the current study, the new compounds 16a-h exhibited potent antiproliferative activity against one or more of the three cancer cell lines at IC 50 values much lower than those of NSAIDs [2,38]. These results could be attributed to the ability of the new compounds to act on other targets besides the COX-2 enzyme.
To identify the potential targets which could contribute to the antiproliferative activity of the new compounds, target prediction of compound 16g, the most active in MTT assay, was evaluated using SwissTargetPrediction [39]. The prediction of the most probable molecular targets by SwissTargetPrediction is done based on 2/3D similarity with a library of 370,000 active compounds. The results in the form of a pie chart showing the top potential targets classified as electrochemical transporters, family A G coupled proteins, enzymes, membrane receptor, oxidoreductase enzyme, protease, and unclassified protein, Figure 15.  OR PEER REVIEW 19 of 39

Target Prediction
In the current study, the new compounds 16a-h exhibited potent antiproliferative activity against one or more of the three cancer cell lines at IC50 values much lower than those of NSAIDs [2,38]. These results could be attributed to the ability of the new compounds to act on other targets besides the COX-2 enzyme.
To identify the potential targets which could contribute to the antiproliferative activity of the new compounds, target prediction of compound 16g, the most active in MTT assay, was evaluated using SwissTargetPrediction [39]. The prediction of the most probable molecular targets by SwissTargetPrediction is done based on 2/3D similarity with a library of 370,000 active compounds. The results in the form of a pie chart showing the top potential targets classified as electrochemical transporters, family A G coupled proteins, enzymes, membrane receptor, oxidoreductase enzyme, protease, and unclassified protein, Figure 15. In addition, the results of SwissTargetPrediction also include a detailed report showing the potential targets that could mediate the activity of compound 16g. The report includes the targets arranged in descending order of their probabilities, their common names, and classification. The results of the target prediction of compound 16g are provided in supplementary data (Table S1).
Among 100 entries, the potential targets of compound 16g were investigated to identify the targets which could contribute to the anti-inflammatory/antiproliferative activities. Among these targets, COX-2, and MAP P38α could contribute to anti-inflammatory and antiproliferative activities of 16g. In addition, four oncogenic kinases (EGFR, CDK2, BRAF, and VEGFR1) are also among the potential targets that may contribute to the antiproliferative activity of compound 16g, Table 4.
Similarly, target prediction was also performed for the remaining compounds (16af, and 16h). The results of this study were analyzed to identify the shared molecular tar- In addition, the results of SwissTargetPrediction also include a detailed report showing the potential targets that could mediate the activity of compound 16g. The report includes the targets arranged in descending order of their probabilities, their common names, and classification. The results of the target prediction of compound 16g are provided in Supplementary Data (Table S1).
Among 100 entries, the potential targets of compound 16g were investigated to identify the targets which could contribute to the anti-inflammatory/antiproliferative activities. Among these targets, COX-2, and MAP P38α could contribute to anti-inflammatory and antiproliferative activities of 16g. In addition, four oncogenic kinases (EGFR, CDK2, BRAF, and VEGFR1) are also among the potential targets that may contribute to the antiproliferative activity of compound 16g, Table 4. Table 4. Potential targets of compounds 16a-h based on the results of SwissTargetPrediction.

Comp.
Molecular Targets   COX-2  P38α  EGFR  CDK2  BRAF  VEGFR1 16a indicated that the enzyme/kinase could be a potential target; (−) indicated that the target is not a potential target for the new compounds.
Similarly, target prediction was also performed for the remaining compounds (16a-f, and 16h). The results of this study were analyzed to identify the shared molecular targets with those of compound 16g, Table 4.
The roles of COX-2 and MAP P38α in different types of solid tumors were discussed in several reports [2,40]. Accordingly, small molecule inhibitors of these two target proteins exhibited antiproliferative activities against several types of cancer cells [2,40,41]. Moreover, many of the small molecule inhibitors that target EGFR, CDK2, BRAF, and VEGFR1 kinases were also reported with potent anticancer activity [42][43][44][45]. Among these six targets (Table 4), COX-2 and CDK2 were identified as potential targets for all the new compounds, while seven of the new compounds were expected to act on EGFR and VEGFR1. On the other hand, MAP P38α and BRAF were identified as the potential targets for five and four of the new compounds, respectively, Table 4.

Molecular Docking and Binding Mode Analysis
In the current study, the new compounds 16a-h exhibited potent antiproliferative activity against one or more of the three cancer cell lines, Table 2. Considering the weak antiproliferative activity of most of the NSAIDs [2,38], the high potency of the new compounds 16a-h against the tested cancer cell lines could be attributed to their multitarget activity. This conclusion was supported with the results of the target prediction, Table 4. These findings are also in concordance with the findings in our previous report [32]. Accordingly, a molecular docking study of compound 16g was performed into the molecular targets identified in the target prediction test, Table 4. The aim of this study was to evaluate the binding affinity, orientation, and interactions of compound 16g against those of the co-crystallized ligands of the six target proteins. The crystal structures of the target proteins were downloaded for the Protein Data Bank. The study was done by AutoDock 4.2 [46], and the results were analyzed and visualized using DSV [28].
Validation of the docking procedures was performed for each of the six targets. The co-crystallized ligands were re-docked into their corresponding proteins and the binding orientation and interactions were compared with those of the co-crystallized ligands. The results of the validation process are provided in the Supplementary Data (Figures S46-S51).

Docking into the Targets Involved in Inflammation
In addition to their role in inflammation, COX-2 and MAP p38α also play important roles in cell proliferation and apoptosis [2,40]. Several small molecule inhibitors of COXs and MAP p38α exhibited potent antiproliferative activity against several types of solid tumors [2,40]. In the current study, COXs and MAP p38α were identified as potential targets for compound 16g. Accordingly, a molecular docking study of compound 16g into the two proteins was performed.

Docking into COX-2
Compound 16g was docked into the active site of COX-2 (pdb: 3LN1). The results revealed a binding free energy of −10.13 kcal/mol compared to −10.27 kcal/mol for celecoxib. Analysis of the binding mode of compound 16g revealed superposition of the phenyl rings over the two phenyl rings of celecoxib, Figure 16. The phenyl ring (A) superposed with the sulfonamide-bearing ring in celecoxib, where the methoxy group extended with the sulfonamide group into the side pocket and formed one carbon-hydrogen bond with Gln178. Moreover, the chlorophenyl moiety in compound 16g superposed partially with the tolyl moiety of celecoxib, where the chloro atom also formed similar hydrophobic interactions with Leu370 and Trp373 like the methyl group in celecoxib, Figure 16.  Investigation of the binding interaction of compound 16g into COX-2 revealed four conventional hydrogen bonds between the cyano and carbonyl groups with Ser516 Arg106, and Tyr341 amino acids in COX-2, Figure 16. In addition, the methoxy group in compound 16g formed two carbon-hydrogen bonds with Gln178 and Leu338, similar to the hydrogen bonds formed by the sulfonamido group in celecoxib, supplementary data ( Figure S46).
Analysis of the binding mode revealed partial superposition of the phenyl rings of 16g over the two phenyl rings in SB2. Moreover, the nitrogen atom of the cyano group in 16g occupied the same position of the pyridinyl nitrogen of SB2 and also formed one conventional hydrogen bond with Met109, Figure 17. Investigation of the binding interaction of compound 16g into COX-2 revealed four conventional hydrogen bonds between the cyano and carbonyl groups with Ser516, Arg106, and Tyr341 amino acids in COX-2, Figure 16. In addition, the methoxy group in compound 16g formed two carbon-hydrogen bonds with Gln178 and Leu338, similar to the hydrogen bonds formed by the sulfonamido group in celecoxib, Supplementary Data ( Figure S46).
Analysis of the binding mode revealed partial superposition of the phenyl rings of 16g over the two phenyl rings in SB2. Moreover, the nitrogen atom of the cyano group in 16g occupied the same position of the pyridinyl nitrogen of SB2 and also formed one conventional hydrogen bond with Met109, Figure 17.
Investigation of the binding interactions of compound 16g into MAP p38 showed four conventional hydrogen bonds with the key amino acids Lys53, Met109, and Leu171, Figure 17 Unlike the co-crystallized ligand SB2, no unfavorable interactions or steric clashes were observed between compound 16g and MAP p38α, Supplementary Data ( Figure S47). Investigation of the binding interactions of compound 16g into MAP p38 showed four conventional hydrogen bonds with the key amino acids Lys53, Met109, and Leu171, Figure 17 Unlike the co-crystallized ligand SB2, no unfavorable interactions or steric clashes were observed between compound 16g and MAP p38α, supplementary data (Figure S47).

Docking Study into Oncogenic Kinases
In addition, a docking study of compound 16g was performed on the four oncogenic kinases (EGFR, CDK2, BRAF, and VEGFR1) identified in the target prediction test, Table  4. The results were also compared with those of the co-crystallized ligands of the four kinases.

Docking into EGFR
Compound 16g was docked into EGFR (pdb: 1M17) [42], and the results revealed a significantly higher binding affinity (ΔGb = −9.52 kcal/mol) compared to the co-crystallized ligand, erlotinib (ΔGb = −7.39 kcal/mol). The higher affinity of compound 16g could be attributed to the higher number of hydrogen bonds and to the electrostatic interaction with Asp831. Analysis of the binding interactions of 16g into EGFR revealed two conventional hydrogen bonds with the key amino acids Cys773 and Asp831 ( Figure 18) compared to one conventional hydrogen bond for erlotinib, supplementary data ( Figure S48).
Investigation of the binding mode of compound 16g revealed partial superposition of the pyrrole ring over the ethynylphenyl moiety of erlotinib, where the two moieties formed similar binding interactions with Thr830 (carbon-hydrogen bond) and Lys721 (pialkyl interaction), Figure 18. The chlorophenyl moiety in 16g also superposed with the phenyl ring of the quinazoline nucleus in erlotinib and extended into the front pockets, forming similar hydrophobic interactions with Leu694, Ala719, and Leu820 like the 7methoxyethoxy moieties of erlotinib.

Docking Study into Oncogenic Kinases
In addition, a docking study of compound 16g was performed on the four oncogenic kinases (EGFR, CDK2, BRAF, and VEGFR1) identified in the target prediction test, Table 4. The results were also compared with those of the co-crystallized ligands of the four kinases.

Docking into EGFR
Compound 16g was docked into EGFR (pdb: 1M17) [42], and the results revealed a significantly higher binding affinity (∆G b = −9.52 kcal/mol) compared to the co-crystallized ligand, erlotinib (∆G b = −7.39 kcal/mol). The higher affinity of compound 16g could be attributed to the higher number of hydrogen bonds and to the electrostatic interaction with Asp831. Analysis of the binding interactions of 16g into EGFR revealed two conventional hydrogen bonds with the key amino acids Cys773 and Asp831 ( Figure 18) compared to one conventional hydrogen bond for erlotinib, Supplementary Data ( Figure S48).

Docking into CDK2
To perform a docking study of compound 16g into CDK2, the pdb: 2VTP [43] was used. The results revealed a significantly higher binding affinity for compound 16g (ΔGb = −10.0 kcal/mol) compared to the co-crystallized ligand, LZ9 (ΔGb = −7.57 kcal/mol). An overlay view of the best ranked pose of compound 16g into CDK2 revealed partial superposition of the pyrrolizine nucleus over the pyrazole ring in LZ9, where the nitrogen of Investigation of the binding mode of compound 16g revealed partial superposition of the pyrrole ring over the ethynylphenyl moiety of erlotinib, where the two moieties formed similar binding interactions with Thr830 (carbon-hydrogen bond) and Lys721 (pialkyl interaction), Figure 18. The chlorophenyl moiety in 16g also superposed with the phenyl ring of the quinazoline nucleus in erlotinib and extended into the front pockets, forming similar hydrophobic interactions with Leu694, Ala719, and Leu820 like the 7methoxyethoxy moieties of erlotinib.

Docking into CDK2
To perform a docking study of compound 16g into CDK2, the pdb: 2VTP [43] was used. The results revealed a significantly higher binding affinity for compound 16g (∆G b = −10.0 kcal/mol) compared to the co-crystallized ligand, LZ9 (∆G b = −7.57 kcal/mol). An overlay view of the best ranked pose of compound 16g into CDK2 revealed partial superposition of the pyrrolizine nucleus over the pyrazole ring in LZ9, where the nitrogen of the cyano group occupied the same position of N 10 of the pyrazole ring. This nitrogen atom also formed one conventional hydrogen bond with Leu83 like LZ9, Figure 19. The higher affinity of compound 16g compared to the co-crystallized ligand (LZ9) is also attributed to the higher amount of hydrogen, supplementary data ( Figure S49). Compound 16g formed three conventional hydrogen bonds with Leu83, Lys129, and Asp145 in CDK2, Figure 19.

Docking into BRAF
Contrary to the docking results for EGFR and CDK2, the results of the docking study of compound 16g for the wild type BRAF (pdb: 4RZV) [44] revealed a lower binding affinity (ΔGb = −11.02 kcal/mol) compared to the co-crystallized ligand, vemurafenib (ΔGb = −12.77 kcal/mol). Investigation of the binding orientation revealed that the pyrrolizine nucleus of compound 16g partially superposed over the pyrrolo [2,3-b]pyridine nucleus of vemurafenib, Figure 20. The higher affinity of compound 16g compared to the co-crystallized ligand (LZ9) is also attributed to the higher amount of hydrogen, Supplementary Data ( Figure S49). Compound 16g formed three conventional hydrogen bonds with Leu83, Lys129, and Asp145 in CDK2, Figure 19.

Docking into BRAF
Contrary to the docking results for EGFR and CDK2, the results of the docking study of compound 16g for the wild type BRAF (pdb: 4RZV) [44] revealed a lower binding affinity (∆G b = −11.02 kcal/mol) compared to the co-crystallized ligand, vemurafenib (∆G b = −12.77 kcal/mol). Investigation of the binding orientation revealed that the pyrrolizine nucleus of compound 16g partially superposed over the pyrrolo [2,3-b]pyridine nucleus of vemurafenib, Figure 20.
In addition, the nitrogen atom of the cyano group in compound 16g was able to form one conventional hydrogen bond with the key amino acid Cys532 like the pyridine nitrogen in vemurafenib. However, 16g displayed only two conventional hydrogens with BRAF compared to six hydrogen bonds for vemurafenib, Supplementary Data ( Figure S50). In addition, the nitrogen atom of the cyano group in compound 16g was able to form one conventional hydrogen bond with the key amino acid Cys532 like the pyridine nitrogen in vemurafenib. However, 16g displayed only two conventional hydrogens with BRAF compared to six hydrogen bonds for vemurafenib, supplementary data ( Figure  S50).

Docking into VEGFR1
The docking study into VEGFR1 (pdb: 3HNG) also revealed a lower binding affinity for compound 16g (ΔGb = −10.67 kcal/mol) compared to -12.06 kcal/mol for the co-crystallized ligand, 8ST. The best fit conformation of compound 16g superposed only partially with the chlorophenyl-carboxamide moiety of 8ST, Figure 21. In addition, compound 16g also showed a fewer number of hydrogen bonds with VEGFR1 than 8ST, supplementary data ( Figure S51). These results were also similar to those obtained from the docking study of compounds 16g into BRAF.

Docking into VEGFR1
The docking study into VEGFR1 (pdb: 3HNG) also revealed a lower binding affinity for compound 16g (∆G b = −10.67 kcal/mol) compared to -12.06 kcal/mol for the co-crystallized ligand, 8ST. The best fit conformation of compound 16g superposed only partially with the chlorophenyl-carboxamide moiety of 8ST, Figure 21. In addition, compound 16g also showed a fewer number of hydrogen bonds with VEGFR1 than 8ST, Supplementary Data ( Figure S51). These results were also similar to those obtained from the docking study of compounds 16g into BRAF. In conclusion, the results of the docking study revealed higher binding affinities for compound 16g toward four (COX-2, MAP P38α, EGFR, and CDK2) of the six targets compared to their co-crystallized ligands. However, additional in silico studies are needed to support these results.

RMSD Analysis and Hydrogen Bond Interaction Estimation
Molecular dynamic simulations (MDS) have proven to have a high value in many computational studies especially in the accurate determination of the binding affinity and the stability of the ligand-protein complexes. Accordingly, six molecular dynamic experi- In conclusion, the results of the docking study revealed higher binding affinities for compound 16g toward four (COX-2, MAP P38α, EGFR, and CDK2) of the six targets compared to their co-crystallized ligands. However, additional in silico studies are needed to support these results.

RMSD Analysis and Hydrogen Bond Interaction Estimation
Molecular dynamic simulations (MDS) have proven to have a high value in many computational studies especially in the accurate determination of the binding affinity and the stability of the ligand-protein complexes. Accordingly, six molecular dynamic experiments were conducted for the synthesized compound 16g in complex with each of the six potential targets as generated from the docking study. Based on the calculated RMSD values for the enzyme Cα atoms as well as compound 16g heavy atoms, the two complexes of CDK2 and EGFR bound to 16g were predicted to be the most stable as depicted from Figure 22. The maximum RMSD values of the two complexes mentioned only reached 1.5 and 1.7 Å for CDK2 and EGFR, respectively. In contrast, the maximum RMSD values of the other four targets in complex with 16g reached 2.6, 2.7, 3 and 3.2 Å for COX-2, MAP p38α, BRAF and VEGFR1, respectively. Furthermore, the stability of the hydrogen bond interactions between the lead compound 16g, and the six potential targets were monitored through the entire MDS. Generally, a hydrogen bond is considered stable and valid when the distance between the hydrogen bond donor and acceptor is kept less than 3.5 Å. This criterion was maintained only in all the formed hydrogen bonds of 16g-CDK2 and 16g-EGFR complexes. Some of the hydrogen bonds in the other four complexes were significantly unstable Table 5. The MDS results highlighted a high potentiality for compound 16g to inhibit both CDK2 and EGFR.  The maximum RMSD values of the two complexes mentioned only reached 1.5 and 1.7 Å for CDK2 and EGFR, respectively. In contrast, the maximum RMSD values of the other four targets in complex with 16g reached 2.6, 2.7, 3 and 3.2 Å for COX-2, MAP p38α, BRAF and VEGFR1, respectively. Furthermore, the stability of the hydrogen bond interactions between the lead compound 16g, and the six potential targets were monitored through the entire MDS. Generally, a hydrogen bond is considered stable and valid when the distance between the hydrogen bond donor and acceptor is kept less than 3.5 Å. This criterion was maintained only in all the formed hydrogen bonds of 16g-CDK2 and 16g-EGFR complexes. Some of the hydrogen bonds in the other four complexes were significantly unstable Table 5. The MDS results highlighted a high potentiality for compound 16g to inhibit both CDK2 and EGFR.

MM-PBSA Calculations
The binding free energy that resulted from the binding of 16g to each of the potential targets were calculated using the MM-PBSA approach. In this approach, the calculations are based on all trajectories extracted from the MDS, in contrast to the docking score that is based on a single conformation. Accordingly, this technique is considered as a more reliable indicator compared to the energy score obtained from the docking studies. As part of these calculations, the free energy of each component (receptor, ligand, and complex) was calculated by summing its free energy of solvation and molecular mechanics' potential energy in vacuum. The free energy of solvation includes both the polar solvation energy and nonpolar solvation energy (non-electrostatic, calculated by the solvent accessible surface area; SASA model). Finally, the free energy of binding was calculated by subtracting the energy of the receptor and ligand from the energy of the complex. The calculated types of energies and binding free energy values for the six complexes are summarized in Table 6. The results indicated a higher stability for the 16g-CDK2 and 16g-EGFR complexes than the other four complexes as revealed by the high negative binding free energy of the two complexes. The average binding free energy of the 16g-CDK2 and 16g-EGFR complexes were −401 and 387 KJ/mol, respectively, which suggests a strong and stable binding for 16g with the two targets. On the other hand, 16g achieved average binding free energy of −360, −354, −337, and −328 KJ/mol with COX-2, p38α, BRAF, and VEGFR-1. We believe that the results from the MD simulations validated our design and supported our hypothesis for CDK2 and EGFR as potential targets for 16g. To calculate physicochemical properties related to drug-likeness, several web-based tools are available for free. Of these tools, SwissADME (http://www.swissadme.ch) [47] was used to calculate physicochemical properties of compounds 10, and 16a-h, Table 7. Meanwhile, the calculation of molecular volumes and drug-likeness score (DLS) of these compounds was done using Molsoft (http://molsoft.com/mprop). The molecular weight of all the new compounds 16a-h was in the range of 386.42-427.50 daltons, which fall within the limits of Lipinski's rule. They also showed comparable or slightly higher molecular volumes compared to compound 10. In addition, the number of the hydrogen bond acceptors (H A )/donors (H D ) were also within the limits of Lipinski's rule. Accordingly, no violations from Lipinski's rule was observed for any of the new compounds 16a-h.
All of the new compounds exhibited similar bioavailability scores with 81.60-84.79% oral absorption compared to compound 10 (84.79%). The new compounds also exhibited drug-likeness score in the range of 0.25-0.95 compared to 0.80 for compound 10. Among the new compounds, compound 16g exhibited the highest DLS.

Metabolic Study
In the current study, Biotransformer (http://biotransformer.ca) [48] was used to predict the metabolic pathways and the expected metabolites of compound 16g. The test compound was submitted to the server. Thereafter, phase I metabolic transformation was selected. The output results in the form of metabolic pathways were obtained describing the expected metabolites and the transforming phase I enzymes which could perform this action in humans. The results were collected in one figure, Figure 23.
The results of the metabolic study of compounds 16g revealed five potential metabolic pathways and eight expected metabolites. The metabolic pathways included epoxidation of the aromatic rings, aromatic hydroxylation, O-demethylation of the methoxy group, Nhydroxylation of the carboxamide nitrogen, and oxidation of one of the secondary carbons of the pyrrolidine ring, Figure 23.
The prediction of phase II metabolites of compound 16g did not show any results. This was due to the absence of the polar functional groups (OH and COOH), which can undergo conjugation metabolism. However, most of the expected phase I metabolites of compound 16g (Figure 23) have alcoholic/phenolic OH groups that can undergo glucuronidation of sulfate conjugation. dict the metabolic pathways and the expected metabolites of compound 16g. The test compound was submitted to the server. Thereafter, phase I metabolic transformation was selected. The output results in the form of metabolic pathways were obtained describing the expected metabolites and the transforming phase I enzymes which could perform this action in humans. The results were collected in one figure, Figure 23.

Pharmacophore Search
The pharmacophore search was done using Pharmit [24]. The compound library was initially uploaded to Pharmit as a compressed file in "sdf.gz." format. The study was started by selecting the crystal structure of COX-2 (pdb: 1CX2) with the co-crystallized ligand (SC-558). The displayed pharmacophore features were edited to include the pharmacophore features in Figure 4. The settings including hits reduction and screening were set to the default values. The search MolPort was used to select the uploaded compound library. The search type was set to pharmacophore search. The results of the pharmacophore search appeared in a tabular form including codes of the hits, RMSD, MW, and number of RBs.

Chemistry
Chemical reagents and solvents were purchased from Sigma Aldrich (Darmstadt, Germany). Melting points (uncorrected) of compounds 16a-h were determined by an IA 9100MK-Digital melting point apparatus (Cole-Parmer, Vernon Hills, IL, USA). Absorption bands in the infrared (IR) spectra were recorded on a RUKER TENSOR 37 FTIR spectrophotometer. Molecular ions and mass spectra (MS) of compounds 16a-h were analyzed using the Shimadzu Qp-2010 Plus mass spectrometer (EI ionization mode). The elemental analyses (C, H, and N) of compounds 16a-h were measured in a Microanalytical Center, Cairo University. 1 H-NMR, 13 C-NMR, and DEPT C 135 spectra were recorded using BRUKER AVANCE III at 500 MHz, 125 and 125 MHz, respectively.
Copies of the spectra data including IR, mass, 1 H-NMR, 13 C-NMR, and DEPT C 135 spectra of compounds 16a-h are provided in Supplementary Data (Figures S3-S45).
The readers are advised to consider that the two phenyl rings are assigned as phenyl ring A and B (Figure 9).

General Procedure (A) for Preparation of Compounds (16a-h)
The new compounds were prepared from compound 15a,b according to the reported procedures [20]. A mixture of the pyrrolizine 15a,b (3.4 mmol), the appropriate aldehyde (4.4 mmol), and glacial acetic acid (0.5 mL) in absolute ethanol (30 mL) was refluxed for 4-6 h. The reaction mixture was then concentrated and set aside to cool. The solid product obtained was recrystallized from acetone-chloroform.

Antiproliferative Activity Cell Culture
In the current study, the three cancer cell lines (MCF7, A2780, and HT29) used were obtained from the ATCC. The cell lines were cultured following our previous report [20]. On the other hand, the normal MRC5 cells used in the current study were maintained in Eagle's minimum essential medium following the previous report [33].

Antiproliferative Activity Assay
The antiproliferative activity was measured using the MTT assay according to the previous report [35]. Briefly, the cancer/normal cell lines were cultured in 96-well (3 × 10 3 /well) separately. The cells were treated by tested compounds 16a-h in final concentrations of 0, 0.1, 1, 10, 25, and 50 µM and incubated at 37 • C for 72 h. The MTT was added and the antiproliferative activities (IC 50 values) of the new compounds were calculated.

Cell Cycle Analysis
The effect of compound 16g on the cell cycle distribution of MCF7 was performed following the previous report [36]. Briefly, the cancer cells were cultured for 72 h with 0.00, 0.01, 0.05, and 0.10 µM final concentration of compound 16g. Next, the cells were washed with PBS x1 and trypsinized. The collected cells were spinned and fixed with 70% ethanol. Ribonuclease A was added (15 min) after suspending the cells in cold PBS x1. Following the addition of propidium iodide, analysis of the ice-cold cells was done by flow cytometry (BC, FC500, Brea, CA, USA).

Annexin V FITC/PI Assay
MCF7 cells were treated by compound 16g at 0.00, 0.01, 0.05 and 0.10 µM final concentrations for 72 h. Briefly, the superanent of the cells was collected in ice-cold tubes, trypsinized and incubated (37 • C). The procedures were completed following the previous report [37]. The samples were analyzed using flow cytometry (BC, FC500, Brea, CA, USA).

Target Prediction
The SwissTargetPrediction (http://www.swisstargetprediction.ch) [39] was used to predict the molecular targets of the new compounds 16a-h. The test compounds were submitted one by one to the server. After running the prediction process, the results were obtained in the form of a pie chart for each compound including the major classes of the potential targets ranging with their percent. In addition, a detailed report for the first 100 entry of these targets was also obtained.

RMSD Analysis and Hydrogen Bond Interaction Estimation
In this section, GROMACS 5.1 software was used to conduct all the molecular dynamics simulations [51]. Six MDS experiments were conducted on the six complexes retrieved from the docking step; each complex contains 16g bound to a different target (CDK2, COX2, EGFR, VEGFR1, BRAF, or MAP p38α). The Automated Topology Builder (ATB) and Repository version 3 [52] were implemented to generate a topology file for 16g under the GROMOS96 force field. The generated ligand topology was joined with each of the six enzymes' topology using the standard published protocol [53]. Solvation of the six complexes was done using single point charge (SPC) water model. A proper number of ions was added to the processed systems using the gmx genion script. The neutralized solvated systems were energy minimized using the steepest descent minimization algorithm with a maximum of 50,000 steps and <10.0 kJ/mol under GROMOS96 43a1 force field [54]. After that, two equilibration ensembles were conducted to ensure proper equilibrations for all the processed systems. At the beginning, an NVT ensemble with a constant number of particles, volume, and temperature (310 K) was done for 1 ns then followed by an NPT ensemble with a constant number of particles, pressure, and temperature for 4 ns. The Particle Mesh Ewald (PME) method with a 12 Å cut-off and 12 Å Fourier spacing were used to get the long range electrostatic [53]. The six equilibrated systems entered the production stage without any restraints for 50 ns with a time step of 2 fs, and the structural coordinates were saved every 5 ps. Both the temperature (310k) and the pressure (1atm) were regulated throughout the simulation V-rescale weak coupling method (modified Berendsen thermostat) and the Parrinello-Rahman method [55,56]. The RMSD of the whole system was calculated from the generated trajectories from the production step as well as the distances of the formed HBs.

MM-PBSA Calculation
Binding free energy calculations were performed using the MM-PBSA which applies the following equation: G (Complex) is the total free energy of the ligand-protein complex. G (Receptor) and G (Ligand) are total free energies of the isolated protein and ligand in solvent, respectively. The total free energy was calculated for all MD trajectories from its molecular mechanics potential energy plus the energy of solvation using the g_mmpbsa package implemented in GROMACS software [57]. Individual energies and the corresponding SD were calculated and then summed together to yield the average total free energy of each.

Physicochemical Properties and Drug-Likeness
The physicochemical properties of compounds 16a-h were calculated by the Swis-sADME webserver (http://www.swissadme.ch/) [47]. Each compound was submitted to the server followed by running the calculations. On the other hand, the Molsoft webserver (http://molsoft.com/mprop/) was used in the calculation of drug-likeness scores (DLSs) of the final compounds.

Metabolic Study
Biotransformer (http://biotransformer.ca) [48] was used to predict the metabolic pathways and the metabolites of compound 16g. The test compound was submitted in the server. Thereafter, the task of metabolic transformation was selected (phase I/II), and the number of reaction steps was set to one. The output results in the form of metabolic pathways/expected metabolites were obtained.

Conclusions
In the current study, an in silico approach based on free software was used to optimize the antiproliferative activity and investigate the potential mechanism of action of a series of pyrrolizine-based Schiff bases. A compound library of 288 Schiff base derivatives was designed based on compound 10. A pharmacophore of the compound library search was performed. Structural analysis of the top-scoring hits was conducted, and a preliminary docking study into COX-2 was performed to select the promising hits for the synthesis. The chemical synthesis and structural elucidation of the new compounds 16a-h were discussed. The MTT assay was used to evaluate the antiproliferative activity of compounds 16a-h against MCF-7, A2780, and HT29 cancer lines (IC 50 = 0.01-40.50 µM). Amongst the new compounds, compound 16g exhibited the highest antiproliferative activity against MCF7 cells (IC 50 = 0.01 µM). To assess the toxicity and selectivity of the new compounds, their growth inhibitory activity was evaluated against normal MRC5 cells (IC 50 = 1.27-24.06 µM). Compound 16c showed the highest selectivity index against MCF7 and A2780 cells, while 16g was the most selective for the HT29 cell line. To investigate the potential mechanism of action of the new compounds, 16g, the most active in the MTT assay was evaluated for its effects on the cell cycle distribution of MCF7 cells. The results revealed cell cycle arrest at the G 2 /M phase. Compound 16g also induced a dose-dependent increase in the apoptotic events in MCF7 cells compared to the control (7.9-17.8%). SwissTargetPrediction was used to predict the potential molecular targets which could mediate the anticancer potential of 16g. The results revealed six potential targets including COX-2, MAP P38α, EGFR, CDK2, BRAF, and VEGFR1. A comparative molecular docking study of compound 16g was performed into the six targets, where the results revealed high binding affinities for 16g toward four of these targets (COX-2, MAP P38α, EGFR, CDK2). The molecular dynamic simulation also revealed favored stability and binding energy for 16g in complex with CDK2 and EGFR, while COX-2 was in the third order. These findings suggested that compound 16g could serve as a potential anticancer agent.
Supplementary Materials: Supplementary material including spectral data (IR, Mass, 1 H-NMR, 13 C-NMR, and DEPT C 135 spectra) of compounds 16a-h are provided with this manuscript (Figures S1-S61).  Data Availability Statement: The supplementary materials were provided with the manuscript "were uploaded with the first and revised submission". You can use the supplementary material file uploaded to the system with the revised manuscript.