Computational investigation of a ternary model of SnoN-SMAD3-SMAD4 complex

Computational investigation of a ternary model of SnoN-SMAD3-SMAD4 complex

Computational Biology and Chemistry 83 (2019) 107159 Contents lists available at ScienceDirect Computational Biology and Chemistry journal homepage:...

2MB Sizes 0 Downloads 16 Views

Computational Biology and Chemistry 83 (2019) 107159

Contents lists available at ScienceDirect

Computational Biology and Chemistry journal homepage: www.elsevier.com/locate/cbac

Research article

Computational investigation of a ternary model of SnoN-SMAD3-SMAD4 complex

T

Mingfei Jia,1, Yelei Dinga,1, Xiaolong Lib,*, Ningfang Maob,*, Jie Chena,* a b

Department of Urology, Changzheng Hospital, Naval Military Medical University, Shanghai, 200003, China Department of Orthopedics, Changhai Hospital, Naval Military Medical University, Shanghai, 200433, China

A R T I C LE I N FO

A B S T R A C T

Keywords: Molecular dynamics simulations Conformational ensemble Free energy landscapes Protein-protein interactions Transcriptional factors

The transforming growth factor β (TGFβ) plays an essential role in the regulation of cellular processes such as cell proliferation, migration, differentiation, and apoptosis by association with SMAD transcriptional factors that are regulated by the transcriptional regulator SnoN. The crystal structure of SnoN-SMAD4 reveals that SnoN can adopt two binding modes, the open and closed forms, at the interfaces of SMAD4 subunits. Accumulating evidence indicates that SnoN can interact with both SMAD3 and SMAD4 to form a ternary SnoN-SMAD3-SMAD4 complex in the TGFβ signaling pathway. However, how the interaction of SnoN with the SMAD3 and SMAD4 remains unclear. Here, molecular dynamics (MD) simulations and molecular modeling methods were performed to figure out this issue. The simulations reveal that SnoNopen exists in two, open and semi-closed, conformations. Molecular modeling and MD simulation studies suggest that the SnoNclosed form interferes with the SMAD3SMAD4 protein; in contract, the SnoNopen can form a stable SnoN-SMAD3-SMAD4 complex. These mechanistic mechanisms may help elucidate the detailed engagement of SnoN with two SMAD3 and SMAD4 transcriptional factors in the regulation of TGFβ signaling pathway.

1. Introduction As a multifunctional cytokine, transforming growth factor β (TGFβ) plays an important role in the regulation of a broad spectrum of cellular processes, encompassing cell proliferation, migration, differentiation, and apoptosis (Tielemans et al., 2018). Dysregulation of TGFβ signaling is closely linked to metastatic cancer (Ahmadi et al., 2019). It is established that malfunction within the TGFβ signaling pathway is responsible for approximatively 40% of cancers. SMADs, a group of transcriptional factors, are downstream proteins of TGFβ signaling that directly transduce TGFβ-initiated signals from the cell surface receptors to the cell nucleus in which the formed multi-protein complex can modulate the expression of diverse genes to control a variety of biological activities (Heldin and Moustakas, 2016; Wu et al., 2007; Zhu et al., 2018). SMADs contains a family of eight structurally similar proteins that can be classified into three different subtypes: receptor-regulated SMADs (R-SMADs, SMAD1, SMAD2, SMAD3, SMAD5, and SMAD8/9),

common partner SMAD (Co-SMAD, SMAD4), and inhibitory SMADs (ISMADs, SMAD6/7) (Heldin and Moustakas, 2016). In the canonical pathway, R-SMADs such as SMAD2 and SMAD3 are activated by dual phosphorylation of Ser423 and Ser425 of SMAD3 as an example at the C-terminus of the protein via an active heteromeric complex with cell surface type II and type I Ser/Thr kinase receptors (Kamato et al., 2013). Two phosphorylated SMAD3 proteins interact with one CoSMAD (i.e. SMAD4) to form a heterotrimer, which, acting as a transcriptional regulatory complex, translocate into the nucleus to regulate gene transcription. In recent years, accumulating evidence suggests that the transcriptional regulator SnoN serves as a corepressor or coactivator of TGFβmediated transcription via the association with the SMAD3-SMAD4 complex (Marrow et al., 2014; Tecalco-Cruz et al., 2018). At high SnoN levels when overexpressed, it inhibits gene transcription regulated by the TGFβ-SMAD signaling pathway by the disruption of the formation of SMAD3-SMAD4 heterotrimer. In contrast, at basal SnoN levels, it also behaves as a transcriptional coactivator and thus promotes TGFβ-

Abbreviations: Co-SMAD, common partner SMAD; GPU, graphic processing unit; MD, molecular dynamics; MM/GBSA, molecular mechanics/generalized born surface area; PPIs, protein-protein interactions; RMSD, root-mean-square deviation; R-SMAD, receptor-regulated SMAD; SASA, solvent-accessible surface area; TGFβ, transforming growth factor β ⁎ Corresponding authors. E-mail addresses: [email protected] (X. Li), [email protected] (N. Mao), [email protected] (J. Chen). 1 These authors contributed equally. https://doi.org/10.1016/j.compbiolchem.2019.107159 Received 17 September 2019; Received in revised form 28 October 2019; Accepted 3 November 2019 Available online 09 November 2019 1476-9271/ © 2019 Elsevier Ltd. All rights reserved.

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

covalent bonds involving hydrogen atoms. The coordinates were saved every 5.0 ps for analysis. To test the convergence of sampling of the conformational space, we performed three replicas of MD simulations with different initial velocities. The analyses of the free energy landscapes of the SnoN conformation in the three replicas revealed the similar distribution (Figure S1), suggesting the convergence of sampling of the conformational space. Thus, in the production result analysis, we focus on the first replica of simulations.

