Dimerization and ligand binding in tyrosylprotein sulfotransferase-2 are influenced by molecular motions

Copyright © and Moral Rights are retained by the author(s) and/ or other copyright owners. A copy can be downloaded for personal non-commercial research or study, without prior permission or charge. This item cannot be reproduced or quoted extensively from without first obtaining permission in writing from the copyright holder(s). The content must not be changed in any way or sold commercially in any format or medium without the formal permission of the copyright holders.


Introduction
Sulfation of tyrosine residues in proteins is a very important post synthetic modification, which is still less explored. The reaction is catalyzed by an enzyme -Tyrosylprotein Sulfotransferase (TPST) 1,2 . Tyrosine sulfation is important for cell signaling and is also target for drug design 3 4 . Sulfated tyrosine residues are important for the interaction of HIV's gp120 protein to CCR5 in order to get entry to CD4+ Tlymphocytes 5 . Sulfated PSGL-1 is involved in the infection process of Enterovirus71 virus, which is responsible for neurological diseases in children 6 . Sulfated tyrosine residues are found in approximately sixty immune system proteins 7 8 . The sulfation of protein tyrosine residues is carried out as enzyme catalyzed transfer of negatively charged sulfate group from the sulfate donor 3'-phosphoadenosine 5'phosphosulfate (PAPS) to the hydroxyl group of a tyrosine leading to a tyrosine O4-sulfate ester and Adenosine 3',5'-diphosphate 2 . Two isoforms TPTST-1 and TPST-2, exist in humans 9 . TPST has a type II transmembrane topology i.e. a short N-terminal cytoplasmic domain, a 17 residue transmembrane domain, and a luminal catalytic domain 10 . The enzyme has two N-glycosylation sites, that are four cysteine residues on the luminal oriented side of enzyme. 4 The crystal structure of the catalytic domain of human TPST-2 (denoted as TPST2∆C18,containing the sequence from G43 to L359) in complex with PAP (3'-phosphoadenosine 5'-phosphate) and the substrate (the peptide C4P5Y3) was solved at a resolution of 1.9 Å 11 ( Figure 1). The substrate peptide contains only one tyrosine residue which is a sulfate acceptor and is surrounded by six acidic residues, which provide the substrate negative charge of minus six 11 . The catalytic domain of TPTS-2, consists ofa α/β motif with a five-stranded parallel β-sheet, surrounded by α helices 11 . The TPST-2 exists as a homodimer, which is the catalytically active state of the enzyme, where two subunits of the dimer are designated as protomer A and protomer B 11 and both participate in the formation of the binding site and the catalytically active complex. The alpha helices α2-α4 play an important role in dimer formation in TPST-2 and residues R118 and R122 harboring α2-α3 also participate in the substrate binding of the other protomer and vice versa ( Figure 1). Residue W113 plays a key role in the dimerization process and W113A mutant appears to exists as monomer in solution based on experimental studies 12 . The crystal structure of  overall fold, atomistic interactions and contacts between its  both protomers, however it shows no information about how  conformational flexibility affects key interactions responsible  for the formation of the catalytically active dimer and the  structural changes, associated with the binding of the ligands  (cofactor and substrate). Molecular Dynamics is one of the most appropriate and broadly implemented method for studying protein flexibility and dynamics which provides atomistic insights that cannot be obtained experimentally 13 14 15 16 17 18 . In a recently accepted paper 16 we comprehensively analyzed the conformational flexibility of the wild type full complex of TPST2 (WTFC), which includes the enzyme, substrate and the cofactor in their catalytically active dimer form. We studied the mutations of key residues in the substrate binding site and the catalytic site and analyzed the correlated motions in the WTFC. However, it is essential to understand the atomistic nature of the structural changes which occur upon the binding of each ligand, PAPS and C4P5Y3. MD simulations can provide the needed structural insights, as a fast alternative to the experimentally determined structures, as far as the only available crystal structure so far is the one of the full complex (WTFC) 12 . In addition it is very important to understand atomistic details of the conformational effects associated with the formation of the protein dimer, which is the catalytically active form of TPST-2. Computational analysis of the effect of W113A mutation can provide unique atomistic insights into the dimerization process. As far as there is no available experimentally derived structure of this key mutant, atomistic MD can be reliable and fast alternative to gain the vital missing structural information. In order to complete the missing knowledge, we extended our study 16 and made simulations to answer the following questions: i) how the conformational flexibility influences the protein-protein interactions in the dimer formation; ii) what is the atomistic nature of W113A which is a key mutation preventing the dimerization iii) how does cofactor and substrate binding influence the structure of the apoenzyme (APO), the enzyme with bound cofactor complex (EC) and the enzyme with bound substarte complex (ES); iii) how does the correlated motions in TPST2 change in the different states of the enzyme upon binding of the ligands.

