Abstract
Triple-negative breast cancer is the leading worldwide cause of cancer-related deaths in women. The prospection and development of new substances with antitumoral potential is of great importance for the treatment of this disease. The objective of this work was to identify a commercial drug or ligand that could potentially bind to the FZD7 transmembrane protein and inactivate the Wnt signaling pathway in triple-negative breast cancer cells. We aimed at computationally modeling the FZD7, Wnt3, and Wnt3a proteins, at making them available in protein model databases, and at conducting docking analysis to assess the binding free energy between FZD7 and the selected ligands. The Wnt3 and Wnt3a proteins were modeled by homology modeling, and the FZD7 protein was modeled by homology modeling and ab initio modeling. The ligands were selected based on their similarity to the palmitoleic acid and were gathered from the ZINC database. A total of 30 commercially available ligands were found in the ZINC database. The docking results show that the ligands zinc08221009, zinc13546050, zinc05260769, zinc04529321, and zinc05972969 are good candidates for novel drug development. The created models and conducted analysis by this work will most certainly help in future research on the Wnt signaling pathway and its components.
1. Introduction
Breast cancer is the world's leading cancer in women and the second most common cancer type in the world, representing more than 20% of new cases every year. Although of good prognosis when diagnosed early, with cure rates that can reach 90%, in Brazil, the mortality rates are still high, mostly due to late diagnosis. In developed countries, the 5-year survival rate after treatment is 85%, whereas in developing countries this rate drops to 70% (World Health Organization, 2017).
Triple-negative breast cancer cells (TNBC) do not have a therapeutic molecular marker. In this type of tumor there is the superexpression of the canonical and noncanonical Wnt signaling pathway, which contributes to cell transformation and tumor progression. This pathway is activated by ligands of the transmembrane FZD7 protein, which activate the Wnt signaling cascade with subsequent accumulation of β-catenin in the nucleus of the cell. This β-catenin increase activates transcription factors that contribute to the development of tumors (Dey et al., 2013).
Bioinformatics, while still a somewhat new science, has been proving to be of great value to research since it allows for rapid and lower cost analysis. With the advance of technology and high-throughput techniques, there is massive generation of data that are stored in public databases, such as the Gene Expression Omnibus (GEO) database, Protein Database (PDB), and ZINC database (Berman et al., 2000; Irwin et al., 2012; Barrett et al., 2013). This extensive amount of data requires bioinformatics techniques for their handling, interpretation, organization, and analysis. This work aims at identifying, using bioinformatic tools, potential commercially available ligands to the FZD7 protein. It is expected that these ligands hinder the binding of Wnt3 and Wnt3a to FZD7, thus inactivating the Wnt signaling pathway.
2. Methodology
2.1. Triple-negative breast cancer cells superexpressed protein
An initial search for superexpressed proteins in TNBC was done in the GEO database (Barrett et al., 2013), a public database for gene expression. The results were filtered, and the protein Frizzled 7 (FZD7) was chosen (Yang et al., 2011), registered under GDS4069 in the database.
2.2. Molecular modeling
The primary protein sequence of the proteins FZD7, Wnt3, and Wnt3a, obtained from the UniProt database (Uniprot, 2007), with the codes O75084, P56703, and P56704, respectively, were used to search for those proteins' tridimensional structure in the PDB (Berman et al., 2000). There are no available experimental structures of the complete sequences of the proteins searched. The primary sequences were then used in a PSI-BLAST (Altschul et al., 1997) search against the PDB to find homologue proteins to them. The global alignment between the target proteins' sequences and their respective protein templates was done by using the EMBOSS Needle online tool (Li et al., 2015).
The homologue proteins found through the PSI-BLAST search were used as templates for homology modeling in the Modeller software (Sali and Blundell, 1993). Using a single template for Wnt3 and another single template for Wnt3a, 100 models were generated for each. With the use of multiple templates, another 100 models were generated for both Wnt3 and Wnt3a. As for FZD7, 100 models were generated by using multiple templates. Due to low quality of the models generated for FZD7, the online tool I-TASSER (Yang et al., 2015) was used for ab initio modeling. PyMOL software (Delano, 2015) was used in the editing and visualization of templates and generated models and in calculating the root mean square deviation (RMSD) between them.
2.3. Model quality assessment
The models created by homology were ranked based on their DOPE (Discrete Optimized Protein Energy) statistical potential (Shen and Sali, 2006) generated by the Modeller software, the best model being the one presenting the lowest value. DOPE is a statistical potential dependent on atomic distances; it is calculated by using a sample of the native structure of a protein and is based on the probabilistic theory. ERRAT (Macarthur et al., 1994) and VERIFY3D (Eisenberg et al., 1997) analyses were done to assess the atom distribution and 3D profiles of the models. For ERRAT analyses, an overall quality ≥80% is considered a good result. For VERIFY3D analysis, a good result is when 80% or more of the residues present a mean 3D-1D score ≥0.2.
The best models had their stereochemical quality evaluated and their Ramachandran plot generated by using the PROCHECK online tool (Laskowski et al., 1993) through the PDBsum database (De Beer et al., 2014). A Ramachandran plot, also known as Ramachandran diagram, is a way to visualize dihedral angles φ versus ψ of amino-acid residues of a protein, which allows visualization of energetically allowed regions for said angles per residue. In a Ramachandran plot, atoms are considered as rigid spheres with dimensions equal to their van der Waals radius. The φ e ψ torsion angles that cause collision between said spheres are represented as disallowed regions in the plot (Ramachandran et al., 1963). The best models were chosen as the ones with the lowest percentages of residues in the disallowed regions of the Ramachandran plot and with the lowest RMSD values calculated by PyMOL (Delano, 2015).
2.4. Ligand search
For the virtual screening of ligands, the ZINC database was used (Irwin et al., 2012). A search was done for commercially available compounds that had 99% similarity to the palmitoleic acid using the Tanimoto coefficient (Willett, 2006).
2.5. Molecular docking
A docking analysis was done between the ligands found in the virtual screening step and the cysteine rich domain (CRD) of the FZD7 protein, which was represented by the protein pdb5t44 with a resolution of 1.9 Å (Nile et al., 2017). The rigid-flexible docking methodology was applied by using the Autodock Vina software (Trott and Olson, 2010) and for graphical visualization and file editing AutoDockTools (Morris et al., 2009) was used. The docking box was visually established based on the binding region between the protein XWnt8 and the CRD of the protein FZD8 from the pdb4f0a template (Janda et al., 2012) used during the homology modeling of the Wnt3a protein. The dimensions were 22 × 24 × 26, with the x center at −12.207, the y center at −14,689, and the z center at 15.32. The palmitoleic acid (ZINC08221009) was considered as the control ligand in the docking analysis due to its role in the Wnt pathway activation (Doubravska et al., 2011).
The files of the ligands obtained from the ZINC database were in mol2 format and were then converted to the pdbqt format by using OpenBabel (O'Boyle et al., 2011), which is the format used by AutoDockTools. The X-Score software (Wang et al., 2002) was used to validate the binding free energy results estimated by AutodockVina. The software PyMOL (Delano, 2015) was used to add hydrogen atoms to the ligands after the docking process and to convert the file format from pdb to mol2 for those to be read by X-Score.
3. Results
3.1. Molecular docking and quality assessment
Due to the lack of available tridimensional structures for the Wnt3 and Wnt3a proteins, a BLAST (NCBI) search against the PDB for similar proteins was done. For such, the FASTA sequences of the Wnt3 and Wnt3a proteins were used and redundant results were discarded (Altschul et al., 1990). The pdb4f0a_B protein was the only significant alignment found for the proteins Wnt3 and Wnt3a, with a sequence identity of 41% for both target sequences and a 3.25 Å resolution (Janda et al., 2012). Another protein, pdb4krr_A with a 2.1 Å (Chu et al., 2013), was also a result from the PSI-BLAST search; however, it only presented 24% of sequence identity with the target sequences.
For the Wnt3 protein, there was no equivalence for residues 1–59 and residue 355 with the structure of the aligned protein, therefore those residues were not considered during its modeling step. As for Wnt3a, there was no equivalence for residues 1–56 and residue 352 with the structure of the aligned protein, therefore those were not considered during its modeling step. None of the aforementioned residues has an active role in the binding of those proteins tp FZD7. The primary sequences of proteins Wnt3 and Wnt3a were aligned with the sequence of the template pdb4f0a_B (Fig. 1) by using the EMBOSS Needle online tool (Li et al., 2015).

Sequence alignment between the pdb4f0a template and proteins Wnt3
FZD7 has three defined experimental structures of fractions of its structure: pdb4z33, with a 2.45 Å resolution, representing residues 569 to 574; pdb5t44, with a 1.99 Å resolution, representing residues 31 to 168; and pdb5urv, with a 2.20 Å resolution, representing residues 30 to 168. A BLAST (NCBI) search against the PDB was done by using the primary sequence of FZD7 to find proteins with similar sequences while excluding redundant results. Only two alignments resulted from this search, which were used as templates for the modeling step (pdb4jkv_A e 5t44_A); however, the FZD7 residues 1–29, 169–229, and 566–574 did not have matching residues in the templates.
Of the 100 models created for Wnt3 using multiple templates (pdb4f0a_B and pdb4krr), model number 91 was chosen as the best since it showed the lowest DOPE value, −26226.7. The model was then submitted to stereochemical analysis by using PROCHECK through the PDBsum online tool, which generated its Ramachandran plot. Its plot showed 90.3% of residues in the most favorable regions, 8.1% in additionally allowed regions, 0% in generously allowed regions, and 1.6% in disallowed regions (Fig. 2). The RMSD analysis, which was done by using PyMOL, between the created Wnt3 model and its templates resulted in 0.341 Å. The model's overall quality factor was calculated at 71.58% by using the ERRAT tool. The 3D profile analysis, using the VERIFY3D software, resulted in 79.66% of residues with a mean 3D-1D score ≥0.2.

Homology modeling results of protein Wnt3.
As for the 100 Wnt3a models generated by using multiple templates (pdf4f0a_B e pdb4krr), model number 1 was chosen since it showed the lowest DOPE value, −25848.1. The model's Ramachandran plot shows 86.8% of residues in the most favorable region, 8.9% of residues in additionally allowed regions, 2.7% in generously allowed regions, and 1.6% in disallowed regions. The ERRAT analysis resulted in an overall quality of 66.19% for the model, whereas its 3D profile analysis using VERIFY3D showed 78.31% of residues with a mean 3D-1D score ≥0.2. The RMSD analysis between the best model and its templates resulted in 0.430 Å (Fig. 3).

Homology modeling results of protein Wnt3a.
FZD7 model number 49 was chosen from the 100 models created by using multiple templates (pdb4jkv_a and pdb5t44_a), since it showed the lowest DOPE value, −61017.5. The Ramachandran plot created by using PROCHECK showed 88.1% of residues in the most favorable region, 9.2% of residues in additionally allowed regions, 1.2% of residues in generously allowed regions, and 1.5% of residues in disallowed regions. The model was then submitted to ERRAT analysis, which resulted in a score of 43.2% (Fig. 4). The model's 3D profile was evaluated by using VERIFY3D, which resulted in 52.44% of residues with a mean 3D-1D score ≥0.2. Due to the low quality of the homology generated model, it was decided to use the ab initio modeling approach for FZD7 using the I-TASSER software.

FZD7 model 49 ERRAT analysis.
The best models created by the ab initio approach were ranked by I-TASSER using their C-score. A variation of three templates (pdb517dA, pdb4jkvA, pdb4qinA) was used to generate the models. The best model showed an RMSD value of 9.3 ± 4.6 Å and a C-score of −0.74. The ERRAT analysis of the best model had a score of 83.6%, whereas its VERIFY3D analysis resulted in 64.29% of residues showing a mean 3D-1D score ≥0.2. A Ramachandran plot was generated and it showed 75.4% of residues in the most favorable region, 19.4% of residues in additionally allowed regions, 2.5% of residues in generously allowed regions, and 2.7% of residues in disallowed regions (Fig. 5).

FZD7 ab initio modeling results.
To increase the bioinformatics analysis reliability, it was decided to use an experimentally defined tridimensional structure (pdb5t44) for the docking analysis. This structure has a 1.99 Å resolution and it corresponds to the CRD of the FZD7 protein covering residues 31 to 168. The CRD is the FZD7 binding site to the proteins Wnt3 and Wnt3a (Janda et al., 2012).
3.2. Ligand search and virtual screening
The ZINC database search resulted in 29 compounds, and their structure files were downloaded in mol2 and sdf formats. The palmitoleic acid structure file was also downloaded. FAF-Drugs 4 (Alland et al., 2005; Lagorce et al., 2015), an online server, was used to check the compounds for the Lipinski's rule of five (Table 1). The Lipinski's rule of five is a set of rules that describe molecular properties to determine whether a chemical compound has characteristics similar to those of known drugs. For a molecule to be rejected as a drug candidate, it must violate two or more of the following rules: molecular weight ≤500 g/mol; LogP ≤5; number of hydrogen bond donors (HBD) ≤5; and number of hydrogen bond acceptors (HBA) ≤10 (Lipinski et al., 2001). None of the evaluated ligands violated more than one rule, which indicates that these compounds have drug-like properties.
Lipinski's Values of 30 Compounds from the ZINC Database
HBA, hydrogen bond acceptors; HBD, hydrogen bond donors. Boldface is palmitoleic acid.
The literature search to find a protein with similar activity to that of Wnt3a resulted in only one antibody, anti-FZD7 (CRD) named vantictumab (Anti-Fzd or OMP-18R5). However, the researchers responsible for its creation do not have its tridimensional structure available (Gurney et al., 2012). There are no antibodies available in public databases that are able to bind to the CRD of the FZD7 protein.
3.3. Molecular docking
The molecular docking results between the 30 ligands, and FZD7 are shown in Table 2 along with their X-Score values. The best five ligands are shown in Figure 6. The OpenBable software (O'Boyle et al., 2011) was used to convert the ligand file formats from mol2 to pdbqt.

Chemical structure of the five best ligands shown in Table 2.
Binding Free Energy Between 30 Compounds from the ZINC Database and the Cysteine Rich Domain of the FZD7 Protein Using AutoDock Vina and X-Score
Boldface is palmitoleic acid.
The palmitoleic acid (zinc08221009), which was here considered a control ligand, had a calculated binding free energy of −8.0 kcal/mol. Two ligands, zinc13546050 and zinc05260769, had the same binding free energy as the control molecule. Two other ligands, zinc04529321 and zinc05972969, had binding free energies of −8.1 kcal/mol. The LigPlot software (Wallace et al., 1995) was used to visually demonstrate the interaction between the ligands and FZD7. Images of the best interactions are shown in Figure 7 along with a visualization of the docking created by PyMOL.

Interaction between FZD7's CRD (pdb5t44) and the best ligands from the molecular docking analysis. In the right-hand panel, the semicircles indicate van der Waals interactions between FZD7 residues and the ligand. Doted lines represent hydrogen bonds. CRD, cysteine rich domain.
The LigPlot representations show that the interaction between the palmitoleic acid (zinc08221009) and FZD7's CRD is mediated by hydrophobic contacts. The same is true to ligand zinc13546050. Ligands zinc05260769 and zinc04529321 show a hydrogen bond with residue Lys61, besides the hydrophobic contacts. Ligand zinc05972969 shows a hydrogen bond with residue Lys61 and another with Gln55 besides the hydrophobic contacts.
4. Discussion
Triple-negative breast cancer is the breast malignancy with the worst prognosis to date. This is partially due to the lack of molecular markers in these cells and no specific and efficient treatment for them. Since there are currently no target receptors for drugs in these cells, chemotherapy is still the best course of treatment despite its adverse effects. Thus, the prospect of new and more effective, specific, and less cytotoxic drugs is of great importance. Development of new treatments is, however, slowed by the aforementioned lack of well-defined molecular markers in this type of cancer. A better understanding of the underlying molecular mechanisms involved in the development and progression of this cancer type is necessary to allow the creation of new and better drugs and therapies. This type of study can also lead to the discovery of new therapeutic targets.
The relationship between the canonical Wnt signaling pathway and the development and progression of tumors has been studied for more than 50 years (Nusse and Varmus, 1982). The discovery and development of new therapeutic agents is hindered by the lack of a detailed understanding of receptor
There are currently no studies, to the best of our knowledge, that identify lead molecules that target the CRD of the FZD7 protein. There are studies that identify molecules that block the binding between Wnt3a and the co-receptor LRP6 (Lu et al., 2011, 2012; Zhang et al., 2013). However, as shown by Kuhnert et al. (2004), interrupting the Wnt pathway by blocking interactions to LRP6 results in the loss of gastrointestinal homeostasis in rats, which makes the CRD of the protein FZD7 a more attractive target for the development of treatments.
The models of the proteins Wnt3, Wnt3a, and FZD7 created by this work can be used for docking analysis and molecular dynamics in future studies. In the light of that, this study can serve as a basis for the development of new models and further molecular studies targeting the Wnt signaling pathway for new therapies. The best models of the proteins Wnt3 and Wnt3a were submitted to the Protein Model Database (PMDB) under the identification PM0081550 and PM0081549, respectively. The PMDB requires that submitted models have a maximum distance between consecutive carbon alpha (CA) atoms in the range of 3.6 and 4.0 Å, and a distance higher than 3.6 Å between nonconsecutive CA atoms. The databank automatically checks for those requirements and for major stereochemical problems in submitted models (Castrignanò et al., 2006). Molecular dynamics analysis can be done in the future to optimize the quality of the generated models.
The results here presented reveal as potential lead molecules the ligands zinc13546050, zinc05972969, zinc05260769, zinc0452932, and the palmitoleic acid by their binding free energies against the CRD of the FZD7 in the docking analysis. The binding free energy of ligands zinc05972969 and zinc04529321 was superior to that of the control molecules. The four ligands and the palmitoleic acid have drug-like properties as shown by their compliance to the Lipinski's rule of five. Future in vitro and in vivo studies are required to assess the effects of these ligands in TNBC. The results here presented should serve as a starting point for future studies involving the Wnt signaling pathway, the molecular interaction between its components, and the development of new drugs that target this pathway.
Footnotes
Author Disclosure Statement
The authors declare they have no competing financial interests.
Funding Information
No funding was received for this article.