induced transcription by means of stabilization of the SMAD3-SMAD4 heterotrimer. As such, SnoN exerts dual roles in the regulation of gene transcription in the TGFβ signaling pathway (Kamato et al., 2013). Because of SnoN’s ability to activate or inhibit the TGFβ signaling pathway, it is thus hypothesized that SnoN can stabilize or disrupt the association of SMAD3-SMAD4 complex. To rationalize this hypothesis, SnoN should have two different binding poses with SMAD4. In fact, the recently structural determination of SnoN-SMAD4 complex have revealed that SnoN binds to SMAD4 with two distinct forms in either an open or a closed binding mode. According to the crystal structure of SnoN-SMAD4 complex, the open form of SnoN makes less contacts with SMAD4 than the closed form. As such, whether the open form of SnoN is stable in the SnoN-SMAD4 complex is unclear and how the interaction of SnoN with the SMAD3-SMAD4 is unresolved. To figure out the aforementioned questions, we used molecular dynamics (MD) simulations and molecular modeling methods to investigate the binary SnoN-SMAD4 and SMAD3-SMAD4 complexes and the ternary SnoN-SMAD3-SMAD4 complex. MD simulations have widely utilized to explore the conformational dynamics and proteins and protein-protein interactions (PPIs) (Ho and Hamelberg, 2018; Li et al., 2019a; Lu et al., 2019a, b; Ni et al., 2019a; Rodriguez-Bussey et al., 2018; van Gunsteren et al., 2018). Based on extensive MD simulations of SnoN-SMAD4 and SMAD3-SMAD4 complexes, we can obtain the different interactions of SnoN with SMAD4 and the stable features of SnoN within SMAD4 in distinct binding modes. Furthermore, molecular modeling of SnoN-SMAD3-SMAD4 interaction using the stable SnoN-SMAD4 and SMAD3-SMAD4 complexes from MD simulations can uncover the dynamic properties of the ternary SnoNSMAD3-SMAD4 complex in the regulation of TGFβ signaling pathway.

2.3. Binding free energy calculations The binding free energy between SnoN and SMAD4 was assessed using the molecular mechanics/generalized born surface area (MM/ GBSA) method (Hou et al., 2011; Kong et al., 2015; Liu et al., 2018a, b; Momin et al., 2018; Shi et al., 2018; Xie et al., 2019) with the MMPBSA.py module in AMBER 14. ΔGbind = Gcom – (GSnoN + GSMAD4)

(1)

ΔGbind = ΔH – T⋅ΔS ≈ ΔEMM + ΔGsolv ― T⋅ΔS

(2)

ΔEMM = ΔEint + ΔEvdw + ΔEele

(3)

ΔGsolv = ΔGGB + ΔGSA

(4)

where the intramolecular interactions (ΔEint) that consist of bond, angle, and dihedral energies in Eq. 3 did not compute because the single trajectory method was used in the MM/GBSA calculations. The nonpolar contribution to the solvation free energy (ΔGSA) was assessed via the solvent-accessible surface area (SASA) by means of the LCPO method: ΔGSA = γ × SASA + β, where the surface tension constants γ and β were set to 0.0072 and 0, respectively. The polar contribution to the solvation free energy (ΔGGB) was calculated by the Generalized Born (GB) model (Onufriev’s GB, IGB = 2). The interior and exterior dielectric constants were set to 4 and 78.5, respectively. Owing to the nonconvergence of normal-mode analysis of large PPIs (Du et al., 2019; Wang et al., 2019a, b), the conformational entropy (T⋅ΔS) did not calculate between the SnoN and the SMAD4.

2. Materials and methods 2.1. Simulated system preparations The X-ray crystal structures of SnoN-SMAD4 (PDB ID: 5C4V) (Walldén et al., 2017) and SMAD3-SMAD4 (PDB ID: 1U7F) (Chacko et al., 2004) complexes were downloaded from the Protein Data Bank. The Amber ff14ipq force field was used for protein residues. Force field parameters from the AMBER Parameter Database were used for phosphoserine. Zinc ions were treated using the cationic dummy atom approach as previously reported (Li and Merz, 2017). The two systems were solvated in a truncated octahedral box of TIP3P (Jorgensen et al., 1983) water with a 10 Å buffer. Both systems were neutralized using counterions and subsequently a physiological saline concentration equal to 0.15 mol/L was added to both systems to simulate proteins under the physiological condition.

2.4. Free energy landscapes Free energy landscapes (Wang et al., 2019a, b; Zhang et al., 2019) were calculated for both systems by plotting the 2D histogram of p(D1, D2) according to F(D1, D2) = ― kBTln p(D1, D2)

(5)

where p(D1, D2) is the occupancy probability of visiting each grid point on the two-dimensional distance plane. D1 and D2 are the distances between the Cα atoms of SnoNH317 and SMAD4Y430 and between the Cα atoms of SnoNE341 and SMAD4R378, respectively.

2.2. MD simulations AMBER 14 (Case et al., 2005) was used for energy minimization, heating, equilibration, and production steps by means of graphic processing unit (GPU) code. Both systems were minimized using 10,000 steps of solvent minimization, 10,000 steps of hydrogen-only minimization, 10,000 steps of side chain minimization, and 50,000 steps of all-atom minimization (Li et al., 2019b; Lu et al., 2015). Both systems were then heated from 0 to 300 K over 1000 ps with 2 fs time-steps and 10.0 kcal mol-1 Å−2 position restraints on protein atoms were used. This was followed by constant temperature equilibration at 300 K for 2000 ps, with a 10.0 kcal mol-1 Å−2 position restraints on protein atoms using a canonical ensemble (NVT). In the production of 500 ns simulations using an isothermal isobaric ensemble (NPT) with periodic boundary conditions, the Langevin thermostat was used to maintain temperature (Wu and Brooks, 2003). An 10 Å cut-off was used for shortrange nonbonded interactions and particle mesh Ewald (Darden et al., 1993) was used to treat long-range electrostatic interactions. The SHAKE (Ryckaert et al., 1977) method was applied to constrain all

