Exploring the Molecular Basis of dsRNA Recognition by Mss116p

Article pubs.acs.org/Langmuir

Exploring the Molecular Basis of dsRNA Recognition by Mss116p Using Molecular Dynamics Simulations and Free-Energy Calculations Qiao Xue, Ji-Long Zhang, Qing-Chuan Zheng,* Ying-Lu Cui, Lin Chen, Wen-Ting Chu, and Hong-Xing Zhang* State Key Laboratory of Theoretical and Computational Chemistry, Institute of Theoretical Chemistry, Jilin University, Changchun 130023, China S Supporting Information *

ABSTRACT: DEAD-box proteins are the largest family of helicase that are important in nearly all aspects of RNA metabolism. However, it is unclear how these proteins recognize and bind RNA. Here, we present a detailed analysis of the related DEAD-box protein Mss116p-RNA interaction, using molecular dynamics simulations with MM-GBSA calculations. The energetic analysis indicates that the two strands of double strands RNA (dsRNA) are recognized asymmetrically by Mss116p. The strand 1 of dsRNA provides the main binding affinity. Meanwhile, the nonpolar interaction provides the main driving force for the binding process. Although the contribution of polar interaction is small, it is vital in stabilizing the protein−RNA interaction. Compared with the wild type Mss116p, two studied mutants Q412A and D441A have obviously reduced binding free energies with dsRNA because of the decreasing of polar interaction. Three important residues Lys409, Arg415 and Arg438 lose their binding affinity significantly in mutants. In conclusion, these results complement previous experiments to advance comprehensive understanding of Mss116p-dsRNA interaction. The results also would provide support for the application of similar approaches to the understanding of other DEAD-box protein-RNA complexes.

■

short N-terminal extension (NTE), followed by an α-helical Cterminal extension (CTE) and a basic tail (Figure 1).12 In previous research, a series of biochemical experiments have proposed the preconditions and the mechanism of the unwinding process.13−15 ATP binding is required, which would favor the closed conformation of helicase cores.16 Then D1 would wedge dsRNA and utilize the energy of ATP hydrolysis to separate it.13 However, most of these studies focus on the combined RNA-ATP-Mss116p complex. At this point, the major unanswered questions are the following: How does Mss116p recognize and bind its native substrates? What are the functions of the domains and conserved motifs in the binding mechanism? Mallam and co-workers have proposed that D1 acts as an ATP-binding domain and D2 functions as a RNAduplex recognition domain.17 Mohr and co-workers have also found that the CTE domain functions to stabilize the RNA, which is important in recognition.12 In addition, Mallam and co-workers also solve the crystallographic structure of the Mss116p (contains D2 and CTE domains)-dsRNA complex (Figure 1). However, the static crystal structure and the limited resolution are insufficient to provide enough atomic-level

INTRODUCTION DEAD-box proteins are ATP-dependent RNA helicases that modulate RNA and ribonucleoprotein (RNP) structures in nearly all aspects of RNA metabolism, such as translation, RNA splicing, ribosome assembly, RNA degradation, and nuclear transport.1,2 Each DEAD-box protein contains a conserved helicase core that consists of two tandem RecA-like domains (D1 and D2) connected by a short linker. Many proteins also contain N-/C-terminal extensions that may target proteins to specific substrates by RNA−protein or protein−protein interactions.3,4 DEAD-box proteins unwind RNA duplexes by locally nonprocessive strand separation.5,6 They are different from processive RNA helicases that unwind RNA duplexes by translocation through the duplexes. This unwinding mechanism is helpful in amending the local structures of RNA and the RNA−protein complex without globally disrupting the whole structure. Related DEAD-box protein Mss116p of Saccharomyces cerevisiae is often used as an important model system for studying the function and mechanism of DEAD-box protein.7−9 Mss116p functions as an RNA chaperone that facilitates the splicing of mitochondrial group I and group II introns.10 This protein can bind RNA nonspecifically and unwind doublestranded RNA (dsRNA) to disrupt the stable but inactive RNA structures.11 In Mss116p, the helicase core is preceded by a © 2013 American Chemical Society

Received: June 21, 2013 Revised: July 27, 2013 Published: July 29, 2013 11135

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

■

Article

METHODS Molecular Systems Preparation. The Mss116p-dsRNA structure for MD simulations was taken from the Protein Data Bank (PDB ID 4DB2). The structures of Mss116p with singlestranded RNA as well as mutants Q412A and D441A were generated from the crystal structure using Discovery Studio 2.5 software.25 All crystal water molecules were considered during the MD simulations. MD Simulations. Explicit solvent MD simulations were performed using the AMBER 11 program,26 with the ff10 force field (RNA)27,28 and ff99SBildnnmr force field (protein).29,30 force field ff10 collects a variety of updates and modifications that lead to a better representation of RNA. For protein, ff99SBildnnmr has two recent modifications compared to ff10. The first, ff99SBildn, changes side-chain torsions for several amino acids. The second, ff99SBnmr, suggests modifications to backbone torsions based on NMR data. ff99SBildnnmr has high accuracy in MD simulations compared to experimental data.31 Each system was solvated in a box of a TIP3P water model with at least a 10 Å distance around the complex. To keep system neutral, Na+ counterions were added. Whole systems were subjected to 2000 steps of steepest decent minimization followed by 3000 steps of conjugate gradient minimization. Subsequently, the systems were heated from 0 to 300 K in 300 ps by Langevin dynamics with a collision frequency of 1 ps−1.32 Then, the systems went through 500 ps equilibrium MD simulations. Finally, a total of 20 ns was simulated for each system under NPT ensemble conditions using periodic boundary conditions and particle-mesh Ewald (PME) for long-range electrostatics.33,34 Short-range interactions were cut off at 12 Å, and bonds involving hydrogen were held fixed using SHAKE.35,36 The time step was set to 2 fs. All trajectories were recorded every 2 ps.37 All of the figures were created with Pymol and Chimera.38,39 The hydrogen bond is defined if the distance between the donor and acceptor is shorter than 3.5 Å and the angle is less than 60°. The ionic bond is defined if the distance is less than 5 Å. The ionic bond is considered to be formed between ribonucleotides and positively charged amino acids. MM-GBSA Calculations. A total of 5000 snapshots were extracted from each equilibrated trajectory to perform MMGBSA calculations.40,41 In this method, the binding free energy is computed as follows:

