J Mol Model (2015) 21: 190 DOI 10.1007/s00894-015-2739-5
ORIGINAL PAPER
The impact of ligands on the structure and flexibility of sulfotransferases: a molecular dynamics simulation study Li Zhao 1 & Pupu Zhang 1,4 & Shiyang Long 1 & Linlin Wang 1,2 & Pu Tian 1,3
Received: 17 March 2015 / Accepted: 15 June 2015 / Published online: 8 July 2015 # Springer-Verlag Berlin Heidelberg 2015
Abstract Sulfotransferases catalyze transfer of the sulfurylgroup (-SO3) from 3′-phosphoadenosine 5′-phosphosulfate (PAPS) to a large number of substrates. They play an important role in phase II metabolic process. The impact of the cofactor (PAPS) on the structure and flexibility of the enzyme has been studied extensively, and the response of the activecap region to cofactor binding was proposed as the molecular basis for substrate selectivity. In this study, individual and cooperative effects of the cofactor and substrate on the structure and flexibility of the enzyme were investigated. Molecular dynamics simulations were performed for four systems, including free enzyme, binary complexes (cofactor or substrate bound enzyme) and ternary complex (both cofactor and substrate bound enzyme). The influence of ligands (the cofactor and the substrate) on the structure and flexibility of the enzyme, especially that of the active-site cap region, was analyzed. Moreover, mutual structural impact of the ligands was examined as well. The results show that the impact of both the cofactor and the substrate was significant. Our study indicated that the substrate, such as lithocholic acid (LCA), participated in regulating the structure and flexibility of the enzyme actively rather than merely being selected passively. Additionally, the observed synergistic effects of the cofactor * Pu Tian [email protected] 1
School of Life Sciences, Jilin University, Changchun, China
2
Ultrasound Department, China-Japan Union Hospital of Jilin University, Changchun, China
3
Key Laboratory for Molecular Enzymology and Engineering, the Ministry of Education, Jilin University, Changchun, China
4
Present address: Zhongshan Ophthalmic Center, Sun Yat-Sen University, Zhongshan, China
and the substrate demonstrated the importance of examining both ligands in understanding enzymes. Keywords Cofactor and substrate . Human cytosolic sulfotransferase 2A1 . Molecular dynamics simulation . Structure and flexibility
Introduction Cytosolic sulfotransferases (SULTs) catalyze the transfer of a sulfonate group from the unique cofactor 3′phosphoadenosine 5′-phosphosulfate (PAPS) to a wide variety of endogenous and xenobiotic substrates [1–3], including drugs, toxic compounds, steroid hormones, and neurotransmitters. They are phase II metabolic enzymes that participate in chemical defense. The catalytic reaction are broadly involved in both modulating biological activity and increasing solubility of substrates [4]. However, there are reports that link sulfonation imbalance to induction of carcinogenic responses [5–7]. Human cytosolic sulfotransferases (hSULTs) could be divided into four families, including hSULT1, hSULT2, hSULT4, and hSULT6 [8, 9]. Many of them have very broad and overlapping substrate specificities as demonstrated by recent progress in structural biology [8]. Allali-Hassani et al. found evidence for structural Bpriming^ of hSULTs active site by cofactor binding, which may further influence the spectrum of substrates [8]. Ilana Berger et al. found high active site plasticity and structural flexibility of hSULT1A1, and suggested that such character enables binding of different acceptors and plays a dominant role in controlling the specificity and activity [10]. Martiny et al. explored the flexibility of hSULT1A1, hSULT1A3, and hSULT1E1 and identified multiple receptor conformations, which were proposed to be the
190 Page 2 of 9
J Mol Model (2015) 21: 190
most suitable ones for ligand binding prediction [11]. Besides substrate selectivity, substrate inhibition of hSULTs was found to be related to the characteristic structure and flexibility of the enzyme [12]. hSULTs possess a stretch of about 30 amino acids long loop, which serves as an active cap, interacts with and encapsulates both the cofactor and substrates [13–15]. Cook and coworkers observed the response of this active-cap to cofactor binding and suggested its important role in substrate selectivity and recognition [12]. Upon binding of the cofactor PAPS, the active-site cap closes and the binding probability of large substrates (e.g., Raloxifene) decreased [13]. When interactions between the cap and their adjacent residues were weakened via mutagenesis, the enzyme no longer discriminated between large and small substrates and gained increased catalytic efficiency toward large substrates [16]. However, up till now, little is known about the influence of substrates on the structure and dynamics of the enzyme, which are fundamental to understand both the catalytic and substrate inhibition mechanisms. To date SULT2A1 has the broadest known substrate spectrum among SULTs [8, 14]. hSULT2A1 is found in numerous tissues and is concentrated in liver. Its substrates include small planar molecules, steroid, cyclic amines, numerous drugs, and small peptides. The crystal structures of SULT2A1 in complex with 3′-phosphoadenosine 5′-phosphate (PAP) [17], dehydroepiandrosterone (DHEA) [18], androsterone (ADT) [19], and lithocholic acid (LCA) are available in protein data bank with PDB codes 1EFH, 1J99, 1OV4, and 3F3Y respectively. Previously, we investigated the impact of ligands on the thermal stability of SULT2A1 with MD simulations [20], and observed that while both ligands increased the thermal
stability of the enzyme, the extent of stabilization is different for different ligands binding scenarios (cofactor, substrate and both). The thermo-denaturation order of studied complexes was consistent with experimental observations [8]. Here extensive molecular dynamics (MD) simulations were used to investigate the impact of ligands (cofactor and substrate) on the structure and flexibility of human SULT2A1 at physiological temperature. Four systems, including the free enzyme, binary complexes (cofactor or substrate LCA bound enzyme), and the ternary complex (enzyme with both cofactor and substrate LCA bound), were explored. Firstly the structure and flexibility variation among four systems were analyzed based on the root-mean-square fluctuation (RMSF) and evolution of secondary structures with time. Secondly, based on rootmean-square deviation (RMSD) of the active-site cap region and inter-residue distances between the active-site cap region and the main part of the enzyme, we characterized the structural variation of the active-site cap region induced by the ligands. Finally, the mutual structural impact of the ligands was analyzed. RMSD of each ligand as well as the interactions between ligands and the enzyme were examined. The results demonstrated the significant influence of ligands, including both the cofactor and the substrate, on the structure and flexibility of the enzyme.
Fig. 1 Schematic representation of the procedures for building the topology and force field parameters of cofactor PAP. The phosphate moiety in 3PO3 and TH5P were joined onto ADE to obtain topology
and force field parameters of PAP. Topology and parameters of ADE 3PO3 and TH5P are available in CHARMM27 force field [23]
System setup and simulation methods Four models (free enzyme, PAP-bound enzyme, LCA-bound enzyme and both PAP and LCA bound enzyme) were derived from the crystal structure of SULT2A1 complex with both PAP and LCA (PDB code: 3F3Y). Missing atoms were added
J Mol Model (2015) 21: 190
Page 3 of 9 190
Fig. 2 Secondary structure for each residue in a free enzyme, b cofactor bound enzyme, c substrate bound enzyme, and d ternary complex. Secondary structure identities were assigned using DSSP [27]. T: bend and hydrogen bonded turn, E: extended strand participated in beta ladder, B: isolated betabridge, H: alpha helix, G: 3-helix, I: 5-helix, C: coil
using MODELLER [21]. All four models were solvated to obtain a cube large enough to allow at least 12 Å of water between protein and cube surfaces, requiring about 12,500 water molecules. The net charge of the systems was neutralized with Na+, and NaCl was added to the box to 0.15 mM. The resulting four systems are denoted 2A1, 2A1 +PAP, 2A1+LCA, and 2A1+PAP+LCA respectively. All MD simulations were performed with NAMD 2.7 software package [22] using CHARMM27 force fields [23] and TIP3P water model [24]. Topology and parameters of PAP were obtained by joining three segments ADE, 3PO3, and TH5P according to standard CHAR MM procedures as shown in Fig. 1. Topology and parameters of LCA were taken from CHARMM CGenFF files. Assuming a neutral PH as reported in literature [8], the cofactor PAP in this study was not protonated and it had four net charges while the substrate LCA was protonated. PAP rather than PAPS was chosen to study
the influence of PAP and PAPS on the structure and flexibility of hSULTs based on the experimental convenience and the belief that their effects are nearly identical [8], and we followed the same practice in our simulations. Non-covalent van der Waals interactions were cut off at 12 Å. Particle mesh Ewald (PME) [25] was used to calculate long range electrostatic interactions. Bond length involving hydrogen atoms was constrained using SHAKE algorithm [26] with the integration time step of 2 fs. The starting configurations were energy-minimized and equilibrated as previously reported [20]. Then 1-μs simulations at 310 K under NVT condition for 2A1, 2A1 + PAP, 2A1 + LCA, and 2A1 + PAP + LCA systems were performed and collectively 4 μs trajectories were generated. Due to the fact that 2A1, 2A1 + PAP and 2A1 + LCA systems are built from the ternary crystal structure, the first 500-ns of each trajectory is discarded,
Table 1 The probabilities of listed segments in α-helical state for 2A1, 2A1+PAP, 2A1+LCA, and 2A1+PAP+LCA. These segments have a complementary probability of existing in the coil state. The uncertainty is the standard error (SE). A sample of 2000 snapshots were randomly
selected from 10,000 snapshots and the average probabilities of segments being in α-helical state was calculated. The SE was the standard deviations of 1000 samples. (With a sample size of 500, SE was smaller than 0.07)
Residues
2A1
2A1+PAP
2A1+LCA
2A1+PAP+LCA
LYS64-SER68 ILE71-ARG74 ILE102-LEU104 LYS107-SER111 GLY169-TRP170 MET171-MET173
0.65±0.03 0.64±0.02 0.54±0.02 0.76±0.03 0.43±0.02 0.83±0.01
0.46±0.03 0.96±0.00 0.35±0.04 0.78±0.02 0.74±0.03 0.89±0.01
0.46±0.03 0.68±0.03 0.35±0.03 0.77±0.02 0.65±0.02 0.85±0.01
0.39±0.04 0.92±0.00 0.37±0.04 0.69±0.03 0.88±0.02 0.93±0.00
190 Page 4 of 9
J Mol Model (2015) 21: 190
and snapshots of the last 500 ns were taken with an interval of 100 ps to analyze the impact of the cofactor and the substrate on the structure and dynamics of the enzyme.
Results and discussion Structure and flexibility of the enzyme To compare the structure and flexibility of the enzyme in the four systems, secondary structure assignment was carried out using DSSP [27] and its evolution with time was shown in Fig. 2. The notations of secondary structure elements are the same as reported in literature [8, 9]. Overall, there is no obvious difference of secondary structure alignments among four systems except that some peptide segments interconverting between α-helical and coil states exhibited different weight allocation in these two forms. These segments include Lys64-Ser68 in α4, Ile71-Arg74 in α5, Lys107-Ser111 in α7, and Trp163-Met173 in α11. Their probabilities in αhelical state are shown in Table 1. Lys64-Ser68 in α4 and Ile102-Leu104 in α7 in free enzyme tend to exist in αhelical state while they tend to be in coil state in binary or ternary complexes. The probability of Ile71-Arg74 in α5 as αhelix in free enzyme and substrate-bound systems is slightly higher than that in the cofactor-bound and the ternary system. The N-terminal two residues of α11 tend to be in α-helical state in ternary complex, while they exhibit coil state in the free enzyme; the C-terminal three residues possess higher tendency to form α-helix in ternary complex than the other three systems. RMSF shown in Fig. 3 was calculated after minimizing rotation and translation of all backbone heavy atoms. At first glance, free enzyme exhibits larger flexibility as reflected by larger RMSF, especially in the segment Asn226 to Gly252, which is the active-site cap region discussed by Ian Cook [14, 15], while RMSF of the enzyme in ternary system is the smallest. Binding of the cofactor was observed to decrease RMSF of some segments (Trp72-Glu79, Lys188-Asp190, Asn217-Lys248, Phe258-Leu269, and Glu280-Glu285) and increase that of the segment (Ser88-Arg94). The effect of substrate binding on the enzyme RMSF is less extensive, with decrease for the segment Phe220-Ser259 and increase for segments Val69-Ser92, Gly156-Pro172, and Ala275-Phe282 respectively. It is worth noting that cofactor exhibits significant influence on the substrate bound enzyme almost throughout the whole sequence, and particularly for residues Val69Ala86, Gln189-Asp190, and Glu207-Asp267, while the impact of substrate on cofactor bound enzyme is limited to decreasing the flexibility of the active-cap region. In summary, both the cofactor and the substrate decrease flexibility of the
Fig. 3 RMSF of residues in free enzyme (red), cofactor bound enzyme (green), substrate bound enzyme (blue), and enzyme in ternary complex (pink)
enzyme to different extents based on reduction of RMSF, they together show significant synergistic effect.
Dynamics of the active-cap hSULT2A1 possesses a dynamic active-site cap comprising 27 residues (Asn226 to Gly252) [13]. This region is the largest loop of the enzyme and is often referred to as loop3 in literature [8, 9]. Much attention has been paid to the impact of cofactor on the structure and dynamics of the enzyme, especially on this loop region [12]. The primary influence of the cofactor on the enzyme is reflected by the isomerization of the active-site cap between open and close forms [13–16, 28]. This isomerization was suggested as the molecular basis for substrate selectivity of the enzyme [14, 15]. Here root-meansquare deviation (RMSD) of this region was calculated to describe its conformational variation induced by the cofactor, the substrate and both. Again, rotation and translation of backbone heavy atoms from the initial snapshot were minimized and the RMSD values shown in Fig. 4 were calculated using backbone heavy atoms of the active-site cap. RMSD and its fluctuation for the free enzyme is the largest while this value in ternary complex is the smallest. The cofactor and the substrate
Fig. 4 RMSD evolution for the active-cap region (Asn226 to Gly252) in free enzyme (red), cofactor bound enzyme (green), substrate bound enzyme (blue), and enzyme in ternary complex (pink)
J Mol Model (2015) 21: 190
Page 5 of 9 190
Table 2 Inter-residue distances and corresponding standard deviations. In each residue pair, one belongs to the active-cap region and the other belongs to the main part of the enzyme Residues
2A1(Å)
2A1+PAP (Å)
2A1+LCA(Å)
2A1+PAP+LCA(Å)
Tyr231(Oβ)-Ala243(Cβ) Leu233(Cδ)-Pro70(Cγ) Leu234(Cδ)-Pro70(Cγ) Asp237(Oδ)-Lys138(Nε)
7.054±1.531 8.319±1.440 5.097±1.828 8.381±3.987
8.643±2.430 4.998±1.048 9.606±1.071 6.442±1.399
9.937±1.604 5.501±1.261 9.455±1.090 6.803±1.259
5.212±1.044 5.560±1.187 9.258±1.019 7.011±1.031
Tyr238(Oβ)-Met16(Cε) Val239(Cγ)-Met137(Cε) Val239(Cγ)-Lys138(Nε) Val240(Cγ)-Lys138(Nε) Asp241(Oδ)-Asn136(Nδ) Arg247(Nη)-Asp125(Oδ) Arg247(Nη)-Asp253(Oδ) Lys248(Nε)-Asp253(Oδ)
10.608±4.391 7.955±2.296 11.495±2.327 10.218±3.115 7.038±1.816 15.545±2.947 11.928±2.517 8.257±4.171
6.562±2.433 6.349±1.604 9.011±1.838 7.016±1.500 9.584±1.384 7.470±1.663 6.000±1.775 8.070±5.015
6.084±2.335 5.042±1.467 9.296±1.459 7.029±1.448 6.095±1.022 6.953±3.779 3.886±0.831 11.485±1.560
6.089±1.942 6.097±1.194 8.518±1.424 6.776±1.339 8.155±1.272 7.982±1.168 6.446±0.772 3.082±0.736
show significant cooperative effect in decreasing conformational fluctuations of this region. The impact of ligands on the structure and dynamics of the active-cap region was further monitored by inter-residue distances. Inter-residue distances between all active-cap residues and all other residues of the enzyme were calculated and averaged over time. There are a total of ten residue pairs that have average distances smaller than 10.0 Å and exhibit significant differences among free enzyme, binary, and ternary complexes. Atoms used to calculate inter-residue distance together with the corresponding residues discussed here were listed in Table 2. The listed inter-residue distances in cofactor and substrate bound enzymes are very different from that of the free enzyme (see Table 2), the observation indicates the significant impact of cofactor or substrate on the structure of the enzyme. Both the substrate and the cofactor significantly impact relative positions of the active cap with respect to the main part of the enzyme. Furthermore the inter-residue distances in ternary complex are also very different from that of the binary complex (see Table 2), especially for residue pairs Tyr231/Ala243, Asp241/Lys242, and Lys248/Asp253. This again demonstrates the synergistic effect of cofactor and substrate on the structure of the enzyme.
Fig. 5 a RMSD evolution for cofactor in binary complex (green) and ternary complex (pink). b RMSD evolution for substrate in binary (blue) and ternary complex (pink)
Flexibility of the ligands and their interface with the enzyme To describe the stability variation of ligand in binary and ternary complexes, RMSD of ligands was calculated. Time dependence of RMSD was shown in Fig. 5. The results demonstrate the mutual stabilizing effect of cofactor and substrate. Contacts between ligands and the enzyme were calculated to characterize the interface between them. A contact between a ligand and the enzyme is defined as occurring when one atom of the ligand and any one atom of a residue in the enzyme come within 3.5 Å. After evaluating the range from 2.0 Å to 5.0 Å with the interval of 0.5 Å, 3.5 Å was chosen as the optimal cutoff distance, which enabled us to include all conserved binding sites referred to in literature [17–19]. Residues forming contacts between the cofactor and the enzyme, as well as existence probabilities of the contacts were shown in Table 3a. Residues Lys44, Ser45, Gly46, Thr47, Asn48, Trp49, Arg121, Ser129, Tyr184, Ser218, Phe220, Met223, and Arg247 form stable contacts with the cofactor in both binary and ternary complexes with probabilities larger than 80 %. Among these residues, Lys44 to Trp49 belong to the 5′-phosphosphate-binding (PSB) loop, Arg247 belongs to
190 Page 6 of 9
J Mol Model (2015) 21: 190
Table 3 Residues forming contacts with (a) cofactor and (b) substrate in binary and ternary complexes. Existence probabilities of the contacts were also shown. All the standard errors of the existence probabilities are smaller than 0.01 and not shown here. The SE was calculated as the process mentioned in Table 1. In addition, if 500 snapshots were selected each time, SE was smaller than 0.02 (a) Residues
2A1+PAP
2A1+PAP+LCA
LYS44
1.00
1.00
SER45 GLY46
0.88 1.00
0.99 1.00
THR47 ASN48
0.97 0.97
1.00 1.00
TRP49
0.89
0.91
ARG121 SER129
1.00 1.00
1.00 1.00
TYR184 LYS188
0.99 0.33
1.00 0.42
SER218 SER219
0.97 0.14
0.96 0.09
PHE220 MET223 TYR231 LEU245
0.98 0.95 0.01 0.32
0.95 0.91 0.00 1.00
LEU246 ARG247 LYS248 GLY249 VAL250
0.59 0.99 0.48 0.49 0.02
0.99 1.00 1.00 1.00 0.04
SER251 (b) Residues PRO14 MET16
0.04
0.00
2A1+LCA 0.47 0.80
2A1+PAP+LCA 0.80 0.11
GLY17 PHE18 PRO43
0.47 0.53 0.00
0.99 0.55 0.01
THR47 ASN48 ILE71 TRP72 GLU73 PRO76 TRP77 SER80 GLU81 ILE82 GLY83 ALA86 HSE99 PHE133 TRP134
0.00 0.11 0.06 0.66 0.12 0.01 0.98 0.21 0.01 0.67 0.01 0.00 0.35 0.95 0.66
0.02 0.95 0.00 0.91 0.01 0.01 0.99 0.73 0.00 0.99 0.01 0.01 0.97 0.78 0.98
Table 3 (continued) MET137 PHE139 ILE140 LEU159 TYR160 TYR231 LEU234 SER235 TYR238 VAL239 LEU245
0.82 0.91 0.03 0.00 0.29 0.58 0.96 0.03 0.93 0.05 0.02
0.16 0.96 0.01 0.05 0.58 0.98 0.88 0.01 0.96 0.00 0.00
the 3′-phosphate-binding (PB) loop, and other residues belong to the adenine-ring-binding (SLH3) loop [3]. All of these residues are family-wide conserved residues and the results here show that they play an important role in binding of cofactor. In addition, for residues Leu245, Leu246, Lys248, and Gly249 (they belonging to PB loop), the probabilities of contacts between them and the cofactor were increased significantly after binding of substrate. The substrate thus significantly increases the binding of the cofactor with the enzyme, mainly through the PB loop. Table 3b showed the residues forming contacts between substrate and the enzyme, together with their existence probabilities. For residues Trp77, Phe139, Leu234, and Tyr238, they form stable contacts with the substrates in both binary and ternary complexes. All of these residues are critical for binding of many substrates in addition to LCA [18, 19], and Tyr238 is relevant to the substrate inhibition [29]. For residues Pro14, Gly17, Asn48, Trp72, Ser80, Ile82, His99, Trp134, Tyr160, and Tyr231, their contact probabilities with the substrate were increased significantly in the presence of the cofactor, and His99 is a key catalytic residue. The cofactor could significantly increase the binding of His99 with the substrate LCA as indicated by the increased existence probability from 35 to 97 %. It is worth noting that Asn48 is the only residue that is observed in this study to bind with both cofactor and substrate. For residues Met16, Phe133, and Met137, however, their contact probabilities between the substrate and the enzyme were decreased by the addition of cofactor. On the whole, the impact of the cofactor on the binding of substrate with the enzyme is more complex, in contrast to the sole stabilization effect observed for the substrate on the binding of cofactor. Interactions between the ligands and hSULT2A1 Interaction energies between two ligands and the enzyme were calculated to quantify direct interactions between them. In general, the interaction of both the cofactor and the substrate with the enzyme in ternary complex is stronger than that in corresponding binary complexes. The interaction energies between the cofactor and the enzyme were −1046.73±51.00 kcal mol−1 in binary complex and −1063.18±28.27 kcal mol−1 in
J Mol Model (2015) 21: 190
Page 7 of 9 190
Fig. 6 The interaction energy between (a) cofactor and (b) substrate and interacting enzyme residues. Upper panels: data for binary complexes; Middle panel: data for the ternary complex; Lower panel: their differences (E(ternary) - E(binary))
ternary complex. The interaction energies between the substrate and the enzyme were −43.68±27.38 kcal mol−1 and −65.55±15.87 kcal mol−1 in binary and ternary complexes respectively. Because of the four net charges of the cofactor, the electrostatic interaction energy between the cofactor and the enzyme was much larger than that between the substrate and the enzyme. More than 97 % (97.53 % for binary complex and 97.78 % for ternary complex) of the interaction energy between the enzyme and the cofactor was the electrostatic energy. This ratio was less than 50 % in complexes (24.61 % for binary complex and 48.21 % for ternary complex) between the enzyme and the substrate. It was worth noting that the electrostatic interaction between the enzyme and the substrate became more important in ternary complex. To estimate the effect of each residue on the interactions between the enzyme and the ligands, the interaction energy of each interacting residue with the ligands were calculated and the results were shown in Fig. 6. There were three key residues (Lys44, Arg121, and Arg247, see Fig. 6a) that formed very stable interactions between the cofactor and the enzyme both in binary and ternary complexes, among which Lys44 is a key residue related to the catalytic reaction [3]. All these residues were located along the interface between the enzyme and the cofactor identified above. Figure 6a (lower panel) showed residues that experienced significant substrate induced change of interaction energy between the cofactor and the enzyme, they were mainly located along the active-cap region. For the enzyme complex with the substrate, residues Arg19, Lys44,
Trp72, Lys138, and Tyr238 formed very stable interactions between the enzyme and the substrate both in binary and ternary complexes (see Fig. 6b). Among these residues only Tyr238 was located along the interface between the enzyme and the substrate. As shown in Fig. 6b (lower panel) the interaction of more comprehensive residues with the substrate was varied by the cofactor.
Table 4 Residues forming hydrogen bonds with (a) cofactor and (b) substrate in binary and ternary complexes. Existence probabilities of the hydrogen bonds were also shown. All the standard errors of the existence probabilities are smaller than 0.01 and not shown here. The SE was calculated as the process mentioned in Table 1. In addition, if 500 snapshots were selected each time, SE was smaller than 0.02 (a) Residues LYS44 SER45 GLY46 THR47 ASN48 (b) Residues ASN48 TYR231 TYR238
2A1+PAP 1.43 0.36 0.56 1.21 1.50
2A1+PAP+LCA 1.19 0.49 1.00 1.80 1.85
2A1+LCA 0.00 0.11 0.30
2A1+PAP+LCA 0.48 0.49 0.71
190 Page 8 of 9
Hydrogen bonds between ligands and the enzyme were also calculated using VMD plugin [30]. A hydrogen bond is defined by a 3.0 Å or shorter distance between a donor and an acceptor atom, and by a 120° or larger donor-hydrogenacceptor angle. The residues forming hydrogen bond with ligands, as well as the existence probabilities of these hydrogen bonds were presented in Table 4. For hydrogen bonds with existence probabilities larger than 100 %, the corresponding residues could form more than one hydrogen bond with the ligand at the same time. It is the residues Lys44 to Asn48 (PSB loop) that form hydrogen bonds with the cofactor as shown in Table 4a. This was consistent with that in literature [17]. Comparing existence probabilities of hydrogen bonds formed by the cofactor and the enzyme in binary and ternary complexes, it was found that the substrate could significantly increase the hydrogen bond interactions of the cofactor and the enzyme. Residues Tyr231 and Tyr238 form hydrogen bonds with the substrate in binary and ternary complexes. The cofactor also increases the hydrogen bonds stability. It is worth noting that Asn48 forms hydrogen bond with the substrate in the presence of the cofactor. Asn48 is the only residue forming hydrogen bond with both cofactor and substrate in ternary complex.
J Mol Model (2015) 21: 190
References 1. 2.
3.
4. 5. 6.
7.
8.
9.
10.
Conclusions In this study, extensive MD simulations have been used to investigate interactions between hSULT2A1 and its ligands, including both the cofactor and the substrate LCA. Through analysis of both intra- and intermolecular dynamical structural and energetic fluctuations, we provided an atomically detailed description of how the cofactor impact dynamical behavior of the enzyme, especially the active cap region of the enzyme. Taking LCA as an example, we demonstrated that the substrate may affect the structure and flexibility of the enzyme solely, and more significantly when in cooperation with the cofactor. These observations suggested that the substrate, such as LCA, participates in regulating the structure and flexibility of the enzyme actively rather than merely being selected passively. Since hSULT2A1 has a rather broad substrate spectrum, the impact of other substrates needs to be explored to clarify the role of substrate more clearly.
11.
12. 13.
14.
15.
16. 17.
18.
19.
20. Acknowledgments This work is supported by the National Science Foundation of China (Nos. 31270758, 21403085), Jilin Province Science and Technology Development Plan (20130522009JH) and by High Performance Computing Center of Jilin University, China.
21.
Falany C N (1 997 ) En zymolo gy o f h uman cyto soli c sulfotransferases. FASEB J 11:206–216 Negishi M, Pedersen LG, Petrotchenko E, Shevtsov S, Gorokhov A, Kakuta Y, Pedersen LC (2001) Structure and function of sulfotransferases. Arch Biochem Biophys 390:149–157 Chapman E, Best MD, Hanson SR, Wong CH (2004) Sulfotransferases: structure, mechanism, biological activity, inhibition, and synthetic utility. Angew Chem Int Ed Engl 43:3526–3548 Di L (2014) The role of drug metabolizing enzymes in clearance. Expert Opin Drug Metab Toxicol 10:379–393 Glatt H (2000) Sulfotransferases in the bioactivation of xenobiotics. Chem Biol Interact 129:141–170 Yoshinari K, Petrotchenko EV, Pedersen LC, Negishi M (2001) Crystal structure-based studies of cytosolic sulfotransferase. J Biochem Mol Toxicol 15:67–75 Falany JL, Macrina N, Falany CN (2002) Regulation of MCF-7 breast cancer cell growth by beta-estradiol sulfation. Breast Cancer Res Treat 74:167–176 Allali-Hassani A, Pan PW, Dombrovski L, Najmanovich R, Tempel W, Dong A, Loppnau P, Martin F, Thornton J, Edwards AM, Bochkarev A, Plotnikov AN, Vedadi M, Arrowsmith CH (2007) Structural and chemical profiling of the human cytosolic sulfotransferases. PLoS Biol 5:e97 Dong D, Ako R, Wu B (2012) Crystal structures of human sulfotransferases: insights into the mechanisms of action and substrate selectivity. Expert Opin Drug Metab Toxicol 8:635–646 Berger I, Guttman C, Amar D, Zarivach R, Aharoni A (2011) The molecular basis for the broad substrate specificity of human sulfotransferase 1A1. PLoS ONE 6:e26794 Martiny VY, Carbonell P, Lagorce D, Villoutreix BO, Moroy G, Miteva MA (2013) In silico mechanistic profiling to probe small molecule binding to sulfotransferases. PLoS One 8:e73587 Leyh TS, Cook I, Wang T (2013) Structure, dynamics and selectivity in the sulfotransferase family. Drug Metab Rev 45:423–430 Cook IT, Wang T, Falany CN and Leyh T (2013) Active site gating controls substrate selectivity in cytosolic sulfotransferases. FASEB J 27:270.6 Cook IT, Wang T, Almo SC, Kim J, Falany CN, Leyh TS (2013) The gate that governs sulfotransferase selectivity. Biochemistry 52: 415–424 Cook I, Wang T, Falany CN, Leyh TS (2012) A nucleotide-gated molecular pore selects sulfotransferase substrates. Biochemistry 51: 5674–5683 Cook I, Wang T, Falany CN, Leyh TS (2013) High accuracy in silico sulfotransferase models. J Biol Chem 288:34494–34501 Pedersen LC, Petrotchenko EV, Negishi M (2000) Crystal structure of SULT2A3, human hydroxysteroid sulfotransferase. FEBS Lett 475:61–64 Rehse PH, Zhou M, Lin SX (2002) Crystal structure of human dehydroepiandrosterone sulphotransferase in complex with substrate. Biochem J 364:165–171 Chang HJ, Shi R, Rehse P, Lin SX (2004) Identifying androsterone (ADT) as a cognate substrate for human dehydroepiandrosterone sulfotransferase (DHEA-ST) important for steroid homeostasis: structure of the enzyme-ADT complex. J Biol Chem 279:2689– 2696 Zhang P-p, Zhao L, Long S-y, Tian P (2015) The effect of ligands on the thermal stability of sulfotransferases: a molecular dynamics simulation study. J Mol Model 21:1–7 Fiser A, Sali A (2003) MODELLER: generation and refinement of homology-based protein structure models. Macromol Crystallogr D 374:461–491
J Mol Model (2015) 21: 190 22.
Phillips JC, Braun R, Wang W, Gumbart J, Tajkhorshid E, Villa E, Chipot C, Skeel RD, Kale L, Schulten K (2005) Scalable molecular dynamics with NAMD. J Comput Chem 26:1781–1802 23. MacKerell AD, Bashford D, Bellott M, Dunbrack RL, Evanseck JD, Field MJ, Fischer S, Gao J, Guo H, Ha S, Joseph-McCarthy D, Kuchnir L, Kuczera K, Lau FTK, Mattos C, Michnick S, Ngo T, Nguyen DT, Prodhom B, Reiher WE, Roux B, Schlenkrich M, Smith JC, Stote R, Straub J, Watanabe M, Wiórkiewicz-Kuczera J, Yin D, Karplus M (1998) All-atom empirical potential for molecular modeling and dynamics studies of proteins. J Phys Chem B 102:3586–3616 24. Jorgensen WL, Chandrasekhar J, Madura JD, Impey RW, Klein ML (1983) Comparison of simple potential functions for simulating liquid water. J Chem Phys 79:926–935 25. 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
Page 9 of 9 190 26.
27.
28.
29.
30.
Ryckaert J-P, Ciccotti G, Berendsen HJC (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 Kabsch W, Sander C (1983) Dictionary of protein secondary structure: pattern recognition of hydrogen-bonded and geometrical features. Biopolymers 22:2577–2637 Cook I, Wang T, Almo SC, Kim J, Falany CN, Leyh TS (2013) Testing the sulfotransferase molecular pore hypothesis. J Biol Chem 288:8619–8626 Lu LY, Hsieh YC, Liu MY, Lin YH, Chen CJ, Yang YS (2008) Identification and characterization of two amino acids critical for the substrate inhibition of human dehydroepiandrosterone sulfotransferase (SULT2A1). Mol Pharmacol 73:660–668 Humphrey W, Dalke A, Schulten K (1996) VMD: visual molecular dynamics. J Mol Graph 14(33-8):27–28
ORIGINAL PAPER
The impact of ligands on the structure and flexibility of sulfotransferases: a molecular dynamics simulation study Li Zhao 1 & Pupu Zhang 1,4 & Shiyang Long 1 & Linlin Wang 1,2 & Pu Tian 1,3
Received: 17 March 2015 / Accepted: 15 June 2015 / Published online: 8 July 2015 # Springer-Verlag Berlin Heidelberg 2015
Abstract Sulfotransferases catalyze transfer of the sulfurylgroup (-SO3) from 3′-phosphoadenosine 5′-phosphosulfate (PAPS) to a large number of substrates. They play an important role in phase II metabolic process. The impact of the cofactor (PAPS) on the structure and flexibility of the enzyme has been studied extensively, and the response of the activecap region to cofactor binding was proposed as the molecular basis for substrate selectivity. In this study, individual and cooperative effects of the cofactor and substrate on the structure and flexibility of the enzyme were investigated. Molecular dynamics simulations were performed for four systems, including free enzyme, binary complexes (cofactor or substrate bound enzyme) and ternary complex (both cofactor and substrate bound enzyme). The influence of ligands (the cofactor and the substrate) on the structure and flexibility of the enzyme, especially that of the active-site cap region, was analyzed. Moreover, mutual structural impact of the ligands was examined as well. The results show that the impact of both the cofactor and the substrate was significant. Our study indicated that the substrate, such as lithocholic acid (LCA), participated in regulating the structure and flexibility of the enzyme actively rather than merely being selected passively. Additionally, the observed synergistic effects of the cofactor * Pu Tian [email protected] 1
School of Life Sciences, Jilin University, Changchun, China
2
Ultrasound Department, China-Japan Union Hospital of Jilin University, Changchun, China
3
Key Laboratory for Molecular Enzymology and Engineering, the Ministry of Education, Jilin University, Changchun, China
4
Present address: Zhongshan Ophthalmic Center, Sun Yat-Sen University, Zhongshan, China
and the substrate demonstrated the importance of examining both ligands in understanding enzymes. Keywords Cofactor and substrate . Human cytosolic sulfotransferase 2A1 . Molecular dynamics simulation . Structure and flexibility
Introduction Cytosolic sulfotransferases (SULTs) catalyze the transfer of a sulfonate group from the unique cofactor 3′phosphoadenosine 5′-phosphosulfate (PAPS) to a wide variety of endogenous and xenobiotic substrates [1–3], including drugs, toxic compounds, steroid hormones, and neurotransmitters. They are phase II metabolic enzymes that participate in chemical defense. The catalytic reaction are broadly involved in both modulating biological activity and increasing solubility of substrates [4]. However, there are reports that link sulfonation imbalance to induction of carcinogenic responses [5–7]. Human cytosolic sulfotransferases (hSULTs) could be divided into four families, including hSULT1, hSULT2, hSULT4, and hSULT6 [8, 9]. Many of them have very broad and overlapping substrate specificities as demonstrated by recent progress in structural biology [8]. Allali-Hassani et al. found evidence for structural Bpriming^ of hSULTs active site by cofactor binding, which may further influence the spectrum of substrates [8]. Ilana Berger et al. found high active site plasticity and structural flexibility of hSULT1A1, and suggested that such character enables binding of different acceptors and plays a dominant role in controlling the specificity and activity [10]. Martiny et al. explored the flexibility of hSULT1A1, hSULT1A3, and hSULT1E1 and identified multiple receptor conformations, which were proposed to be the
190 Page 2 of 9
J Mol Model (2015) 21: 190
most suitable ones for ligand binding prediction [11]. Besides substrate selectivity, substrate inhibition of hSULTs was found to be related to the characteristic structure and flexibility of the enzyme [12]. hSULTs possess a stretch of about 30 amino acids long loop, which serves as an active cap, interacts with and encapsulates both the cofactor and substrates [13–15]. Cook and coworkers observed the response of this active-cap to cofactor binding and suggested its important role in substrate selectivity and recognition [12]. Upon binding of the cofactor PAPS, the active-site cap closes and the binding probability of large substrates (e.g., Raloxifene) decreased [13]. When interactions between the cap and their adjacent residues were weakened via mutagenesis, the enzyme no longer discriminated between large and small substrates and gained increased catalytic efficiency toward large substrates [16]. However, up till now, little is known about the influence of substrates on the structure and dynamics of the enzyme, which are fundamental to understand both the catalytic and substrate inhibition mechanisms. To date SULT2A1 has the broadest known substrate spectrum among SULTs [8, 14]. hSULT2A1 is found in numerous tissues and is concentrated in liver. Its substrates include small planar molecules, steroid, cyclic amines, numerous drugs, and small peptides. The crystal structures of SULT2A1 in complex with 3′-phosphoadenosine 5′-phosphate (PAP) [17], dehydroepiandrosterone (DHEA) [18], androsterone (ADT) [19], and lithocholic acid (LCA) are available in protein data bank with PDB codes 1EFH, 1J99, 1OV4, and 3F3Y respectively. Previously, we investigated the impact of ligands on the thermal stability of SULT2A1 with MD simulations [20], and observed that while both ligands increased the thermal
stability of the enzyme, the extent of stabilization is different for different ligands binding scenarios (cofactor, substrate and both). The thermo-denaturation order of studied complexes was consistent with experimental observations [8]. Here extensive molecular dynamics (MD) simulations were used to investigate the impact of ligands (cofactor and substrate) on the structure and flexibility of human SULT2A1 at physiological temperature. Four systems, including the free enzyme, binary complexes (cofactor or substrate LCA bound enzyme), and the ternary complex (enzyme with both cofactor and substrate LCA bound), were explored. Firstly the structure and flexibility variation among four systems were analyzed based on the root-mean-square fluctuation (RMSF) and evolution of secondary structures with time. Secondly, based on rootmean-square deviation (RMSD) of the active-site cap region and inter-residue distances between the active-site cap region and the main part of the enzyme, we characterized the structural variation of the active-site cap region induced by the ligands. Finally, the mutual structural impact of the ligands was analyzed. RMSD of each ligand as well as the interactions between ligands and the enzyme were examined. The results demonstrated the significant influence of ligands, including both the cofactor and the substrate, on the structure and flexibility of the enzyme.
Fig. 1 Schematic representation of the procedures for building the topology and force field parameters of cofactor PAP. The phosphate moiety in 3PO3 and TH5P were joined onto ADE to obtain topology
and force field parameters of PAP. Topology and parameters of ADE 3PO3 and TH5P are available in CHARMM27 force field [23]
System setup and simulation methods Four models (free enzyme, PAP-bound enzyme, LCA-bound enzyme and both PAP and LCA bound enzyme) were derived from the crystal structure of SULT2A1 complex with both PAP and LCA (PDB code: 3F3Y). Missing atoms were added
J Mol Model (2015) 21: 190
Page 3 of 9 190
Fig. 2 Secondary structure for each residue in a free enzyme, b cofactor bound enzyme, c substrate bound enzyme, and d ternary complex. Secondary structure identities were assigned using DSSP [27]. T: bend and hydrogen bonded turn, E: extended strand participated in beta ladder, B: isolated betabridge, H: alpha helix, G: 3-helix, I: 5-helix, C: coil
using MODELLER [21]. All four models were solvated to obtain a cube large enough to allow at least 12 Å of water between protein and cube surfaces, requiring about 12,500 water molecules. The net charge of the systems was neutralized with Na+, and NaCl was added to the box to 0.15 mM. The resulting four systems are denoted 2A1, 2A1 +PAP, 2A1+LCA, and 2A1+PAP+LCA respectively. All MD simulations were performed with NAMD 2.7 software package [22] using CHARMM27 force fields [23] and TIP3P water model [24]. Topology and parameters of PAP were obtained by joining three segments ADE, 3PO3, and TH5P according to standard CHAR MM procedures as shown in Fig. 1. Topology and parameters of LCA were taken from CHARMM CGenFF files. Assuming a neutral PH as reported in literature [8], the cofactor PAP in this study was not protonated and it had four net charges while the substrate LCA was protonated. PAP rather than PAPS was chosen to study
the influence of PAP and PAPS on the structure and flexibility of hSULTs based on the experimental convenience and the belief that their effects are nearly identical [8], and we followed the same practice in our simulations. Non-covalent van der Waals interactions were cut off at 12 Å. Particle mesh Ewald (PME) [25] was used to calculate long range electrostatic interactions. Bond length involving hydrogen atoms was constrained using SHAKE algorithm [26] with the integration time step of 2 fs. The starting configurations were energy-minimized and equilibrated as previously reported [20]. Then 1-μs simulations at 310 K under NVT condition for 2A1, 2A1 + PAP, 2A1 + LCA, and 2A1 + PAP + LCA systems were performed and collectively 4 μs trajectories were generated. Due to the fact that 2A1, 2A1 + PAP and 2A1 + LCA systems are built from the ternary crystal structure, the first 500-ns of each trajectory is discarded,
Table 1 The probabilities of listed segments in α-helical state for 2A1, 2A1+PAP, 2A1+LCA, and 2A1+PAP+LCA. These segments have a complementary probability of existing in the coil state. The uncertainty is the standard error (SE). A sample of 2000 snapshots were randomly
selected from 10,000 snapshots and the average probabilities of segments being in α-helical state was calculated. The SE was the standard deviations of 1000 samples. (With a sample size of 500, SE was smaller than 0.07)
Residues
2A1
2A1+PAP
2A1+LCA
2A1+PAP+LCA
LYS64-SER68 ILE71-ARG74 ILE102-LEU104 LYS107-SER111 GLY169-TRP170 MET171-MET173
0.65±0.03 0.64±0.02 0.54±0.02 0.76±0.03 0.43±0.02 0.83±0.01
0.46±0.03 0.96±0.00 0.35±0.04 0.78±0.02 0.74±0.03 0.89±0.01
0.46±0.03 0.68±0.03 0.35±0.03 0.77±0.02 0.65±0.02 0.85±0.01
0.39±0.04 0.92±0.00 0.37±0.04 0.69±0.03 0.88±0.02 0.93±0.00
190 Page 4 of 9
J Mol Model (2015) 21: 190
and snapshots of the last 500 ns were taken with an interval of 100 ps to analyze the impact of the cofactor and the substrate on the structure and dynamics of the enzyme.
Results and discussion Structure and flexibility of the enzyme To compare the structure and flexibility of the enzyme in the four systems, secondary structure assignment was carried out using DSSP [27] and its evolution with time was shown in Fig. 2. The notations of secondary structure elements are the same as reported in literature [8, 9]. Overall, there is no obvious difference of secondary structure alignments among four systems except that some peptide segments interconverting between α-helical and coil states exhibited different weight allocation in these two forms. These segments include Lys64-Ser68 in α4, Ile71-Arg74 in α5, Lys107-Ser111 in α7, and Trp163-Met173 in α11. Their probabilities in αhelical state are shown in Table 1. Lys64-Ser68 in α4 and Ile102-Leu104 in α7 in free enzyme tend to exist in αhelical state while they tend to be in coil state in binary or ternary complexes. The probability of Ile71-Arg74 in α5 as αhelix in free enzyme and substrate-bound systems is slightly higher than that in the cofactor-bound and the ternary system. The N-terminal two residues of α11 tend to be in α-helical state in ternary complex, while they exhibit coil state in the free enzyme; the C-terminal three residues possess higher tendency to form α-helix in ternary complex than the other three systems. RMSF shown in Fig. 3 was calculated after minimizing rotation and translation of all backbone heavy atoms. At first glance, free enzyme exhibits larger flexibility as reflected by larger RMSF, especially in the segment Asn226 to Gly252, which is the active-site cap region discussed by Ian Cook [14, 15], while RMSF of the enzyme in ternary system is the smallest. Binding of the cofactor was observed to decrease RMSF of some segments (Trp72-Glu79, Lys188-Asp190, Asn217-Lys248, Phe258-Leu269, and Glu280-Glu285) and increase that of the segment (Ser88-Arg94). The effect of substrate binding on the enzyme RMSF is less extensive, with decrease for the segment Phe220-Ser259 and increase for segments Val69-Ser92, Gly156-Pro172, and Ala275-Phe282 respectively. It is worth noting that cofactor exhibits significant influence on the substrate bound enzyme almost throughout the whole sequence, and particularly for residues Val69Ala86, Gln189-Asp190, and Glu207-Asp267, while the impact of substrate on cofactor bound enzyme is limited to decreasing the flexibility of the active-cap region. In summary, both the cofactor and the substrate decrease flexibility of the
Fig. 3 RMSF of residues in free enzyme (red), cofactor bound enzyme (green), substrate bound enzyme (blue), and enzyme in ternary complex (pink)
enzyme to different extents based on reduction of RMSF, they together show significant synergistic effect.
Dynamics of the active-cap hSULT2A1 possesses a dynamic active-site cap comprising 27 residues (Asn226 to Gly252) [13]. This region is the largest loop of the enzyme and is often referred to as loop3 in literature [8, 9]. Much attention has been paid to the impact of cofactor on the structure and dynamics of the enzyme, especially on this loop region [12]. The primary influence of the cofactor on the enzyme is reflected by the isomerization of the active-site cap between open and close forms [13–16, 28]. This isomerization was suggested as the molecular basis for substrate selectivity of the enzyme [14, 15]. Here root-meansquare deviation (RMSD) of this region was calculated to describe its conformational variation induced by the cofactor, the substrate and both. Again, rotation and translation of backbone heavy atoms from the initial snapshot were minimized and the RMSD values shown in Fig. 4 were calculated using backbone heavy atoms of the active-site cap. RMSD and its fluctuation for the free enzyme is the largest while this value in ternary complex is the smallest. The cofactor and the substrate
Fig. 4 RMSD evolution for the active-cap region (Asn226 to Gly252) in free enzyme (red), cofactor bound enzyme (green), substrate bound enzyme (blue), and enzyme in ternary complex (pink)
J Mol Model (2015) 21: 190
Page 5 of 9 190
Table 2 Inter-residue distances and corresponding standard deviations. In each residue pair, one belongs to the active-cap region and the other belongs to the main part of the enzyme Residues
2A1(Å)
2A1+PAP (Å)
2A1+LCA(Å)
2A1+PAP+LCA(Å)
Tyr231(Oβ)-Ala243(Cβ) Leu233(Cδ)-Pro70(Cγ) Leu234(Cδ)-Pro70(Cγ) Asp237(Oδ)-Lys138(Nε)
7.054±1.531 8.319±1.440 5.097±1.828 8.381±3.987
8.643±2.430 4.998±1.048 9.606±1.071 6.442±1.399
9.937±1.604 5.501±1.261 9.455±1.090 6.803±1.259
5.212±1.044 5.560±1.187 9.258±1.019 7.011±1.031
Tyr238(Oβ)-Met16(Cε) Val239(Cγ)-Met137(Cε) Val239(Cγ)-Lys138(Nε) Val240(Cγ)-Lys138(Nε) Asp241(Oδ)-Asn136(Nδ) Arg247(Nη)-Asp125(Oδ) Arg247(Nη)-Asp253(Oδ) Lys248(Nε)-Asp253(Oδ)
10.608±4.391 7.955±2.296 11.495±2.327 10.218±3.115 7.038±1.816 15.545±2.947 11.928±2.517 8.257±4.171
6.562±2.433 6.349±1.604 9.011±1.838 7.016±1.500 9.584±1.384 7.470±1.663 6.000±1.775 8.070±5.015
6.084±2.335 5.042±1.467 9.296±1.459 7.029±1.448 6.095±1.022 6.953±3.779 3.886±0.831 11.485±1.560
6.089±1.942 6.097±1.194 8.518±1.424 6.776±1.339 8.155±1.272 7.982±1.168 6.446±0.772 3.082±0.736
show significant cooperative effect in decreasing conformational fluctuations of this region. The impact of ligands on the structure and dynamics of the active-cap region was further monitored by inter-residue distances. Inter-residue distances between all active-cap residues and all other residues of the enzyme were calculated and averaged over time. There are a total of ten residue pairs that have average distances smaller than 10.0 Å and exhibit significant differences among free enzyme, binary, and ternary complexes. Atoms used to calculate inter-residue distance together with the corresponding residues discussed here were listed in Table 2. The listed inter-residue distances in cofactor and substrate bound enzymes are very different from that of the free enzyme (see Table 2), the observation indicates the significant impact of cofactor or substrate on the structure of the enzyme. Both the substrate and the cofactor significantly impact relative positions of the active cap with respect to the main part of the enzyme. Furthermore the inter-residue distances in ternary complex are also very different from that of the binary complex (see Table 2), especially for residue pairs Tyr231/Ala243, Asp241/Lys242, and Lys248/Asp253. This again demonstrates the synergistic effect of cofactor and substrate on the structure of the enzyme.
Fig. 5 a RMSD evolution for cofactor in binary complex (green) and ternary complex (pink). b RMSD evolution for substrate in binary (blue) and ternary complex (pink)
Flexibility of the ligands and their interface with the enzyme To describe the stability variation of ligand in binary and ternary complexes, RMSD of ligands was calculated. Time dependence of RMSD was shown in Fig. 5. The results demonstrate the mutual stabilizing effect of cofactor and substrate. Contacts between ligands and the enzyme were calculated to characterize the interface between them. A contact between a ligand and the enzyme is defined as occurring when one atom of the ligand and any one atom of a residue in the enzyme come within 3.5 Å. After evaluating the range from 2.0 Å to 5.0 Å with the interval of 0.5 Å, 3.5 Å was chosen as the optimal cutoff distance, which enabled us to include all conserved binding sites referred to in literature [17–19]. Residues forming contacts between the cofactor and the enzyme, as well as existence probabilities of the contacts were shown in Table 3a. Residues Lys44, Ser45, Gly46, Thr47, Asn48, Trp49, Arg121, Ser129, Tyr184, Ser218, Phe220, Met223, and Arg247 form stable contacts with the cofactor in both binary and ternary complexes with probabilities larger than 80 %. Among these residues, Lys44 to Trp49 belong to the 5′-phosphosphate-binding (PSB) loop, Arg247 belongs to
190 Page 6 of 9
J Mol Model (2015) 21: 190
Table 3 Residues forming contacts with (a) cofactor and (b) substrate in binary and ternary complexes. Existence probabilities of the contacts were also shown. All the standard errors of the existence probabilities are smaller than 0.01 and not shown here. The SE was calculated as the process mentioned in Table 1. In addition, if 500 snapshots were selected each time, SE was smaller than 0.02 (a) Residues
2A1+PAP
2A1+PAP+LCA
LYS44
1.00
1.00
SER45 GLY46
0.88 1.00
0.99 1.00
THR47 ASN48
0.97 0.97
1.00 1.00
TRP49
0.89
0.91
ARG121 SER129
1.00 1.00
1.00 1.00
TYR184 LYS188
0.99 0.33
1.00 0.42
SER218 SER219
0.97 0.14
0.96 0.09
PHE220 MET223 TYR231 LEU245
0.98 0.95 0.01 0.32
0.95 0.91 0.00 1.00
LEU246 ARG247 LYS248 GLY249 VAL250
0.59 0.99 0.48 0.49 0.02
0.99 1.00 1.00 1.00 0.04
SER251 (b) Residues PRO14 MET16
0.04
0.00
2A1+LCA 0.47 0.80
2A1+PAP+LCA 0.80 0.11
GLY17 PHE18 PRO43
0.47 0.53 0.00
0.99 0.55 0.01
THR47 ASN48 ILE71 TRP72 GLU73 PRO76 TRP77 SER80 GLU81 ILE82 GLY83 ALA86 HSE99 PHE133 TRP134
0.00 0.11 0.06 0.66 0.12 0.01 0.98 0.21 0.01 0.67 0.01 0.00 0.35 0.95 0.66
0.02 0.95 0.00 0.91 0.01 0.01 0.99 0.73 0.00 0.99 0.01 0.01 0.97 0.78 0.98
Table 3 (continued) MET137 PHE139 ILE140 LEU159 TYR160 TYR231 LEU234 SER235 TYR238 VAL239 LEU245
0.82 0.91 0.03 0.00 0.29 0.58 0.96 0.03 0.93 0.05 0.02
0.16 0.96 0.01 0.05 0.58 0.98 0.88 0.01 0.96 0.00 0.00
the 3′-phosphate-binding (PB) loop, and other residues belong to the adenine-ring-binding (SLH3) loop [3]. All of these residues are family-wide conserved residues and the results here show that they play an important role in binding of cofactor. In addition, for residues Leu245, Leu246, Lys248, and Gly249 (they belonging to PB loop), the probabilities of contacts between them and the cofactor were increased significantly after binding of substrate. The substrate thus significantly increases the binding of the cofactor with the enzyme, mainly through the PB loop. Table 3b showed the residues forming contacts between substrate and the enzyme, together with their existence probabilities. For residues Trp77, Phe139, Leu234, and Tyr238, they form stable contacts with the substrates in both binary and ternary complexes. All of these residues are critical for binding of many substrates in addition to LCA [18, 19], and Tyr238 is relevant to the substrate inhibition [29]. For residues Pro14, Gly17, Asn48, Trp72, Ser80, Ile82, His99, Trp134, Tyr160, and Tyr231, their contact probabilities with the substrate were increased significantly in the presence of the cofactor, and His99 is a key catalytic residue. The cofactor could significantly increase the binding of His99 with the substrate LCA as indicated by the increased existence probability from 35 to 97 %. It is worth noting that Asn48 is the only residue that is observed in this study to bind with both cofactor and substrate. For residues Met16, Phe133, and Met137, however, their contact probabilities between the substrate and the enzyme were decreased by the addition of cofactor. On the whole, the impact of the cofactor on the binding of substrate with the enzyme is more complex, in contrast to the sole stabilization effect observed for the substrate on the binding of cofactor. Interactions between the ligands and hSULT2A1 Interaction energies between two ligands and the enzyme were calculated to quantify direct interactions between them. In general, the interaction of both the cofactor and the substrate with the enzyme in ternary complex is stronger than that in corresponding binary complexes. The interaction energies between the cofactor and the enzyme were −1046.73±51.00 kcal mol−1 in binary complex and −1063.18±28.27 kcal mol−1 in
J Mol Model (2015) 21: 190
Page 7 of 9 190
Fig. 6 The interaction energy between (a) cofactor and (b) substrate and interacting enzyme residues. Upper panels: data for binary complexes; Middle panel: data for the ternary complex; Lower panel: their differences (E(ternary) - E(binary))
ternary complex. The interaction energies between the substrate and the enzyme were −43.68±27.38 kcal mol−1 and −65.55±15.87 kcal mol−1 in binary and ternary complexes respectively. Because of the four net charges of the cofactor, the electrostatic interaction energy between the cofactor and the enzyme was much larger than that between the substrate and the enzyme. More than 97 % (97.53 % for binary complex and 97.78 % for ternary complex) of the interaction energy between the enzyme and the cofactor was the electrostatic energy. This ratio was less than 50 % in complexes (24.61 % for binary complex and 48.21 % for ternary complex) between the enzyme and the substrate. It was worth noting that the electrostatic interaction between the enzyme and the substrate became more important in ternary complex. To estimate the effect of each residue on the interactions between the enzyme and the ligands, the interaction energy of each interacting residue with the ligands were calculated and the results were shown in Fig. 6. There were three key residues (Lys44, Arg121, and Arg247, see Fig. 6a) that formed very stable interactions between the cofactor and the enzyme both in binary and ternary complexes, among which Lys44 is a key residue related to the catalytic reaction [3]. All these residues were located along the interface between the enzyme and the cofactor identified above. Figure 6a (lower panel) showed residues that experienced significant substrate induced change of interaction energy between the cofactor and the enzyme, they were mainly located along the active-cap region. For the enzyme complex with the substrate, residues Arg19, Lys44,
Trp72, Lys138, and Tyr238 formed very stable interactions between the enzyme and the substrate both in binary and ternary complexes (see Fig. 6b). Among these residues only Tyr238 was located along the interface between the enzyme and the substrate. As shown in Fig. 6b (lower panel) the interaction of more comprehensive residues with the substrate was varied by the cofactor.
Table 4 Residues forming hydrogen bonds with (a) cofactor and (b) substrate in binary and ternary complexes. Existence probabilities of the hydrogen bonds were also shown. All the standard errors of the existence probabilities are smaller than 0.01 and not shown here. The SE was calculated as the process mentioned in Table 1. In addition, if 500 snapshots were selected each time, SE was smaller than 0.02 (a) Residues LYS44 SER45 GLY46 THR47 ASN48 (b) Residues ASN48 TYR231 TYR238
2A1+PAP 1.43 0.36 0.56 1.21 1.50
2A1+PAP+LCA 1.19 0.49 1.00 1.80 1.85
2A1+LCA 0.00 0.11 0.30
2A1+PAP+LCA 0.48 0.49 0.71
190 Page 8 of 9
Hydrogen bonds between ligands and the enzyme were also calculated using VMD plugin [30]. A hydrogen bond is defined by a 3.0 Å or shorter distance between a donor and an acceptor atom, and by a 120° or larger donor-hydrogenacceptor angle. The residues forming hydrogen bond with ligands, as well as the existence probabilities of these hydrogen bonds were presented in Table 4. For hydrogen bonds with existence probabilities larger than 100 %, the corresponding residues could form more than one hydrogen bond with the ligand at the same time. It is the residues Lys44 to Asn48 (PSB loop) that form hydrogen bonds with the cofactor as shown in Table 4a. This was consistent with that in literature [17]. Comparing existence probabilities of hydrogen bonds formed by the cofactor and the enzyme in binary and ternary complexes, it was found that the substrate could significantly increase the hydrogen bond interactions of the cofactor and the enzyme. Residues Tyr231 and Tyr238 form hydrogen bonds with the substrate in binary and ternary complexes. The cofactor also increases the hydrogen bonds stability. It is worth noting that Asn48 forms hydrogen bond with the substrate in the presence of the cofactor. Asn48 is the only residue forming hydrogen bond with both cofactor and substrate in ternary complex.
J Mol Model (2015) 21: 190
References 1. 2.
3.
4. 5. 6.
7.
8.
9.
10.
Conclusions In this study, extensive MD simulations have been used to investigate interactions between hSULT2A1 and its ligands, including both the cofactor and the substrate LCA. Through analysis of both intra- and intermolecular dynamical structural and energetic fluctuations, we provided an atomically detailed description of how the cofactor impact dynamical behavior of the enzyme, especially the active cap region of the enzyme. Taking LCA as an example, we demonstrated that the substrate may affect the structure and flexibility of the enzyme solely, and more significantly when in cooperation with the cofactor. These observations suggested that the substrate, such as LCA, participates in regulating the structure and flexibility of the enzyme actively rather than merely being selected passively. Since hSULT2A1 has a rather broad substrate spectrum, the impact of other substrates needs to be explored to clarify the role of substrate more clearly.
11.
12. 13.
14.
15.
16. 17.
18.
19.
20. Acknowledgments This work is supported by the National Science Foundation of China (Nos. 31270758, 21403085), Jilin Province Science and Technology Development Plan (20130522009JH) and by High Performance Computing Center of Jilin University, China.
21.
Falany C N (1 997 ) En zymolo gy o f h uman cyto soli c sulfotransferases. FASEB J 11:206–216 Negishi M, Pedersen LG, Petrotchenko E, Shevtsov S, Gorokhov A, Kakuta Y, Pedersen LC (2001) Structure and function of sulfotransferases. Arch Biochem Biophys 390:149–157 Chapman E, Best MD, Hanson SR, Wong CH (2004) Sulfotransferases: structure, mechanism, biological activity, inhibition, and synthetic utility. Angew Chem Int Ed Engl 43:3526–3548 Di L (2014) The role of drug metabolizing enzymes in clearance. Expert Opin Drug Metab Toxicol 10:379–393 Glatt H (2000) Sulfotransferases in the bioactivation of xenobiotics. Chem Biol Interact 129:141–170 Yoshinari K, Petrotchenko EV, Pedersen LC, Negishi M (2001) Crystal structure-based studies of cytosolic sulfotransferase. J Biochem Mol Toxicol 15:67–75 Falany JL, Macrina N, Falany CN (2002) Regulation of MCF-7 breast cancer cell growth by beta-estradiol sulfation. Breast Cancer Res Treat 74:167–176 Allali-Hassani A, Pan PW, Dombrovski L, Najmanovich R, Tempel W, Dong A, Loppnau P, Martin F, Thornton J, Edwards AM, Bochkarev A, Plotnikov AN, Vedadi M, Arrowsmith CH (2007) Structural and chemical profiling of the human cytosolic sulfotransferases. PLoS Biol 5:e97 Dong D, Ako R, Wu B (2012) Crystal structures of human sulfotransferases: insights into the mechanisms of action and substrate selectivity. Expert Opin Drug Metab Toxicol 8:635–646 Berger I, Guttman C, Amar D, Zarivach R, Aharoni A (2011) The molecular basis for the broad substrate specificity of human sulfotransferase 1A1. PLoS ONE 6:e26794 Martiny VY, Carbonell P, Lagorce D, Villoutreix BO, Moroy G, Miteva MA (2013) In silico mechanistic profiling to probe small molecule binding to sulfotransferases. PLoS One 8:e73587 Leyh TS, Cook I, Wang T (2013) Structure, dynamics and selectivity in the sulfotransferase family. Drug Metab Rev 45:423–430 Cook IT, Wang T, Falany CN and Leyh T (2013) Active site gating controls substrate selectivity in cytosolic sulfotransferases. FASEB J 27:270.6 Cook IT, Wang T, Almo SC, Kim J, Falany CN, Leyh TS (2013) The gate that governs sulfotransferase selectivity. Biochemistry 52: 415–424 Cook I, Wang T, Falany CN, Leyh TS (2012) A nucleotide-gated molecular pore selects sulfotransferase substrates. Biochemistry 51: 5674–5683 Cook I, Wang T, Falany CN, Leyh TS (2013) High accuracy in silico sulfotransferase models. J Biol Chem 288:34494–34501 Pedersen LC, Petrotchenko EV, Negishi M (2000) Crystal structure of SULT2A3, human hydroxysteroid sulfotransferase. FEBS Lett 475:61–64 Rehse PH, Zhou M, Lin SX (2002) Crystal structure of human dehydroepiandrosterone sulphotransferase in complex with substrate. Biochem J 364:165–171 Chang HJ, Shi R, Rehse P, Lin SX (2004) Identifying androsterone (ADT) as a cognate substrate for human dehydroepiandrosterone sulfotransferase (DHEA-ST) important for steroid homeostasis: structure of the enzyme-ADT complex. J Biol Chem 279:2689– 2696 Zhang P-p, Zhao L, Long S-y, Tian P (2015) The effect of ligands on the thermal stability of sulfotransferases: a molecular dynamics simulation study. J Mol Model 21:1–7 Fiser A, Sali A (2003) MODELLER: generation and refinement of homology-based protein structure models. Macromol Crystallogr D 374:461–491
J Mol Model (2015) 21: 190 22.
Phillips JC, Braun R, Wang W, Gumbart J, Tajkhorshid E, Villa E, Chipot C, Skeel RD, Kale L, Schulten K (2005) Scalable molecular dynamics with NAMD. J Comput Chem 26:1781–1802 23. MacKerell AD, Bashford D, Bellott M, Dunbrack RL, Evanseck JD, Field MJ, Fischer S, Gao J, Guo H, Ha S, Joseph-McCarthy D, Kuchnir L, Kuczera K, Lau FTK, Mattos C, Michnick S, Ngo T, Nguyen DT, Prodhom B, Reiher WE, Roux B, Schlenkrich M, Smith JC, Stote R, Straub J, Watanabe M, Wiórkiewicz-Kuczera J, Yin D, Karplus M (1998) All-atom empirical potential for molecular modeling and dynamics studies of proteins. J Phys Chem B 102:3586–3616 24. Jorgensen WL, Chandrasekhar J, Madura JD, Impey RW, Klein ML (1983) Comparison of simple potential functions for simulating liquid water. J Chem Phys 79:926–935 25. 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
Page 9 of 9 190 26.
27.
28.
29.
30.
Ryckaert J-P, Ciccotti G, Berendsen HJC (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 Kabsch W, Sander C (1983) Dictionary of protein secondary structure: pattern recognition of hydrogen-bonded and geometrical features. Biopolymers 22:2577–2637 Cook I, Wang T, Almo SC, Kim J, Falany CN, Leyh TS (2013) Testing the sulfotransferase molecular pore hypothesis. J Biol Chem 288:8619–8626 Lu LY, Hsieh YC, Liu MY, Lin YH, Chen CJ, Yang YS (2008) Identification and characterization of two amino acids critical for the substrate inhibition of human dehydroepiandrosterone sulfotransferase (SULT2A1). Mol Pharmacol 73:660–668 Humphrey W, Dalke A, Schulten K (1996) VMD: visual molecular dynamics. J Mol Graph 14(33-8):27–28