2.5. ZDOCK and RDOCK protein-protein docking The molecular modeling of the ternary complex of SnoN-SMAD3SMAD4 was carried out using ZDOCK and RDOCK protein-protein docking methods (Hwang et al., 2014). ZDOCK is a rigid-body proteinprotein method that can probe the rotational and translational space of a protein surface. RODCK, a CHARMM-based energy minimization method, can refine and score the output poses generated by ZDOCK. The initial structures of SnoN and SMAD3-SMAD4 used in the docking study were taken from the representative structures from MD simulations of SnoN-SMAD4 and SMAD3-SMAD4 complexes, respectively. The whole surface of SMAD4 in the SMAD3-SMAD4 complex was set as “receptors” to explore the potential binding sites of SnoN. The detailed docking parameters were used as previous studies (Mou et al., 2019). 2

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

Fig. 1. Surface and cartoon representations of the crystal structure of SnoN-SMAD4 as trimers of heterodimers (PDB ID: 5C4V), with one open (organge) and two closed (lime) SnoN proteins bound to the interfaces of SMAD4 trimer (pink). The tight insert shows the interactions between the open SnoN (SnoNopen) and SMAD4. The left insert shows the interactions between the closed SnoN (SnoNclosed) and the interfaces of the subunits B and C of SMAD4. Secondary structure elements of αhelix (α) and β-strand (β) are depicted with numbers increasing from the N-terminal to the C-terminal tails.

3. Results

SMAD4 subunits (Chakrabarti et al., 2016; Liu et al., 2018a, b; Lu et al., 2019b; Trusova and Gorbenko, 2018). Then, ZDOCK and RDOCK protein-protein docking methods were used to construct the potential ternary complex of SnoN-SMAD3-SMAD4.

3.1. Overview of SnoN-SMAD4 and SMAD3-SMAD4 structural complexes The X-ray crystal structural complex between the SAND domain of SnoN (residues: T238-S356) and the MH2 domain of SMAD4 (residues: I314-Q549) has recently been solved at 2.6 Å resolution (PDB ID: 5C4V) (Walldén et al., 2017). The structure presents as a trimer of SnoNSMAD4 heterodimers (Fig. 1). The secondary structural elements of SnoN include five α-helices, five β-strands, and loops connecting them. Typically, α1, α3, and α4 are irregular helices. Notably, SnoN adopts distinct binding modes at the interfaces of SMAD4 subunits. In the interface of the subunits A and B of SMAD4, SnoN is in the open form, designated as SnoNopen. In this state, SnoNopen forms less contacts with SMAD4, only by the interactions of the β5 of SnoNopen with the β8 of SMAD4A. By sharp contrast, in both the interfaces of the subunits B and C as well as the subunits A and C of SMAD4, SnoN is in the closed form, designated as SnoNclosed. In addition to the β5, the helices α3-α5 of the two SnoNclosed proteins interact with the helix α1 of SMAD4 and the helix α2 of the adjacent SMAD4. Owing to the similar binding modes of the two SnoNclosed at the interfaces of SMAD4 subunits, we only focus on the SnoNclosed at the interface of the subunits B and C of SMAD4 in the following analyses. The crystal structural complex of SMAD3-SMAD4 was solved at 2.5 Å resolution (PDB ID: 1U7F) (Chacko et al., 2004). This complex presents a heterotrimer with a molar ratio of SMAD3:SMAD4 to 2:1 (Fig. 2). Notably, the C-terminal tail of SMAD3 is dually phosphorylated at residues Ser423 and Ser425. This dual phosphorylation of SMAD3 facilitates its binding to SMAD4. As clearly showed in Fig. 1, the SnoNopen form engages in much less contacts with SMAD4 than the SnoNclosed form. This leads us to hypothesize that whether the SnoNopen form would stabilize within the interface of SMAD4 subunits. Furthermore, how in the detailed interaction of SnoN with the SMAD3-SMAD4 protein in the ternary complex. To answer these issues, we used extensive MD simulations to probe the dynamic characterizations of different SnoN proteins at the interface of

3.2. Conformational dynamics of SnoN-SMAD4 and SMAD3-SMAD4 complexes To reveal the conformational dynamics of both SnoNopen and SnoNclosed forms at the interface of SMAD4 subunits, the root-meansquare deviation (RMSD) for all Cα atoms of SnoNopen and SnoNclosed was monitored. As shown in Fig. 3A, overall, both SnoNopen and SnoNclosed forms were stable at the interface of SMAD4 subunits with small fluctuations of RMSD plots. However, the RMSD of SnoNopen was a little larger than that of SnoNclosed, with the values of 2.40 ± 0.26 Å for the SnoNopen and 2.07 ± 0.23 Å for the SnoNclosed. This suggested that the SnoNopen form had the potential to undergo conformational changes along the simulations. Fig. 3B shows the RMSD plot for all Cα atoms of SMAD3-SMAD4 complex. The RMSD plot indicated that the SMAD3-SMAD4 complex was stable during the simulations with the RMSD value of 1.99 ± 0.25 Å. The stable SMAD3-SMAD4 complex obtained from the MD simulations will be used in the construction of ternary SnoNSMAD3-SMAD4 complex. 3.3. Free energy landscapes of SnoN within the interfaces of SMAD4 subunits Inspection of RMSD plots of different SnoN proteins provided a clue that the SnoNopen form had a conformational plasticity. To further unravel the conformational diversity of the SnoNopen form in the simulations, the free energy landscapes of both SnoNopen and SnoNclosed forms were analyzed. Two distance indexes were used to construct the free energy landscapes. One distance (D1) is between the Cα atoms of SnoN His317 (SnoNH317) from the β5 and SMAD4 Tyr430 (SMAD4Y430) from the β8, which monitors the interactions between the β5 of SnoN 3

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

