Sushmita Pradhan, Himakshi Sarma, Barsha Bharadwaz and Venkata Satish Kumar Mattaparthi*
Molecular Modelling and Simulation Lab, Department of Molecular Biology and Biotechnology, Tezpur University, Tezpur-784 028, Assam, India
Received date: June 1, 2017; Accepted date: June 21, 2017; Published date: June 26, 2017
Citation: Pradhan S, Sarma H, Bharadwaz B, Mattaparthi VSK (2017) Comparative Study on the Binding Affinity of Methimazole and Propylthiouracil to Thyroid Peroxidase as an Anti-Thyroid Drug: An In-silico Approach. J Mol Imag Dynamic 7:131. doi: 10.4172/2155-9937.1000131
Copyright: © Pradhan S, et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use,distribution, and reproduction in any medium, provided the original author and source are credited.
Visit for more related articles at Journal of Molecular Imaging & Dynamics
Graves’ disease (GD), an autoimmune disorder, scars majority of women worldwide, causing hyperthyroidism, Graves’s ophthalmopathy and goitre. Thyroid Peroxidase (TPO) is an active target of anti-thyroid drugs, Methimazole and Propylthiouracil, which inhibit the enzyme function of catalysing the thyroid hormones synthesis. Most of the protein-drug interaction studies so far have been focussed mainly at in vivo level, or by using Myeloperoxidase and Lactose peroxidases as TPO surrogates for the same. This makes the molecular interaction of TPO with the drugs crucial to understand. In this study, we used the molecular dynamics (MD) to study the molecular interaction differences between TPO201-500 and both drugs. The binding free energy calculation done using Molecular Mechanics Poisson–Boltzmann (MM-PBSA) and generalized Born and surface area continuum solvation (GBSA) results indicated that both drugs bind strongly to TPO201-500, with Propylthiouracil having slightly higher binding energy than Methimazole. We found that both drugs interacted with the residues- Asp238, His239, Phe243, Thr487 and His494 through hydrophobic interactions and formed stable hydrogen bonds with residue Arg491 of TPO201-500. Since both drugs engage residues- Asp238, His239 and His494 which falls within the proximal heme binding site and catalytic site of TPO201-500, we can conclude that these drugs may be conducive in inhibiting the enzymatic activity of TPO201-500.
Graves’ disease; Thyroid peroxidase; Methimazole; Propylthiouracil; Molecular dynamics; MM-PBSA/GBSA
Graves' disease (GD) is an autoimmune Thyroid disorder (AITD), typically characterized by the presence of thyrotoxicosis, goiter, and ophthalmopathy (Grave’s ophthalmopathy) . It is one of the major causes of hyperthyroidism in the geographical areas with Iodine abundance . GD has reported a frequency of 20-30 cases per 100,000 individuals each year, with the possibility 3% of women and 0.5% of men acquiring GD during their lifetime . The main ground cause of GD is the binding of circulating IgG antibodies to G-protein– coupled thyrotropin receptor (TSHR), resulting in the activation of Thyroglobulin gene which codes for thyroid hormonestriiodothyronine (T3) and thyroxine (T4) [4,5].
These antibodies are thyroid-stimulating antibodies (TSAbs) , - often called TSH receptor autoantibody (TSHRAb or TRAb)  that replicates the effect of the TSHR’s natural ligand Thyrotropin (TSH) by stimulating cyclic adenosine monophosphate (cAMP)–dependent signal transduction, thereby activating the Thyroglobulin gene present within the thyroid gland [6,8]. Thyroglobulin gene is naturally activated by hypothalamus upon the binding of TSH to TSHR. Since the activation of Thyroglobulin gene by the binding of TSIs to THSR is not natural, there is no feedback system leading to the excess production of thyroid hormones: triiodothyronine (T3) and thyroxine (T4) that is catalyzed by Thyroid Peroxidase (TPO) [6,9]. This activation also triggers follicular hypertrophy and hyperplasia, thyroid enlargement .
The synthesis of triiodothyronine T3 and thyroxine T4 hormones from Thyroglobulin gene is catalyzed by Thyroid Peroxidase (TPO) upon the iodination and coupling of tyrosyl residue in Thyroglobulin  and is hence the target site of anti-thyroid drugs for the treatment of GD. TPO201-500 is a 933 residues transmembrane homodimer [11-13], affixed to the apical membrane of thyroid follicle cells by its TMD region on its C-terminal end (residues 847–871) [14,15]. The Xray crystallised structure of TPO was determined 16 years ago, but of low resolution [16,17]. TPO shares 47%, homology with eosinophil peroxidase (EPO), 48% with Lactose Peroxidase (LPO), and 47% with myeloperoxidase (MPO) [18,19]. Based on these homology studies, it was found that Asp238, His239, Glu399, and His494 were important for proximal heme linkage. The proximal histidine residues present in heme peroxidases are crucial in the redox properties of heme iron for catalysis [20,21]. The homology modeling of TPO was done recently, reporting that its cis and trans conformations are involved equally in its antigenicity .
The antithyroid drug (ATD), Propyl- thiouracil (6-propyl-2- thiouracil) (PTU) and Methimazole (1-methyl-2-mercaptoimidazole, Tapazole) (MMZ), are prescribed exclusively for controlling the hyperthyroidism in GD [23-26]. Both ATDs obstructs the iodination of tyrosyl residues in Thyroglobulin catalyzed by Thyroid Peroxidase, thereby arresting the production of thyroid hormones [27-29]. MMZ is preferred to PTU, due to its better adherence and bioavailability [23,30-32]; however, PTU is also recommended for pregnant women [33,34]. Irradiation studies too have been done to treat hyperthyroidism .
The drug binding mechanism of PTU/MMZ to TPO have only been done in-vivo [36-39], however the molecular insights are yet to be understood. Due to the close homology, LPO has been used as surrogates to study interaction of TPO with anti-Thyroid drug MMZ. MMZ-LPO catalyzed oxidation of 2,2- azio-bis-3-ethylbenthiazoline- 6-sulfonic acid (ABTS) and H2O2 was carried out where the turnover was estimated to be ~7.0 lM. . The enzyme inhibition study between LPO- MMZ showed that inhibition occurred at the heme group and side chains Arg255, Glu258, and Leu262 of LPO. Gln105 and His109 of the LPO too were involved in the enzyme inhibition .
In this study, we have taken up PTU and MMZ as our drug of choice considering the disease. The protein structure of TPO201-500 was constructed first, and its stability was asserted. The complex structures of TPO201-500 with both the drugs were constructed using PatchDock server [42,43] which gave the results based on the maximum surface area and minimum atomic contact energy. LigPlot+ software v.1.4.5  was employed to investigate the various bonding and non-bonding interactions between the drugs and TPO201-500. Since the residues- Asp238, His239 and His494 belonging to the heme linkage region of the protein TPO were seen to be involved in the bond interaction with both drugs, these sites will now no longer be available for the iodination and coupling of tyrosyl residues by TPO for the synthesis of thyroid hormones. R491 was also found to play crucial role in the drug-protein communication by forming hydrogen bond. Binding energy calculation further confirmed that both drugs possessed high affinity towards TPO201-500.
The antechamber  module of AMBER12  software package was used for prepping the ATDs prior to the simulation. Each of the TPO201-500-ATD complexes was subjected to Molecular Dynamics simulation to study its stabilizing properties. In addition, Poisson– Boltzmann and generalized Born and surface area continuum solvation (MM/PBSA and MM/GBSA) [47-51] techniques were endorsed to calculate binding free energies of MMZ and PTU to TPO201-500. PBTOT for TPO201-500-MMZ and TPO201-500-PTU was -11.05 and -14.02 kcal mol-1 while GBTOT for TPO201-500- MMZ and TPO201-500-PTU was found -9.43 and -15.47 kcal mol-1 respectively, which put forth the inhibitors to have high affinity to TPO201-500 in the dynamic system. Our findings affirm that both the drugs are equally effective as an anti-thyroid drug. Moreover, our results may further provide guidance to developing more potent inhibitors of TPO and preventing thyroid hormone production in GD.
The amino acid sequences were obtained from UniProt , ID No. P07202 by querying the TPO201-500 protein sequence against RCSB PDB web server [53,54] domain of Thyroid Peroxidase protein sequence (residues 201-500) was submitted to Iterative Threading ASSembly Refinement server (ITASSER) [55-57] for the 3D structure prediction. From I-TASSER we obtained 5 best models ranked according to their C-scores out of which best model is chosen based on C-score, Template Modelling (TM) score  and Root Mean Square Deviation (RMSD) score. Using RAMPAGE  and VERIFY 3D server [60,61] we validated the structure of TPO201-500.
Based upon the 3D structural validation checks, Model 1 showed consistently good results among the rest, hence chosen as an ideal structure of TPO201-500 for further studies. SOPMA (Self-optimized prediction method alignment) server  was adopted in studying the secondary structure of TPO201-500. Anti-thyroid drug: Propylthiouracil (PTU) and Methimazole (MMZ/MMI) were chosen for this study. 3D conformer of both the ATDs were procured from PubChem server  in SDF format, which was later converted to pdb file using open Babel software v.2.3.2 .
The initially constructed 3D structures of TPO201-500 with PTU and MMZ were processed separately for the rigid molecular docking using the PatchDock [42,43] web server. PatchDock works on the principle of geometric surface docking algorithm. The clustering score of 1.5 Å RMSD was given to remove the unwanted models. The ten complexes of TPO201-500-PTU and TPO201-500-MMZ were obtained from the server. Based on highest surface area and lowest atomic contact energy, one complex each of TPO201-500-MMZ and TPO201-500-PTU were selected.
Determining interacting residues and bond types
The best structure of TPO201-500 in complex with PTU and MMZ processed by PatchDock was selected among the rest considering the minimum atomic contact energy and maximum surface area. The best structure was visualized in LigPlot+ software v.1.4.5  for evaluating the interacting residues of TPO201-500 interacting with the drugs. Interface residues are defined as the residues with a contact distance less than 6 Å from the interacting partner.
Molecular dynamics simulations of thyroid peroxidase (TPO201-500) with anti-thyroid drugs (ATD) - methimazole and propylthiouracil
The molecular dynamics study was performed using AMBER12  software package where ff99SB force field  parameters were used for protein (TPO201-500) and generalized AMBER force field (GAFF)  parameters of the anti-thyroid drugs−Methimazole and Propylthiouracil were established with the Antechamber  program, wherein AM1-BCC [67,68] was used to fix the proportionate partial atomic charges. Using xLEaP and antechamber module in AMBER12, we prepared the initial coordinate and the topology files for two the complex systems.
Xleap also added all the missing parameters to the complexes. All the canonical systems were then neutralized with Na+ ions, followed by the hydration with TIP3P  water molecules in a box with a buer size of 10 from the solute in x, y, and z axis. The initial minimization of the solvated structures was done imposing the constraints over the solute for the first 500 steps using the steepest descent algorithm followed by another 500 steps with a conjugate gradient method. Second minimization was carried out freely for 500 steps. 8 Å was applied as the cut off for both the minimization steps. Gradual heating was done to both the complexes with the successive increase in the temperature from 0-300K. The generated ensemble was then equilibrated for 100ps using the NPT ensemble at 300 K and 1 atmospheric pressure. Full MD production run was carried out till 15ns.
Periodic boundary conditions along with all the electrostatic interactions were maintained throughout the simulations using the particle mesh Ewald (PME) method [70,71] under isothermal and isobaric conditions respectively. The temperature was governed using Berendsen thermostat . The production run was attained with the NPT ensemble for all the formerly mentioned Thyroid Peroxidase-Anti thyroid drug complexes for 15 ns with 2 fs time step. Furthermore, the shake algorithm  was used to restrain all the bonds including the hydrogen atoms.
To check the stability of the obtained MD trajectories, each trajectory were clustered to acquire comparative Root Mean Square Deviation (RMSD), Radius of Gyration (Rg), B-factor and Solvent Accessible Surface Area (SASA) using the CPPTRAJ algorithm  of Amber12 package. We also proceeded further to analyze hydrogen bond properties based on all possible hydrogen donors (HD) and acceptors (HA) in the system, the percentage occupancy of hydrogen during each step of MD trajectory and considering bond distances and bond angle made by HA–H–HD atoms. Here, the hydrogen bonds having consistency for more than 5% of the simulation period were taken as an output result. Molecular graphics and analyses of the trajectories obtained were performed using UCSF Chimera package  and VMD v.1.9.3 .
Binding Free Energy (BFE) calculations
The binding free energy between the protein and the drugs were assessed by considering the thermodynamics pattern of Figure 1. The binding free energies of two anti-thyroid drugs- Methimazole and Propylthiouracil to Thyroid peroxidase were computed using MMGBSA (Molecular Mechanics Generalized Born Surface Area) and MM-PBSA (Molecular Mechanics Poisson–Boltzmann Surface Area approach [47-51] in AMBER program. We took 100 conformational snapshots extracted from the last 1 ns of MD trajectories. Independent MM-GBSA/PBSA analysis was done for (i) Thyroid peroxidase (TPO201-500) (ii) anti-thyroid drugs (ATD) - MMZ/PTU and (iii) the complex- TPO201-500-MMZ/PTU. The binding free energies of TPO201-500 with both the drugs were calculated using equation 1.
ΔGbinding = ΔGcomplex- [ΔGTPO+ ΔGATD] --------------->(1)
Where G-binding corresponds to the total binding free energy, Gcomplex, GTPO and GATD are the solvation free energies of the complex, thyroid peroxidase and anti-thyroid drugs, respectively. BFE as stated by second law of thermodynamics can also be summed as the difference between enthalpy ΔH and the conformational entropy -TΔS [48,51] of the system which can be represented by the equation 2.
ΔGbinding = ΔH - TΔS ---------------> (2)
ΔH can be explained as sum of the changes of the gas-phase molecular mechanics energies (ΔEMM), polar and non-polar solvation energy (ΔGsolv); see equation 3. Also, ΔEMM is the sum total, of the differences in the bond, angle and dihedral energies (ΔEinternal), van der Waals (ΔEvdW), and electrostatic energies (ΔEele), as shown in the equation 4. Accordingly, equation 5 divides ΔGsolv into polar (ΔGPB/GB) and non-polar solvation free energy (ΔGsurf);
ΔH = ΔEMM + ΔGsolv ---------------> (3)
ΔEMM = ΔEinternal + ΔEvdW + ΔEele ---------------> (4)
ΔGsolv = ΔGPB/GB + ΔGsurf ---------------> (5)
The terms Eint, EvdW, and Eele are internal energy, van der Waals interaction, and electrostatic interaction of the systems respectively. The terms ΔGPB/GB and ΔGsurf are polar and non-polar solvation free energies. The modified GB model of Onufriev  was considered in the determination of the ΔGGB. A value of 80 was considered for the internal dielectric constant of the solvent (water) while, the external dielectric constant of the solute was set to 1. The non-polar solvation free energy (ΔGsurf) for both the systems was estimated equation 6.
ΔGsurf = γ × SASA +β ---------------> (6)
Where SASA is the abbreviation of solvent-accessible surface-area. In case of PB method, γ was to 0.00542 kcal (mol-1 Å-2) and β = 0.92 kcal mol-1, while for the GB method, γ was customized to 0.0072 kcal (mol-1Å-2) and β = 0 kcal mol-1 respectively . 1.4 Å was the probe radius of the solvent used for these calculations.
Structural modelling and validation of thyroid peroxidase
We modelled the structure of Thyroid Peroxidase (TPO201-500) from I-TASSER (Figure 2). The server provided us with five models using the threading algorithm. The ideal model was selected based upon the C-score computed from the comparative structural significance. The C-score for TPO201-500 was 1.20. The overall quality of model was achieved by various bioinformatics tools. TM scored TPO201-500 in the rage of 0.88 ± 0.07 and RMSD estimated was 3.8 ± 2.6 Å. The structural justification of the ITASSER yielded structure was done by RAMPAGE server. We observed 81.5% of residues were within most favored region, 11.4% in allowed region and 7.0% in disallowed region. Verify-3D server ascertains the accuracy of the predicted secondary structure (3D) model with its respective residues (1D) by assigning a structural class based on its location and environment. 81.6% of residues had an averaged 3D-1D score >= 0.2 which confirms the accuracy of the predicted model to near correct to the actual structure. SOPMA server identified 39.33% α-helices, 47.67% random coil, 7.33% or extended strand and 5.67% of β-strands respectively in the modeled structure of Thyroid Peroxidase.
Protein-drug interaction study
It is an evident fact that Methimazole and Propylthiouracil interferes in the synthesis of thyroid hormones by regulating the enzyme inhibition of thyroid peroxidase. The in-vivo/in-vitro binding mechanisms of PTU/MMZ to TPO201-500 have also been done exclusively, yet the interactions at the molecular level are not known. So as to understand the drug-protein interaction at the molecular stage, we docked Thyroid peroxidase with both the drugs using PatchDock, server. This particular server works on the principle of geometric hashing which is implied as a best molecular shape complementarity between the docked structures considering the wide interface area and lesser steric hindrances. The docking gave us array of docked structures for TPO201-500-MMZ/PTU complex. These complexes were graded basis of geometric shape, interface area and atomic contact energy as shown in Figures 3 and 4. The best ranked structure in case of TPO201-500-MMZ complex had a geometric shape complementarity score of 2300, approximate interface area of 260.2 Å2 and atomic contact energy of (ACE) of -152.3 kcal mol-1 (Solution Structure 1 from Figure 3). Likewise, the top scored structure for TPO201-500-PTU was found to have the geometric score of 2912, interface area of about of 328.8 Å2 and ACE of -146.6kcal mol-1 (Solution Structure 1 from Figure 4).
Next, we examined the hydrogen and hydrophobic interaction occurring in protein-drug complex using the Ligplot+ tool as shown in Figures 5A and 5B. Figure 5A details the interaction between Thyroid Peroxidase and Methimazole. Arg491 formed a hydrogen bond with the thiol atom i.e., Sulphur at the distance of 3.25 Å. Hydrophobic bonds are formed between the protein and the drug by residues Asp238, His239, Phe243, Thr487 and His494. Figure 5B showcases the molecular interaction between TPO201-500 and Propylthiouracil. We can clearly observe the hydrogen bond formed between residue Arg491 and Sulphur atom of the drug PTU at the distance of 3.04 Å. Amino acids Asp238, Gln235, His239, Phe243, Glu399, Thr487 and His494 formed the hydrophobic bonds with the drug. In both the scenarios, drug is seen to have bonded with proximal heme binding residues of the protein- Asp238, His239 and His494.
Figure 5: LIGPLOT representation of Anti-Thyroid drugs-Methimazole (MMZ) and Propylthiouracil (PTU) binding to Thyroid Peroxidase (TPO). Ligand, protein, and hydrogen bonds are in thick purple, thick brown and green dash lines, respectively. Hydrophobic contacts forming residues are shown in sky blue colour as an arc with spikes towards ligand atoms (A) Thyroid Peroxidase (TPO) bound to Methimazole (MMZ). (B) Thyroid Peroxidase (TPO) bound to Propylthiouracil (PTU).
Since these residues are bound to the drug, which renders the active site of the protein unavailable for its ritual catalyzing activity thereby inhibiting the production of the thyroid hormones. Another homolog study of TPO was done in LPO by Singh et al., where they found that the heme moiety residues held MMZ tightly on one side while hydrophobic amino acids - Arg255, Glu258, and Leu262 held it on the other end. Many extensive studies on LPO as TPO surrogate has supported the importance of heme binding and catalytic site residues mainly Hisidines and Aspartic acid in mammalian peroxidases [79-81]. A review on MPO highlighted its ligand binding site to be present within the distal heme pocket which upon the binding of PTU, isoniazid could lead to irreversible inactivation of the enzyme . The homology studies done by Taurog in 1999 and Furtmüller et al. have pointed out the essentiality of the proximal histidine residues for the prime functioning of TPO201-500 to catalyse the iodination of tyroysl residues in thyroglobulin. And we have noted these histidine residues His239 and His494 to have been engaged by anti-thyroid drugs due to hydrophobic linkage, also relieving the protein of its catalytic role. Comparative study on the effective dosage of Methimazole and Propylthiouracil in Patients with GD by Nakamura et al.  revealed that 15 mg/d of MMI is prescribed for mild and moderate GD, 30 mg/d of PTU is recommended for severe forms of GD. A team of scientists from Spain and Argentina  determined the pharmacological activities of a copper based propylthiouracil ([Cu(PTU)2]2, and concluded that this combination of PTU would be good candidate due to its antioxidant, antimicrobial and alkaline phosphatase activity.
MD simulation study on the TPO201-500-MMZ/PTU complex
The molecular dynamics simulations were administered to probe the efficacy of both the anti-thyroid drugs- Methimazole and Proplythiouracil as an inhibitor for thyroid peroxidase. Extracting the trajectory output, we analyzed the stability of TPO201-500-MMZ/PTU complexes in an explicit environment. The conformational dynamics of TPO201-500 complexed with the drugs-MMZ and PTU were estimated by plotting the Root Mean Square Deviation (RMSD) by taking their respective equilibrated structures as a reference. The comparative RMSD graph for the Cα atoms of the TPO201-500- MMZ/PTU complexes as a function of time can been seen in Figure 6A. MMZ bound TPO201-500 started to settle at the RMSD value of around 5 Å while MMZ bound XPA was exhibiting fluctuations even at the RMSD value of ~ 7.5 Å. The stability of PTU over MMZ in the biological system can also be explained from the study done by Yoshihara et al., where they found that over exposure of the in utero foetus to MMZ increased the risks of congenital malformations and least with PTU prescriptions.
Figure 6: Comparative MD analyses of TPO-MMZ and TPO-PTU. (A) RMSD values for TPO-MMZ and TPO-PTU complexes relative to the starting structure during MD simulation. (B) Solvent accessible surface area (SASA) values for TPO-MMZ and TPOPTU complexes. (C) Radius of Gyration (Rg) values of TPO-MMZ and TPO-PTU. (D) B-factor distribution for TPO-MMZ and TPOPTU.
Figure 6B shows the solvent accessible surface area (SASA) for the both the protein-drug complexes. Both the systems, be it TPO201-500- MMZ or TPO201-500-PTU had proper access to the solvent provided in the water box throughout the simulation. TPO201-500- MMZ presented the SASA value of within 20000-22000 Å2 while PTUshowed SASA value in the range of 21000-22500 Å2. TPO201-500- PTU complex required much larger surface area to access the solvent in the system.
We determined the comparative radius of gyration (Rg) for both the drugs complexed with TPO201-500. Rg quantifies the compactness and the spatial dispersal of the atoms, of a particular component molecule, around an axis. Figure 6C shows the radius of gyration for TPO201-500 bound to Methimzole and Propylthiouracil as a function of time. Methimazole bound TPO201-500 showed oscillations from 6 ns onwards, while TPO201-500-PTU was settled till 9 ns, then it started oscillating.
Next, we analyzed B-factor graph for the backbone Cα-atoms for the TPO201-500 in complex with both the drugs to observe the variations in their structural dynamics from their average position. Figure 6D shows that both the systems presented more or less coherent atomic fluctuations, and the only difference between them being the higher peaks as seen in the TPO201-500-PTU complex, whereas TPO201-500-MMZ had shorter peaks meaning less atomic fluctuations.
Hydrogen bond analysis as seen in Figure 7 showed intermolecular hydrogen bonds formed between the enzyme and the drug during the simulation period. As per Figure 7A, we saw two hydrogen bonds formation in TPO201-500-PTU complex form, while only single hydrogen bond formed in case of TPO201-500-MMZ (Figure 7B) during simulation. We further calculated the inter-molecular hydrogen bond stability of both the systems.
Table 1 shows hydrogen bond occupancy of the inter-molecular hydrogen bonds formation as the function of time period. We found the important amino acids Arg491 had high fraction of occupancy in the complex of TPO201-500-MMZ while Asp238 and Phe485 had greater fractions for hydrogen bond formation in case of TPO201-500- PTU. Higher fraction of hydrogen bond definitely suggested the stronger interaction between the protein TPO201-500 and the ATDs. These occupancies align exactly with the Figures 7A and 7B for the total number of the hydrogen bonds formed within the systems. We simultaneously took comparative snapshots of both the TPO201-500- MMZ and TPO201-500-PTU systems at different intervals of simulation period and are presented in Figure 8.
|Anti-Thyroid drugs||Atoms of Hydrogens bonds||Fractions|
|Methimazole (MMZ)||N1:LIG- HG2:ARG491||0.6100|
|N1:LIG- OE1: GLN316||0.0021|
|Propylthiouracil (PTU)||N:LIG- OD1:ASP238||0.5346|
Table 1: H-bond occupancy among all MD simulations.
Both the systems presented the β-sheets formations. PTU bound TPO201-500 had loss of the α-helices over the simulation period.
Calculations of binding free energies
To measure the efficacy of the anti-thyroid drugs- Methimazole and Propylthiouracil in inhibiting the enzyme activity of Thyroid peroxidase in Graves’ disease, the binding free energies were evaluated using MM-PBSA/GBSA method as provided in AMBER program. The calculated values are summarized in Table 2. According to Table 2, the binding free energy of the TPO201-500-MMZ complex using MM-PB/ GBSA methods is found to be -11.05 (PBTOT) and -9.43 (GBTOT) kcal mol-1. TPO201-500-PTU complex had the binding free energy of -14.02 (PBTOT) and -15.47 (GBTOT) kcal mol-1. As per the data obtained, the drug-protein binding was considerably enhanced by the contributions of the van der Waals interactions (ΔEvdW), non-polar solvation energy (ΔGsurf) and electrostatic interactions (ΔGele).
|(kcal mol-1)||(kcal mol-1)|
|MM||ELE||-16.90 ± 4.46||-13.45 ± 3.61|
|VDW||-12.48 ± 1.78||-21.62 ± 1.25|
|GAS||-29.38 ± 4.30||-35.08 ± 3.09|
|PBSOR||-1.53 ± 0.05||-2.13 ± 0.11|
|PBSA||PBCAL||19.86 ± 3.67||23.18 ± 2.19|
|PBSOL||18.33 ± 3.69||21.05 ± 2.22|
|PBELE||2.96 ± 1.62||9.73 ± 2.56|
|PBTOT||-11.05 ± 1.81||-14.02 ± 2.09|
|GBSURF||-1.53 ± 0.05||-2.21 ± 0.11|
|GBSA||GB||21.48 ± 3.62||21.73 ± 2.73|
|GBSOL||19.95 ± 3.64||19.60 ± 2.80|
|GBELE||4.59 ± 1.30||8.28 ± 1.36|
|GBTOT||-9.43 ± 1.26||-15.47 ± 1.24|
ELE:Electrostatic energy as calculated by the MM force field,
VDW = Van der Waals contribution from MM,
INT = Internal energy arising from bond, angle, and dihedral terms in the MM force field. (This term always amounts to zero in the single trajectory approach),
GAS = Total gas phase energy (sum of ELE, VDW, and INT),
PBSUR/GBSUR = Non-polar contribution to the solvation free energy calculated by an empirical model.
PBCAL/GB = The electrostatic contribution to the solvation free energy calculated by PB or GB, respectively.
PBSOL/GBSOL = Sum of non-polar and polar contributions to solvation,
PBELE/GBELE = Sum of the electrostatic solvation free energy and MM electrostatic energy
PBTOT/GBTOT = Final estimated binding free energy in kcal mol-1calculated from the terms above.
Table 2: Binding free energies (kcal mol-1) between Thyroid peroxidase (TPO) and Anti-Thyroid drugs- Methimazole (MMZ) and Propylthiouracil (PTU).
The PBTOT or GBTOT values of both the complexes showed not much difference in their overall binding energy. For PBTOT, the differences were 14.02-11.05 = 2.97 kcal mol-1, while their GBTOT difference was however quite large i.e., 15.47-9.43 = 6.04 kcal mol-1. Based on these results, Propylthiouracil showed better results in enzyme inhibition than Methimazole. Nonetheless, both the drugs are equally qualified in target binding of Thyroid peroxidase and thereby preventing the production of thyroid hormones. A population study conducted by Azizi and his team in 2014 for the management of hyperthyroidism in Asia too showed the preference of PTU over MMZ by the clinician’s due to its better outcome. A density functional theory method was opted by Tappa and his group  to evaluate the inhibitory effect of the MMZ tautomers-1-methyl-1H-imidazole-2 (3H)-thione (M1) and 1-methyl-1H-imidazole-2-thiol (M2) where the effect of M1 was superior to M2 tautomer. Gupta and Kumar  did metabolic modelliing and tested the drug activity of PTU against TPO, TSHR and sodium iodide symporter (NIS) using systems biology approach and deduced that PTU is effective on TPO. They also observed TSHR to be potent target for therapeutic study. A study done on PTU, KClO4, or TSH  revealed that these compounds reduced the oxidative stress due to high concentrations of iodide in the thyroids of both Metallothioneins(MTs) MT-I/II KO and WT mice. Bhabak and Mugesh in 2010 , and Manna, Roy and Mugesh in 2013 , reported Selenium analogues of ATDs to have better inhibitory effect compared to Sulphur based thinamide drugs on TPO suggesting that this selenium based drugs may be the next potential ATDs.
In this computational study, the protein structure of thyroid peroxidase was modelled first, followed by the independent docking of the protein with MMZ and PTU. We have identified the residues that are involved in the drug-protein interaction. Both the drugs engaged the heme linkage active site residues- Asp238, His239 and His494, depleting the protein of its function to catalyze the iodination of tyrosyl residues on thyroglobulin for T3 and T4 hormone production. The hydrophobic interaction with residues Asp238, His239, Phe243, Thr487 and His494 was seen in both the cases. Proximal histidine residues His239 and His494 were also found to be interacting with the drug making the enzyme unfit to carry out its biological role. Arg491 was also found to be binding with the drug by hydrogen bonding. RMSD showed the stability of MMZ to be much higher than PTU. TPO201-500-PTU had higher values for SASA and B-factor graphs. Hydrogen bond analyses too substantiated the stronger and stable interaction shared by these anti-thyroid drugs to the enzyme. MMPBSA/ GBSA binding free energy results showed both drugs to be equipped in binding to TPO201-500 and in inhibiting the catalytic activity of TPO201-500 in thyroid hormone synthesis as seen in the literatures. Propylthiouracil showed comparatively higher binding efficiency than methimazole as per GBTOT analysis, yet there was not much difference between their binding energies considering the PBTOT results. After studying all the obtained evidences, we deduce that both the anti-thyroid drugs are competent enough to inhibit the enzyme activity of Thyroid peroxidase, and in process keeping the unwanted production of thyroid hormone at bay.
The authors extend their deepest gratitude towards Tezpur University for the start-up grant. This research was also supported by Department of Biotechnology (DBT) funded Bioinformatics Infrastructure facility (BIF) in the Department of Molecular Biology and Biotechnology at Tezpur University.