Figure 1. Schematic of the Mss116p (contains D2 and CTE)-dsRNA complex. D2, yellow; CTE, blue; dsRNA, green; conserved motifs, orange.

information about the molecular basis of RNA recognition by Mss116p. Besides, some studies also indicate that mutations of residues such as Q412A and D441A would influence the unwinding activity.18,19 Why and how do the mutants reduce the dsRNA binding? So far, little is known about the essence of dsRNA binding, which is the initial step in the whole unwinding process. To understand the mechanism comprehensively, studies about the origin of the binding affinity in the Mss116p-dsRNA complex are necessary. In recent years, molecular dynamics (MD) simulations combined with binding free-energy calculations have been proven to be successful in researching protein−protein, protein−DNA, and protein−RNA interactions.20−23 This method can provide not only plentiful dynamic structural information about the protein-RNA complex but also abundant energy information. Detailed information would contribute to the further understanding of the essence of protein−RNA interactions and demonstrate the origin of the binding affinity.24 To characterize the recognition and binding mechanism of the Mss116p-dsRNA interaction and identify the key residues in the interaction, 20 ns MD simulations followed by the molecular mechanics generalized Born surface area (MM-GBSA) calculation were performed on the Mss116pdsRNA complex. Furthermore, MD simulations of mutants Q412A and D441A with dsRNA complexes were executed to illustrate the origin of the reduced binding affinity. These simulations and MM-GBSA calculations would complement the experiments for a better understanding of the binding mechanism of the Mss116p-dsRNA complex by providing atomic details and experimentally inaccessible dynamic conformations. The application of computational methods in explaining the structural basis of Mss116p recognition would provide strong motivation for the use of similar approaches to address this kind of question regarding other DEAD-box proteins.

ΔG bind = Gcomplex − Gprotein − G RNA

(1)

Gcomplex, Gprotein, and GRNA are the free energy of the complex, protein, and RNA, respectively. The free energy, G, can be calculated using the scheme as follows: G = EMM + Gsol − TS

(2)

EMM = E int + Eele + Evdw

(3)

Gsol = GGB + GSA

(4)

EMM, Gsol, and TS represent the molecular mechanics component in the gas phase, the stabilization energy due to solvation, and a vibrational entropy term, respectively. EMM is given as the sum of Eint, Eele, and Evdw, which are the internal, Coulomb, and van der Waals interaction terms, respectively. Gsol is the solvation energy and can be separated into an electrostatic solvation free energy (GGB) and a nonpolar solvation free energy (GSA). GGB can be calculated with the generalized Born (GB) method.42 GSA is considered to be 11136

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

proportional to the molecular solvent-accessible surface area (SASA) buried on binding.43 Normal-mode analysis was performed to estimate the change in conformational entropy upon RNA binding (TS) using the nmode module of AMBER 11. Because the normal-mode analysis is computationally expensive, only 32 snapshots of the 5000 snapshots were extracted equally for the entropy calculation.44,45 To obtain a detailed view of Mss116p and RNA interactions, the MM-GBSA method was also used to calculate the binding free energy of each residue at the interface between protein and RNA. Each residue can be partitioned into two groups: backbone and side chain. The snapshots used in the binding free-energy decomposition are the same as those used in the binding free-energy calculation.

■

RESULTS Stability of the Binding Pocket and Fluctuations of Two Conserved Motifs. The systems quickly achieve equilibrium after 1 ns (Figure S1) The RMSD value of Mss116p in the Mss116p-dsRNA complex is 1.22 Å (standard deviation 0.13). The results of root-mean-square fluctuations (RMSF) are displayed in Figure 2. It is easy to notice that the

Figure 3. Molecular surface of Mss116p with bound dsRNA colored by the electrostatic potential (blue, positive; red, negative). dsRNA is shown in red.

ssRNA in which Mss116p was in a closed state.13 We superimpose this structure on the representative structure of our simulation (Figure 4). The flexibilities of these two motifs

Figure 2. RMSF relative to the initial structure for backbone atoms of the Mss116p-dsRNA complex and Mss116p.

RMSF values of most residues are less than 1 Å. For the conserved motifs, IV (residues 378−384), IVa (408−415), V (429−434), and Vb (448−454) have little flexibility whereas Va (437−443) and VI (462−469) motifs are elastic. Figure 3 presents the electrostatic potential at the Mss116p-dsRNA complex surface, which shows an evident positively charged binding pocket. This binding pocket is primarily composed of the three conserved motifs IV, IVa, and V. Similar RMSF values between the Mss116p protein and Mss116p-dsRNA complex indicate that the binding pocket is conserved whether RNA binds Mss116p. We have also analyzed the secondary structures during the simulation (Figure S2.) It is obvious that all of the secondary structures were stable during the simulation. On the basis of the above structural analysis, the global protein structure and the overall topology for Mss116p, especially the binding pocket, are reserved. From Figure 2, it is obvious that motifs Va and VI fluctuate significantly during the simulations. Va is located at the interface of D1, D2, and dsRNA. VI is located at the interface of D1 and D2. Recent research solved the crystal structure of the helicase cores of Mss116p bound to adenosine nucleotide and

Figure 4. Superposition of the representative structure of the Mss116p-dsRNA complex (Mss116p, yellow) to the closed-state Mss116p−ssRNA−ADP complex (PDB ID 3I62) (Mss116p, blue; ADP, orange).

may be to avoid the clash with D1 when Mss116p underwent a conformational transition. In addition, the high flexibility would help these two motifs make up the ATP binding site in “closed” Mss116p. Arg438 in motif Va should be noticed because it has been deemed to interact strongly with dsRNA. The RMSF value of Arg438 (2.32 Å in the Mss116p-dsRNA complex and 1.62 Å in the Mss116p protein) also shows that it has a different fluctuation in the Mss116-dsRNA complex and the Mss116p single protein. Figure 5 monitors the dihedral changes of Arg438 in the Mss116p-dsRNA complex and Mss116p. A comparison of the representative structures of the Mss116pdsRNA complex to Mss116p indicates that the side chain of Arg438 is flexible in Mss116p whereas it is fixed in the Mss116p-dsRNA complex. 11137

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