Methods and Computational Details
The conformational flexibility of different states of TPST-2 enzyme (APO, EC, ES and W113A mutant) was explored using atomistic molecular dynamics (MD) simulations for 100 ns, by means of Gromacs 4.5.5 package 19 20 21 with GROMOS96 43a1 22 force field (Table S1). The system setup and preparation was done as in reference 16 . The crystal structure of TPST-2 12 with PDB code 3AP1 was used. Energy minimization in gas phase was performed by using the steepest descent 23 and conjugate gradient 24 until the maximum force was smaller than 100 KJ/mol -1 /nm -1 . Subsequently a periodic box was defined and periodic boundary conditions were applied. The box size was set to 1.0 nm by using Single Point Charge (SPC) 25 , water model was used to solvate the system and Na + and Cl ions were used to neutralize its total charge (Table S1). New energy minimization was performed for the whole system using again the steepest descent and the conjugate gradient algorithms. After that position restrain molecular dynamics was carried out for 50 ps in NVT ensemble 26 at constant temperature of 300K with time step of 0.002 ps. The productive 100ns MD simulations were performed in NPT ensemble at constant temperature of 300K with initial velocities taken from Maxwell velocity distribution at 300K and an integration time step of 0.002 ps. The system was kept at 300K by Berendsen temperature coupling and Parrinello-Rahman pressure coupling were with constant (τ T ) of 0.1 ps and 1 bar pressure, time constant (τ P ) of 0.5 ps. Electrostatic interactions were treated by the Particle Mesh Ewald (PME) 27 . The cut off for the Coulomb interactions was 1.0 nm, the Fourier spacing was 0.135 nm. The van der Waals interactions were modelled using Lennard-Jones potential with cut off distance set to 1.4 nm (rvdw) and using switching function. The LINCS algorithm 28 was used to maintain all the covalent bonds involving hydrogen atoms rigid. The coordinates were saved after every 20 ps.

Analysis of Molecular dynamics simulations
The analyses of the trajectories obtained from the simulations such as root mean square deviation (RMSD) of Cα atoms of the protein with respect to minimized crystal structure, root mean square fluctuations (RMSF), radius of gyration (Rg), electrostatic interactions, hydrogen bonds, solvent accessible surface area (SASA), and cluster analysis were carried out with the Gromacs tools. The visualization of molecular structures and trajectories was done using VMD 29 software. The domain cross-correlation analysis (DCCA) was performed using the Bio3D package 30 . The Cij matrix were used to represent a Please do not adjust margins Please do not adjust margins cross correlation between the ith and jth atoms and ranges from -1 to +1. To create the cross correlation Cij matrix the Cα atoms of the protein, substrate and all atoms of the cofactor were included. The positive values in DCCA map reflect the correlated motions and the negative values reflect the anticorrelated motions .