Fig. 2. Cartoon representation of the crystal structure of SMAD3-SMAD4 heterotrimer (PDB ID: 1U7F) with one SMAD4 (pink) molecule bound to two SMAD3 molecules (light cyan). The two phosphoserine residues Ser423 (pSer423) and Ser425 (pSer425) at the C-terminal tail of SMAD3 are depicted by stick models.

and the β8 of SMAD4, and the other (D2) is between the Cα atoms of SnoNE341 from the α5 and SMAD4R378 from the α1, which describes the conformational dynamics of the α5 of SnoN. Fig. 4 showed the free energy landscapes projected onto the two distance indexes for the SnoNopen (A) and SnoNclosed (B) forms. Dramatically, the free energy landscapes are quite different for these two

SnoN proteins. For example, the two free-energy basins (C1 and C2) for the SnoNopen form were observed compared with the SnoNclosed form; the latter mainly had only one major state. C1 and C2 are approximately located in regions of the corresponding D1 and D2 values of 7.5–7.7 Å and 26.0–30.0 Å as well as 7.3–7.8 Å and 8.0–13.0 Å, respectively. Interestingly, the C2 state of the SnoNopen from is close to the major state of the SnoNclosed from. Taken together, such remarkable

Fig. 3. (A) Time dependence of root-mean-square deviation (RMSD) for all Cα atoms of SnoNopen (red line) and SnoNclosed (black line) along the 500 ns MD simulations. (B) Time dependence of RMSD for all Cα atoms of SMAD3-SMAD4 complex along the 500 ns MD simulations. 4

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

Fig. 4. Constructed 2D free-energy contour of SnoNopen (A) and SnoNclosed (B) with respect to D1 (the Cα atoms between SnoNH317 and SMAD4Y430) and D2 (the Cα atoms between SnoNE341 and SMAD4R378). The SnoNopen features two major substates, designated as C1 (the open form) and C2 (the semi-closed form). The solid blue square represents the position of the initial structure of SnoNopen and SnoNclosed conformations. The unit of free energy is in kcal/mol.

differences between the conformational distributions of the SnoNopen and SnoNclosed forms suggested that the SnoNopen form can adopt different binding modes at the interfaces of SMAD4 subunits along the simulations. To confirm that our MD simulations are enough to construct a reliable free energy landscape, we assessed the convergence of the free energy landscape by choosing different subsets from the complete simulation dataset. Six subsets of the complete dataset were extracted with different accumulated simulation times by truncated each MD simulation into varied length, namely 400 ns, 420 ns, 440 ns, 460 ns, 480 ns, and 500 ns. To check whether the sampling is converged or not, we projected the MD snapshots from each subset onto the same two order parameter pairs. As shown in Figures S1 (SnoNopen) and S2 (SnoNclosed), all datasets demonstrate similar free energy landscapes and no extra metastable state appears while increasing the simulation time, suggesting the convergence of current MD simulations.

subsequently. As shown in Fig. 5A, the most representative structure of the C1 was superimposed to the crystal structure of SnoN-SMAD4 complex. One can see that the conformation of the SnoN in the C1 is similar to that of the crystal SnoNopen from in the SnoN-SMAD4 complex, with the RMSD between the two states less than 1 Å. This comparison indicated that the SnoN protein in the C1 still adopted the open form. In a similar vein, the most representative structure of the C2 was superimposed to the crystal structure of SnoN-SMAD4 complex. As shown in Fig. 5B, the conformation of the SnoN in the C2 is significantly different from that of the crystal SnoNopen from in the SnoN-SMAD4 complex. The RMSD between the two states is 2.47 Å. This comparison suggested that the SnoN protein in the C2 has changed its binding mode rather than in the open form. Inspection of the SnoN in the C2 found that the two irregular helices α3 and α4 are in the loop conformation that make additional interactions with the helices α1 and α2 of SMAD4, a similar situation observed in the SnoNclosed from at the interface of SMAD4 subunits. However, the SnoN in the C2 is not equivalent to the SnoNclosed from in the crystal structure. As such, we defined the SnoN in the C2 as a semi-closed (SnoNsemi-closed) from. To further detect whether the existence of new ligand binding sites in the SnoNsemi-closed from, the Fpocket program (Schmidtke et al., 2010) was used to analyze the potential binding sites in the most

3.4. The SnoNopen from has two states: open and semi-closed forms In order to elucidate the distinct binding modes of the SnoNopen from with the SMAD4, the cluster analyses of C1 and C2 were performed and the most representative structures of the C1 and C2 were

Fig. 5. (A) Structural superimposition of the most representative structure of the C1 (pink) to the crystal structure of SnoN-SMAD4 complex (light cyan). (B) Structural superimposition of the most representative structure of the C2 (green) to the crystal structure of SnoN-SMAD4 complex (light cyan). 5

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