Arg415 to C9, and Lys409 to G7, should be noticed. These ionic bonds were not shown in the crystal structure although they may be important to the Mss116p-dsRNA interaction. In addition, most of the ionic bonds formed between the residues belong to the binding pocket and Va with strand 1. Notably, Lys414 formed one ionic bond with G13′ and Lys542 formed two ionic bonds with C2′ and C3′ in strand 2 about half of the time. Though these ionic bonds are not as stable as above, they reflect the function of the CTE domain in binding RNA, which has been proven by experiment.12 The hydrogen bonds in the crystal structure were also observed to be stable during the MD simulation, and most of them formed between D2 and strand 1. Remarkably, Arg438 also formed a hydrogen bond with strand 1, which is not observed in the crystal structure. This hydrogen bond may induce the different RMSF fluctuations and dihedral changes in Arg438 in the Mss116p protein and the Mss116p-dsRNA complex. Calculation of Binding Free Energy. To investigate the source of the recognition and binding process, simulations of Mss116p-dsRNA, Mss116p-strand 1, and Mss116p-strand 2 were performed for 20 ns, and binding free energies were computed for 5000 snapshots using the MM-GBSA method. To minimize potential bias from the initial structure, the first 10 ns of each simulation was excluded from further analysis. Overall binding free energies were broken down into contributions from electrostatics, van der Waals interactions, nonpolar solvation interactions, and electrostatic solvation interactions. The contribution of entropy is also considered. For all computed energies, the results are presented as ensemble-averaged values over all 5000 snapshots. Decomposing the binding free energies would help us identify the key residues in the binding process. Generally, the MM-GBSA method is not able to reproduce the absolute value of the experimental results exactly.46 However, there is a good correlation between the calculated binding free energy by this method and experimental results, some of which are almost linear. In other words, this method can give a good rank of binding free energy for different systems.47,48 In this study, our aim is to use MM-GBSA calculations to investigate the sources dominating binding rather than to calculate the absolute binding free energies. Thus, we are more concerned about the contributions of different energy items and key residues than the absolute binding free energy. Because the MM-GBSA method can give relatively accurate binding free energies as well as the decomposing energy items and individual contributions, it is

Figure 5. Dihedral changes (psi) during MD simulations in the Mss116p-dsRNA complex and the Mss116p protein.

Intermolecular Interactions. To investigate the interaction characteristics between Mss116p and dsRNA, we analyzed the ionic bonds and hydrogen bonds formed at the interface of Mss116p and dsRNA (Figure 6). It is clear that

Figure 6. Ionic bond and hydrogen bond occupancies between Mss116p and dsRNA in the Mss116p-dsRNA complex.

several ionic bonds existed for a long time during the MD simulations. Arg415 formed an ionic bond with C8 during the whole 20 ns, and this ionic bond has been confirmed in the crystal structure. Three other strong ionic bonds, Arg438 to C9,

Table 1. Binding Free Energy and Its Components of Mss116p-dsRNA, Mss116p-strand1, Mss116p-strand 2, Q412A-dsRNA, and D441A-dsRNA Complexes (kcal/mol) ΔEvdw ΔEele ΔGGB ΔGSA ΔGpola ΔGnonpolb ΔGc TΔS ΔΔGd a

dsRNA

strand 1

strand 2

Q412A

D441A

−109.47 ± 7.14 −2820.53 ± 114.56 2818.59 ± 104.96 −17.58 ± 0.79 −1.94 −127.05 −128.99 ± 17.57 −59.6 ± 6.62 −69.39 ± 18.94

−107.04 ± 7.73 −1979.67 ± 74.34 1978.17 ± 71.17 −14.02 ± 0.93 −1.5 −121.06 −122.56 ± 12.80 −51.57 ± 6.58 −70.96 ± 14.40

−70.08 ± 7.00 −1353.95 ± 81.85 1385.68 ± 80.70 −10.28 ± 0.96 31.73 −80.36 −48.63 ± 9.97 −48.69 ±4.79 0.06 ± 11.02

−104.07 ± 7.06 −2489.38 ± 110.01 2510.82 ± 105.35 −14.87 ± 0.71 21.44 −118.94 −97.51± 10.35 −47.72± 5.66 −49.78 ± 11.80

−113.54 ± 6.73 −3133.95 ± 130.58 3152.73 ± 125.21 −17.39 ± 0.85 18.78 −130.93 −112.15 ± 13.39 −60.28 ± 4.86 −51.87 ± 14.25

ΔGpol = ΔEele + ΔGGB bΔGnonpol = ΔEvdw + ΔGSA cΔG = ΔEele + ΔGGB + ΔEvdw + ΔGSA dΔΔG = ΔG − TΔS 11138

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

Table 2. Ribonucleotides Have More Than a 2 kcal/mol Binding Affinity (kcal/mol) ribonucleotide G5 G6 G7 C8 C9 C10 C3′ a

ΔEvdw −6.70 −6.27 −5.50 −3.73 −1.83 −2.60 −4.40

± ± ± ± ± ± ±

2.74 2.66 2.79 2.99 2.73 2.69 2.90

ΔEele −62.78 −88.49 −104.67 −104.95 −106.03 −83.32 −51.46

± ± ± ± ± ± ±

ΔGGB 11.46 9.93 10.11 6.73 8.71 7.73 12.67

67.73 91.51 103.78 105.31 106.33 83.94 50.43

± ± ± ± ± ± ±

ΔGSA 9.92 8.29 7.62 5.34 6.77 6.11 9.46

−0.79 −0.86 −0.86 −0.68 −0.50 −0.54 −0.67

± ± ± ± ± ± ±

0.08 0.06 0.05 0.07 0.08 0.06 0.09

ΔGpola

ΔGnonpolb

ΔGc

4.94 3.01 −0.88 0.36 0.30 0.61 −1.03

