Ling Yang
Xiaoting Sun
Ying Ye3†
Shun Zhu- 1Department of Cellular and Genetic Medicine, School of Basic Medical Sciences, Fudan University, Shanghai, China
- 2Department of Medical Oncology, Shuguang Hospital, Shanghai University of Traditional Chinese Medicine, Shanghai, China
- 3Department of Oral Implantology, Shanghai Engineering Research Center of Tooth Restoration and Regeneration, School and Hospital of Stomatology, Tongji University, Shanghai, China
- 4Department of ENT, Second People's Hospital of Shenzhen, First Affiliated Hospital of Shenzhen University, Shenzhen, China
- 5Department of Biochemistry, University of Iowa, Iowa City, IA, United States
p38 mitogen-activated protein kinases are signaling molecules with major involvement in cancer. A detailed mechanistic understanding of how p38 MAPK family members function is urgently warranted for cancer targeted therapy. The conformational dynamics of the most common member of p38 MAPK family, p38α, are crucial for its function but poorly understood. Here we found that, unlike in other cancer types, p38α is significantly activated in pancreatic adenocarcinoma samples, suggesting its potential for anti-pancreatic cancer therapy. Using a state of the art supercomputer, Anton, long-timescale (39 μs) unbiased molecular dynamics simulations of p38α show that apo p38α has high structural flexibility in six regions, and reveal potential catalysis mechanism involving a “butterfly” motion. Moreover, in vitro studies show the low-selectivity of the current p38α inhibitors in both human and mouse pancreatic cancer cell lines, while computational solvent mapping identified 17 novel pockets for drug design. Taken together, our study reveals the conformational dynamics and potentially druggable pockets of p38α, which may potentiate p38α-targeting drug development and benefit pancreatic cancer patients.
Introduction
p38 mitogen-activated protein kinases (MAPKs) play critical roles in cellular responses, proliferation, survival, cell cycle, and migration in cancer. p38 MAPK family includes p38α (MAPK14), p38β (MAPK11), p38γ (MAPK12), and p38δ (MAPK13). The four p38 MAPK family members have different tissue expression patterns, with p38α being ubiquitously expressed at significant levels in most cell types. The p38 MAPKs function in a cell context-dependent manner (1–4). However, a detailed mechanistic understanding of how p38 MAPK family members function is still not well-understood. A major challenge will be to determine when and how to specifically target p38 MAPK for disease therapy.
To gain insights into its molecular mechanisms and design potential therapeutics, the structure of MAP kinase p38α has been extensively studied over the last two decades (5–8). X-ray crystallography and nuclear magnetic resonance (NMR) showed that p38α consists of two domains, a 135-residue N-terminal domain mainly composed of β-sheets and a 225-residue C-terminal domain mainly composed of α-helices, in between of which lies the catalytic site, i.e., the ATP-binding pocket (9, 10). Despite its seemingly rigid crystal structures, p38α is highly dynamic, which is supported by various evidence. Firstly, in the majority (~78%) of the p38α crystal structures, e.g., 3L8X (11), 3OCG (12), and 2NPQ (13), the glycine-rich loop and/or the activation loop is invisible in the electron density map, indicating the structural flexibility of these loops. Secondly, hydrogen-exchange mass spectrometry (HX/MS) study showed that phosphorylation of the activation loop induced conformational changes in p38α (14); While X-ray crystallographic study showed that phosphorylation brings the N-terminal domain and C-terminal domain closer (15). Thirdly, NMR experiments of p38α in apo and ligand-bound forms suggested that the ATP-binding pocket is highly flexible even after ligand binding (16, 17). These studies reveal the important roles of conformational dynamics in the activation and catalysis of p38α.
Pharmacologically, p38α is considered as a potential target for various diseases such as inflammatory diseases and cancer. The potent small-molecule inhibitor SB203580 competitively binds the ATP-binding pocket and has been widely used to study p38 MAPK functions (5). Nevertheless, given the similarity of the ATP-binding site in different kinases, SB203580 also targets non-p38 protein kinases, usually at higher concentrations (18). Over the years, many p38α inhibitors have been developed, for example BIRB796, which leads to a conformational reorganization that prevents ATP binding and activation (19). The C-terminal domain has also been predicted to have the flexibility for potentially binding differently shaped compounds (20). For p38α inhibitor development, allosteric small-molecule inhibitors that target other regions of the kinase are warranted and they might reduce the off-target effects in drug applications.
Clinically, p38 MAPK inhibitors show significant effects in pre-clinical animal models but repeatedly failed in clinical trials (21). In cancer study, p38α and p38β increase cell proliferation and invasion of colon cancer, follicular lymphoma, ovarian cancer, and more recently, pancreatic cancer (22–25). It is indicated that targeting p38 MAPKs, especially p38α, should be explored for cancer therapy. However, several phase I clinical studies with p38 MAPK small-molecule direct or indirect inhibitors show hepatic, neurological, gastrointestinal, and cardiovascular toxicities, indicating that the commonly studied inhibitors are not highly selective (26). Currently, one of the major challenges in pancreatic cancer drug development is to overcome drug off-target effects and drug resistance (26, 27). It is rational to speculate that next-generation highly selective p38 MAPK inhibitors may exhibit less adverse effects. However, it still requires further investigation in cancer patients.
In this work, by screening pancreatic adenocarcinoma (PDAC) patient samples (n = 40) and TCGA database, we demonstrated that p38 MAPKs, especially p38α are highly expressed and activated. p38α blockades show significant anti-tumor effect in PDAC cells, but are not selective enough. To pave the way for highly selective inhibitor development, conformational dynamics of p38α was examined. Anton supercomputer long-timescale (39 μs) MD simulations using both AMBER and OPLS force fields show the p38α flexibility in six regions. Mechanistically, p38α MD simulations reveal a “butterfly” motion that might be important for p38α catalytic function. In addition, computational solvent mapping reveals 17 novel pockets that are potentially druggable for cancer therapy. To our knowledge, this is the first comprehensive study of both the conformational dynamics and potentially druggable pockets in p38α. This study provides insights into understanding the molecular mechanism of p38α function and into developing potential drugs with high specificity and selectivity against PDAC.
Results
p38 MAPKs Expression Correlates With Poor Prognosis in PDAC Patients
To investigate the role of p38 MAPKs in cancer, we screened a panel of human tumor tissues in The Cancer Genome Atlas (TCGA) that spontaneously express p38 MAPKs. We have found that the majority of the tumor tissues express similar level of p38 MAPKs to its adjacent healthy tissues (black labeled, Figure 1A, Figures S1A–C), while uterine carcinosarcoma (UCS), uterine corpus endometrial carcinoma (UCEC), and chromophobe renal cell carcinoma (KICH) show decreased p38α (MAPK14), p38β (MAPK11), and p38γ (MAPK12) expression compared with control pancreas (green labeled, Figure 1A, Figures S1A,B). Surprisingly, between more than 30 types of cancer, only pancreatic adenocarcinoma (PAAD or PDAC) shows a significant increase of p38α, p38β, and p38γ compared with healthy control pancreas (red labeled, Figure 1A, Figures S1A,B), suggesting that p38 MAPKs may be involved in PDAC development. To further validate these findings, we detected the mRNA level for all four members of the p38 MAPK family. Indeed, p38α, p38β, and p38γ are all increased markedly, while p38δ (MAPK13) shows a trend of increase (Figure 1B). p38α shows a significant dominance among these genes (Figure 1C), indicating its importance among four members. We next investigated the relationship of p38 with PDAC clinical outcomes. Surprisingly, p38α showed a trend of increasing with advanced stage, and strongly correlated with poor overall survival in PDAC patients (Figures 1D,E). These results suggest that p38α plays a role in PDAC progression and correlates with poor prognosis in PDAC patients.
Figure 1. p38α expression correlates with poor prognosis and adipose markers in PDAC patients. (A) Transcriptomic expression levels of MAPK14 across multiple cancer types and paired normal samples, with each dot representing a distinct tumor or normal sample. Red dot, tumor sample; Green dot, control sample; Red group name, significantly upregulated; Green group name, significantly downregulated; Black group name, not significant. (B) Transcriptomic expression levels of MAPK14, MAPK11, MAPK12, and MAPK13 in PDAC and adjacent pancreas samples [n(Control) = 179 samples; n (PAAD) = 171 samples]. (C) Transcriptomic expression levels of MAPK14, MAPK11 and MAPK12 in PDAC samples [n (PAAD) = 171 samples]. (D) Violin plots of MAPK14 based on PDAC patient pathological stage [n (PAAD) = 171 samples]. (E) Overall survival (OS) analysis of PDAC patients based on MAPK14 expression. (F) Correlation analysis of MAPK14 and adipose markers (PLIN1, PLIN2, PLIN3, PLIN4, PLIN5) expression in human PDAC tissues and adjacent pancreas. Control group = 179 samples; PAAD group = 171 samples. *p < 0.05. NS = not significant. Data presented as mean ± s.e.m.
p38α Expression Correlates With Adipose Markers in PDAC Tissues
It was reported that proliferating cancer cells may take up exogenous lipids and activating endogenous lipid biosynthesis (28), and tumor implanted in adipose environment show significant lipid metabolic reprogramming (29). Considering that PDAC is one of the tumors that adjacent to the adipose environment, we tested the correlation of p38 MAPKs and lipid droplet marker perilipin (PLIN) family in the PDAC database. Surprisingly, p38α strongly correlated with PLIN 2 and 3, two small lipid droplets markers, but not PLIN1, 4, and 5, which are the big lipid droplets marker (Figure 1F). It was plausible that small lipid droplets exist inside pancreatic tumor cells or in mobilized adipocytes in the tumor microenvironment. These findings support that p38α correlates with adipose-rich PDAC and may be involved in cancer lipid metabolism.
p38α Is Activated in PDAC Patient Samples
To validate our findings, we measured p38α levels in human PDAC patients. In this study, the patient PDAC sample and adjacent pancreas were compared. In accordance with previous report (30), PDAC tissues show significant infiltration of inflammatory cells, stromal cellular components, and activated fibroblasts (Figure 2A). Next, we detected p38α in these samples. As expected, p38α is highly expressed in human PDAC tissues (Figures 2B,D; Figure S2). Consistent with this result, western blot and immunohistochemistry staining show highly activated phospho-p38α in human PDAC tissues (Figures 2B,D). Furthermore, pathological analysis shows that both p38α and phospho-p38α are mainly located in epithelial cancer cells (Figure 2C). Western blot of various cell lines demonstrates that PDAC cancer cells and endothelial primary cells expressed high levels of p38α protein, whereas human THP-1 macrophage-like cell line, and human stromal fibroblasts lacked p38α expression (Figure 2E). These findings further support our notion that p38α is activated in PDAC cells. These pilot clinical findings validate our TCGA data.
Figure 2. p38α is activated in cancer cells in PDAC patient samples. (A) Micrographs of H&E and immunohistochemistry staining with CD163, PDGFRβ, and αSMA. Quantification of CD163+, PDGFRβ+, and αSMA+ signals (n = 8 random fields per group). (B) Micrographs of immunohistochemistry staining with p38 and Phospho-p38 in PDAC and adjacent pancreas samples. Quantification of p38+, Phospho-p38+ signals (n = 8 random fields per group). (C) Pathological analysis of p38 and Phospho-p38 in PDAC and adjacent pancreas samples (n = 20 samples per group). (D) Protein expression levels of p38 and Phospho-p38 in human PDAC tissues and adjacent pancreas (n = 20 samples per group). (E) protein expression levels of p38 and Phospho-p38 in various human cell lines (n = 3 samples per group). (F) Cell viabilities of Pan02 cell lines treated with 3.125–200 μM SB203580 for 24 h (n = 6 samples per group). *p < 0.05; **p < 0.01; ***p < 0.001. NS = not significant. Data presented as mean ± s.e.m.
p38α Blockade Enhances Apoptosis of Human and Mouse PDAC Cells
To gain further mechanistic insights of p38α blockade on cancer, human pancreatic adenocarcinoma cells MiaPaca-2, Panc-1, and mouse PDAC cell line Panc02 were projected for p38α blockade. The canonical inhibitor SB203580, which can block both p38α and p38β (31), and ralimetinib, a small molecule inhibitor specifically for p38α under phase 2 clinical trial (32), were tested. As expected, both drugs significantly inhibit PDAC cell viability, with the IC50 around 56.89μM (Figure 2F). These findings suggest the potential of p38α as a therapeutic target for PDAC.
p38α Is Highly Dynamic Locally and Globally
Despite its potential for PDAC therapy, our understanding of the structural and functional basis of p38α is still poor. To investigate the structural dynamics of p38α, we applied the state of the art supercomputer, Anton for simulation of our potential targets (33). Three simulation runs of apo p38α and one simulation run of ligand-bound p38α (Table 1) for the AMBERff99SB-ILDN and OPLS-AA/L force fields (abbreviated as AMBER and OPLS) were performed (Table 2) and the convergence was evaluated (Figure S3). We plot the root mean square deviations (RMSDs) of the protein backbone from their crystallographic positions for all simulation runs. In AMBER simulations, RMSDs of apo p38α range from 2 to 5 Å (Figure 3A, left panel, red, orange, and green lines); in OPLS simulations, RMSDs of the apo p38α run go up to 8 Å (Figure 3A, right panel, red line). In contrast, ligand-bound p38α is below 4 Å in both AMBER and OPLS simulations (Figure 3A, blue lines). The RMSD discrepancy between AMBER and OPLS simulations is possibly due to their different preference for secondary structure propensities. The RMSD discrepancy between apo and ligand-bound p38α suggests that the conformational flexibility of p38α is dampened upon ligand binding. It should be noted that as the N terminal and C terminal ends of p38α can be quite floppy (Figure S4), these ends were truncated in MD simulations. These results show that apo p38α exhibits higher structural flexibility than ligand-bound p38α.
Table 1. Equilibrium angles for the eight restrained dihedral angles of the ligand SB203580.
Table 2. Simulation parameters for MD simulations of apo and SB203580-bound p38α.
Figure 3. Dynamics of p38α kinase in MD simulations. (A) RMSDs of p38α in AMBER and OPLS simulations. Three runs of apo p38α (apo1, apo2, and apo3) and one run of ligand-bound p38α are shown in red, orange, green, and blue, respectively. The probability density functions are shown in the right panel. (B) Representative snapshots (pink) aligned to the initial structure (light blue). Regions of interest are marked by magenta circle. (C) RMSFs mapped onto p38α structure. Sausage representation is used for the protein. The thickness of the tube is normalized within each simulation run with thicker tubes corresponding to higher RMSFs. The color of the tube is normalized using all simulation runs using the “rainbow” gradient with warmer color corresponding to higher RMSFs. Color normalization was done by setting the color range to 0–8.3 Å (the largest RMSF value in all simulations). The N-terminal end (residue 4–13) and C-terminal end (residue 345–354) are not shown and excluded from the color normalization for better visualization. (D) RMSFs of p38α backbone in AMBER and OPLS simulations. The color scheme is same as in RMSD.
To show global structural excursions of apo p38α, we visually inspected trajectories of simulation runs and found several representative conformations (Figure 3B, pink ribbons) that were aligned to the crystal structure (Figure 3B, light blue ribbons). In AMBER apo3, we see a spike in the RMSD plot at ~1,790 ns (Figure 3A, left panel, green line), which attributes to the loss of helicity in αD helix and “closing” motion of the glycine-rich loop (Figure 3B, AMBER apo3, magenta circles; Figure S5). In OPLS apo1, RMSD rises rapidly to ~8 Å (Figure 3A, right panel, red line), which can be mainly attributable to the “closing” motion of the glycine-rich loop that is accompanied by large conformational changes in both the N-terminal and C-terminal domains (Figure 3B, OPLS apo1 magenta circle). In OPLS apo2, RMSD climbs up to ~4.5 Å at ~1,000 ns and descends to 3.5 Å at ~2,000 ns (Figure 3A, right panel, orange line), which identify the “closing” and “opening” of the glycine-rich loop (Figure 3B, OPLS apo2, magenta circle). In OPLS apo3, RMSD rises steadily to ~4 Å (Figure 3A, right panel, green line), which classifies the “closing” of the glycine-rich loop (Figure 3B, OPLS apo3, magenta circle). Among these structural excursions, OPLS apo1 seems to have the highest RMSD at 7.0 Å (Figure 3B, OPLS apo1, magenta circle), probably due to the global conformational change of both the N-terminal and C-terminal domains. It is also intriguing to see that in both OPLS apo1 and OPLS apo2, the glycine-rich loop is deformed significantly compared with the β-hairpin structure in the crystal structure (Figure 3B, OPLS apo1 and OPLS apo2).
Furthermore, to understand the local structural dynamics of p38α, residue-based root mean square fluctuations (RMSFs) were calculated and mapped onto the p38α structure (Figure 3C). In general, we found that local structural dynamics is highly dependent on the secondary structure of the local region. In α-helices- or β-sheets-rich local regions, low structural flexibility is observed, such as β-sheets 1–5 in the N-terminal domain (Figure S6) (34), and C-terminal domain that is made up of α-helices E, F, H, and G (Figure 3C, blue and cyan tubes). In local regions mainly consisting of disordered secondary structures such as loops and turns, high structural flexibility is observed (Figure 3C, orange and red tubes). In AMBER and OPLS simulations of p38α, high structural flexibility is observed mainly in 6 local regions (Figures 3C,D): (1) the glycine-rich loop (residue 30–38); (2) the L6 loop (residue 93–99); (3) the αD helix (residue 113–119); (4) the activation loop (residue 169–183); (5) the MAP kinase insert (residue 243–261); (6) the L16 loop (residue 305–330). Taken together, these data suggested p38α is highly dynamic in the entire kinase structure and in the specific secondary structures of local regions.
“Butterfly” Motions Potentially Contribute to p38α Enzymatic Catalysis
To gain further insights on the catalysis mechanism of p38α, we performed a principal component analysis (PCA) on p38α trajectories. For both AMBER and OPLS simulations, the trajectories of three apo p38α runs and one ligand-bound p38α run were combined to build a structure ensemble for PCA. The two most dominant PCs, PC1 and PC2, account for a large fraction of the correlated motions (S7A). Both AMBER and OPLS force fields produce considerable structural excursions (as big as 100 Å) from the crystal structures (Figures 4A,B). AMBER simulations show a more restricted structural excursion pattern (Figure 4A) while OPLS simulations have a more expanded pattern (Figure 4B). These results are in accordance with our previous findings (Figure 3A).
Figure 4. Principal component analysis. PC1-PC2 plots for AMBER simulations (A) and OPLS simulations (B). Data for apo1, apo2, apo3 and bound simulations are shown as red, orange, green, and blue dots. Data for the crystal structures of apo p38α (PDB code: 1P38) and SB203580-bound p38α (PDB code: 1A9U) are shown as white circle and triangle, respectively on the PC1-PC2 plot. Probability density functions are shown in upper and right panels. (C–F) show the two extreme structures along the PC1 or PC2 axis as blue and pink ribbons. Blue and pink indicate the – and + axis, respectively.
Structural morphing between the two extreme structures along the PC1 and PC2 axis was performed (Movies S1–S4). For PC1 and PC2, we overlaid the two extreme structures onto each other (Figures 4C–F). Note that the pink and blue structures represent the + and – axis. In AMBER simulations, PC1 represents the “butterfly” motion between the N-terminal and C-terminal domains (Figure 4C; Movie S1), where the N-terminal domain and the C-terminal domain resembling the wings of “butterfly”; PC2 represents the “twist” motion, where the N-terminal domain and the C-terminal domain rotates relative to each other (Figure 4D; Movie S2). In OPLS simulations, PC1 represents a large conformational change where both the N-terminal domain and the C-terminal domain move toward to each other, completely occluding the ATP-binding pocket (Figure 4E; Movie S3); PC2 represents a “butterfly” motion (Figure 4F; Movie S4). Notably, although the similarity of PCs in AMBER and OPLS simulations are limited (Figure S7B), it appears that both force fields have a reasonable agreement with the experimental crystal structures (Figures S7C,D).
Simulations Agree With Experimental NMR Observables
Chemical shifts have been important indicators of local backbone conformations (35). For verification of our results, we compared simulated and experimental chemical shifts of apo p38α. Linear regression of chemical shifts for five atom types in two representative simulation runs is performed (Figures 5A,B). The correlation coefficient (r2) suggest that there is excellent agreement for atom CB (r2: 0.98, 0.99) and atom CA (r2: 0.88, 0.89); reasonable agreement for atom C and N (r2: 0.42, 0.64); and bad agreement for atom H (r2: 0.19, 0.28). The correlation coefficient of the five atoms ranks similarly in other simulation runs (Figure S8). The top six outlier amino acid residues were highlighted (Figures 5A,B, red dots). To confirm the PPM results, we also performed chemical shift predictions using SHIFTX2. As expected, similar results were observed, indicating the accuracy of our simulation data (Figures S9–S11). Despite different structural flexibility in AMBER and OPLS simulations, there seems to be no significant difference between the accuracy of the chemical shifts. Furthermore, experimental residual dipolar couplings (RDCs) of apo p38α were applied for accuracy validation (Table 3). We compare simulated and experimental RDCs for 39 residues (Figure S12A) in apo p38α (36). The probability density function of Pearson's R (between simulated and experimental RDCs) indicated that MD simulations using both force fields agree generally well with RDC experiments (Figure S12B). Linear regressions of simulated and experimental RDCs show the same trend except OPLS apo1 (r2 is 0.00), possibly due to the surface charge distribution difference for all 85 charged residues (Table 4, Figure S13). The other five simulation runs agree generally well (r2 is 0.47, 0.34, 0.38, 0.35, 0.19 for AMBER apo1, AMBER apo2, AMBER apo3, OPLS apo2, OPLS apo3; Figures 5C,D). Taken together, these data suggest the high quality and accuracy of our simulations.
Figure 5. Comparison between simulated and experimental NMR observables. (A,B) Comparison between simulated and experimental chemical shifts of apo p38α for simulation runs AMBER apo1 and OPLS apo1. Chemical shifts were calculated from MD snapshots at an interval of 1 ns using PPM. The x-axis and y-axis of the correlation plot are the simulated and experimental chemical shifts, respectively. The linear regression lines are shown as gray, with slope a, intersection b, and correlation coefficient r2 shown in top left corner. The top six outlier residues with the largest RMS error are shown as red dots. (C,D) Comparison between simulated and experimental RDCs of apo p38α. Results are shown for all AMBER simulations and OPLS simulations of apo p38α. The x-axis and y-axis of the correlation plot are the simulated and experimental RDCs, respectively. The linear regression lines are shown as red, with slope a, intersection b, and correlation coefficient r2 shown in top left corner.
Table 3. Experimental 1H–15N residual dipolar couplings of the 39 residues in apo p38α.
Table 4. List of charged residues used in distance matrix calculations.
Dynamics of p38α Inhibitor in the ATP-Binding Pocket
The commonly used p38α inhibitor, SB203580, was suggested as a potential drug for cancer therapy (37). To study the drug-kinase dynamics of p38α, we measured the distances of two major interactions: (1) the hydrogen bond between the backbone of MET109 and SB203580; (2) the π-π stacking interactions between the phenyl sidechain of TYR35 and SB203580. In AMBER simulations, the hydrogen bond is well-maintained; the π-π stacking interaction is frequently disrupted with ring-ring distance up to 13 Å (Figure 6A). In OPLS simulations, although with minor fluctuations, both the hydrogen bond and the π-π stacking interaction are well-maintained (Figure 6A). SB203580 is less stable in AMBER simulations than in OPLS simulations, which is consistent with the wider range of RMSD fluctuations in AMBER simulations (Figure 3A, blue lines), suggesting that the ligand stability might affect the protein flexibility and vice versa.
Figure 6. p38α-SB203580 interaction and potential binding pockets in p38α. (A) Distance plots of hydrogen bond and stacking interactions between p38α and the inhibitor SB203580. The hydrogen bond distance is measured between backbone amide of MET109 and atom NB1 in the inhibitor (red lines). The stacking interaction distance is measured between the center of mass of the six-member ring in TYR35 and the inhibitor (blue lines). Probability density functions are shown in the side panel. (B) Cell viabilities of Pan02, RAW, MS5 cell lines treated with 3.125–50 μM SB203580 for 24 h (n = 6 samples per group). Cell viabilities of Pan02, RAW, MS5 cell lines treated with 1.25–10 μM ralimetinib for 24 h (n = 6 samples per group). (C) SB203580 and ralimetinib inhibit phosphorylation of ERK and AKT at 4 h and 30 min in Pan02 cells. GAPDH indicates the loading level in each lane. (n = 3 samples per group) (D) Violin plot of consensus cluster strength for potential ligand-binding pockets identified from AMBER (blue) and OPLS (red) simulations of apo p38α. Results with consensus cluster strength S ≥ 16 are plotted with quartiles shown as dashed lines inside the violins. Note that for pockets with very few occurrences in simulations, only single lines are displayed. (E) Representative snapshots showing the location of potential ligand-binding pockets identified in FTMap analysis of AMBER and OPLS simulations. p38α protein is shown in white ribbon representation and pockets are indicated by corresponding solvent probes shown in green surface representation. Note that the apo p38α crystal structure (PDB code: 1P38) is used here for consistency in visualization.
Current p38α Inhibitor Is Not Highly Selective
Although a variety of p38α inhibitors have been developed with enhanced specificity (38), most of these inhibitors are ATP competitors. Due to the similarity of the ATP-binding site of different kinases, off-target effects remain one of the biggest obstacles for the clinical application of p38α inhibitors. We tested the p38α inhibitors in PDAC cancer cells and various host cell lines. Indeed, both first and second generation of p38α inhibitors impedes cell viability of healthy host fibroblasts, and monocytes in a similar pattern as of cancer cells (Figure 6B). Surprisingly, in PDAC cells treated with two classical inhibitors, the phosphorylation of other kinases such as AKT, ERK is significantly altered under a relatively modest concentration (Figure 6C). Interestingly, the total amount of AKT was reduced under 100 μM SB203580 and 1 μM ralimetinib treatment, suggesting that an indirect regulatory pathway might be involved (Figure 6C). These findings suggest current p38α inhibitor is not selective enough. New approaches to inhibit the p38α needs to be explored for PDAC treatment.
Potential Drug Binding Sites in p38α
To explore potential ligand-binding pockets in p38α, computational solvent mapping was performed using FTMap (39, 40). Pockets with consensus cluster strength S ≥ 16 are considered druggable. Using these criteria, 19 druggable pockets were found (Table 5): αD-L13 pocket, L4-L7 pocket, P-αC pocket, αG-L14 pocket, αE-L16 pocket, αE-L16b pocket, P-L12 pocket, L12 pocket, αE-β7 pocket, αC-αL16 pocket, β5-αL16 pocket, P-L16 pocket, MKI pocket, ATP pocket, αE-αF pocket, αH-MKI pocket, αF-αG pocket, β2-L4 pocket, αE-L4 pocket (Figure 6D). These pockets have a varied frequency of occurrences in simulations (Table 6) and are scattered in the entire p38α structure (Figure 6E). For comparison, FTMap analysis performed on crystal structure ensemble consisting of 196 crystal structures of p38α identified 6 druggable pockets: P-L16 pocket, ATP pocket, L4–L7 pocket, αG-L14 pocket, β5-αL16 pocket, MKI pocket, which partially verified our simulation results. Besides the well-studied ATP pocket and MKI pocket, the remaining 17 pockets are novel pockets that may be explored for further cancer drug development. Interestingly, a comparison between apo and ligand-bound simulations indicates that the presence of the ligand seems to prevent the exposure of a few binding pockets in AMBER (Figure S14A) and OPLS simulations (Figure S14B). Taken together, we provide potential new drug-binding sites that may pave the way for new generation drug design targeting p38α.
Table 5. Details of pocket lining residues of potential ligand-binding pockets identified from FTMap analysis.
Table 6. Number of occurrences of potential ligand-binding pockets in MD simulations.
Discussion
The role of p38 MAPKs in cancer is still under intense investigation. What is the role of p38 MAPK family members in cancer? The one-word answer is “context-dependent.” Some human tumors, such as HCCs, have lower p38 MAPK activity than non-tumorigenic tissues (41). Several reports showed that p38α is a tumor suppressor. In support of that view, negative regulators of p38α, such as the phosphatases PPM1D and DUSP26, are overexpressed in human tumors (42, 43). However, bearing the evidence of p38α signaling in tumor suppression, mutations of p38α have not been consistently identified in human tumors. It is suspected that cancer cells may benefit from the p38α signaling pathway. In support of that view, p38α blockades show significant tumor suppressing effects in vivo in various cancer types (44, 45), suggesting the dual roles of p38α signaling in cancer. Here in our work, a thorough screening for more than 30 types of tumors and their counterparts shows that p38α expression is significantly increased in PDAC tissues, which drives our curiosity. We originally hypothesized that p38α may also be activated in PDAC cells and p38α blockades may benefit PDAC patients. Our experimental data support this hypothesis, and SB203580 significantly inhibit human and mouse PDAC cell growth. Moreover, in our recently unpublished data, SB203580 reduces the activity of pancreatic tumor derived-macrophages in in vivo models. Our data suggest that p38α blockades may block various cell components in PDAC tumor microenvironment, and may serve as a potential target for PDAC therapy.
Clinically, obesity is associated with cancer risk (46), and the adipose tissue microenvironment supports cancer development, metastasis (47), and drug resistance (29). For tumors predominantly occur in the adipocyte-rich microenvironments such as breast, prostate, ovarian, colon, and pancreatic cancers, the degree of adipose tissue involvement is often correlated with poor prognosis (48). Our data shows that MAPK14 is correlated with PLIN2 and PLIN3, but not PLIN1, 4, and 5. Notably, PLIN2 and 3 is the marker for small lipid droplets. These results may reflect that cancer cells induce lipolysis of surrounding adipocyte, which is previously reported in breast cancer (49).
After establishing the role of p38α in PDAC cells, we turn to computational biology tools for more detailed information on the atomistic level for this interesting target. We performed molecular dynamics study on p38α conformational dynamics on the state of the art, highly parallel supercomputer, Anton. Supercomputer-powered biomolecular simulation provides us abundant information related to p38α function. Initially, a powerful and direct metrics, all-too-all RMSDs were calculated to examine the convergence of MD simulations. We observed the frequent appearance of low RMSD stripes in both AMBER and OPLS simulations (Figure S3A, dark blue areas), indicating a reasonable convergence. The results show that across different simulation runs, the conformational similarity is quite low, especially for OPLS apo1 (6.3 μs), which is significantly different from all other simulation runs (Figure S3B, mid panel, lower panel, red areas), possibly due to large conformational change of the glycine-rich loop in OPLS apo1 simulation (Figure S5). We found 11 p38α crystal structures with glycine-rich loop conformations that at least partially resemble the unfolded conformation found in MD simulations, suggesting that the novel conformation of the glycine-rich loop might be realistic.
The high flexibility we found in local regions of apo p38α is generally in line with previous simulation studies. Kuzmanic et al., performed metadynamics simulations of p38α to find that the glycine-rich loop, the activation loop, and the L16 loop are highly flexible (50). McClendon et al. carried out multiple microsecond-scale MD simulations of protein kinase A (PKA) and observed the opening and closing motions of the glycine-rich loop (51). Shan et al. performed long MD simulations of EGFR kinase to find that it is intrinsically disordered in the glycine-rich loop and the activation loop (52). Kumar et al. performed long MD simulations of Aurora-A kinase and identified several highly flexible regions including the glycine-rich loop and the activation loop (53). These data support that the conformational flexibility of local regions might be a common feature in protein kinases, presumably conferring adaptability in enzyme catalysis and protein-protein interactions.
The “butterfly” motion identified in PCA analysis is a quite interesting protein dynamics motion. Similar motions were seen in simulations of other enzymes as well such as the adenylate kinase (54). Many two-lobe enzymes may show similar type of motions during catalysis. Meanwhile, it is possible that both “conformational selection” and “induced fit” (55–57) contribute to the ligand-binding because the PC1-PC2 subspace of apo and ligand-bound p38α has both overlapping and non-overlapping regions. Nevertheless, the PCA result is complicated by two factors. Firstly, the initial structure of apo and SB203580-bound p38α simulations are similar (with RMSD of only 0.5 Å). Therefore, it is highly likely to have structural overlaps. Secondly, the apo and SB203580-bound p38α crystal structures come from mouse and human, respectively. The two structures differ at residue 48 and 263 (HIS and ALA for the mouse form; LEU and THR for the human form), which may be a potential source of error.
Despite the different conformational changes, linear regression of simulated and experimental chemical shifts shows similarly good agreement in AMBER and OPLS simulations, indicating that chemical shift is not sensitive to conformational changes. In contrast, linear regression of simulated and experimental RDCs shows a lower correlation coefficient in OPLS apo1 than other simulation runs, suggesting that RDC is sensitive to conformational changes. Large conformational changes of both N-terminal and C-terminal domains in OPLS apo1 may also contribute to the relatively worse RDC agreement. Since accurate RDC predictions are dependent on the precise distribution of the surfaces charges (58), we calculated all-to-all distance matrix for all 85 charged residues (Table 4) in p38α including ARG, LYS, ASP, and GLU for the simulation run AMBER apo1, AMBER apo2, and OPLS apo1. The charged residue distance matrices have a much bigger discrepancy between AMBER apo1 and AMBER apo2 than between AMBER apo1 and OPLS apo1 (Figure S13), suggesting that surface charge distribution can be an important contributing factor in RDC predictions.
Newly discovered potential druggable pockets in this study may benefit millions of PDAC patients. Notably, ATP pocket, MAP kinase insert (MKI) pocket, and αC-αL16 pocket have been experimentally verified (13, 59). Several novel pockets are in proximity to the binding site of upstream activators or downstream substrates of p38α and can potentially be explored to modulate its activity. There are at least four pockets that can be potentially used for this purpose: (1) αD-L13 pocket, which can be used to block docking of kinase substrates; (2) αC-αL16 pocket, which can be used to lock p38α in inactivate conformation; (3) αG-L14 pocket, which can be used to block p38α-TAB1 interaction (60); (4) αE-β7 pocket, which can be used to block p38α-MK2 interaction (61) and p38α-MKP5 interaction (62). Interestingly, the binding of the inhibitor seems to veil many potential binding pockets in AMBER and OPLS simulations (Figure S14). Further, virtual screenings and experimental validations are needed to confirm the druggability of these novel pockets.
Taken together, our study provides an interesting potential target for combating PDAC and offers detailed conformational dynamics information on p38α protein. Our findings have also paved potential avenues for developing the new classes of p38α-targeting drug, which may overcome the side effects of current therapy and benefit PDAC patients.
Materials and Methods
Cell Culture
Murine Panc02 pancreatic ductal adenocarcinoma cell line was kindly provided by Dr. Maximilian Schnurr from Munich University, Germany. Murine monocyte/macrophage-like cell line RAW264.7 and human monocyte/macrophage-like cell line THP-1 were kindly provided by Dr. Dapeng Yan at the Fudan University, China. Human TERT-immortalized fibroblasts Whirly, human PDAC cell lines MiaPaCa-2, Murine PDAC cell Pan02, and Human umbilical vein endothelial cell HUVEC were kindly provided by Dr. Yihai Cao from the Karolinska Institutet, Sweden. Murine fibroblast cell line MS5 were purchased from the ATCC. Panc02, RAW264.7, MS5, Whirly, MiaPaCa-2 were cultured in 10% FBS-DMEM (Cat. No. TBD10569, TBD, China) containing 100 U/ml penicillin, 100 μg/ml streptomycin (Cat. No. MA0110, Meilunbio, China). THP-1 cells were cultured in 10 % FBS-RPMI1640 (Cat. No. HY1640, TBD, China) containing 100 U/ml penicillin, 100 μg/ml streptomycin. HUVEC cells were cultured in 10% FBS-Medium 199/EBSS (Cat. No. SH30253.01, HyClone) containing 100 U/ml penicillin, 100 μg/ml streptomycin. All cell lines used in our study were negative for mycoplasma (Cat. No. LT07-318; Lonsa).
Human Patient Samples
All studies related to clinical human samples were approved by the Ethical Review Committee in Shuguang Hospital, Shanghai University of Traditional Chinese Medicine, Shanghai, China. Pancreatic tumor samples and adjacent pancreatic tissues were collected from cancer patients with written informed permission.
Database Analysis
Transcriptome data from patient samples of pancreatic cancer were analyzed using the online database, The Cancer Genome Atlas (TCGA, https://www.cancer.gov/about-nci/organization/ccg/research/structural-genomics/tcga) to investigate whether the expression of interesting markers is altered in tumor tissue. RNA sequencing analysis and visualization platform Gene Expression Profiling Interactive Analysis (GEPIA, http://gepia.cancer-pku.cn/) (63) were used to perform the correlation analysis and survival analysis in pancreatic cancer cohort. A gene expression profile across various cancer types and paired normal samples was generated from GEPIA.
Immunoblot
Fresh tumor tissues and cultured cells were lysed in RIPA buffer and the proteinase and phosphatase inhibitor cocktail (Cat. No. MA0151, Meilunbio, China; Cat. No. MB2678, Meilunbio, China; 1:100). For immunoblot, each protein sample and a standard molecular weight marker (Cat. No. WJ102, EpiZyme, China) were loaded onto a 10% SDS-PAGE gel (Cat. No. PG112, EpiZyme, China), followed by transferring onto a polyvinylidene difluoride (PVDF) membrane (Cat. No. IPVH00010, Millipore). The membranes were blocked with 5% skimmed milk for 2 h, and were probed overnight at 4°C with a rabbit anti-phospho-p38 antibody (Cat. No. 4631S, Cell Signaling; 1:1,000), a rabbit anti-p38 antibody (Cat. No. 9212S, Cell Signaling; 1:1,000), a rabbit anti-phospho-Erk1/2 antibody (Cat. No. 9101S, Cell Signaling; 1:1,000), a rabbit anti-Erk1/2 antibody (Cat. No. 4695S, Cell Signaling; 1:1,000), a rabbit anti-phospho-AKT antibody (Cat. No. 2118-1, Epitomics; 1:1,000), a rabbit anti-AKT antibody (Cat. No. 1085-1, Epitomics; 1:1,000) and a mouse anti-GAPDH antibody (Cat. No. A01020, Abbkine; 1:1,000) in 5% skimmed milk. After rigorous washing with PBS containing 0.1% Tween-20 (Cat. No. T8220, Solarbio, China), membranes were incubated with a goat anti-mouse HRP-conjugated IgG antibody (Cat. No. AS003, ABclonal; 1:5,000) and a goat anti-rabbit HRP-conjugated IgG antibody (Cat. No. AS014, ABclonal; 1:5,000). Target proteins were visualized using a super-sensitive ECL luminescence reagent (Cat. No. MA0186, Meilunbio, China) with a Molecular Imager ChemiDoc XRS System (Bio-Rad).
RNA Extraction and Quantitative Real-Time PCR
Total RNAs were extracted from tumor tissues and cultured cells using RNAsimple Total RNA kit (Cat. No. DP419, TIANGEN, China). Total RNA from each sample was reversely transcribed using an All-in-One cDNA Synthesis SuperMix (Cat. No. B24408, Bimake, China). Reverse transcription was performed at 42°C for 60 min, followed by 70°C for 5 min to inactivate the enzyme activity. cDNA samples were stored at −20°C and subjected to qPCR using a StepOnePlus Real-Time PCR System (Applied Biosystems). Each qPCR sample was performed in a triplicate and 10 μl reaction containing 2 × SYBR Green qPCR Master Mix (Cat. No. B21202, Bimake, China), 50 nM forward and reverse primers, and 4 μl cDNA. The qPCR protocol was executed for 45 cycles and each cycle consisted of denaturation at 95°C for 15 s, annealing at 60°C for 1 min, and extension at 72°C for 1 min. The primer pairs used in our experiments included: mouse Gapdh forward: 5′-CCAGCAAGGACACTGAGCAA-3′; mouse Gapdh reverse: 5′-GGGATGGAAATTGTGAGGGA-3′; mouse Mapk14 forward: 5′- GGGACACCCCCTGCTTATCT-3′; mouse Mapk14 reverse: 5′- TCCCTGCTTTCAAAGGACTGG-3′; mouse Mapk11 forward: 5′- GGACCTGAACAGGATCGTGAA-3′; mouse Mapk11 reverse: 5′- CTCAGCCATGAAGCCTCCC-3′; mouse Mapk12 forward: 5′- CGCCGTGTACCAAGACCTG-3′; mouse Mapk12 reverse: 5′- GAGGCGCAACTCTCTGTAGG-3′; mouse Mapk13 forward: 5′- GAGGCGCAACTCTCTGTAGG-3′; mouse Mapk13 reverse: 5′- CACTCAGGGTCTCATGCTTCA-3′; human GAPDH forward: 5′-AGGGCTGCTTTTAACTCTGGT-3′; human GAPDH reverse: 5′-CCCCACTTGATTTTGGAGGGA-3′; human MAPK14 forward: 5′- CCCGAGCGTTACCAGAACC-3′; human MAPK14 reverse: 5′- TCGCATGAATGATGGACTGAAAT-3′; human MAPK11 forward: 5′-AAGCACGAGAACGTCATCGG-3′; human MAPK11 reverse: 5′- TCACCAAGTACACTTCGCTGA-3′; human MAPK12 forward: 5′- CCCAGACATCAGGGAGTAATGG-3′; human MAPK12 reverse: 5′- TCTATCGGATACTTCAGCGTCA-3′; human MAPK13 forward: 5′- CACTCAGGGTCTCATGCTTCA-3′; human MAPK13 reverse: 5′- GCTTGCGTTGGTCAGGACA-3′.
Histology and Immunohistochemistry
Paraffin-embedded tissues were cut in the thickness of 5 μm, mounted onto glass slides, baked for 1 h at 60°C, deparaffinized in Xylene (Cat. No. 10023418, SCR, China), and sequentially rehydrated in 99, 95, and 70% ethanol (Cat. No. 10009218, SCR, China). Tissue slides were counterstained with Haematoxylin (Mayer's) (Cat. No. MB9897, Meilunbio, China) and Eosin (Cat. No. MA0164, Meilunbio, China) before dehydration with 95 and 99% ethanol, and were mounted with neutral balsam (Cat. No. 1004160, SCR, China). Stained tissues were analyzed under a light microscope (Leica DM IL LED). For immunohistochemical staining, paraffin-embedded tissue sections were stained with a rabbit anti-PDGFRβ antibody (Cat. No. ab32570, Abcam; 1:200); a mouse anti-CD163 antibody (Cat. No. ab156769, Abcam; 1:200); a rabbit anti-αSMA antibody (Cat. No. ab32575, Abcam; 1:200); a rabbit anti-phospho-p38 antibody (Cat. No. 4631S, Cell Signaling; 1:100), and a rabbit anti-p38 antibody (Cat. No. 9212S, Cell Signaling; 1:100). A ready to use HP IHC detection kit (Cat. No. abs957, Absin, China) was used for visualization. Captured images were further analyzed using Adobe Photoshop CS software.
Cell Viability
Panc02, RAW, and MS5 cells at the density of 1 × 104 cells per well were seeded in a 96-well plate. After 24 h of incubation, cells were treated with 3.125–50 μM SB203580 SB203580 (Cat.No.S1076, selleck,USA) and 1.25–10 μM ralimetinib (Cat. No.S1494, selleck,USA). Cell viability was measured using a Cell Counting Kit-8 (Cat. No. MA0218-5, Meilunbio, China). In brief, 10 μl of the solution containing 2-(2-methoxy-4-nitrophenyl)-3-(4-nitrophenyl)-5-(2,4-disulfophenyl)-2H-tetrazolium (WST-8) was added to each well of a 96-well plate, followed by incubation for 4 h. Densitometry was measured at the wavelength of 450 nm by a Synergy 2 Multi-Mode Microplate Reader (BioTek).
Simulation System Preparation
p38α protein structures were prepared using Maestro 9.0 (Schrödinger, LLC). The crystal structures of apo p38α (PDB code: 1P38) (9) and SB203580-bound (PDB code: 1A9U) p38α (10) were used as initial structures. The protein was capped with acetyl group (ACE) at the N-terminus and N-methyl group (NME) at the C-terminus to prevent artificial charge-charge interactions. AMBERff99SB-ILDN (64) and OPLS-AA/L (65) force fields were used for the protein. A box size of 85.0 × 85.0 × 85.0 Å3 was used to ensure the enclosure of the entire protein. Apo p38α or SB203580-bound p38α system was solvated in water described by the TIP3P model (66) to mimic physiological condition and neutralized using nine sodium ions (67).
Ligand Parameterization
To derive parameters for the SB203580, we used different approaches. For AMBER simulation, the bonded parameters and non-bonded parameters for SB203580 were parameterized using General Amber Force Field (GAFF) (68). The partial charges were derived using the restrained electrostatic potential (RESP) method on the RED server (69, 70). For OPLS simulations, the bonded parameters, non-bonded parameters, and partial charges were derived by analogy method as previously described (71). In both simulations of ligand-bound p38α, dihedral restraints were applied to eight dihedral angles to keep the ligand conformations close to the crystal structures. The force constant for dihedral restraints is 4184 kJ/mol.rad. Each equilibrium angle (Table 1) was obtained from four SB203580-bound p38α crystal structures [PDB codes: 1A9U (10), 1PME (72), 2EWA (17), and 3GCP (73)].
Molecular Dynamics Simulations
Before the production simulation, energy minimization, and equilibration were performed on our in-house clusters using Desmond 2.4 (74). Energy minimization was initially carried out using the steepest descent algorithm for 1,000 steps and the L-BFGS algorithm for 1,000 steps. The energy-minimized structure was incrementally heated up from 50 to 298 K (50, 100, 150, 200, 250, 298 K) with 100 ps at each temperature in the NVT ensemble using Berendsen thermostat (75). Next, the protein was equilibrated at 298 K for an additional 10 ns in the NVT ensemble using the Nosé-Hoover thermostat (76, 77). For van der Waals interactions and short-range electrostatic interactions, a 9 Å cutoff was used; for long-range electrostatic interactions, Particle Mesh Ewald (PME) was used (78). The “bonded,” “near,” and “far” time steps are 1, 1, and 3 fs, respectively.
Anton Supercomputer Molecular Dynamics Simulations
For production simulations, we performed simulations on the Anton supercomputer (Pittsburgh supercomputing center, Pittsburgh, PA) (33). Production simulations were performed in the NVT ensemble with the temperature maintained at 298 K using Nosé-Hoover thermostat (76, 77). For van der Waals interactions and short-range electrostatic interactions, the cutoff was automatically determined by the “guess_chem” utility on Anton for optimized performance. Gaussian Split Ewald (GSE) (79) method was used for long-range electrostatic interactions. M-SHAKE algorithm (80) was used to constrain all bonds, enabling a 2 fs time step. Snapshots of the simulations were saved at an interval of 0.1 ns for analysis. Visualizations of simulation snapshots and trajectories were done using both PyMOL 1.8 (PyMOL 2015) and VMD 1.9 (81). For AMBER and OPLS simulations, three independent simulations of apo p38α and one simulation of SB203580-bound p38α were performed with detailed parameters shown in Table 2.
Calculations of RMSD and RMSF
To obtain root mean square deviations (RMSDs) in p38α simulations, the trajectory (non-terminal residues, residue 14–344) at an interval of 0.1 ns was fitted to the initial structure on protein backbone; RMSDs were calculated for each snapshot over simulation time. To obtain all-to-all RMSDs, the trajectory at an interval of 1 ns was fitted to the initial structure on the protein backbone before all-to-all RMSD matrix was generated using pytraj (82), a Python package binding to cpptraj program (83). To obtain root mean square fluctuations (RMSFs) in p38α simulations, the trajectory at an interval of 0.1 ns was fitted to the initial structure on protein backbone; RMSFs were then calculated for each residue.
Principal Component Analysis
Principal components analysis (PCA) was performed for a combined trajectory of apo and ligand-bound p38α for each force field. The input trajectory was initially aligned on Cα atoms of p38α and used to extract eigenvectors and eigenvalues. Then, the trajectory of Cα atoms was projected onto the two most dominant eigenvectors to obtain principal components PC1 and PC2. For experimental comparison, PCA of 44 full-length crystal structures of p38α was performed with PDB entries as follows: 5LAR, 5ETI, 5ETC, 4LOO, 3U8W, 3S3I, 3RIN, 3PY3, 3NNW, 3NNV, 3NNU, 3KQ7, 3ITZ, 3GFE, 3GC7, 3DT1, 3D83, 3D7Z, 2ZB0, 2YIX, 2I0H, 1ZZL, 1YQJ, 1WBW, 1WBV, 1WBT, 1WBS, 1WBO, 1WBN, 1W84, 1W83, 1W82, 1W7H, 1R3C, 1P38, 1OVE, 1OUY, 1OUK, 1M7Q, 1DI9, 1BMK, 1BL7, 1BL6, 1A9U. All PCA calculations were performed using the analysis tools g_covar and g_anaeig in Gromacs 5.0 (84, 85).
Calculations of Chemical Shifts and Residual Dipolar Couplings
To obtain chemical shifts of the backbone atoms in p38α from MD simulations, we performed ensemble-based chemical shift calculations on MD trajectory at an interval of 1 ns using the PPM program (86). PPM predicts the chemical shifts of CA, CB, C, N, and H based solely on the physical and chemical properties of the protein structure. For comparison, we also performed chemical shift calculations on MD trajectory at an interval of 1 ns using the SHIFTX2 program (87). SHIFTX2 predicts the chemical shifts of CA, CB, C, N, and H using both structure- and sequence- based criteria. Experimental p38α chemical shifts were obtained from Biological Magnetic Resonance data bank (BMRB entry number: 6468) (88); the entry includes 219 1H chemical shifts, 684 13C (including CA, CB, and C) chemical shifts, and 219 15N chemical shifts for apo p38α. For comparison of simulated and experimental chemical shifts of atom CA, CB, C, N, and H, different residues were used depending on the availability of experimental chemical shifts for specific atoms in BMRB entry 6,468. To obtain residual dipolar couplings (RDCs) of each residue in p38α, we performed RDC calculations on structural snapshots at an interval of 1 ns using the PALES program (58, 89). Experimental RDC data for 39 residues of apo p38α and the parameters in RDC (Table 3) were obtained from previous report (36). Specifically, a concentration of bacteriophage Pf1 at 20 mg/mL and a pH of 6.0 was used. The sodium chloride concentration was optimized using a range from 0.005 to 0.5 M and 0.2 M was used for final RDC calculations.
FTMap Analysis
A standalone version of FTMap (39, 40) was used to perform calculations on MD snapshots at an interval of 1 ns. FTMap analysis mainly consists of three steps: (1) independent docking of 16 organic probes on the protein using a fast Fourier transform correlation approach; (2) refinement of the probe positions using energy minimizations; (3) clustering and ranking of the resulting poses to identify consensus sites (CSs). The 16 solvent probes used in FTMap algorithm are: acetamide, acetonitrile, acetone, acetaldehyde, methylamine, benzaldehyde, benzene, isobutanol, cyclohexane, N, N-dimethylformamide, dimethyl ether, ethanol, ethane, phenol, isopropanol, urea. Consensus cluster strength (S) is defined as the number of probe clusters in a specific binding pocket on the protein. For experimental comparison, 196 crystal structures were also used for FTMap analysis. Their PDB codes are: 1A9U, 1BL6, 1BL7, 1BMK, 1DI9, 1IAN, 1KV1, 1KV2, 1LEW, 1LEZ, 1M7Q, 1OUK, 1OUY, 1OVE, 1OZ1, 1P38, 1R39, 1R3C, 1W7H, 1W82, 1W83, 1W84, 1WBN, 1WBO, 1WBS, 1WBT, 1WBV, 1WBW, 1WFC, 1YQJ, 1YW2, 1YWR, 1ZYJ, 1ZZ2, 1ZZL, 2BAK, 2FSL, 2FSM, 2FSO, 2FST, 2GFS, 2GHL, 2GHM, 2GTM, 2GTN, 2I0H, 2NPQ, 2PUU, 2QD9, 2RG5, 2RG6, 2Y8O, 2YIS, 2YIW, 2YIX, 2ZAZ, 2ZB0, 2ZB1, 3BV2, 3BV3, 3BX5, 3CTQ, 3D7Z, 3DS6, 3DT1, 3FC1, 3FI4, 3FKL, 3FKN, 3FKO, 3FLN, 3FLQ, 3FLS, 3FLW, 3FLY, 3FLZ, 3FMH, 3FMJ, 3FMK, 3FML, 3FMM, 3FMN, 3FSF, 3FSK, 3GC7, 3GCP, 3GCQ, 3GCS, 3GCV, 3GFE, 3GI3, 3HA8, 3HEC, 3HEG, 3HL7, 3HLL, 3HP2, 3HP5, 3HRB, 3HUB, 3HUC, 3HV3, 3HV4, 3HV5, 3HV6, 3HV7, 3HVC, 3ITZ, 3IW5, 3IW6, 3IW7, 3IW8, 3K3I, 3K3J, 3KF7, 3KQ7, 3L8S, 3LFA, 3LFB, 3LFC, 3LFE, 3LFF, 3MGY, 3MH0, 3MH1, 3MH2, 3MH3, 3MPA, 3MPT, 3MVL, 3MVM, 3MW1, 3NEW, 3NNU, 3NNV, 3NNW, 3NNX, 3NWW, 3O8P, 3O8T, 3O8U, 3OBG, 3OBJ, 3OC1, 3OCG, 3OD6, 3ODZ, 3P4K, 3P5K, 3P78, 3P79, 3P7A, 3P7B, 3P7C, 3PG3, 3PY3, 3QUD, 3QUE, 3RIN, 3ROC, 3S3I, 3U8W, 3UVP, 3UVQ, 3ZS5, 3ZSG, 3ZSH, 3ZYA, 4A9Y, 4AA0, 4AA4, 4AAC, 4DLI, 4DLJ, 4E5A, 4E5B, 4E6A, 4E6C, 4E8A, 4EH2, 4EH3, 4EH4, 4EH5, 4EH6, 4EH7, 4EH8, 4EH9, 4EHV, 4EWQ, 4F9W, 4F9Y, 4FA2, 4GEO, 4KIN, 4KIP, 4KIQ.
Classification of Binding Pockets
Binding pockets are named and classified based on its adjacent secondary structures, which harbor pocket lining residues (Table 5). To determine the identity of a pocket, the intersecting volume of the bounding sphere for the probe clusters and the pocket lining residues were calculated for each pocket using the equation as described previously (90):
Vo is the intersecting volume; r1, r2 are the radii of the bounding sphere of the probe cluster and of the pocket lining residues, respectively; d is the distance between the two sphere centers. A specific probe cluster falls in a pocket category if it has the maximum intersecting volume with that pocket. A scaling factor of 0.75 on r2 was used to avoid the overestimation of the bound site volume. All pockets and subpockets in the ATP-binding site are classified as ATP pocket.
Statistical Analysis
Statistical analyses of biological experiments were performed using the standard two-tailed Student t-test, and p < 0.05 was considered statistically significant.
Data Availability Statement
The datasets analyzed in this study can be found in The Cancer Genome Atlas (https://portal.gdc.cancer.gov/).
Ethics Statement
The studies involving human participants were reviewed and approved by Ethics committee, Shuguang Hospital, Shanghai University of Traditional Chinese Medicine. The patients/participants provided their written informed consent to participate in this study.
Author Contributions
SZ generated the ideas, designed experiments, and wrote the manuscript. LY, XS, YY, and SZ performed most experiments and organized all figures. YL, JZ, WL, and AE participated in discussions.
Funding
This work was supported by a grant from the National Natural Science Foundation of China (Grant No. 81501096) awarded to SZ. This work was also supported by a grant from the National Institutes of Health (R01GM068648) awarded to AE. YY was supported by National Natural Science Foundation of China (81600839) and Shanghai Municipal Health Commission (20164Y0222).
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Acknowledgments
Anton computer was generously made available by D. E. Shaw Research hosted by the National Institutes of Health and the Pittsburgh Supercomputing Center (PSC). We thank David Hall for kindly providing a standalone version of FTMap.
Supplementary Material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fonc.2019.01294/full#supplementary-material
References
1. Nebreda AR, Porras A. p38 MAP kinases: beyond the stress response. Trends Biochem Sci. (2000) 25:257–60. doi: 10.1016/S0968-0004(00)01595-4
2. Kyriakis JM, Avruch J. Mammalian mitogen-activated protein kinase signal transduction pathways activated by stress and inflammation. Physiol Rev. (2001) 81:807–69. doi: 10.1152/physrev.2001.81.2.807
3. Cuenda A, Rousseau S. p38 MAP-kinases pathway regulation, function and role in human diseases. Biochim Biophys Acta. (2007) 1773:1358–75. doi: 10.1016/j.bbamcr.2007.03.010
4. Dhillon AS, Hagan S, Rath O, Kolch W. MAP kinase signalling pathways in cancer. Oncogene. (2007) 26:3279–90. doi: 10.1038/sj.onc.1210421
5. Kumar S, Boehm J, Lee JC. p38 MAP kinases: key signalling molecules as therapeutic targets for inflammatory diseases. Nat Rev Drug Discov. (2003) 2:717–26. doi: 10.1038/nrd1177
6. Zarubin T, Han J. Activation and signaling of the p38 MAP kinase pathway. Cell Res. (2005) 15:11–8. doi: 10.1038/sj.cr.7290257
7. Olson JM, Hallahan AR. p38 MAP kinase: a convergence point in cancer therapy. Trends Mol Med. (2004) 10:125–9. doi: 10.1016/j.molmed.2004.01.007
8. Kaminska B. MAPK signalling pathways as molecular targets for anti-inflammatory therapy–from molecular mechanisms to therapeutic benefits. Biochim Biophys Acta. (2005) 1754:253–62. doi: 10.1016/j.bbapap.2005.08.017
9. Wang Z, Harkins PC, Ulevitch RJ, Han J, Cobb MH, Goldsmith EJ. The structure of mitogen-activated protein kinase p38 at 2.1-A resolution. Structure. (1997) 94:2327–32. doi: 10.1073/pnas.94.6.2327
10. Wang Z, Canagarajah BJ, Boehm JC, Kassisà S, Cobb MH, Young PR, et al. Structural basis of inhibitor selectivity in MAP kinases. Structure. (1998) 6:1117–28. doi: 10.1016/S0969-2126(98)00113-0
11. Das J, Moquin RV, Pitt S, Zhang R, Shen DR, McIntyre KW, et al. Pyrazolo-pyrimidines: a novel heterocyclic scaffold for potent and selective p38 alpha inhibitors. Bioorg Med Chem Lett. (2008) 18:2652–7. doi: 10.1016/j.bmcl.2008.03.019
12. Das J, Moquin RV, Dyckman AJ, Li T, Pitt S, Zhang R, et al. 5-amino-pyrazoles as potent and selective p38alpha inhibitors. Bioorg Med Chem Lett. (2010) 20:6886–9. doi: 10.1016/j.bmcl.2010.10.034
13. Diskin R, Engelberg D, Livnah O. A novel lipid binding site formed by the MAP kinase insert in p38 alpha. J Mol Biol. (2008) 375:70–9. doi: 10.1016/j.jmb.2007.09.002
14. Sours KM, Kwok SC, Rachidi T, Lee T, Ring A, Hoofnagle AN, et al. Hydrogen-exchange mass spectrometry reveals activation-induced changes in the conformational mobility of p38alpha MAP kinase. J Mol Biol. (2008) 379:1075–93. doi: 10.1016/j.jmb.2008.04.044
15. Akella R, Min X, Wu Q, Gardner KH, Goldsmith EJ. The third conformation of p38alpha MAP kinase observed in phosphorylated p38alpha and in solution. Structure. (2010) 18:1571–8. doi: 10.1016/j.str.2010.09.015
16. Honndorf VS, Coudevylle N, Laufer S, Becker S, Griesinger C. Dynamics in the p38alpha MAP kinase-SB203580 complex observed by liquid-state NMR spectroscopy. Angew Chem Int Ed Engl. (2008) 47:3548–51. doi: 10.1002/anie.200705614
17. Vogtherr M, Saxena K, Hoelder S, Grimme S, Betz M, Schieborr U, et al. NMR characterization of kinase p38 dynamics in free and ligand-bound forms. Angew Chem Int Ed Engl. (2006) 45:993–7. doi: 10.1002/anie.200502770
18. Fabian MA, Biggs WH, Treiber DK, Atteridge CE, Azimioara MD, Benedetti MG, et al. A small molecule-kinase interaction map for clinical kinase inhibitors. Nat Biotechnol. (2005) 23:329–36. doi: 10.1038/nbt1068
19. Pargellis C, Tong L, Churchill L, Cirillo PF, Gilmore T, Graham AG, et al. Inhibition of p38 MAP kinase by utilizing a novel allosteric binding site. Nat Struct Biol. (2002) 9:268–72. doi: 10.1038/nsb770
20. Perry JJ, Harris RM, Moiani D, Olson AJ, Tainer JA. p38alpha MAP kinase C-terminal domain binding pocket characterized by crystallographic and computational analyses. J Mol Biol. (2009) 391:1–11. doi: 10.1016/j.jmb.2009.06.005
21. Yong HY, Koh MS, Moon A. The p38 MAPK inhibitors for the treatment of inflammatory diseases and cancer. Expert Opin Investig Drugs. (2009) 18:1893–905. doi: 10.1517/13543780903321490
22. Nemoto T, Kubota S, Ishida H, Murata N, Hashimoto D. Ornithine decarboxylase, mitogen-activated protein kinase and matrix metalloproteinase-2 expressions in human colon tumors. World J Gastroenterol. (2005) 11:3065–9. doi: 10.3748/wjg.v11.i20.3065
23. Elenitoba-Johnson KS, Jenson SD, Abbott RT, Palais RA, Bohling SD, Lin Z, et al. Involvement of multiple signaling pathways in follicular lymphoma transformation: p38-mitogen-activated protein kinase as a target for therapy. Proc Natl Acad Sci USA. (2003) 100:7259–64. doi: 10.1073/pnas.1137463100
24. Zhou HY, Pon YL, Wong AS. Synergistic effects of epidermal growth factor and hepatocyte growth factor on human ovarian cancer cell invasion and migration: role of extracellular signal-regulated kinase 1/2 and p38 mitogen-activated protein kinase. Endocrinology. (2007) 148:5195–208. doi: 10.1210/en.2007-0361
25. Wang L, Bai YY, Yang Y, Hu F, Wang Y, Yu Z, et al. Diabetes mellitus stimulates pancreatic cancer growth and epithelial-mesenchymal transition-mediated metastasis via a p38 MAPK pathway. Oncotarget. (2016) 7:38539–50. doi: 10.18632/oncotarget.9533
26. Noel JK, Crean S, Claflin JE, Ranganathan G, Linz H, Lahn M. Systematic review to establish the safety profiles for direct and indirect inhibitors of p38 Mitogen-activated protein kinases for treatment of cancer. A systematic review of the literature. Med Oncol. (2008) 25:323–30. doi: 10.1007/s12032-008-9039-1
27. Long J, Zhang Y, Yu X, Yang J, LeBrun DG, Chen C, et al. Overcoming drug resistance in pancreatic cancer. Expert Opin Ther Targets. (2011) 15:817–28. doi: 10.1517/14728222.2011.566216
28. Beloribi-Djefaflia S, Vasseur S, Guillaumond F. Lipid metabolic reprogramming in cancer cells. Oncogenesis. (2016) 5:e189. doi: 10.1038/oncsis.2015.49
29. Iwamoto H, Abe M, Yang Y, Cui D, Seki T, Nakamura M, et al. Cancer lipid metabolism confers antiangiogenic drug resistance. Cell Metab. (2018) 28:104–17.e5. doi: 10.1016/j.cmet.2018.05.005
30. Bardeesy N, DePinho RA. Pancreatic cancer biology and genetics. Nat Rev Cancer. (2002) 2:897–909. doi: 10.1038/nrc949
31. Cohen P. The search for physiological substrates of MAP and SAP kinases in mammalian cells. Trends Cell Biol. (1997) 7:353–61. doi: 10.1016/S0962-8924(97)01105-7
32. Patnaik A, Haluska P, Tolcher AW, Erlichman C, Papadopoulos KP, Lensing JL, et al. A first-in-human phase I study of the oral p38 MAPK inhibitor, ralimetinib (LY2228820 Dimesylate), in patients with advanced cancer. Clin Cancer Res. (2016) 22:1095–102. doi: 10.1158/1078-0432.CCR-15-1718
33. Shaw DE, Deneroff MM, Dror RO, Kuskin JS, Larson RH, Salmon JK, et al. Anton, a special-purpose machine for molecular dynamics simulation. Commun Acm. (2008) 51:91–7. doi: 10.1145/1364782.1364802
34. Laskowski RA, Macarthur MW, Mossn DS, Thornton J. Procheck - a program to check the stereochemical quality of protein structures. J Appl Crystallogr. (1993) 26:283–91. doi: 10.1107/S0021889892009944
35. Wishart DS, Sykes BD, Richards FM. Relationship between nuclear magnetic resonance chemical shift and protein secondary structure. J Mol Biol. (1991) 222:311–33. doi: 10.1016/0022-2836(91)90214-Q
36. Honndorf VS, Coudevylle N, Laufer S, Becker S, Griesinger C, Habeck M. Inferential NMR/X-ray-based structure determination of a dibenzo[a,d]cycloheptenone inhibitor-p38alpha MAP kinase complex in solution. Angew Chem Int Ed Engl. (2012) 51:2359–62. doi: 10.1002/anie.201105241
37. Lim SJ, Lee YJ, Lee E. p38MAPK inhibitor SB203580 sensitizes human SNU-C4 colon cancer cells to exisulind-induced apoptosis. Oncol Rep. (2006) 16:1131–5. doi: 10.3892/or.16.5.1131
38. Xing L, Shieh HS, Selness SR, Devraj RV, Walker JK, Devadas B, et al. Structural bioinformatics-based prediction of exceptional selectivity of p38 MAP kinase inhibitor PH-797804. Biochemistry. (2009) 48:6402–11. doi: 10.1021/bi900655f
39. Brenke R, Kozakov D, Chuang GY, Beglov D, Hall D, Landon MR, et al. Fragment-based identification of druggable 'hot spots' of proteins using Fourier domain correlation techniques. Bioinformatics. (2009) 25:621–7. doi: 10.1093/bioinformatics/btp036
40. Kozakov D, Grove LE, Hall DR, Bohnuud T, Mottarella SE, Luo L, et al. The FTMap family of web servers for determining and characterizing ligand-binding hot spots of proteins. Nat Protoc. (2015) 10:733–55. doi: 10.1038/nprot.2015.043
41. Iyoda K, Sasaki Y, Horimoto M, Toyama T, Yakushijin T, Sakakibara M, et al. Involvement of the p38 mitogen-activated protein kinase cascade in hepatocellular carcinoma. Cancer. (2003) 97:3017–26. doi: 10.1002/cncr.11425
42. Bulavin DV, Demidov ON, Saito S, Kauraniemi P, Phillips C, Amundson SA, et al. Amplification of PPM1D in human tumors abrogates p53 tumor-suppressor activity. Nat Genet. (2002) 31:210–5. doi: 10.1038/ng894
43. Yu W, Imoto I, Inoue J, Onda M, Emi M, Inazawa J. A novel amplification target, DUSP26, promotes anaplastic thyroid cancer cell growth by inhibiting p38 MAPK activity. Oncogene. (2007) 26:1178–87. doi: 10.1038/sj.onc.1209899
44. Campbell RM, Anderson BD, Brooks NA, Brooks HB, Chan EM, De Dios A, et al. Characterization of LY2228820 dimesylate, a potent and selective inhibitor of p38 MAPK with antitumor activity. Mol Cancer Ther. (2014) 13:364–74. doi: 10.1158/1535-7163.MCT-13-0513
45. Gupta J, del Barco Barrantes I, Igea A, Sakellariou S, Pateras IS, Gorgoulis VG, et al. Dual function of p38alpha MAPK in colon cancer: suppression of colitis-associated tumor initiation but requirement for cancer cell survival. Cancer Cell. (2014) 25:484–500. doi: 10.1016/j.ccr.2014.02.019
46. Lauby-Secretan B, Scoccianti C, Loomis D, Grosse Y, Bianchini F, Straif K, et al. Body fatness and cancer–viewpoint of the iarc working group. N Engl J Med. (2016) 375:794–8. doi: 10.1056/NEJMsr1606602
47. Zhang M, Di Martino JS, Bowman RL, Campbell NR, Baksh SC, Simon-Vermot T, et al. Adipocyte-derived lipids mediate melanoma progression via FATP proteins. Cancer Discov. (2018) 8:1006–25. doi: 10.1158/2159-8290.CD-17-1371
48. Park J, Morley TS, Kim M, Clegg DJ, Scherer PE. Obesity and cancer–mechanisms underlying tumour progression and recurrence. Nat Rev Endocrinol. (2014) 10:455–65. doi: 10.1038/nrendo.2014.94
49. Dirat B, Bochet L, Dabek M, Daviaud D, Dauvillier S, Majed B, et al. Cancer-associated adipocytes exhibit an activated phenotype and contribute to breast cancer invasion. Cancer Res. (2011) 71:2455–65. doi: 10.1158/0008-5472.CAN-10-3323
50. Kuzmanic A, Sutto L, Saladino G, Nebreda AR, Gervasio FL, Orozco M. Changes in the free-energy landscape of p38α MAP kinase through its canonical activation and binding events as studied by enhanced molecular dynamics simulations. Elife. (2017) 6:e22175. doi: 10.7554/eLife.22175
51. McClendon CL, Kornev AP, Gilson MK, Taylor SS. Dynamic architecture of a protein kinase. Proc Natl Acad Sci USA. (2014) 111:E4623–31. doi: 10.1073/pnas.1418402111
52. Shan Y, Eastwood MP, Zhang X, Kim ET, Arkhipov A, Dror RO, et al. Oncogenic mutations counteract intrinsic disorder in the EGFR kinase and promote receptor dimerization. Cell. (2012) 149:860–70. doi: 10.1016/j.cell.2012.02.063
53. Kumar A, Purohit R. Use of long term molecular dynamics simulation in predicting cancer associated SNPs. PLoS Comp Biol. (2014) 10:e1003318. doi: 10.1371/journal.pcbi.1003318
54. Daily MD, Makowski L, Phillips GN, Cui Q. Large-scale motions in the adenylate kinase solution ensemble: coarse-grained simulations and comparison with solution X-ray scattering. Chem Phys. (2012) 396:84–91. doi: 10.1016/j.chemphys.2011.08.015
55. Formaneck MS, Ma L, Cui Q. Reconciling the “old” and “new” views of protein allostery: a molecular simulation study of chemotaxis Y protein (CheY). Proteins. (2006) 63:846–867. doi: 10.1002/prot.20893
56. Boehr DD, Nussinov R, Wright PE. The role of dynamic conformational ensembles in biomolecular recognition. Nat Chem Biol. (2009) 5:789–96. doi: 10.1038/nchembio.232
57. Guo JJ, Zhou HX. Protein allostery and conformational dynamics. Chem Rev. (2016) 116:6503–15. doi: 10.1021/acs.chemrev.5b00590
58. Zweckstetter M, Bax A. Prediction of sterically induced alignment in a dilute liquid crystalline phase: aid to protein structure determination by NMR. J Am Chem Soc. (2000) 122:3791–2. doi: 10.1021/ja0000908
59. Rettenmaier TJ, Fan H, Karpiak J, Doak A, Sali A, Shoichet BK, et al. Small-molecule allosteric modulators of the protein kinase PDK1 from structure-based docking. J Med Chem. (2015) 58:8285–91. doi: 10.1021/acs.jmedchem.5b01216
60. De Nicola GF, Martin ED, Chaikuad A, Bassi R, Clark J, Martino L, et al. Mechanism and consequence of the autoactivation of p38α mitogen-activated protein kinase promoted by TAB1. Nat Struct Mol Biol. (2013) 20:1182–90. doi: 10.1038/nsmb.2668
61. ter Haar E, Prabhakar P, Liu X, Lepre C. Crystal structure of the P38 alpha-MAPKAP kinase 2 heterodimer. J Biol Chem. (2007) 282:9733–9. doi: 10.1074/jbc.M611165200
62. Zhang YY, Wu JW, Wang ZX. A distinct interaction mode revealed by the crystal structure of the kinase p38α with the MAPK binding domain of the phosphatase MKP5. Sci Signal. (2011) 4:ra88. doi: 10.1126/scisignal.2002241
63. Tang ZF, Li C, Kang B, Gao G, Li C, Zhang Z, et al. GEPIA: a web server for cancer and normal gene expression profiling and interactive analyses. Nucleic Acids Res. (2017) 45:W98–102. doi: 10.1093/nar/gkx247
64. Lindorff-Larsen K, Piana S, Palmo K, Maragakis P, Klepeis JL, Dror RO, et al. Improved side-chain torsion potentials for the Amber ff99SB protein force field. Proteins. (2010) 78:1950–8. doi: 10.1002/prot.22711
65. Kaminski GA, Friesner RA, Tirado-Rives J, Jorgensen ΩΛ. Evaluation and reparametrization of the OPLS-AA force field for proteins via comparison with accurate quantum chemical calculations on peptides. J Phys Chem B. (2001) 105:6474–87. doi: 10.1021/jp003919d
66. Jorgensen WL, Chandrasekhar J, Madura JD, Impey RW, Klein ML. Comparison of simple potential functions for simulating liquid water. J Chem Phys. (1983) 79:926–35. doi: 10.1063/1.445869
67. Aqvist J. Ion water interaction potentials derived from free-energy perturbation simulations. J Phys Chem. (1990) 94:8021–4. doi: 10.1021/j100384a009
68. Wang J, Wolf RM, Caldwell JW, Kollman PA, Case DA. Development and testing of a general amber force field. J Comp Chem. (2004) 25:1157–74. doi: 10.1002/jcc.20035
69. Dupradeau FY, Pigache A, Zaffran T, Savineau C, Lelong R, Grivel N, et al. The R.ED. tools: advances in RESP and ESP charge derivation and force field library building. Phys Chem Chem Phys. (2010) 12:7821–39. doi: 10.1039/c0cp00111b
70. Vanquelef E, et al. RED server: a web service for deriving RESP and ESP charges and building force field libraries for new molecules and molecular fragments. Nucleic Acids Res. (2011) 39:W511–7. doi: 10.1093/nar/gkr288
71. Rockey WM, Elcock AH. Rapid computational identification of the targets of protein kinase inhibitors. J Med Chem. (2005) 48:4138–52. doi: 10.1021/jm049461b
72. Fox T, Coll JT, Xie X, Ford PJ, Germann UA, Porter MD, et al. A single amino acid substitution makes ERK2 susceptible to pyridinyl imidazole inhibitors of p38 MAP kinase. Protein Sci. (1998) 7:2249–55. doi: 10.1002/pro.5560071102
73. Simard JR, Getlik M, Grütter C, Pawar V, Wulfert S, Rabiller M, et al. Development of a fluorescent-tagged kinase assay system for the detection and characterization of allosteric kinase inhibitors. J Am Chem Soc. (2009) 131:13286–96. doi: 10.1021/ja902010p
74. Bowers KJ, Chow E, Xu H, Dror RO. Scalable algorithms for molecular dynamics simulations on commodity clusters. In: Proceedings of the 2006 ACM/IEEE Conference on Supercomputing %@ 0-7695-2700-0. Tampa, FL ACM (2006), p. 84.
75. Berendsen HJC, van Postma JPM, van Gunsteren WF, DiNola AD, Haak JR. Molecular-dynamics with coupling to an external bath. J Chem Phys. (1984) 81:3684–90. doi: 10.1063/1.448118
76. Nose S. A unified formulation of the constant temperature molecular-dynamics methods. J Chem Phys. (1984) 81:511–9. doi: 10.1063/1.447334
77. Hoover WG. Canonical dynamics - equilibrium phase-space distributions. Phys Rev A. (1985) 31:1695–7. doi: 10.1103/PhysRevA.31.1695
78. Essmann U, Perera L, Berkowitz M, Darden T, Lee H, Pedersen LG. A smooth particle mesh ewald method. J Chem Phys. (1995) 103:8577–93. doi: 10.1063/1.470117
79. Shan Y, Klepeis JL, Eastwood MP, Dror RO, Shaw DE. Gaussian split Ewald: a fast Ewald mesh method for molecular simulation. J Chem Phys. (2005) 122:54101. doi: 10.1063/1.1839571
80. Krautler V, Van Gunsteren WF, Hunenberger PH. A fast SHAKE: algorithm to solve distance constraint equations for small molecules in molecular dynamics simulations. J Comp Chem. (2001) 22:501–8. doi: 10.1002/1096-987X(20010415)22:5<501::AID-JCC1021>3.0.CO;2-V
81. Humphrey W, Dalke A, Schulten K. VMD: Visual molecular dynamics. J Mol Graphics Model. (1996) 14:33–8. doi: 10.1016/0263-7855(96)00018-5
82. Nguyen H, Roe DR, Swails J, Case DA. PYTRAJ: Interactive Data Analysis for Molecular Dynamics Simulations. New Brunswick, NJ: Rutgers University (2016).
83. Roe DR, Cheatham TE. PTRAJ and CPPTRAJ: software for processing and analysis of molecular dynamics trajectory data. J Chem Theory Comput. (2013) 9:3084–95. doi: 10.1021/ct400341p
84. Van Der Spoel D, Lindahl E, Hess B, Groenhof G, Mark AE, Berendsen HJ. GROMACS: fast, flexible, and free. J Comput Chem. (2005) 26:1701–18. doi: 10.1002/jcc.20291
85. Abraham MJ, Murtola T, Schulz R, Páll S, Smith JC, Hess B, et al. GROMACS: high performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX. (2015) 1–2:19–25. doi: 10.1016/j.softx.2015.06.001
86. Li DW, Bruschweiler R. NMR-based protein potentials. Angew Chem Int Edn. (2010) 49:6778–80. doi: 10.1002/anie.201001898
87. Han B, Liu Y, Ginzinger SW, Wishart DS. SHIFTX2: significantly improved protein chemical shift prediction. J Biomol Nmr. (2011) 50:43–57. doi: 10.1007/s10858-011-9478-4
88. Vogtherr M, Saxena K, Grimme S, Betz M, Schieborr U, Pescatore B, et al. NMR backbone assignment of the mitogen-activated protein (MAP) kinase p38. J Biomol Nmr. (2005) 32:175. doi: 10.1007/s10858-005-2449-x
89. Zweckstetter M. NMR: prediction of molecular alignment from structure using the PALES software. Nat Protoc. (2008) 3:679–90. doi: 10.1038/nprot.2008.36
Keywords: p38α, molecular dynamics, tumor targeted therapy, conformational dynamics, pancreatic cancer
Citation: Yang L, Sun X, Ye Y, Lu Y, Zuo J, Liu W, Elcock A and Zhu S (2019) p38α Mitogen-Activated Protein Kinase Is a Druggable Target in Pancreatic Adenocarcinoma. Front. Oncol. 9:1294. doi: 10.3389/fonc.2019.01294
Received: 31 October 2019; Accepted: 07 November 2019;
Published: 26 November 2019.
Edited by:
Lasse Jensen, Linköping University, SwedenReviewed by:
Xingxin Wu, Nanjing University, ChinaJianfeng Cai, University of South Florida, United States
Copyright © 2019 Yang, Sun, Ye, Lu, Zuo, Liu, Elcock and Zhu. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Shun Zhu, c2h1bl96aHUmI3gwMDA0MDtmdWRhbi5lZHUuY24=
†These authors have contributed equally to this work