structure of SMAD3-SMAD4 complex. After overlapping, the CHARMM force field-based energy optimization (Brooks et al., 2009) was performed to potentially optimize the superimposed complexes. Fig. 6A shows the superimposition between SnoNopen-SMAD4 and SMAD3-SMAD4. Similar with the crystal structure of SnoN-SMAD4, the SnoNopen form interacts with the SMAD4 by the association of the β5 of SnoN and the β8 of SMAD4. Notably, the dually phosphorylated Ser423 and Ser425 at the C-terminal tail of SMAD3 do not interfere with the SnoNopen form. This modelled complex indicated that the SnoNopen form can interact with the SMAD3-SMAD4 complex to form a stable SnoN-SMAD3-SMAD4 complex. Fig. 6B shows the superimposition between SnoNsemi-closed-SMAD4 and SMAD3-SMAD4. One can see that the dually phosphorylated Ser423 and Ser425 at the C-terminal tail of SMAD3 can also interact with the SnoNsemi-closed form. Fig. 7 shows the superimposition between SnoNclosed-SMAD4 and SMAD3-SMAD4. Compared with the SnoNopen form, the dually phosphorylated Ser423 and Ser425 at the C-terminal tail of SMAD3 have severe steric interference with the SMAD4, leading to the instability of the SnoNclosed form with the SMAD3-SMAD4 complex.

Table 1 Binding free energy calculations (kcal/mol) between SMAD4 and the open, the semi-closed, as well as the closed forms of SnoN. Energy items

Open

Semi-closed

Closed

ΔEvdw ΔEele ΔGgas ΔGsolv ΔGbind

−30.31 ± 9.12 −1241.16 ± 52.42 −1271.47 ± 54.68 953.63 ± 47.74 −317.84 ± 14.91

−42.54 ± 12.09 −1596.65 ± 49.45 −1639.19 ± 52.16 1215.09 ± 40.52 −424.11 ± 30.11

−58.57 ± 14.70 −2077.84 ± 50.67 −2136.42 ± 45.34 1557.41 ± 43.26 −579.01 ± 37.08

representative structure of SnoNsemi-closed from and the results were compared with the SnoNopen from. As shown in Figure S4, Fpocket results revealed that there is a ligand binding site formed by the H5 helix, the H4 helix, and the loop connecting the H1 helix and the β1 strand observed in the SnoNsemi-closed from, which was not observed in the SnoNopen from. This site may represent a hidden allosteric site that can be used for allosteric drug design (Lu et al., 2018; Ni et al., 2019b). 3.5. Interactions of the open, semi-closed, and closed forms of SnoN with SMAD4

3.7. Molecular docking study of SnoN-SMAD3-SMAD4 complex

To further reveal the different binding abilities of the open, the semi-closed, and the closed forms of SnoN to the SMAD4, the MM_GBSA binding free energy (ΔGbind) calculations were carried out, which have previously widely used to assess the binding affinities of PPIs (Kato et al., 2018; Liu et al., 2018a, b; Ni et al., 2018; Zhao et al., 2018). As listed in Table 1, the most favorable contribution to the ΔGbind is derived from the electrostatic interaction (ΔEele). This reflects that in all complexes the electrostatic interactions between SnoN and SMAD4 are mainly responsible for binding. Overall, the calculated ΔGbind values for the SnoNopen, SnoNsemi-closed, and SnoNclosed forms were -317.84 ± 14.91, -424.11 ± 90.11, and -579.01 ± 37.08 kcal/mol. As a result, the strength of ΔGbind decreased in the order SnoNclosed > SnoNsemi-closed > SnoNopen. As expected, the SnoNclosed form has the strongest interaction with the interfaces of SMAD4 subunits, which is consistent with the pronounced contacts between the SnoNclosed form and SMAD4 in the crystal structure.

The ternary complex of SnoN-SMAD3-SMAD4 generated by molecular superimposition method might yield steric conflicts in some local regions. Thus, the stability of the ternary SnoN-SMAD3-SMAD4 complex is unclear. To confirm the dynamic stability of the three transcriptional factors in the ternary complex and to further validate the preferred SnoN conformation in the SnoN-SMAD3-SMAD4 complex, we used ZDOCK and RDOCK protein-protein docking methods to construct the potential SnoN-SMAD3-SMAD4 complexes. A total of 8 potential SnoN-SMAD3-SMAD4 complexes (models 0–7) were generated by the ZDOCK and RDOCK method. Among these ternary complexes, the docked SnoN protein in the model 0 of the docked SnoN-SMAD3SMAD4 complex was similar with the conformation of SnoN in the crystal structure of SnoN-SMAD4. Each 500 ns MD simulation for the 8 SnoN-SMAD3-SMAD4 complexes was then performed. After simulations, the RMSD for all Cα atoms of the SnoN, SMAD3, and SMAD4 in each system was analyzed. As shown in Figure S5, the RMSD plots suggested that SnoN can stabilize in each system, which can also be validated through structural inspection of the SnoN in the SnoN-SMAD3-SMAD4 complex (Figure S6). Furthermore, the MM_GBSA binding free energy (ΔGbind) calculations between the SnoN and SMAD3-SMAD4 were performed to reveal which model was the most preferred binding mode. As shown in Table S1, the ΔGbind decreased in the order models 0 > 3 > 7 > 1 > 5 > 2 > 6 > 4. This result suggested that the model 0 with the docked SnoN conformation similar with the crystal structure had the lowest binding free energy, consistent with the structural superimposition results and further reflecting the strongest

3.6. Superimposition study of SnoN-SMAD3-SMAD4 complex Up to date, the crystal structure of ternary SnoN-SMAD3-SMAD4 complex has remained unavailable. Therefore, we first used the molecular superimposition method to construct the SnoN-SMAD3-SMAD4 complex, as the same procedure previously performed by Walldén et al. (2017). To accomplish this goal, the most representative structural complex of SMAD3-SMAD4 was extracted from the MD simulations. Then, the structural complexes of SnoNopen-SMAD4, SnoNsemi-closedSMAD4, and SnoNclosed-SMAD4 were respectively overlapped with the