−7.49 −7.13 −6.37 −4.42 −2.34 −3.15 −5.08

−2.55 ±1.59 −4.11 ± 2.54 −7.26 ± 2.34 −4.06 ± 1.91 −2.04 ± 1.80 −2.53 ± 1.89 −6.11 ± 2.13

ΔGpol = ΔEele + ΔGGB. bΔGnonpol = ΔEvdw + ΔGSA. cΔG = ΔEele + ΔGGB + ΔEvdw + ΔGSA.

Figure 7. Contribution of each residue and ribonucleotide to wild-type Mss116p-dsRNA binding.

that the core ribonucleotides that interact intermolecularly with the binding pocket make notable contributions to the binding affinity. From Table 2 and Figure 7, we can see there are seven ribonucleotides in dsRNA that make more than a 2 kcal/mol free-energy contribution to binding. It is easy to notice that G6, G7, and C8 with more than a 4 kcal/mol contribution to binding make strong ionic bonds or hydrogen bonds with the binding pocket. G7 makes the largest contribution to binding, and it forms a strong ionic bond with Lys409 and two hydrogen bonds with Gly408 and Thr433. G6 forms an ionic bond with Lys409 and a hydrogen bond with Val383 whereas C8 forms an ionic bond with Arg415 and hydrogen bonds with Arg415 and Gly436 (Figure 6). From the electrostatic potential figure (Figure 3), we can see that these residues act as a clamp to fix the core ribonucleotides. In addition, C3′ also makes a large contribution to the binding affinity. This ribonucleotide makes an ionic bond with Lys542 in the CTE domain. This reflects the fact that the CTE domain also has an important function in dsRNA binding. Dominant Residues and the Contribution of Conserved Motifs. To identify the hotspot residues of the protein−dsRNA interaction interface, decomposing energy analysis was employed for the simulation trajectory and the corresponding result is given in Table S1. It is clear that 15 residues of Mss116p make more than a 2 kcal/mol free-energy contribution to binding. Figure 8 presents the visualized representation for the distribution and contributions of the 15 residues. From the distribution, 12 of these residues are located in the D2 domain and the other 3 residues belong to the CTE domain. Notably, the residues in different domains undergo different interactions with dsRNA. The 12 residues in the D2 domain all undergo strong favorable electrostatic interactions

suitable to investigate the source of the binding process in the protein-dsRNA complex.49,50 Different Binding Affinities of Strand 1 and Strand 2. Table 1 details the binding free energies for dsRNA, strand 1, and strand 2 to Mss116p. When the contribution of entropy is considered, the binding free energy of strand 2 to Mss116p is almost zero, which means that Mss116p almost has no binding affinity to strand 2. However, the binding free energy of strand 1 is much lower than that of strand 2. Obviously, strand 1 is the main RNA chain that binds with Mss116p. By decomposing the binding free energy into different items, it is clear that the differences mostly come from the contribution of van der Waals and hydrophobic interactions to binding free energies. The nonpolar or hydrophobic interactions (ΔGnonpol = ΔEvdw + ΔGSA) with −127.05 kcal/mol provide the main driving force for binding. In comparison to nonpolar contributions, the electrostatic interaction (ΔGpol = ΔEele + ΔGGB) has a much smaller but favorable contribution. Actually, the direct intermolecular electrostatic interactions are highly favorable to binding, but their contributions are counteracted by the large desolvation penalties associated with the binding process. From the results of the individual energy contribution for dsRNA and the single strand, we could determine the origin of the binding affinity difference. Both the polar and nonpolar contributions of binding energies are largely reduced for strand 2 compared to that for strand 1. The nonpolar contribution decreased by 40.7 kcal/mol, and the polar interaction contribution decreased by 33.23 kcal/mol. Importance of Core Ribonucleotides. The result above shows that nearly all of the binding affinity for dsRNA can be attributed to strand 1. When the binding free energies are decomposed into individual ribonucleotides, the result indicates 11139

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

