Proteome-wide prediction of targets for aspirin: new insight into the molecular mechanism of aspirin

Overview of pipeline for proteome-wide prediction of aspirin targets

The pipeline for proteome-wide prediction of aspirin targets is outlined in Fig. 1. Firstly, we constructed the structural database of human proteins (17,425 non-redundant structures), and the binding sites of aspirin (BSiteAs) were used to search against this database using the program CMASA. Secondly, the binding abilities of aspirin to these putative off-targets were estimated using molecular docking. If aspirin docked unsuccessfully to the predicted binding pocket of a particular protein, this protein was removed from the target list of aspirin. Thirdly, the remaining putative targets were subject to molecular mechanics Poisson–Boltzmann surface area (MM-PBSA) free energy calculation. Finally, we performed the pathway enrichment analysis for these putative targets.

The structure database of human proteins

The query database was built by integrating three available sources of protein structure. The first source is the Protein Data Bank managed by Research Collaboratory for Structural Bioinformatics (RCSB PDB) (http://www.rcsb.org) (Rose et al., 2011), which contains 6,983 structures of human proteins. The other two repositories are structure models built by homology modeling techniques from Structure Atlas of Human Genome (Motono et al., 2011) (http://bird.cbrc.jp/sahg) and GPCR Research Database (GPCR-RD) (Zhang & Zhang, 2010), which have 16,529 and 1,028 structures, respectively. This procedure yielded a total of 24,540 human protein structures. The three available sources of protein structure were integrated and removed the redundancy using CD-HIT Suite (Huang et al., 2010). The identity cut-off was set to 0.95. There are totally 17,425 structures in our self-built ‘non-redundant’ structural database of human proteins. The structural database represents a relatively complete and accurate library of human proteins.

To detect and compare the aspirin binding sites in the structure database

We used all representative aspirin-protein complexes in our study. There are a total of six aspirin-protein complex structures currently available (1OXR, 1TGM, 2QQT, 3GCL, 4NSB and 3IAZ). Both 1OXR and 1TGM are structures of phospholipase A2, we only used 1OXR in this study (Fig. 2A). 2QQT and 3GCL are bovine lactoperoxidase, and we used 3GCL in this study. 4NSB and 3IAZ are structures of buffalo chitinase-3-like protein 1 and bovine lactotransferrin, respectively. The two structures were also used in our study. In addition, 1PTH, the complex structure of salicylic acid and COX-1, was also used in this study (Fig. 2A). Therefore, a total of five complex structures (1OXR, 3GCL, 4NSB, 3IAZ and 1PTH) were used in this study. The binding sites of aspirin (BSiteAs) were defined as amino acid residues lying within 5 Å from atoms of the aspirin. The BSiteAs were used to search against 17,425 non-redundant structures using the program CMASA (http://bsb.kiz.ac.cn/CMASA_flex) which was developed in our lab. The accuracy and sensitivity of this tool were 0.96 and 0.86, respectively. The workflow of CMASA was as following. First, the CMASA parsed the query binding site and the database to be searched. Second, for each structure in database, based on the query binding site, the CMASA enumerated all possible combination of residues in the structure using amino acid substitute. The blocks substitution matrix 62 (BLOSUM62, cutoff = 1) was used. These residue combinations are similar to that in the query binding site. Each residue combination forms a possible local structure. Third, the CMASA used contact matrix average deviation (CMAD) between the query binding site and these local structures to filter the candidate matches. Then, the CMASA calculated the RMSD and the RMSD based P-value if the CMAD < cutoff. At last, the CMASA ranked all of the hits and add their information. Hits are considered significant if the CMASA E-value <0.01.

Molecular docking

The hit proteins have similar local structures with BSiteAs and potential to bind to aspirin. However, it does not mean that aspirin can certainly bind to these proteins. To assess whether aspirin can bind to these proteins, flexible ligand-rigid protein docking was performed using CHARMm-based DOCKER (CDOCKER) (Wu et al., 2003) in the Discovery Studio v3.1. The following steps are included in the CDOCKER protocol. (1) A set of ligand conformations are generated using high-temperature molecular dynamics with different random seeds. (2) Random orientations of the conformations are produced by translating the center of the ligand to a specified location within the receptor active site, and performing a series of random rotations. A softened energy is calculated and the orientation is kept if the energy is less than a specified threshold. This process continues until either the desired number of low-energy orientations is found, or the maximum number of bad orientations have been tried. (3) Each orientation is subjected to simulated annealing molecular dynamics. The temperature is heated up to a high temperature then cooled to the target temperature. (4) A final minimization of the ligand in the rigid receptor using non-softened potential is performed. (5) For each final pose, the CHARMm energy (interaction energy plus ligand strain) and the interaction energy alone are calculated. The poses are sorted by CHARMm energy and the top scoring (most negative, thus favorable to binding) poses are retained. In this study, we generated 10 random conformations for each ligand. The parameters of the dynamic steps, the dynamic target temperature and include electrostatics are set to 1,000 steps, 1,000 K and True, respectively, which is the default setting in CDOCKER. The binding sphere of CDOCKER is defined around the local structure detected by CMASA. The center of binding sphere is set as the center of the local structure. And the radius of binding sphere is set as 10 Å, which allows the free rotation of aspirin. If aspirin docked successfully to a particular protein, the binding poses of aspirin showing the lowest energy were retained and used for MM-PBSA (Kollman et al., 2000; Kuhn et al., 2005) free energy calculation. The 3-D structures of docked complexes were visualized using PyMol v1.5.

MM-PBSA free energy calculation and entropy change estimation

The complex structures of aspirin and putative targets were further used to obtain more accurate estimate of binding free energy. Binding free energies (ΔG_bind) of aspirin at the binding site on these target proteins were calculated by using the MM-PBSA free energy calculation protocol in Pipeline pilot v8.5 (http://accelrys.com/products/collaborative-science/biovia-pipeline-pilot/) as follows: (1)${\displaystyle \Delta G_{\mathrm{binding}}^{\mathrm{PB}}=G_{\mathrm{complex}}-G_{\mathrm{receptor}}-G_{\mathrm{ligand}}}$ where G_complex, G_receptor, G_ligand are the free energies of the complex, receptor and ligand respectively. The free energy of each molecule on the right hand side can be considered as the sum of molecular mechanics energy in gas phase(E_MM) and solvation free energy(Δsol).

As the entropy change (TΔS) is the most time-intensive part, it was not included in the above calculation. We estimate the entropy change in another way. In the study of Chang & Gilson (2004), they computed the entropy changes for different receptors upon same ligand binding. The results showed the entropy changes for different receptors were similar. Therefore, we can estimate an approximate entropy change. The process of entropy change estimation is as following: (1) Based on the formula:—RTlnK = ΔG = ΔH − TΔS (where R is the gas constant, T is the absolute temperature), the entropy changes (TΔS) is equal to the sum of ΔH and RTlnK; (2) The experiment study have shown that aspirin binds to PLA2 enzyme (PDBID: 1OXR) specifically with a binding constant (K) of 1.56 × 10⁵M⁻¹ (Singh et al., 2005). For this complex, the enthalpy change (ΔH) was calculated as −2.327 kcal/mol using the MM-PBSA and RTln K was calculated as 7.055 kcal/mol. So the entropy change (TΔS) is 4.728 kcal/mol (equal to −2.327 kcal/mol + 7.055 kcal/mol) for this complex. 3) Based on above assumption and computation, we set the entropy change in the process of aspirin binding to various putative targets as approximate value 4.728 kcal/mol.

Pathway enrichment analysis and interaction network construction

To analyze the significance of KEGG pathways associated with our predicted targets, we collected UniProtKB accession number (AC) of these targets and performed KEGG pathway annotation using the DAVID tool (http://david.abcc.ncifcrf.gov/). The significantly over-represented KEGG pathways were identified based on the Bonferroni-adjusted P value (P < 0.01) (Huang, Sherman & Lempicki, 2009). In addition, based on these pathways, an integrated targets-cellular effect interaction network was constructed.

Identification of putative targets of aspirin in human proteome

We presented a proteome-wide prediction of aspirin targets using structural bioinformatics and system biology approaches. We used comparison of BSiteAs to recognize putative targets and further refined by docking and MM-PBSA in structural bioinformatics part whereas pathway enrichment analysis and interaction network construction were performed in system biology section. The steps in our pipeline for proteome-wide prediction of aspirin-binding proteins are shown in Fig. 1. Firstly, the binding sites of aspirin (BSiteAs) were used as queries to search against 17,425 non-redundant structures of human proteins in our self-build structure database using the program CMASA. Totally, 79 proteins with putative BSiteAs were identified (Table S1). Of these proteins, the top 10 ranked proteins are members of the phospholipase A2. cyclooxygenase, lactoperoxidase and Chitotriosidase families, which are the primary targets of aspirin. The remaining 69 proteins have different structural folds from the primary targets.

The hit proteins have similar local structures with BSiteAs and potential to bind to aspirin. However, it does not mean that aspirin can certainly bind to these proteins. In the second step, molecular docking was used to assess whether aspirin can bind to these proteins. CDOCKER in the Discovery Studio v3.1 was used to dock aspirin to the predicted binding site on these proteins. Proteins that failed to dock aspirin were removed from the target list. Only 26 proteins were considered for further analysis after filtering by molecular docking, 10 proteins of which are the primary targets of aspirin (Table S1).

Finally, MM-PBSA free energy calculation was performed for the lowest-energy protein-aspirin complex obtained in the docking step. In total, 23 proteins bind to aspirin with binding free energies (ΔG_bind) < 0 listed in Table 1 and selected as the putative targets of aspirin. Next, we analyzed the binding modes and affinities of aspirin to these targets, the structural similarity and pathway enrichment of these targets and clinical outcomes of aspirin.

Structural diversity of the putative targets

In order to analyze the structural similarity of the 23 putative targets, each pair of these targets were structurally aligned using the program Combinatorial Extension (CE) (Shindyalov & Bourne, 2001). The targets were considered to be similar if the root-mean-square deviation (RMSD) between two structures is less than 4 Å. A network of targets was generated by linking structurally similar targets (Fig. 2A). The primary targets of aspirin are colored with green, and the newly identified targets are colored with red. Only two proteins HLA class I histocompatibility antigen, A-2 alpha chain (HLA-A) and axin interactor, dorsalization-associated protein (AIDA) clustered with the phospholipase A2 family. Proteins cyclooxygenase-1 and 2 are linked together but not similar with the other proteins. Therefore, the overall structures of the newly identified targets and the primary targets were not similar. It indicates the structural diversity of the 23 putative targets. However, the 23 targets have similar local structures. For example, proteins group IB phospholipase A2 (PLA2G1B) and cyclin-dependent kinase 13 (CDK13) which belong to different protein families have very similar binding sites of aspirin (BSiteAs) (Fig. 2B).

Diverse binding modes of aspirin to the putative targets

In our study, aspirin was docked to the predicted binding sites on putative targets. Results of these docking experiments reveal diverse binding modes of aspirin to these targets (Fig. 3). Some examples are given follows. As a first example, aspirin binds to protein CDK13 with a novel mode compared to ATP-analog inhibitor of the kinase (Fig. 3A). ATP-analog inhibitors exhibit inhibitory activity of kinase by competitive binding to its ATP binding site. In contrast, aspirin binds to CDK13 in the vicinity of the ATP binding site and interact with loops (L1 and L2) which are important for ATP binding. Another example is that aspirin binds to protein ras-related C3 botulinum toxin substrate 1 (RAC1), which is very different from GTP binding (Fig. 3B). Aspirin binds to the other side of RAC1 and interact with the N-terminal part of switch II (sequence ⁵⁶ WDTAG), which is crucial for interaction between RAC1 and protein Arfaptin (Tarricone et al., 2001). The Arfaptin mediates cross-talk between Rac and Arf signaling pathways. The last example is that aspirin binds to protein integrin alpha-L (ITGAL) in the binding site of its inhibitor BQM (Guckian et al., 2008) (Fig. 3C). The analysis of the binding mode of aspirin to the other 19 putative targets is shown in the Figs. S1–S23. The statistic of 23 putative binding pockets of aspirin is shown in the Fig. S24. The coordinates of 23 putative binding pockets are shown in Table S2. The size of the 23 binding pockets is ranging from 10 to 26 residues, if we define the residues lying within 5 Å from aspirin as binding pockets. The pocket size of aspirin binding to phospholipase A2 (PDBID: 1OXR) and cyclooxygenase-1 (PDBID: 1PTH) is fall in this range. The most used amino acids in pockets are leucine (L), cysteine (C), valine (V), glycine (G), tyrosine (Y), isoleucine (I), alanine (A) and phenylalanine (F). The eight amino acids appeared in the pockets more than 20 times. All the eight amino acids are with hydrophobic side chain except cysteine (C). It suggests that hydrophobic interaction is important for aspirin binding. Additionally, the H-bond is also important for aspirin binding because 78% (18/23) binding pockets formed H-bonds with aspirin. The most used amino acids involved in H-bond formation are aspartic acid (D), lysine (K) and histidine (H). Therefore, the combinations of these charged amino acids (D, K and H) and hydrophobic amino acids (L, V, G, Y, I, A and F) may form the pocket that aspirin prefers to bind to.

The binding affinities of aspirin to putative targets

To obtain accurate estimate of the binding energy of different putative targets with aspirin, MM-PBSA free energy calculation protocol was used. In the combination of the experimental study (Singh et al., 2005) and MM-PBSA free energy calculation, we obtained the estimated entropy change (TΔS = 4.728 kcal/mol) upon aspirin binding. The entropy changes do not have large fluctuations when the same ligand binds to a different acceptor based on the study of Chang & Gilson (2004). Therefore, the entropy changes when aspirin binds to various putative targets was assumed as 4.728 kcal/mol to compare free energies associated with different aspirin binding putative targets. The binding free energies including entropy change for the 23 proteins binding to aspirin were calculated and listed in Table 1. The binding free energies of the 23 proteins with aspirin are varied from −6.0 (group IID secretory phospholipase A2, PLA2G2D) to −33.0 (exosome component 3, EXOSC3) kcal/mol. Overall, the binding free energies for newly identified targets (the average −18.4 kcal/mol) are comparable to that for the primary targets (the average −15.3 kcal/mol).

Pathway enrichment and interaction network of putative targets

Using the DAVID tool, we find that our predicted targets are significantly overrepresented for several pathways (p < 0.01) (Table 2). Some of these pathways are strongly involved in inflammation, cardiovascular disease and cancer, such as VEGF signaling, Fc epsilon RI signaling, arachidonic acid metabolism, gonadotropin-releasing hormone (GnRH) signaling and MAPK signaling. To illustrate the relationship between the putative targets and their cellular effect, an integrated interaction network of targets-cellular effect based on their associated pathways was constructed (Fig. 4). The interactions between predicted targets and the major effects involved in cancer development, inflammation and cardiovascular disease were present in this network. Represented by green circles in the network, the predicted targets regulate VEGF, epsilon RI signaling, arachidonic acid metabolism, and MAPK pathways through interactions with other proteins (gray circles) connecting the pathways. Inhibition of predicted targets is expected to down-regulate these pathways, and then prevent inflammation and decrease the risk of cardiovascular disease and cancer.

We use the protein RAC1 as an example to show how inhibition by aspirin leads to anti-inflammation and anti-cancer effect. RAC1, which belongs to the Rho GTPase family, can regulate the machinery that controls the assembly and disassembly of cytoskeletal elements (Bid et al., 2013). RAC1 is required for the VEGF-induced increase in vascular permeability (Monaghan-Benson & Burridge, 2009). Increased vascular permeability is often observed in inflammation (Wilhelm, 1973). Therefore, RAC1 can promote the development of inflammation through VEGF and Fc epsilon RI signaling pathway. Furthermore, RAC1 is involved in gonadotropins expression, cytokine transcription and actin reorganization, which play a critical role for proliferation, migration and evading apoptosis of cancer cell. The inhibition of RAC1 by aspirin will attenuate these processes, which may provide an alternative explanation for the anti-inflammation and cancer effect of aspirin.