Imidazothiazole-based potent inhibitors of V600E-B-RAF kinase with promising anti-melanoma activity: biological and computational studies

Abstract A series of imidazothiazole derivatives possessing potential activity against melanoma cells were investigated for molecular mechanism of action. The target compounds were tested against V600E-B-RAF and RAF1 kinases. Compound 1zb is the most potent against both kinases with IC50 values 0.978 and 8.2 nM, respectively. It showed relative selectivity against V600E mutant B-RAF kinase. Compound 1zb was also tested against four melanoma cell lines and exerted superior potency (IC50 0.18-0.59 µM) compared to the reference standard drug, sorafenib (IC50 1.95-5.45 µM). Compound 1zb demonstrated also prominent selectivity towards melanoma cells than normal skin cells. It was further tested in whole-cell kinase assay and showed in-cell V600E-B-RAF kinase inhibition with IC50 of 0.19 µM. Compound 1zb induces apoptosis not necrosis in the most sensitive melanoma cell line, UACC-62. Furthermore, molecular dynamic and 3D-QSAR studies were done to investigate the binding mode and understand the pharmacophoric features of this series of compounds.


Introduction
Melanoma is malignant tumour of melanocytes that are present in the skin and the eye. People with fair skin are more prone to melanoma than dark-coloured skin individuals. In addition, blue, green, or hazel-colored eyes are at higher risk to melanoma than dark-coloured eyes 1 . Excessive exposure to sun light and family history of melanoma are also among the risk factors of melanoma. The World Health Organisation (WHO) indicates about 132,000 new melanoma skin cancer cases are diagnosed annually worldwide. Moreover, the American Cancer Society estimates diagnosis of more than 100,000 melanoma cases in USA only in 2020 2 . Much research efforts are required to develop more efficient and less toxic anti-melanoma drugs.
RAS-RAF-MEK-ERK pathway is among the potential targets for melanoma treatment. V600E mutated B-RAF is over-activated and over-expressed in several melanoma patients 3,4 . Inhibitors of V600E-B-RAF and RAF1 (C-RAF) kinases have been reported as antiproliferative agents against melanoma. Among which, vemurafenib 5 , encorafenib 6 , and dabrafenib 7,8 (Figure 1) have been marketed as RAF-inhibitory anti-melanoma therapies. Several clinical trials are currently running for these drugs either alone or in combination with other drugs. In addition, sorafenib ( Figure 1) is one of the potent inhibitors of RAF kinases that was utilised as a reference standard in our cell-based assay 9 . Moreover, several imidazothiazole derivatives were reported as V600E-B-RAF or RAF1 kinase inhibitors [10][11][12] .
We have previously reported a series imidazothiazole derivatives as potential antiproliferative agents against melanoma [13][14][15] . In our previous reports, we did not investigate their molecular mechanism of action. In the present study, we tested their molecular mechanism of action against kinases in cell-free and whole-cell assays. The most potent compound was further tested against melanoma and normal skin cells to investigate the selectivity indexes. The ability of the most potent compound to induce apoptosis and necrosis were also tested. Moreover, different computational studies such as molecular dynamic simulation and QSAR were carried out.

Synthesis of the target compounds
The methods are reported in our previously published articles [13][14][15] .

Cell-free kinase profiling
In a final volume of 25 mL, the tested kinase (5-10 mU) was incubated with 25 mM Tris pH 7.5, 0.02 mM EDTA, 0.66 mg/mL myelin basic protein, 10 mM magnesium acetate and [c 33 P-ATP] (specific activity approx. 500 cpm/pmol, concentration as required). The reaction was initiated by the addition of the Mg-ATP mix. After incubation for 40 min at room temperature, the reaction was stopped by the addition of 5 mL of a 3% phosphoric acid solution. 10 mL of the reaction was then spotted onto a P30 filtermat and washed three times for 5 min in 75 mM phosphoric acid and once in methanol before drying and scintillation counting.