Q412A and D441A were performed for 20 ns in order to explore the origin of the binding affinity loss. As expected, the mutated systems have a smaller binding free energy than wild-type protein. The binding free energies for wild type, Q412A, and D441A are −69.39, −49.78, and −51.87 kcal/mol, respectively (Table 1). By comparing the individual energy item for wild type and mutants, we could identify the origin of the binding affinity loss for the mutants. From Table 1, it can be seen that the nonpolar energies of the three systems are similar; that is, the differences between mutant systems and wild type are both less than 9 kcal/mol. But for polar energies including the direct electrostatic interaction and polar solvation interactions, Q412A and D441A decrease by 23.38 and 20.72 kcal/mol, respectively. More importantly, the polar interaction is favorable to binding in wild type but unfavorable in mutants. To explore which residues are responsible for the reduced binding free energy for the mutants, we further decomposed the binding free energy to a per residue value (Figure S3). For mutants, a large reduction in the binding affinity can be observed in Arg438, Arg415, and Lys409. Table 4 lists the detailed energy information for the three residues. It is obvious that side chains of the three residues dominate the binding affinity. For both of the mutants, the loss of the polar contribution of the three residues is the main reason for the decreased binding affinity. To investigate the structural detail changes between mutant and wild type, clustering analysis was applied to extract the representative conformation of the whole ensemble (Figure 9. It is clear that the position of the binding pocket and motif Va are different between mutants and wild type. The pockets of mutants are smaller than those of wild type. The corresponding ribonucleotides of mutants are farther outside than those of wild type.

Figure 8. Distribution of important residues and ribonucleotides and their contributions (kcal/mol).

with dsRNA whereas the contributions of nonpolar interactions are relatively small. However, the residues in the CTE domain are different from those in the D2 domain in that they undergo stronger favorable nonpolar interactions with dsRNA. Remarkably, three positively charged residues Lys409, Arg415, and Arg438 make the largest contributions to binding. These three residues have strong intermolecular interactions with strand 1 (Figure 6). We also calculate the contributions of the conserved motifs. Motifs IV, IVa, V, Va, Vb, and VI give contributions of −14.82, −36.67, −3.98, −14.38, 0.16, and −0.77 kcal/mol, respectively. In all of the conserved motifs, IVa makes the largest contribution. The sum of motifs IV, IVa, V, and Va produces a binding affinity of −69.85 kcal/mol, which includes a large proportion (69.95%) of the total −99.85 kcal/mol. Conserved Side-Chain Interactions as the Primary Source of Affinity. A large fraction of the binding affinity of residues can be attributed to a common set of interactions primarily involving side-chain contacts. The contributions of these groups are detailed in Table 3. Overall, the side chains account for roughly −71.48 kcal/mol of the affinity whereas the backbones account for only about −28.36 kcal/mol. In the contribution of side chains, about a −32.02 kcal/mol contribution comes from polar interactions and −39.46 kcal/ mol comes from nonpolar interactions. In terms of the polar interaction contribution, this is a more favorable contribution than the total overall polar interaction in the Mss116p-dsRNA interaction. Contrast with Mutants. Previous results cannot give the exact reasons for the reduction of the binding affinity of several reported mutants. Accordingly, another two systems including

■

DISCUSSION Here, we used MD simulations, along with binding free-energy calculations, to investigate the structural basis of the Mss116pdsRNA interaction. The study of this system has revealed many important insights into the recognition mechanism of the Mss116p-dsRNA interaction. Previous experiments mostly focus on the entire change in the Mss116p unwinding pathway.11,51 However, it is still unclear how DEAD-box proteins bind dsRNA and what is the driving force of affinity. Unlike the crystal structure, which can give only static information, we take this structure as the initial structure on which to perform MD simulations and present detailed and dynamic atomic-level information of dsRNA binding to the Mss116p protein. Besides, the result of MM-GBSA also presents energy information about the source of the binding affinity, which is hard to obtain from experiments. In addition, relevant MD simulations of mutants binding to dsRNA are performed to identify the origin of reduced binding affinity. During the MD simulations, one of the most significant results is the stability of the whole Mss116p structure. No matter if the protein-RNA complex or a single protein is considered, the protein structure, especially the binding pocket,

Table 3. Energy Determinants of General Affinity (kcal/mol) side chain backbone a

Evdw

Eele

GGB

GSA

Gpola

Gnonpolb

ΔGc

−32.22 −22.51

−3748.91 2338.66

3716.89 −2342.81

−7.23 −1.69

−32.02 −4.15

−39.46 −24.20

−71.48 −28.36

ΔGpol = ΔEele + ΔGGB. bΔGnonpol = ΔEvdw + ΔGSA. cΔG = ΔEele + ΔGGB + ΔEvdw + ΔGSA. 11140

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

Table 4. Binding Free Energies and Individual Energy Terms of Lys409, Arg415, and Arg 438 in Wild Type and Mutants (kcal/ mol) Lys409

Arg415

Arg438

a

Mss116p Q412A D441A Mss116p Q412A D441A Mss116p Q412A D441A

ΔGbackbone

ΔGside chain

ΔEvdw

ΔEele

ΔGGB

ΔGSA

ΔGpola

ΔGnonpolb

0.35 0.27 0.41 0.14 0.04 0.13 −0.04 0.18 0.03

−9.23 −4.96 −3.52 −16.48 −11.54 −12.99 −13.40 −1.78 −5.79

−1.71 −1.38 −1.54 −0.56 0.33 0.17 −0.53 −0.29 −1.78

−342.77 −321.27 −328.96 −313.97 −296.68 −303.89 −275.47 −208.53 −241.25

336.01 318.36 327.69 298.45 285.01 291.00 263.08 207.32 237.64

−0.41 −0.39 −0.3 −0.26 −0.17 −0.15 −0.51 −0.11 −0.37

−6.76 −2.91 −1.28 −15.52 −11.67 −12.89 −12.39 −1.21 −3.61

−2.12 −1.77 −1.84 −0.82 0.16 0.03 −1.04 −0.39 −2.15

ΔGc −8.88 −4.69 −3.12 −16.34 −11.51 −12.86 −13.43 −1.60 −5.76

± ± ± ± ± ± ± ± ±

3.78 2.97 1.08 2.62 3.04 2.19 2.58 1.12 2.45

ΔGpol = ΔEele + ΔGGB. bΔGnonpol = ΔEvdw + ΔGSA. cΔG = ΔEele + ΔGGB + ΔEvdw + ΔGSA.

Figure 9. (a) Comparison of binding pocket and motif Va between wild-type Mss116p (yellow) and Q412A (cyan). (b) Comparison of binding pocket and motif Va between wild-type Mss116p (yellow) and D441A (magenta).

hydrogen bonds between two RNA strands are cut off by D1, strand 2 is easy to separate from the dsRNA whereas strand 1 maintains its position because of the strong binding affinity. This result is in agreement with the experiment in which Mss116p unwinds dsRNA by local strand separation.5 From the comparison between Mss116p-dsRNA and Mss116p, it is easy to determine the structural change in motif Va. The flexibility of this motif is helpful to the binding process because it must allow dsRNA to step aside to bind to the binding pocket. The changing Arg438 dihedral proves this. After the binding of dsRNA, the side chain of Arg438 changes its orientation to form a strong hydrogen bond and ionic bond with strand 1. The CTE domain also functions to help protein bind dsRNA. During the MD simulations, we can see an ionic bond that is formed between the CTE domain and dsRNA. Interestingly, this ionic bond formed with strand 2 not strand 1, which is far from the interface of D1 and D2. This would help to stabilize the dsRNA and not influence the subsequent unwinding. This result is consistent with the experiment in which the CTE domain contributes to the binding and unwinding process.9 By the calculation of the binding free energy using the MMGBSA method, our results present the source of recognition

is always very stable. This stable structure is favorable in the binding process for several reasons. First, the pocket is mostly composed of positively charged residues that can bind with ribonucleotides by strong electrostatic interactions. Second, the entire unwinding mechanism of the DEAD-box protein relies on the transformation of the protein that D1 undergoes from the open state to the closed state to unwind the dsRNA. Then one strand of dsRNA would depart from the protein and leave a single strand of RNA.17 The rigid binding pocket would ensure that Mss116p undergoes a stable interaction with strand 1 during the whole unwinding process. In addition, the binding free-energy results are consistent with the unwinding mechanism in which strand 1 is the primary source of binding with Mss116p compared to strand 2. The total binding free energies of dsRNA, strand 1, and strand 2 are −69.39, −70.96, and 0.06 kcal/mol, respectively. This result reflects that strand 1 has much stronger binding affinity for Mss116p than does strand 2. By decomposing the energy to individual ribonucleotides (Figure 7), we can see that the ribonucleotides that make large contributions to binding are mostly concentrated in strand 1. Strand 2 has little binding affinity for the D2 domain. The only interactions between two strands are hydrogen bonds. Therefore, when the protein is in the closed state and the 11141

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

strand in the binding process. On the basis of the result of MMGBSA, the contribution of the nonpolar interaction is the main source of the binding affinity. Although the contribution of the polar interaction is small in wild-type Mss116p, it plays an important role in stabilizing the protein−RNA interaction and is the origin of a different binding affinity for mutants. For Q412A and D441A mutants, the polar contribution is largely reduced whereas the nonpolar contribution remains the same as that in wild type. From the decomposition of free energy to a per residue value, conserved motifs IV, IVa, and V in the binding pocket and motif Va contribute about 70% of the total binding free energy. Several key residues are identified to be important for the Mss116p-dsRNA interaction. In these residues, Lys409, Arg415, and Arg438 are considered to play key roles in binding. In mutants, these three residues lose their binding affinity because of the shrinkage of the binding pocket in the mutants. In conclusion, the rigid structure of the binding pocket and several important residues would fix strand 1 with strong intermolecular interactions. Strand 1 is the pivotal strand in the Mss116p-dsRNA binding process. The binding affinity of strand 2 comes from the interaction of the CTE domain, which is very weak compared to that of strand 1. The diversity of the binding affinity in strand 1 and strand 2 is the essence of recognition and would be beneficial to the next unwinding process. Detailed information about the Mss116p-dsRNA interaction obtained from the present simulations is meaningful to the further understanding of the unwinding mechanism.

and the binding process. We also obtain the energy contributions from individual functional components: van der Waals energy, electrostatic interaction energy, and solvation energy. First, a large portion of the affinity is derived from nonpolar energy in which the van der Waals contribution is the most dominant term. The electrostatic energy makes a much smaller contribution whereas the intermolecular electrostatic energy is extremely favorable. Besides, most of the affinity comes from the side chains of residues. Second, when the binding free energies are decomposed to per residue values, the residues with important functions in binding are identified. By dividing these residues into the conserved motifs, we find that IV, IVa, V, and Va make the most significant contributions to the binding affinity. The binding affinity of these four motifs is about 70% of the total binding affinity whereas the contribution of motifs IVa and Va is about 51% of the total binding affinity. Only three residuesLys409, Arg415, and Arg438make greater than 8 kcal/mol contributions to binding. These three residues, which belong to motifs IVa and Va, are all situated in the primary binding domain that includes the binding pocket and motif Va and make strong hydrogen bonds or ionic bonds with the core ribonucleotides of strand 1. From Table 4, we can see that the polar interactions make the largest contributions to their binding. In conclusion, nonpolar energy is the main driving force for the entire binding and polar energy also has an important function in stabilizing the Mss116p-dsRNA interaction. By performing the MD simulations for two mutants Q412A and D441A, we find that the mutations partially decrease the binding affinity. dsRNA can bind to mutants with a lower binding affinity compared to that of wild type. By decomposing the binding free energy to per residue values, the contributions of Lys409, Arg415, and Arg438 are largely reduced. The binding affinity of Arg438 almost disappeared in both mutants. Through the comparison of the representative conformations of wild type and mutants, we can infer that the pocket between IVa and Va shrinks in the mutants. This shrinkage would increase the distance between the protein and RNA and then induce the loss of electrostatic interaction of motifs IVa and Va, especially Lys409, Arg415, and Arg438. Compared to wild-type Mss116p, motifs Va of both mutants are less flexible (Figure S4). Va of both mutants is unable to bind dsRNA as well as wild-type Mss116p. This result also proves that dsRNA is unable to bind closely with a mutant. Although the electrostatic interaction makes a much smaller contribution than the nonpolar interaction in the Mss116p binding to dsRNA, it is responsible for the loss of binding affinity of the mutants. The smaller binding pocket of mutant induces the loss of binding affinity. This result may further result in an insufficient unwinding activity that is consistent with the experiment in which Q412A and D441A do not have the ability to unwind dsRNA.19,52

■

ASSOCIATED CONTENT

S Supporting Information *

Residues giving more than a 2 kcal/mol binding affinity. Timedependent RMS deviation values from the starting structures along the simulation trajectory. Second structure variation along MD simulations. Contribution of each residue and ribonucleotide to Q412A-dsRNA and D441AdsRNA complexes. RMSF relative to the initial structure for backbone atoms of Mss116p-dsRNA, Q412AdsRNA, and D441A-dsRNA complexes. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]; [email protected]. Fax: (+86) 431-8849-8962. Tel: (+86) 431-8849-8962 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work is supported by the Natural Science Foundation of China. (grant nos. 21273095, 20903045, and 21203072)