Fig. 6. The modelled SnoNopen-SMAD3-SMAD4 (A) and SnoNsemi-closed-SMAD3-SMAD4 (B) complexes after energy minimizations. The SnoNopen and SnoNsemi-closed were superimposed based on the SMAD4 with the most representative structure of SMAD3-SMAD4 complex extracted from the MD simulations, respectively. The two phosphoserine residues Ser423 (pSer423) and Ser425 (pSer425) at the C-terminal tail of SMAD3 (magenta) are depicted by stick models.

6

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

Fig. 7. The modelled SnoNclosed-SMAD3-SMAD4 complex. The SnoNclosed was superimposed based on the SMAD4 with the most representative structure of SMAD3-SMAD4 complex extracted from the MD simulations. The two phosphoserine residues Ser423 (pSer423) and Ser425 (pSer425) at the Cterminal tail of SMAD3 (magenta) are depicted by stick models. The C-terminal tail of SMAD3 has severe steric conflicts with the SnoNclosed form.

Declaration of competing Interest

binding affinity between the SnoN and the SMAD3-SMAD4.

The authors declare no competing financial interests. 4. Discussion Acknowledgements Previous biological experiments showed that SnoN has dual roles in the regulation of TGFβ signaling pathway through interaction with transcriptional regulators SMAD3 and SMAD4. Overexpression of SnoN in cells suppresses TGFβ signaling, blocking TGFβ-induced gene expression and cell cycle arrest (Edmiston et al., 2005). This suggests that at high levels, SnoN promotes cell proliferation and transformation. In contrast, at low levels, SnoN acts as a positive mediator to facilitate TGFβ-dependent inhibition of proliferation (Sarker et al., 2005). However, the mechanistic insights into the dual roles of SnoN involvement in the TGFβ signaling pathway remains to be understood. Recently, Walldén et al. (2017) have determined the crystal structure of binary SnoN-SMAD4 complex and found that SnoN can adopt two different states in the SMAD4 active site, the open and the close states. Further molecular modeling of ternary SnoN-SMAD3-SMAD4 complex, based on the crystal structures of SnoN-SMAD4 and SMAD3SMAD4 complexes, implied that the open state of SnoN can interact with the SMAD3-SMAD4 complex, while the close state is unable to form a ternary complex. However, the stability and the dynamic properties of the ternary SnoN-SMAD3-SMAD4 complex are unknown. In the present study, we performed MD simulations and molecular modeling methods to elucidate the associations of SnoN with the SMAD3-SMAD4 complex. Free energy landscape calculations revealed that the SnoNopen form at the interface of SMAD4 subunits can exist two conformational states: the SnoNopen and the SnoNsemi-closed states. The MM_GBSA binding free energy calculations revealed that the binding abilities of SnoN-SMAD4 complexes decreased in the order SnoNclosed > SnoNsemi-closed > SnoNopen, with the SnoNclosed-SMAD4 complex showing the strongest binding. Molecular modeling of the ternary SnoN-SMAD3-SMAD4 complex showed that SnoNclosed can interface with the SMAD3-SMAD4 complex, whereas the SnoNopen can form a stable SnoN-SMAD3-SMAD4 complex. This result can be explained the dual roles of SnoN in the regulation of TGFβ signaling pathway where the low levels of SnoN promote TGFβ-dependent transcription via the formation of a stable ternary SnoN-SMAD3-SMAD4 complex, while the high levels antagonize it owing to the inability of the formation of ternary SnoN-SMAD3-SMAD4 complex.

The study was supported by the grants from the National Natural Science Foundation of China (No.81372005 and No. 81901255). Appendix A. Supplementary data Supplementary material related to this article can be found, in the online version, at doi:https://doi.org/10.1016/j.compbiolchem.2019. 107159. References Ahmadi, A., Najafi, M., Farhood, B., Mortezaee, K., 2019. Transforming growth factor-β signaling: tumorigenesis and targeting for cancer therapy. J. Cell. Physiol. 234, 12173–12187. Brooks, B.R., Brooks, C., Mackerell, A.D., Nilsson, L., Petrella, R.J., Roux, B., Won, Y., Archontis, G., Bartels, C., Boresch, S., Caflisch, A., Caves, L., Cui, Q., Dinner, A.R., Feig, M., Fischer, S., Gao, J., Hodoscek, M., Im, W., Kuczera, K., Lazaridis, T., Ma, J., Ovchinnikov, V., Paci, E., Pastor, R.W., Post, C.B., Pu, J.Z., Schaefer, M., Tidor, B., Venable, R.M., Woodcock, H.L., Wu, X., Yang, W., D.M, Y, Karplus, M., 2009. CHARMM: molecular dynamics simulation package. J. Comput. Chem. 30, 1545–1614. Case, D.A., Cheatham, T.E., Darden, T., Gohlke, H., Luo, R., Merz, K.M., Onufriev, A., Simmerling, C., Wang, B., Woods, R.J., 2005. The Amber biomolecular simulation programs. J. Comput. Chem. 26, 1668–1688. Chacko, B.M., Qin, B.Y., Tiwari, A., Shi, G., Lam, S., Hayward, L.J., De Caestecker, M., Lin, K., 2004. Structural basis of heteromeric smad protein assembly in TGF-beta signaling. Mol. Cell 15, 813–823. Chakrabarti, M., Jang, H., Nussinov, R., 2016. Comparison of the conformations of KRAS isoforms, K-Ras4A and K-Ras4B, points to similarities and significant differences. J. Phys. Chem. B 120, 667–679. Darden, T., York, D., Pedersen, L., 1993. Particle mesh Ewald: an N.log(N) method for Ewald sums in large systems. J. Chem. Phys. 98, 10089–10092. Du, J., Qiu, M., Guo, L., Yao, X., 2019. Computational study of the binding mechanism between farnesoid X receptor α and antagonist N-benzyl-N-(3-(tertbutyl)-4-hydroxyphenyl)-2,6-dichloro-4-(dimethylamino) benzamide. J. Biomol. Struct. Dyn. 37, 1628–1640. Edmiston, J.S., Yeudall, W.A., Chung, T.D., Lebman, D.A., 2005. Inability of transforming growth factor-beta to cuase SnoN degradation leads to resistance to transforming growth factor-beta-induced growth arrest in esophageal cancer cells. Cancer Res. 65, 4782–4788. Heldin, C.H., Moustakas, A., 2016. Signaling receptors for TGF-β family members. Cold Spring Harb. Perspect. Biol. 8, 1–33. Ho, K.C., Hamelberg, D., 2018. Combinatorial coarse-graining of molecular dynamics simulations for detecting relationships between local configurations and overall conformations. J. Chem. Theory Comput. 14, 6026–6034.