Interactions between the Protomers and Dimer Formation
The equilibration of the MD simulations done for the present analysis has been discussed extensively 16 . The TPST-2 enzyme exists as a dimer both in the crystal structure and in vivo 9 . The formation of the substrate binding site in WTFC is a result of the dimerization with the participation of α2, α3 and α4 helices. During the simulations, the aliphatic parts of the side chains of R110 and R105 make hydrophobic interactions with the aromatic ring of W113 ( Figure S1) with average distance of 3.8 Å and 4.1 Å respectively for 100 ns in WTFC. These interactions are a result from conformational flexibility and are not seen in the crystal structure. The side chains of M109, L133 and I102 also make hydrophobic interactions ( Figure S31) with average distance of 3.8 Å, 4.0 Å and 4.4 Å respectively with the aromatic ring of W113 ( Figure 2). In order to understand the structural effect of a key mutation for the dimerization -W113A, we performed 100ns MD simulation of in silico generated structure of the mutant. The structure equilibrated (in two separate runs) in similar pattern as the wild type ( Figure S2, Table S3). In the trajectory of the W113A the residues in the vicinity of the W113A mutation site showed higher RMSFs than in the WT FC ( Figure S3). Importantly RMSFs of W113A also showed sensitive increase in fluctuations of α2, α3 and α4 (residue 105-145) α-helices, which are involved in the dimerization of WTFC (Figure 3). In the W113A mutant the distance between the A113 side chain and the side chains of R110, R105, M109 and L133 is higher than that in the WTFC during the MD trajectory ( Figure  S5 Table S2). This difference indicates about the roles of the above interactions in the formation of the active dimer. The dimerization is very important for the formation of the substrate binding region in WTFC 12 . The RMSFs of the substrate itself and the substrate binding region in the W113A mutant are higher than in the WTFC (Figure 4). These findings indicate the importance of the dimer for the proper formation and the stability of the substrate binding site and are in agreement with the experiments showing that the W113A mutant has reduced the enzyme activity 75 %, as compare to the wild type 11 . Figure 3 The Root Mean Square Fluctuations (RMSF) of mutant W113A (cyan) with respect to wild type TPST-2 (black) for 100 ns trajectory using C α atoms. The region α2-α4 (marked with red circle) involved in dimerization shows increase fluctuation in W113A.
A single mutation might influence not only the very local structure, surrounding the site of the mutation, but also could have a more global effect on the protein structure and interactions. In the context of dynamics a single mutation might influence the local flexibility around the site, but would also have sensitive impact on the correlated motions between different parts of the protein molecule. Unique insight about the correlations between movements of different parts of the protein can be provided by the dynamic cross-correlation Analysis (DCCA). In order to evaluate the effect of the W113A mutation on the correlated motions in TPST-2, we performed DCCA of the W113A mutant and compared it to the pattern of motions in WTFC. Importantly the DCCA reveals that in the mutant, W113A, there are more antcorrelated regions than in the WT FC ( Figure 5). In addition in W113A mutant, the region (α2, α3 and α4) shows reduced correlation with respect to α2, α3 and α4 of WTFC. These residues are involved in dimerization ( Figure 4B). Moreover, residues 180-210 (3'-hosphate binding otif and βe i.e. substrate binding region) of protomer A show anticorrelation towards residues (α2, α3 and α4) of protomer B, which explains the role of W113 in the dimerization and is in agreement with the experimental data about this mutant. The MD analysis of the TPST-2 dimer complements the experimental crystallographic and kinetic studies, providing information about new hydrophobic interactions which are not observed in the crystal structure. In addition the simulations of the W113A dimerization, provide structural explanation of its functional effect, namely the mutation leads to increase of distances to important residues, involved in the dimer stabilization, increase in the flexibility of the substrate and the substrate binding site and enhanced flexibility of the alphahelices (α2-4) involved in the dimer stabilization.

Conformational Flexibility and Correlated Motions upon binding of the Cofactor and the Substrate
The MD study of the whole complex, containing the protein dimer, the substrate and the cofactor (WTFC) 16 has revealed the key importance of the conformational flexibility which cannot be revealed by experimental methods. However, it is crucially important to understand the atomistic nature of the structural changes resulting from the binding of each ligand -PAPS and C4P5Y3. As far as the only available crystal structure so far is the one of the full complex (WTFC) 12 , MD simulations can provide the needed structural insight as fast alternative of the experimentally determined structures 31 . More importantly, an MD analysis of each of the APO, EC and ES in reference to the WTFC would provide unique knowledge about important changes in the correlated motions upon the binding of the substrate and the cofactor which cannot be obtained experimentaly.. In order to provide these insights into the conformational changes associated with binding of the cofactor PAPS and the substrate peptide C4P5Y3, we simulated for 100nsec: i) the apoenzyme (APO) built from protomers A and B; ii) the apoenzyme plus the cofactor PAPS complex (EC); iii) the apoenzyme plus the substrate C4P5Y3 complex (ES) and iv) the substrate in water (S) and compared them with the MD simulations of WT FC 12 .