■

■

CONCLUSIONS Through explicit-solvent MD simulations combined MMGBSA calculations, we have been able to explore the essence of Mss116p-dsRNA interaction and explain a number of important features of the recognition process. The primary binding domain is composed of the stable binding pocket, which is composed of positively charged residues, and a flexible motif Va, which is close to the binding pocket. This binding domain undergoes strong intermolecular interactions with dsRNA, especially strand 1. Strand 1 is the main functional

REFERENCES

(1) Linder, P. Dead-box proteins: a family affairactive and passive players in RNP-remodeling. Nucleic Acids Res. 2006, 34, 4168−4180. (2) Cordin, O.; Banroques, J.; Tanner, N. K.; Linder, P. The DEADbox protein family of RNA helicases. Gene 2006, 367, 17−37. (3) Jankowsky, E.; Fairman, M. E. RNA helicasesone fold for many functions. Curr. Opin. Struct. Biol. 2007, 17, 316−324. (4) Fairman-Williams, M. E.; Guenther, U.-P.; Jankowsky, E. SF1 and SF2 helicases: family matters. Curr. Opin. Struct. Biol. 2010, 20, 313− 324.

11142

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

(5) Yang, Q.; Del Campo, M.; Lambowitz, A. M.; Jankowsky, E. DEAD-box proteins unwind duplexes by local strand separation. Mol. Cell 2007, 28, 253−263. (6) Yang, Q.; Jankowsky, E. The DEAD-box protein Ded1 unwinds RNA duplexes by a mode distinct from translocating helicases. Nat. Struct. Mol. Biol. 2006, 13, 981−986. (7) Mohr, S.; Stryker, J. M.; Lambowitz, A. M. A DEAD-box protein functions as an ATP-dependent RNA chaperone in group I intron splicing. Cell 2002, 109, 769−779. (8) Del Campo, M.; Tijerina, P.; Bhaskaran, H.; Mohr, S.; Yang, Q.; Jankowsky, E.; Russell, R.; Lambowitz, A. M. Do DEAD-box proteins promote group II intron splicing without unwinding RNA? Mol. Cell 2007, 28, 159−166. (9) Del Campo, M.; Mohr, S.; Jiang, Y.; Jia, H.; Jankowsky, E.; Lambowitz, A. M. Unwinding by local strand separation is critical for the function of DEAD-box proteins as RNA chaperones. J. Mol. Biol. 2009, 389, 674−693. (10) Huang, H. R.; Rowe, C. E.; Mohr, S.; Jiang, Y.; Lambowitz, A. M.; Perlman, P. S. The splicing of yeast mitochondrial group I and group II introns requires a DEAD-box protein with RNA chaperone function. Proc. Natl. Acad. Sci. U.S.A. 2005, 102, 163−168. (11) Potratz, J. P.; Del Campo, M.; Wolf, R. Z.; Lambowitz, A. M.; Russell, R. ATP-dependent roles of the DEAD-box protein Mss116p in group II intron splicing in vitro and in vivo. J. Mol. Biol. 2011, 411, 661−679. (12) Mohr, G.; Del Campo, M.; Mohr, S.; Yang, Q.; Jia, H.; Jankowsky, E.; Lambowitz, A. M. Function of the C-terminal domain of the DEAD-box protein Mss116p analyzed in vivo and in vitro. J. Mol. Biol. 2008, 375, 1344. (13) Del Campo, M.; Lambowitz, A. M. Structure of the yeast DEAD box protein Mss116p reveals two wedges that crimp RNA. Mol. Cell 2009, 35, 598−609. (14) Mallam, A. L.; Jarmoskaite, I.; Tijerina, P.; Del Campo, M.; Seifert, S.; Guo, L.; Russell, R.; Lambowitz, A. M. Solution structures of DEAD-box RNA chaperones reveal conformational changes and nucleic acid tethering by a basic tail. Proc. Natl. Acad. Sci. U.S.A. 2011, 108, 12254−12259. (15) Karow, A. R.; Klostermeier, D. A structural model for the DEAD box helicase YxiN in solution: localization of the RNA binding domain. J. Mol. Biol. 2010, 402, 629−637. (16) Liu, F.; Putnam, A.; Jankowsky, E. ATP hydrolysis is required for DEAD-box protein recycling but not for duplex unwinding. Proc. Natl. Acad. Sci. U.S.A. 2008, 105, 20209−20214. (17) Mallam, A. L.; Del Campo, M.; Gilman, B.; Sidote, D. J.; Lambowitz, A. M. Structural basis for RNA-duplex recognition and unwinding by the DEAD-box helicase Mss116p. Nature 2012, 490, 121−125. (18) Mohr, G.; Del Campo, M.; Turner, K. G.; Gilman, B.; Wolf, R. Z.; Lambowitz, A. M. High-throughput genetic identification of functionally important regions of the yeast DEAD-box protein Mss116p. J. Mol. Biol. 2011, 413, 952−972. (19) Bifano, A. L.; Turk, E. M.; Caprara, M. G. Structure-guided mutational analysis of a yeast DEAD-box protein involved in mitochondrial RNA splicing. J. Mol. Biol. 2010, 398, 429−443. (20) Pan, D.; Sun, H.; Shen, Y.; Liu, H.; Yao, X. Exploring the molecular basis of dsRNA recognition by NS1 protein of influenza A virus using molecular dynamics simulation and free energy calculation. Antiviral Res. 2011, 92, 424−433. (21) Xue, W.; Liu, H.; Yao, X. Molecular mechanism of HIV-1 integrase−vDNA interactions and strand transfer inhibitor action: a molecular modeling perspective. J. Comput. Chem. 2012, 33, 527−536. (22) Chen, L.; Zhang, J.-L.; Yu, L.-Y.; Zheng, Q.-C.; Chu, W.-T.; Xue, Q.; Zhang, H.-X.; Sun, C.-C. Influence of hyperthermophilic protein cren7 on the stability and conformation of DNA: insights from molecular dynamics simulation and free energy analysis. J. Phys. Chem. B 2012, 116, 12415−12425. (23) Qin, F.; Ye, W.; Chen, Y.; Chen, X.; Li, Y.; Zhang, J.; Chen, H.F. Specific recognition between intrinsically disordered LEF and DNA. Phys. Chem. Chem. Phys. 2012, 14, 538−545.