7

Computational Biology and Chemistry 83 (2019) 107159

M. Ji, et al.

of protein‐protein interactions (PPIs): A new paradigm for PPI drug discovery. Med. Res. Rev. 39, 2314–2342. Ni, D., Li, X., He, X., Zhang, H., Zhang, J., Lu, S., 2019b. Drugging K-RasG12C through covalent inhibitors: mission possible? Pharma. Thera. 202, 1–17. Rodriguez-Bussey, I., Yao, X.Q., Shouaib, A.D., Lopez, J., Hamelberg, D., 2018. Decoding allosteric communication pathways in Cyclophilin a with a comparative analysis of perturbed conformational ensembles. J. Phys. Chem. B 122, 6528–6535. Ryckaert, J.P., Ciccotti, G., Berendsen, H.J.C., 1977. Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys. 23, 327–341. Sarker, K.P., Wilson, S.M., Bonni, S., 2005. SnoN is a cell type-specific mediator of transforming growth factor-β responses. J. Biol. Chem. 280, 13037–13046. Schmidtke, P., Le Guilloux, V., Maupetit, J., Tufféry, P., 2010. Fpocket: online tools for protein ensemble pocket detection and tracking. Nucleic Acids Res. 38, 582–589. Shi, D., Zhou, S., Liu, X., Zhao, C., Liu, H., Yao, X., 2018. Understanding the structural and energetic basis of PD-1 and monoclonal antibodies bound to PD-L1: a molecular modeling perspective. Biochim. Biophys. Acta - Gen. Subj. 1862, 576–588. Tecalco-Cruz, A.C., Ríos-López, D.G., Vázquez-Victorio, G., Rosales-Alvarez, R.E., MacíasSilva, M., 2018. Transcriptional cofactors Ski and SnoN are major regulators of the TGF-β/Smad signaling pathway in health and disease. Signal Transduct. Target. Ther. 3, 1–15. Tielemans, B., Delcroix, M., Belge, C., Quarck, R., 2018. TGFβ and BMPRII signalling pathways in the pathogenesis of pulmonary arterial hypertension. Drug Discov. Today 24, 703–716. Trusova, V.M., Gorbenko, G.P., 2018. Molecular dynamics simulations of lysozyme–lipid systems: probing the early steps of protein aggregation. J. Biomol. Struct. Dyn. 36, 2249–2260. van Gunsteren, W.F., Daura, X., Hansen, N., Mark, A.E., Oostenbrink, C., Riniker, S., Smith, L.J., 2018. Validation of molecular simulation: an overview of issues. Angew. Chem. Int. Ed. Engl. 57, 884–902. Walldén, K., Nyman, T., Hällberg, B.M., 2017. SnoN stabilizes the SMAD3/SMAD4 protein complex. Sci. Rep. 7, 1–10. Wang, E., Sun, H., Wang, J., Wang, Z., Liu, H., Zhang, J.Z.H., Hou, T., 2019a. End-Point Binding Free Energy Calculation with MM/PBSA and MM/GBSA: Strategies and Applications in Drug Design. Chem. Rev. 119, 9457–9508. Wang, L., Zheng, G., Liu, X., Ni, D., He, X., Cheng, J., Lu, S., 2019b. Molecular dynamics simulations provide insights into the origin of gleevec’s selectivity toward human tyrosine kinases. J. Biomol. Struct. Dyn. 37, 2733–2744. Wu, J., Zhou, Y., Xiang, G., Mitchelson, K., Cheng, J., Wang, Y., Chen, A., Yu, J., Zhang, L., 2007. Identification of the gene transcription and apoptosis mediated by TGF-βSmad2/3-Smad4 signaling. J. Cell. Physiol. 215, 422–433. Wu, X., Brooks, B.R., 2003. Self-guided Langevin dynamics simulation method. Chem. Phys. Lett. 381, 512–518. Xie, T., Yu, J., Fu, W., Wang, Z., Xu, L., Chang, S., Wang, E., Zhu, F., Zeng, S., Kang, Y., Hou, T., 2019. Insight into the selective binding mechanism of DNMT1 and DNMT3A inhibitors: a molecular simulation study. Phys. Chem. Chem. Phys. 21, 12931–12947. Zhang, H., He, X., Ni, D., Mou, L., Chen, X., Lu, S., 2019. How does the novel T315L mutation of breakpoint cluster region-abelson (BCR-ABL) kinase confer resistance to ponatinib: a comparative molecular dynamics simulation study. J. Biomol. Struct. Dyn. https://doi.org/10.1080/07391102.2019.1567390. Zhao, M.L., Wang, W., Nie, H., Cao, S.S., Du, L.F., 2018. In silico structure prediction and inhibition mechanism studies of AtHDA14 as revealed by homology modeling, docking, molecular dynamics simulation. Comput. Biol. Chem. 75, 120–130. Zhu, Q., Chang, A., Xu, A., Luo, K., 2018. The regulatory protein SnoN antagonizes activin/Smad2 protein signaling and thereby promotes adipocyte differentiation and obesity in mice. J. Biol. Chem. 293, 14100–14111.