Stability and Flexibility of APO, EC, ES and S
The RMSDs of APO, EC, ES, and substrate in water are represented in Figure S5 and  The ES complex showed increased fluctuations in the 5'PBS (5'phosphosulfate binding motif) region, which indicates that the binding of PAPS has stabilizing effect on the 5'PBS. The EC shows an increase in fluctuation of R78, S285 and some part of the 3'PBS (3'-phosphate binding motif) region thus suggesting that substrate binding has a stabilizing effect for the cofactor and substrate binding sites in WTFC ( Figure S6). The average value of radius of gyration (Rg) of WTFC is 28.73 Å , the Rg of EC increased to 30.2 Å and for ES is 28.8 Å (Table S4, Figure S7) demonstrating the effect of the substrate binding on compression of the structure of the full complex. These results show that the binding of both the substrate and the cofactor contribute for the stability of the full active complex; however, the effect of the substrate binding is much more prominent. In the EC complex R78 makes a stable interaction with the 5'PBS region of the cofactor as in WT FC ( Figure S8). The S285 hydroxyl group is hydrogen bonded to OAF residue of 5'PBS site during only 54% of the simulation. The hydroxyl group of S285, makes interactions with the OAD oxygen of the 3' phosphate group in contrast to the WTFC where it interacts with the OAF of 5' PBS region. The side chain of R183 in WTFC stabilizes the 3'PBS phosphate OAB oxygen atom (average distance 3.63 Å) however in EC complex it makes a weaker interaction (average distance 5.1 Å). The RMSF of PAPS in the EC complex shows higher fluctuations compared to WTFC ( Figure S9) which indicates the stabilizing role of the substrate for PAPS cofactor. The changes in the above atomistic interactions illustrate the key structural changes which occurs upon binding of the substrate and the cofactor and reinforce the idea of accounting for the conformational flexibility whilst modelling protein-ligand interactions (e.g. in the docking studies). Figure S6 The RMSD of substrate peptide of the wild type TPST-2 (black) compared against substrate peptide only in water (blue) for 100ns using C α atoms.

Correlated Motions
An insight into the correlated motions between different parts of the enzyme molecule can be provided by dynamic cross correlation analysis (DCCA) 32,33 . The DCCA of APO, ES and EC are shown in (Figure S10-S12). The positive regions in the map are represented with red and yellow and are associated with correlated motions whereas anti-correlated motions (negative regions) are shown in light green and blue. In the APO trajectory there are less correlated and more anti correlated motions as compared to WTFC. This result suggests that the binding of the PAPS and C4P5Y3 (substrate peptide) influences not only the local structure and interactions, but also affects more globally the system of protein motions. In the ES complex, the correlated region in ( Figure S11, region a) showed fewer residues involved in correlated motion with respect to the WTFC (region a) and there is no strong correlation between the substrate and the enzyme residues as in the WTFC indicating about the apparent role of the cofactor for good stabilization and orientation of the substrate. In the EC, the extent of both positive and negative correlation is also lower than WTFC again demonstrating the effect of the ligand binding more globally on the correlated motions in the molecule. The above differences with respect to the WTFC are another indication that the protein ligand binding process is dynamical in nature, influences not only the local structure but also have more more global effect on the protein.
The current MD study reveals the importance of conformational flexibility for the dimerization of TPST-2.
Enzymes often exist as dimers in physiological conditions and carry out their biochemical functions in such state; therefore computational modelling of the dimer state provides important insight about their structure-function relationships which can be used for further modelling of their reaction mechanisms, e.g. as in the case of the multiscale modelling of the COX-1 enzyme 34 ).
Our results support the idea that the conformational changes associated with interactions with ligands (substrates, cofactors, inhibitors) can lead to changes of the system of interactions within the entire enzyme molecule, rather than only in the local binding site. In similar way for example, the network of interactions between the aromatic chromophores, which generate the near-UV circular dichroism in class-A βlactamases and human carbonic anhydrase are sensitively influenced by the conformational changes, related to the binding of ligands, 13,[35][36][37] . The MD studies presented here assert about the importance of accounting for the conformational flexibility in computational docking studies of protein-ligand interactions.

Conclusions
The present MD study reveals the key effect of the protein flexibility on the dimerization of TPST-2 and the interactions between the apoenzyme, the cofactor PAPS and the substrate peptide C4P5Y3. New hydrophobic interactions, important for the protomer-protomer interactions are identified and the structural effect of the key dimerization mutation W113A is explained in the absence of experimentally-derived structure of the mutant. Our results show that the binding of the substrate and cofactor influence not only the local interactions in the binding site, but also affect of the correlated motions in the entire molecule. The study asserts about the importance of accounting for conformational dynamics in studying proteinligand and protein-protein interactions.