(24) Singh, A.; Snyder, S.; Lee, L.; Johnston, A. P.; Caruso, F.; Yingling, Y. G. Effect of oligonucleotide length on the assembly of DNA materials: molecular dynamics simulations of layer-by-layer DNA films. Langmuir 2010, 26, 17339−17347. (25) Studio D, version 2.5.; Accelrys Inc.: San Diego, CA, 2009. (26) Case, D.; Darden, T.; Cheatham , T., III; Simmerling, C.; Wang, J.; Duke, R.; Luo, R.; Walker, R.; Zhang, W.; Merz, K. AMBER 11; University of California: San Francisco, 2010; p 142. (27) Banás,̌ P.; Hollas, D.; Zgarbová, M.; Jurečka, P.; Orozco, M.; Cheatham, T. E., III; Šponer, J. i.; Otyepka, M. Performance of molecular mechanics force fields for RNA simulations: stability of UUCG and GNRA hairpins. J. Chem. Theory Comput. 2010, 6, 3836− 3849. (28) Swadling, J. B.; Suter, J. L.; Greenwell, H. C.; Coveney, P. V. Influence of surface chemistry and charge on mineral−RNA interactions. Langmuir 2013, 29, 1573−1583. (29) Lindorff-Larsen, K.; Piana, S.; Palmo, K.; Maragakis, P.; Klepeis, J. L.; Dror, R. O.; Shaw, D. E. Improved side-chain torsion potentials for the Amber ff99SB protein force field. Proteins 2010, 78, 1950− 1958. (30) Li, D. W.; Brüschweiler, R. NMR-based protein potentials. Angew. Chem. 2010, 122, 6930−6932. (31) Beauchamp, K. A.; Lin, Y.-S.; Das, R.; Pande, V. S. Are protein force fields getting better? A systematic benchmark on 524 diverse NMR measurements. J. Chem. Theory Comput. 2012, 8, 1409−1414. (32) Pal, S.; Bandyopadhyay, S. Importance of protein conformational motions and electrostatic anchoring sites on the dynamics and hydrogen bond properties of hydration water. Langmuir 2013, 29, 1162−1173. (33) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577−8593. (34) Tosaka, R.; Yamamoto, H.; Ohdomari, I.; Watanabe, T. Adsorption mechanism of ribosomal protein L2 onto a silica surface: a molecular dynamics simulation study. Langmuir 2010, 26, 9950− 9955. (35) Johnson, S.; Bronowska, A.; Chan, J.; Evans, D.; Davies, A. G.; Wälti, C. Redox-induced conformational change in mercaptoalkanoic acid multilayer films. Langmuir 2012, 28, 6632−6637. (36) Feng, F.; Lee, S. H.; Cho, S. W.; Kömürlü, S.; McCarley, T. D.; Roitberg, A.; Kleiman, V. D.; Schanze, K. S. Conjugated polyelectrolyte dendrimers: aggregation, photophysics, and amplified quenching. Langmuir 2012, 28, 16679−16691. (37) Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H. J. Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. J. Comput. Phys 1977, 23, 327−341. (38) DeLano,W. L. The PyMOL Molecular Graphics System, 2002. (39) Pettersen, E. F.; Goddard, T. D.; Huang, C. C.; Couch, G. S.; Greenblatt, D. M.; Meng, E. C.; Ferrin, T. E. UCSF Chimeraa visualization system for exploratory research and analysis. J. Comput. Chem. 2004, 25, 1605−1612. (40) Kollman, P. A.; Massova, I.; Reyes, C.; Kuhn, B.; Huo, S.; Chong, L.; Lee, M.; Lee, T.; Duan, Y.; Wang, W. Calculating structures and free energies of complex molecules: combining molecular mechanics and continuum models. Acc. Chem. Res. 2000, 33, 889−897. (41) Yao, J.; Chu, Y.; An, R.; Guo, H. Understanding product specificity of protein lysine methyltransferases from QM/MM molecular dynamics and free energy simulations: the effects of mutation on SET7/9 beyond the Tyr/Phe switch. J. Chem. Inf. Model 2012, 52, 449−456. (42) Hawkins, G. D.; Cramer, C. J.; Truhlar, D. G. Parametrized models of aqueous free energies of solvation based on pairwise descreening of solute atomic charges from a dielectric medium. J. Phys. Chem. 1996, 100, 19824−19839. (43) Weiser, J.; Shenkin, P. S.; Still, W. C. Approximate atomic surfaces from linear combinations of pairwise overlaps (LCPO). J. Comput. Chem. 1999, 20, 217−230. 11143

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144