Hou, T., Wang, J., Li, Y., Wang, W., 2011. Assessing the performance of the MM/PBSA and MM/GBSA methods. 1. The accuracy of binding free energy calculations based on molecular dynamics simulations. J. Chem. Inf. Model. 51, 69–82. Hwang, H., Kim, B.-H., Vreven, T., Pierce, B.G., Wiehe, K., Weng, Z., 2014. ZDOCK server: interactive docking prediction of protein-protein complexes and symmetric multimers. Bioinformatics 30, 1771–1773. Jorgensen, W.L., Chandrasekhar, J., Madura, J.D., Impey, R.W., Klein, M.L., 1983. Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 79, 926. Kamato, D., Burch, M.L., Piva, T.J., Rezaei, H.B., Rostam, M.A., Xu, S., Zheng, W., Little, P.J., Osman, N., 2013. Transforming growth factor-β signalling: role and consequences of Smad linker region phosphorylation. Cell. Signal. 25, 2017–2024. Kato, Y., Kihara, H., Fukui, K., Kojima, M., 2018. A ternary complex model of Sirtuin4NAD+-Glutamate dehydrogenase. Comput. Biol. Chem. 74, 94–104. Kong, X., Pan, P., Li, D., Tian, S., Li, Y., Hou, T., 2015. Importance of protein flexibility in ranking inhibitor affinities: modeling the binding mechanisms of piperidine carboxamides as Type I1/2 ALK inhibitors. Phys. Chem. Chem. Phys. 17, 6098–6113. Li, J., Fu, Q., Liang, Y., Cheng, B., Li, X., 2019a. Microsecond molecular dynamics simulations and dynamic network analysis provide understanding of the allosteric inactivation of GSK3β induced by the L343R mutation. J. Mol. Model. 25, 111. Li, J., Wu, G., Fu, Q., Ge, H., Liu, S., Li, X., Cheng, B., 2019b. Exploring the influence of conserved lysine69 on the catalytic activity of the helicobacter pylori shikimate dehydrogenase: A combined QM/MM and MD simulations. Comput. Biol. Chem. 83, 107098. Li, P., Merz, K.M., 2017. Metal ion modeling using classical mechanics. Chem. Rev. 117, 1564–1686. Liu, L., Yu, H., Yang, L., Zhao, X., Huang, X., 2018a. Elucidating proton-mediated conformational changes in an acid-sensing ion channel 1a through molecular dynamics simulation. J. Biomol. Struct. Dyn. 36, 2261–2267. Liu, N., Zhou, W., Guo, Y., Wang, J., Fu, W., Sun, H., Li, D., Duan, M., Hou, T., 2018b. Molecular dynamics simulations revealed the regulation of ligands to the interactions between androgen receptor and its coactivator. J. Chem. Inf. Model. 58, 1652–1661. Lu, S., Banerjee, A., Jang, H., Zhang, J., Gaponenko, V., Nussinov, R., 2015. GTP binding and oncogenic mutations may attenuate hypervariable region (HVR)-Catalytic domain interactions in small GTPase K-Ras4B, exposing the effector binding site. J. Biol. Chem. 290, 28887–28900. Lu, S., Ji, M., Ni, D., Zhang, J., 2018. Discovery of hidden allosteric sites as novel targets for allosteric drug design. Drug Discov. Today 23, 359–365. Lu, S., Ni, D., Wang, C., He, X., Lin, H., Wang, Z., Zhang, J., 2019a. Deactivation pathway of ras GTPase underlies conformational substates as targets for drug design. ACS Catal. 9, 7188–7196. Lu, S., Shen, Q., Zhang, J., 2019b. Allosteric methods and their applications: facilitating the discovery of allosteric drugs and the investigation of allosteric mechanisms. Acc. Chem. Res. 52, 492–500. Marrow, B., Secreted, S., Protect, C., 2014. SnoN signaling in proliferating cells and postmitotic neurons. FEBS Lett. 71, 3831–3840. Momin, M., Yao, X.Q., Thor, W., Hamelberg, D., 2018. Substrate sequence determines catalytic activities, domain-binding preferences, and allosteric mechanisms in Pin1. J. Phys. Chem. B 122, 6521–6527. Mou, L., Dou, W., Meng, G., Sun, K., Chen, X., 2019. The structural basis of the autoinhibition mechanism of glycogen synthase kinase 3β (GSK3β): molecular modeling and molecular dynamics simulation studies. J. Biomol. Struct. Dyn. https://doi.org/ 10.1080/07391102.2019.1615988. Ni, D., Liu, D., Zhang, J., Lu, S., 2018. Computational insights into the interactions between Calmodulin and the c/nSH2 domains of p85α regulatory subunit of PI3Kα: implication for PI3Kα activation by calmodulin. Int. J. Mol. Sci. 19, 151. Ni, D., Lu, S., Zhang, J., 2019a. Emerging roles of allosteric modulators in the regulation

8