In-cell kinase assay
In the cellular V600E-B-RAF phosphorylation assay the murine embryonal fibroblast cell line (MEF) was used, which expresses exogenously introduced constitutively active human V600E-B-RAF fused to the N-terminus of a modified oestrogen receptor (ER). In this construct, V600E-B-RAF is only active in the presence of the oestrogen-analogue hydroxytamoxifen. hydroxytamoxifen treatment results in B-RAF-VE600E kinase activation and subsequent phosphorylation of the direct V600E-B-RAF substrate MEK1 at Ser218/222. MEF-V600E-B-RAF:ER cells were plated in Dulbecco's modified eagle medium (DMEM) supplemented with 10% FCS in multiwell cell culture plates. Compound incubation (90 min at 37 C) was done in serum-free medium. Cells treated with 10 mM sorafenib were used as low control. For stimulation, cells were treated with 1 mM hydroxytamoxifen for 1 h at room temperature. After cell lysis, quantification of MEK1 phosphorylation was assessed in 96-well plates via sandwich-ELISA using a MEK1 specific capture antibody and a phospho-MEK1-Ser218/222 specific detection antibody.

MTT assay
The melanoma and normal skin cell lines were maintained in DMEM including 10% foetal bovine serum (FBS) and 1% penicillin/ streptomycin in a humid atmosphere containing 5% carbon dioxide at 37 C. The cells were taken from culture substrate with 0.05% trypsin-0.02% EDTA and plated at a density of 5 Â 10 3 cells/ well in 96-well plates followed by incubation at 37 C for 24 h in a humid atmosphere with 5% carbon dioxide CO 2 before treatment with various concentrations (3- conventional 3-(4,5-dimethylthiazol-2-yl)-2,5-diphenyltetrazolium bromide (MTT) reduction assay. MTT assays were carried out with CellTiter 96V R (Promega) according to the manufacturer's instructions. The absorbance at 590 nm was recorded using EnVision 2103 (Perkin Elmer; Boston, MA, USA). The IC 50 values were calculated using GraphPad Prism 4.0 software. Triplicate testing was performed.

Caspase-3/7 and LDH release assays
They were performed following the protocols reported in our previously published article 16  2.6.2. Docking and molecular dynamic simulation As a prerequisite for molecular dynamic simulation, docking experiments were carried-out using the program AutoDock Vina 19 . The x-ray crystal structures of V600E-B-RAF kinase (PDB-ID: 3IDP) and RAF1 kinase (PDB-ID: 3OMV) were retrieved from Protein Data Bank (http://www.rcsb.org/pdb/). Bound inhibitors and crystalline water molecules were extracted from the initial structures. Autodock MGL Tools were used for the addition of polar hydrogens and charges. Compound 1zb was treated using the same procedure. Spacing of 1.0 Å between the grid points was used to establish grid boxes covering the active site of the studied macromolecules, centred towards the coordinates of 5.29 (x), 24.13 (y), 32.82 (z) for V600E-B-RAF and towards the coordinates of 28.61 (x), 39.19 (y), 39.09 (z) for RAF1 kinase. Exhaustiveness was set to 12, while number of poses was set to 10.
Molecular dynamic simulations for the compound 1zb with respect to V600E-B-RAF and RAF1 were started from the earlier docked complexes. Desmond software V R (Desmond Molecular Dynamics System, version3.8, D. E. Shaw Research, New York, NY, 2014) embedded within Maestro interface (Schr€ odinger Release 2014-2: Maestro, version 9.8, Schr€ odinger, LLC, New York, NY, 2014) was used to conduct all-atoms molecular dynamic simulations employing OPLS_2005 force field parameters. Each complex was subjected to the same dynamic protocol, in which; TIP3P explicit water molecules as solvent model within an orthorhombic periodic boundary box of the size 10 Å 3 were used to solvate the protein-ligand complex, then, system neutralisation was accomplished by adding appropriate counter-ions followed by adding 0.15 M of salt ions. Prior to the actual dynamic run, system relaxation was achieved by performing a series of short (2000 iterations) restrained and non-restrained solute minimizations steps followed by short 12 ps simulation steps using NVT and NPT ensembles. 50 ns production run was carried out using the NPT ensemble class integrating the equation of motion every 2 fs and setting the temperature and pressure to 300 K and 1 atmosphere, respectively. Short-range interactions cut-off was set to 9 Å and the long-range electrostatic interactions were calculated employing the particle mesh Ewald (PME) method. Results were visualised within maestro environment and analysed using Desmond interaction diagram panel.

FlapV R 3D-QSAR model construction
The 3D-QSAR models for the RAF1 and V600E-B-RAF active compounds were constructed employing the program Flap V R (fingerprint for ligands and proteins) 20,21 . The program calculates fingerprints, which are derived from the GRID Molecular Interaction Fields (MIFs) and are characterised as quadruplets of pharmacophoric features. Model construction was initiated by building a compound database for each enzyme. For the RAF1 enzyme, 9 active compounds with their available pIC 50 values were used for this purpose, whereas 22 active compounds with their corresponding pIC 50 values ( Table 6) were used for the V600E-B-RAF database creation. Next, the protonation states for each molecule at pH 7.4 were calculated using an integrated FlapV R tool called MoKa V R 22 . Subsequently, for each molecule, a maximum of 50 conformers were generated with an RMSD value of 0.3 Å. Ligand-based alignment was achieved employing the most active compound 1zb in its low energy conformation as a template. Later, for each database, a partial least square (PLS) model was constructed employing the GRID-probes; H (shape), DRY (hydrophobic), O (H-bond acceptor) and N1 (H-bond donor). Model validation was accomplished via leave-one-out cross validation and the optimal latent-variables for each model were determined by investigating the R 2 vs Q 2 plot.

Volsurfþ V R QSAR model construction
Two additional QSAR models were constructed; one for the RAF1 enzyme, and another one for the V600E-B-RAF. In each case, the molecular descriptors for the studied compounds were calculated employing the software VolSurfþ V R 23,24 . The software calculates a set of 1D-3D molecular descriptors (ca. 128) of different classes including molecular size/shape, hydrophilic/hydrophobic regions quantification, INTEraction enerGY (INTEGY moments), capacity factors, amphiphilic moments, hydro-lipo balance, molecular diffusivity, LogP, LogD, pH-dependent solubility, molecular flexibility, 3D pharmacophoric descriptors, etc. During model development, the calculated descriptors were recruited as independent variables, whereas the bioactivity values (pIC 50 ) for the compounds under study were recruited as the dependent variable. Development and validation of our QSAR models were achieved employing the QSARINS V R software 25,26 . All subsets were investigated for the first two descriptors selection. Then, optimal combinations of descriptors (greater than two) were reached by means of genetic algorithm (GA). In the GA selection phase, the population size, maximum number of generations and mutation rate were set to 800, 3000 and 0.6, respectively. Multiple linear regression (GA-MLR) method was adapted for the final model building. The robustness for each model was assessed employing the Q 2 LOO (leave one-out) cross validation procedure.

Chemistry
The target compounds 1a-zh were prepared utilising the pathways illustrated in Schemes 1-3 [13][14][15] . The synthetic strategy involved preparation of methyl sulphonyl intermediates 6a,b (Scheme 1) and amine side chain intermediates 13a-c and 15a-f (Scheme 2), followed by interaction of them together to get the aminopyrimidinyl final compounds (Scheme 3). 5a,b. Oxidation of the methylthio group of compounds 5a,b using oxone yielded the corresponding methyl sulphonyl derivatives 6a,b. Scheme 2 illustrates the conversion of 2-aminoethanol into the amino reagents possessing arylamide or arylurea moiety. At the beginning, the amino group of 2-aminoethanol was Cbz-protected using benzyloxycarbonyl chloride in presence of triethylamine (TEA) 28 . The hydroxyl group of compound 7 was converted into methanesulfonate using methanesulfonyl chloride and TEA 29 . Compound 8 was heated with sodium azide to replace the OMs group with azido (compound 9). Reduction of azido compound 9 into the corresponding amino was accomplished using triphenylphosphine in refluxing methanol. The free primary amino group of compound 10 was treated with aryl isocyanates or condensed with carboxylic acids to form the corresponding urea or amide, respectively. The last step of Scheme 2 involves deprotection of Cbz group using hydrogen gas and palladium metal over charcoal.
In Scheme 3, the mesyl intermediates 6a,b were heated with the amino side chain reagents 12a-c and 14a-f in presence of Hunig's base to produce the target compounds possessing fluoro or methoxy on the benzene ring at position 6 of the imidazothiazole nucleus. Boron tribromide was used to demethylate the methoxy compounds 1o-q and 1 u-za to form the corresponding phenolic analogues 1r-t and 1zb-zh, respectively. The structures of the target compounds 1a-zh are drawn in Table 1.

Kinase profiling
Based on our previous knowledge of inhibitory effect of imidazothiazole derivatives against V600E-B-RAF and RAF1 11,12 , we decided to start our biological investigation with testing all the 34 final compounds against both kinases. The target compounds were tested at 1 mM concentration, and the compounds that produced promising inhibitory effect, more than 75% inhibition, were further tested in 10-dose mode to measure their IC 50 values. The inhibition percentage and IC 50 values are demonstrated in Table  1. The target compounds are generally more active against V600E-B-RAF than RAF1. A comparative computational study was carried out to investigate it.
Meta-hydroxyl group on the phenyl ring at position 6 of the imidazothiazole nucleus is more favourable for activity than metamethoxy and para-fluoro substituents against both V600E-B-RAF and RAF1. This is emphasised by the activity of OH compounds 1zb-zh compared with OMe analogues 1 u-za or F derivatives 1e and 1j-n. It is obvious that hydrogen bond donor at the this part of the structure is optimal for stronger affinity with V600E-B-RAF and RAF1 kinases. This was further investigated by computational studies. In addition, amide compounds (1e, 1 u, 1x, 1zb, 1zd, and 1ze) are more active and more potent than the corresponding urea analogues (1a, 1o, and 1q-t). It can be concluded that the length, hydrophobicity, and electronic properties of amide spacer are more optimal for activity than those of urea. Furthermore, terminal benzamido moiety with no substitution (e.g. compounds 1 u and 1zb) or ortho-hydroxyl group (e.g. compounds 1v and 1zc) are more favourable for activity than analogues with bulkier substituents. QSAR study was carried out to further investigate this issue. Among all the final compounds, compound 1zb possessing 3-hydroxyphenyl at position 6 of imidazothiazole and unsubstituted benzamido is the most potent against both V600E-B-RAF and RAF1. It is the only compound that exerted sub-nanomolar IC 50 value against V600E-B-RAF.
The most potent compound, 1zb, against V600E-B-RAF and RAF1 was further tested at 1 mM concentration against a 30-kinase panel belonging to different kinase families. The results are illustrated in Table 2. Apart from V600E-B-RAF, it inhibited 73.14% of FLT3 kinase while its inhibitory effect against all the other kinases was less than 50%. The IC 50 of compound 1zb against FLT3 was further examined and found to be 838 nM (Table 3). So the compound is 856.85 times more selective towards V600E-B-RAF than FLT3. It is also 8.38 folds more selective against V600E-B-RAF than RAF1. Based on these results, compound 1zb can be considered as a relatively selective V600E-B-RAF kinase inhibitor.

Antiproliferative activity against melanoma cell lines
Compound 1zb was tested against four melanoma cell lines and normal skin epithelial cell line. Its potency was compared with sorafenib as a reference standard as it is a potent RAF kinase inhibitor. The results are illustrated in Table 4. Compound 1zb is more potent than sorafenib against all the four tested melanoma cell lines. Its IC 50 values are within sub-micromolar range, while those of sorafenib are in one-digit micromolar scale. The selectivity indexes of compound 1zb against A375, MDA-MB-435, SK-MEL-28, and UACC-62 melanoma cell lines compared to CCD1106K normal skin cells are 20.13, 15.69, 24.37, and 51.44, respectively. Furthermore, the strong potency of compound 1zb against the tested melanoma cell lines that harbour V600E-B-RAF 30 emphasises the postulation that V600E-B-RAF could be at least a part of molecular mechanism(s) of antiproliferative action of compound 1zb.

In-cell kinase assay
The next step of our study was investigation of ability of compound 1zb to inhibit V600E-B-RAF kinase inside the cells after penetrating the cell membrane. It was compared with sorafenib as a reference standard. This assay was done in murine embryonal fibroblast (MEF) cell line and the ability to inhibit V600E-B-RAF was investigated through measuring the concentration of phosphorylated MEK, the downstream protein of V600E-B-RAF. The dose-response curves are illustrated in Figure 2. Compound 1zb could cross the cell membrane and inhibit V600E-B-RAF with IC 50 value of 0.19 mM. This concentration is comparable to its IC 50 value against UACC-62 melanoma cell line (0.18 mM) that has overexpressed V600E-B-RAF 30 . This result can also prove that V600E-B-RAF inhibition could be a molecular mechanism of action against melanoma.

Apoptosis and necrosis assays
As per the results shown in Table 4, UACC-62 is the most sensitive melanoma cell line for compound 1zb. Its IC 50 value is 0.18 mM. We decided to investigate whether this concentration could induce apoptosis and/or necrosis in UACC-62 cells. Caspase-3/7 were measured as proof of apoptosis induction while lactate dehydrogenase (LDH) release was measured as an indicator of necrosis 16 . At 0.18 mM concentration, compound 1zb increased caspase-3/7 release by 53.33% but did not show any significant increase in LDH release (Table 5). This means that compound 1zb at its IC 50 concentration can induce apoptosis of UACC-62 melanoma cells but not necrosis. This finding could help us understand its antiproliferative effect in more details.

Molecular dynamic study
The structure-activity relationships presented in Table 1 were further explored via computational techniques in an attempt to provide plausible explanation for the bioactivity vibrations existed among our compound series against both enzymes (RAF1 and V600E-B-RAF). Thus, we started our endeavour by performing molecular docking procedure for the most active compound 1zb against the two enzymes. Then, best docked-complexes achieved were utilised to carry-on a 50 nanosecond (ns) molecular dynamic simulations under explicit hydration environment. In each case, system stability and trajectories were assessed in term of rootmean-squared deviations (RMSD) of protein Ca-atoms as well as the ligand heavy atoms. Figure 3(A,B) illustrate the RMSD values (compared to the starting frame-0) along the simulation time, where changes of the order of 1-3 Å are perfectly acceptable for small, globular proteins. Accordingly, in the case of RAF1 enzymeligand complex system equilibration was achieved after 5 ns, while for the V600E-B-RAF; system equilibration was achieved after 18 ns and remained stable during the rest of the simulation time.
Interactions of compound 1zb and the RAF1 kinase revealed a number of important findings. The 3-hydroxyphenyl group (OH) was able to establish hydrogen bonding with the Cys-424 residue for about 65% of the simulation time (direct and water bridge). Additional water bridge interaction was observed between the benzamide carbonyl oxygen (C¼O) and the corresponding gatekeeper residue Thr-421 for 58% of the simulation time ( Figure  4(A,B)). Furthermore, an interesting intra-molecular hydrogen bonding was observed between the benzamide nitrogen and the pyrimidine nitrogen last for 41% of time, which is indicative for conformational flexibility restriction. Other significant interaction were the hydrophobic interactions with the residues Phe-475, Leu-406, and Ile-355. On the other hand, an analysis of compound 1zb within the V600E-B-RAF kinase revealed another set of interactions. The 3hydroxyphenyl group (OH) was able to establish two hydrogen bonds (direct and water-bridged); one with the hinge Cys-532 residue and another one with the Gln-530 residues for about 83% and 65% of the simulation time, respectively. Additional direct Hbond interaction was observed between the benzamide carbonyl oxygen (C¼O) and the corresponding residue Ser-536 for about 45% of the simulation time ( Figure 5(A,B)). The residue Gly-534 was able to form H-bond with the amino group (NH) next to pyrimidine ring for about 49% of the simulation time. Other hydrophobic interactions were those with the residues Phe-583 and Phe-595.

Flap V R 3D-QSAR study
At the next level of our investigation, we were interested in gaining more detailed information about the crucial factors affecting selectivity in each enzyme case, and the contribution of structural variations to the observed potency, thus we decided to conduct a 3D-QSAR analysis employing the program Flap V R 20,21 . The program frequently used for building QSAR models based on molecular interaction fields (MIFs) generated through the GRID force-field 31 . Four different GRID-probes were utilised for calculating the final interaction fields including; H-probe for computing the shape constraints on small ligand or protein cavity, DRY-probe for computing favourable hydrophobic interactions, O-probe (carbonyl oxygen) acting as hydrogen bond acceptor (HBA) and N1-probe (amidic NH) acting as hydrogen bond donor (HBD) probe.
In the case of RAF1 compounds, the best developed model was a two-Latent variables model with an R 2 of 0.799 ( Figure 6(A))    Table 5. Caspase-3/7 and LDH release assay results of compound 1zb at 0.18 mM concentration over UACC-62 melanoma cells.

Assay
Result (compared to control cells treated with DMSO) Caspase-3/7 53.33% increase LDH release 0% increase consisted of three interaction fields encoding shape, hydrophobic and HBA regions. Alignment of the most potent compound (1zb, IC 50 ¼ 8.2 nM) to the model (Figure 7(A)) revealed a good match between its structural components and the given fields, where the amino group (NH) next to pyrimidine ring is matching the HBA field (red) and the benzamide ring is matching well with the hydrophobic field (yellow) and there is no shape violation observed at the benzamide ring. On the other hand, the least potent compounds 1ze, 1zg and 1zf exhibited the same structural components match except for the shape field. The trifluoromethyl groups on the benzamide ring of compound 1ze clearly violating the shape constrains (Figure 7(B)). The same phenomenon was observed for the 4-morpholino substitution of 1zg (Figure 7(C)) and the 5-methylimidazole ring substitution of 1zf (Figure 7(D)). Accordingly, it seems that any benzamide substitution violating of the shape constraints would have negative impact on the RAF1 bioactivity and only small compacted substituents could be tolerated.
In the case of V600E-B-RAF compounds, the best developed model was a three-Latent variables model with an R 2 of 0.857 ( Figure 6(B)) consisted of three interaction fields encoding shape, hydrophobic and HBA regions. Again, alignment of the most potent compound (1zb, IC 50 ¼ 0.978 nM) to the model ( Figure  8(A)) revealed a good match between its structural components and the given fields, where the 3-hydroxyphenyl group (OH) is matching the HBA (red) and the shape (blue mesh) fields. The benzamide ring is matching well with the hydrophobic field (yellow) and there is no shape violation observed at the benzamide ring. The least potent compound (1 b, IC 50 ¼ 245 nM) suffers from many drawbacks as observed in Figure 8(B). The 4-fluorophenyl   group deviates from the perfect shape match and provide no complementarity with the HBA field. Additionally, the presence of urea linker appeared inferior leading to misalignment of the phenylurea ring from the perfect hydrophobic ring match. In compound 1z (Figure 8(C)), the presence of 3-methoxyphenyl (HBA) in front of HBA field and the shape constraints violation by the benzamide ring substitutions (3-trifluoromethyl and 4-morpholino) appeared inferior to the bioactivity. The lower potency observed for compound 1t (Figure 8(D)) could be attributed to the presence of urea linker and the phenylurea ring misalignment. In light of this information, it seems that the presence of small HBD group such as (OH) at meta-position of the phenyl ring is highly recommended, and the amide linker is more preferred over the urea linker. Furthermore, small compacted substituents at the ortho and meta positions of the benzamide ring are of more positive impact for the bioactivity.
3.6.3. Volsurfþ V R QSAR study Later, we turned our focus towards developing QSAR models for the RAF1 and V600E-B-RAF enzymes based on multiple linear regressions (MLR) analysis of molecular descriptors (ca. 128) generated by VolSurfþ V R software 23,24 so that to gain more insight of structural trends responsible for activity discrimination. Two statistically significant models were achieved and found to fit well with the experimental biodata (Figures 9(A,B)) as indicated by their  corresponding Fischer's value (F), squared correlation coefficients (R 2 ), and the standard errors of estimate (s). The best developed model capable of describing the inhibitory pattern exerted by our compounds against the C-RAF enzyme is illustrated by the following equation; pIC50 ¼ À35:647 Â IW1 À 1:367 Â FLEX þ 13:505 where in this model; IW1 is the INTEGY (INTEraction enerGY) moment presented as vector pointing from centre of mass to the centre of the hydrophilic region, Low INTEGY descriptors such as IW1 indicates the hydrophilic moieties are either close to the centre of mass or they balance at an opposite ends of the molecule. Flex is the maximum flexibility of the molecule estimated after the random generation of 50 conformers; such descriptor is calculated as the logarithmic averaged differences between maximum and minimum distances of atom "i" in a selected conformer.
The two descriptors found to negatively influence the bioactivity as demonstrated by their corresponding negative coefficients, hence, for more improved RAF1 inhibitory activity, the model suggested lesser molecular flexibility which is in agreement with the intra-molecular HB observed during molecular dynamic simulation. Furthermore, distribution of hydrophilic moieties towards the terminal ends of the structure would greatly improve the activity (Table 6).
On the other hand, the best developed model for the V600E-B-RAF enzyme is illustrated by the following equation.   where, ACDODO is a 3D-pharmacophoric descriptor based on the triplets of pharmacophoric points, where it represents an Acceptor-Donor-Donor triplet. IW2 is the INTEGY (INTEraction enerGY) moment presented as vector pointing from centre of mass to the centre of the hydrophilic region, Low INTEGY descriptors such as IW2 indicates the hydrophilic moieties are either close to the centre of mass or they balance at an opposite ends of the molecule. DD8 is the difference between the maximum hydrophobic volumes (obtained upon variation of ligand conformations) and the hydrophobic volume (D8) of the imported 3D structure calculated at the 8 levels of energy. Furthermore, the hydrophobic volume (D8) is the molecular envelope generating attractive hydrophobic interactions calculated using the DRY-probe at an energy range of À1.6 kcal/mol. PSA is the Polar Surface Area calculated via the sum of polar region contributions.
Among the emerged descriptors, three were found to negatively influence activity (IW2, DD8 and PSA), while one descriptor (ACDODO) found to have positive impact. Thus, the model suggested the needed requirements for improved V600E-B-RAF bioactivity as follow; 1) the presence of three pharmacophoric points (Acceptor-Donor-Donor) is greatly preferred, which is consistent with the dynamic profile discussed earlier, 2) the proper positioning of the small hydrophilic moieties towards the terminal rings is highly advised and 3) the proper hydrophilic-lipophilic balance is highly desirable ( Table 6).

Conclusion
In this study, we studied the molecular mechanism of action of a promising series of imidazothiazole derivatives that are potent against melanoma cells. The target compounds were tested against V600E-B-RAF and RAF1 kinases. The structure-activity relationship study revealed the importance of the following structural moieties for activity; meta-hydroxyphenyl moiety at position 6 of the imidazothiazole nucleus, amide linker, and terminal benzene ring attached to the amide spacer either unsubstituted or substituted with a small group. The most promising compound that possesses all these pharmacophoric features is compound 1zb. It showed extremely high potency and relative selectivity towards V600E-B-RAF. Moreover, compound 1zb exerted high potency against different melanoma cell lines with sub-micromolar IC 50 values. It also showed high selectivity towards melanoma cells than normal skin cells. Furthermore, it could penetrate the cell membrane and inhibit V600E-B-RAF kinase with IC 50 value of 0.19 mM that is comparable to its IC 50 value against UACC-62 melanoma cells (0.18 mM). At 0.18 mM concentration, it induced caspase-3/7 release that lead to apoptosis. These promising results of compound 1zb especially makes it a promising drug candidate for further consideration and optimisation to end up with a new promising anti-melanoma agent.