Langmuir

Article

(44) Hawkins, G. D.; Cramer, C. J.; Truhlar, D. G. Pairwise solute descreening of solute charges from a dielectric medium. Chem. Phys. Lett. 1995, 246, 122−129. (45) Xu, B.; Shen, H.; Zhu, X.; Li, G. Fast and accurate computation schemes for evaluating vibrational entropy of proteins. J. Comput. Chem. 2011, 32, 3188−3193. (46) Hou, T.; Wang, J.; Li, Y.; Wang, W. 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. 2011, 51, 69−82. (47) Hou, T.; Wang, J.; Li, Y.; Wang, W. Assessing the performance of the molecular mechanics/Poisson Boltzmann surface area and molecular mechanics/generalized Born surface area methods. II. The accuracy of ranking poses generated from docking. J. Comput. Chem. 2011, 32, 866−877. (48) Cui, Y. L.; Zhang, J. L.; Zheng, Q. C.; Niu, R. J.; Xu, Y.; Zhang, H. X.; Sun, C. C. Structural and dynamic basis of human cytochrome P450 7B1: a survey of substrate selectivity and major active site access channels. Chem.Eur. J. 2013, 19, 549−557. (49) Hou, T.; Li, N.; Li, Y.; Wang, W. Characterization of domain− peptide interaction interface: prediction of SH3 domain-mediated protein−protein interaction network in yeast by generic structurebased models. J. Proteome Res. 2012, 11, 2982−2995. (50) Hou, T.; Zhang, W. Characterization of domain-peptide interaction interface: a case study on the amphiphysin-1 SH3 domain. J. Mol. Biol. 2008, 376, 1201. (51) Mallam, A. L.; Del Campo, M.; Gilman, B.; Sidote, D. J.; Lambowitz, A. M. Structural basis for RNA-duplex recognition and unwinding by the DEAD-box helicase Mss116p. Nature 2012, 490, 121−125. (52) Bifano, A. L.; Caprara, M. G. A DExH/D-box protein coordinates the two steps of splicing in a group I intron. J. Mol. Biol. 2008, 383, 667−682.

11144

dx.doi.org/10.1021/la402354r | Langmuir 2013, 29, 11135−11144