Homology modeling, simulation and molecular docking studies of catechol-2, 3-Dioxygenase from Burkholderia cepacia: Involved in degradation of Petroleum hydrocarbons

Catechol 2, 3-dioxygenase is present in several types of bacteria and undergoes degradation of environmental pollutants through an important key biochemical pathways. Specifically, this enzyme cleaves aromatic rings of several environmental pollutants such as toluene, xylene, naphthalene and biphenyl derivatives. Hence, the importance of Catechol 2, 3-dioxygenase and its role in the degradation of environmental pollutants made us to predict the three-dimensional structure of Catechol 2, 3-dioxygenase from Burkholderia cepacia. The 10ns molecular dynamics simulation was carried out to check the stability of the modeled Catechol 2, 3- dioxygenase. The results show that the model was energetically stable, and it attains their equilibrium within 2000 ps of production MD run. The docking of various petroleum hydrocarbons into the Catechol 2,3-dioxygenase reveals that the benzene, O-xylene, Toluene, Fluorene, Naphthalene, Carbazol, Pyrene, Dibenzothiophene, Anthracene, Phenanthrene, Biphenyl makes strong hydrogen bond and Van der waals interaction with the active site residues of H150, L152, W198, H206, H220, H252, I254, T255, Y261, E271, L276 and F309. Free energy of binding and estimated inhibition constant of these compounds demonstrates that they are energetically stable in their binding cavity. Chrysene shows positive energy of binding in the active site atom of Fe. Except Pyrene all the substrates made close contact with Fe atom by the distance ranges from 1.67 to 2.43 Å. In addition to that, the above mentioned substrate except pyrene all other made π-π stacking interaction with H252 by the distance ranges from 3.40 to 3.90 Å. All these docking results reveal that, except Chrysene all other substrate has good free energy of binding to hold enough in the active site and makes strong VdW interaction with Catechol-2,3-dioxygenase. These results suggest that, the enzyme is capable of catalyzing the above-mentioned substrate.


Background:
Catechol and its derivatives are key metabolic intermediates in the aerobic degradation of aromatic hydrocarbons by bacteria.Catechol dioxygenases one such enzyme plays a major role in the metabolism of aromatic hydrocarbons [1].Catechol dioxygenases can be divided into two major groups: those that cleave the aromatic ring between the vicinal diols (the intradiol dioxygenase group) and those that cleave the ring to one side of the vicinal diols (the extradiol dioxygenase group) [2].Catechol 2, 3-dioxygenase (C23O), also called meta-cleavage dioxygenase, which belongs to the extradiol dioxygenase group and catalyzes the meta-ring cleavage of (substituted) catechols [3].A wide variety of C23Os was reported from various bacteria in which most of them are able to catalyze catechol more efficiently than the substituted catechols, including methylcatechols, although some rare exceptions.For example, the C23O from Pseudomonas putida strain was capable to grow on toluidine and cleaves 3-methylcatechol at a higher rate than catechol [4].Catechol 2, 3-dioxygenase is present in several types of bacteria that degrade environmental pollutants, and it is one of the key biochemical "pathways" through which the contaminants are degraded.These enzymes metabolize the aromatic rings compounds such as toluene, xylene, naphthalene and biphenyl derivatives, which are the environmental pollutant of serious concern [5].
The Catechol-2, 3-dioxygenase enzyme has potential uses in bioremediation of Petroleum hydrocarbon compounds.However, the lack of characterization of the enzyme has limited its development as a practical biocatalyst.It is desirable to widen the substrate specificity of the enzyme to make it a more generic catalyst for breaking down these pollutants [6].The binding properties of Petroleum hydrocarbon compounds along with identification of active site residues could be used for site directed mutagenesis in order to identify their role in activity and substrate specificity.Consequently, these enzymes were the subject of different studies and have been suitable targets for site directed mutagenesis to improve their activities [7].Selecting the best mutant which had the best activity in oil degradation and use them in the future for cleaning oil spills in the environment [8].Present study utilizes protein-ligand docking as a tool to screen various petroleum hydrocarbons in order predicts the potential target for degradation of the substrates.

Figure 1:
The sequence alignments between Catechol 2, 3dioxygenase from Burkholderia cepacia and Pseudomonas Putida.The symbol asterisks represents single fully conserved residues, colons: fully conserved strong groups, periods: fully conserved and dot: weaker groups.The sequences were aligned with ClustalX 2.0.10.

Structure prediction of Catechol 2, 3-dioxygenase (C23O)
The three dimensional structure of Catechol 2, 3-dioxygenase was generated using the homology modeling method.Protein sequence of C23O (314aa, UniProtKB accession number: Q7M0R6) from Burkholderia cepacia was retrieved from Uniprot (http://www.uniprot.org/).A BLASTp search was made against the PDB database (PDB, http://www.rscb.org/pdb/) to spot the most appropriate template for homology modeling.Accordingly, the structure of Catechol 2, 3-dioxygenase from Pseudomonas Putida [9] was selected as best template based on their sequence percentage identity (41%) and query coverage (99%) (Figure 1).The alignment of both template and target was done using the align2D script in the Modeler9V8.Then the aligned sequence is used for model generation using the model_ligand.pyscript in the Modeler9V8

Molecular dynamics simulation
The molecular dynamics simulation of Catechol-2, 3dioxygenase was carried out using GROMACS 4.5.4[17, 18] molecular dynamics package.The atoms of the molecule was treated with OPLS-AA/L all-atom force field [19] and solvated with TIP4P water model, which extended to 1.50 nm triclinic box from the molecule to edge of the box.In order to neutralize the system 10 Na2+ were replaced with SOL molecule.The 50,000 steps of steepest descent and conjugate gradient energy minimization algorithm was used to minimize the protein with a tolerance of 1000 kJ mol -1 nm -1 .The PME [20] method was used for long-range interaction with a 1.2 nm cut-off and 0.16 nm Fourier spacing also used.Then, the energy minimized solvated system was considered a sensible one in terms of geometry and solvent orientations and used further.The molecular bond angles and geometry of the water were constrained with LINCS [21] and SETTLE [22] algorithm respectively.Consecutively to regulate the temperature (310 K) and pressure (1 atm) of the system, the V-rescale weak coupling and Parrinello-Rahman method [23] were used respectively.Both NVT (constant number of particles, Volume and Temperature) and NPT (constant number of particles, Pressure and Temperature) position restraint method was used for 100 ps.These results show that, the system is well equilibrated in terms of temperature, pressure, density and total energy.In continuation with that, this pre-equilibrated system was subjected to 10000 ps (10 ns) production MDS with a time-step of 2 fs.In every 2 ps the structural coordinates are saved and analyzed.The above-mentioned protocol was referred from our previous work [24].

Molecular Docking
The docking of substrate into the Catechol-2, 3-dioxygenase was carried using Autodock 4.0 software package [25] (Morris, et al. 1998).This software contains AutoTor, AutoGrid and AutoDock, in which AutoTor defines the rotatable bonds and root to detect the rotation.AutoGrid generates grid map for each substrate, here the grid box was generated based on the FE 2+ atom.The box size was set to 50 × 50 × 50 Å in the x, y, and z axes.The spacing between grid points was 0.375 Å. AutoDock menu helps for docking methodology.Lamarckian genetic algorithm (LGA) was used for the conformer searching in the grid box mentioned.This method is used to calculate the binding energy for each molecule as a scoring function.In docking protocol the number of individuals in each population is set to 150 and the maximum number of energy evaluations is set to 2,50,00,000, whereas the maximum number of generations was set to 1000.The maximum number of top individuals, the rate of gene mutation and crossover were set as default value.The docking protocols are used from the work of Mannu et al.

Discussion: Model validation
The modeled structure of C23O was subjected to validation analysis and compared with crystallographic C23O structure using WHATIF server.The modeled structure shows 89.3% favored region and 10.7% in the allowed region of Ramachandran plot and Z-score (-1.635) of Ramachandran plot also shows the quality of the model.In terms of the quality packing, the WHATIF results show that the backbone-backbone contacts (-0.48), backbone-sidechain contacts (-3.27), sidechainbackbone contacts (-2.16) and sidechain-sidechain contacts (-1.95).The overall quality value for all contacts is -2.32, which is a normal range value and the second generation packing score (-2.32) for all contacts also confirms the quality of the packing.In addition, the overall quality factor predicted by ERRAT shows that 81.311 % (Figure 2a) for Catechol-2, 3-dioxygenase from Burkholderia cepacia, which is somewhat lesser than the template (86.288) structure quality and maintained.Moreover, the Verify3D results (Figure 2b) illustrate that the compatibility score of target-template, which gave on clear cut values for the similar compatibility except in the region of 150-200 residues.All these results suggest that the modeled structure is valid one to use for further studies.

Molecular dynamics simulation
In order to evaluate the stability of the Catechol-2, 3dioxygenase from Burkholderia cepacia the molecular dynamics simulation of 10 ns was performed.The minimum potential energy (Figure 4a) of the protein -8.60324e +05 (7038 th ps) prove that the modeled structure was energetically stable during the production MD run.The RMSD profile shows that the first 2000 ps the protein backbone rise to equilibrate and remains quite stable until the 10 th ns of production MD with the maximum RMSD of 0.35 nm.This longest period of equilibration confirms the prolonged stability of the protein (Figure 4b).Moreover, the RMS fluctuation of Catechol-2, 3-dioxygenase was also calculated and described in (Figure 4c

Molecular Docking
The Fe atom containing Catechol-2, 3-dioxygenase from Burkholderia cepacia has the cabability of degrading the hydrocarbons like Benzene, O-xylene, Toluene, Fluorene, Naphthalene, Carbazol, Pyrene, Dibenzothiophene, Anthracene, Phenanthrene, Biphenyl and Chrysene (Figure 5) from crude oil.The present study described the catalytic mechanism of Catechol-2, 3-dioxygenase with its substrate through molecular docking analysis.Accordingly, the above mentioned substrates were prepared and docking with the active site of Catechol-2, 3-dioxygenase and analyzed.From the Table 1 (see supplementary material), each Hydrogen Bond (HB) and Van der Waals interaction made between substrate and enzyme, along with Free energy of binding, estimated inhibition constant (Ki) and substrate distances with Fe atom were shown.The docking results reveal that the benzene, Oxylene, Toluene, Fluorene, Naphthalene, Carbazol, Pyrene, Dibenzothiophene, Anthracene, Phenanthrene, Biphenyl makes HB and VdW interaction with the active site residues of H150, L152, W198, H206, H220, H252, I254, T255, Y261, E271, L276 and F309 and shown the free energy binding of -3.72, -4.00, -4.21, -4.57, -4.63, -4.29, -3.93, -3.50, -3.30, -2.14 and -1.20 respectively.Chrysene shows positive energy of binding in the active site atom of Fe.Except Pyrene, all above-mentioned substrates made close contact with Fe atom by the distance ranges from 1.67 to 2.43 Å.In addition, the above-mentioned substrate except pyrene, all other made π-π stacking interaction with H252 by the distance of ranges from 3.40 to 3.90 Å (Figure 6 & 7).The pyrene molecule is present 7.40 away from the Fe atom and makes T-shaped stacking interaction with F309 with a distance of 3.10 Å.Moreover, the catalytic Fe atom retains their coordinate interaction with H150, H220 and E271 by the distance of 2.45, 2.63 and 2.13 Å respectively.All these docking results reveal that, except Chrysene, all other substrate has good free energy of binding to hold enough in the active site and makes strong VdW interaction with Catechol-2,3-dioxygenase.These results suggest that, the enzyme is capable of catalyzing the above-mentioned substrate.

Conclusion:
In order to understand the mechanism of Catechol 2, 3dioxygenase in the degradation of environmental pollutants, the three-dimensional structure was modeled and analyzed their stability through 10 ns molecular dynamics simulation studies.The validation results ERRAT shows 81.311 % overall quality factor.The compatibility score predicted by Verify3D explains the compatibility nature of modeled enzyme is quite same with the template structure.These two successive validation results suggest that the modeled structure is reliable one for further studies.In addition to that the 10 ns MD simulation also carried out to check the stability of modeled enzyme and shows it has lowest potential energy to sustain its three-dimensional structure.The docking studies of various petroleum hydrocarbons with Catechol 2,3-dioxygenase illustrates that except Chrysene and Pyrene all other are able to form strong interaction with H150, L152, W198, H206, H220, H252, I254, T255, Y261, E271, L276 and F309 with free energy of binding ranges from -4.63 to -1.20.Chrysene shows positive energy of binding with the active site of Fe, whereas the pyrene shows negative energy of binding, though it's docked far from the active site atom of Fe due to its larger molecular nature.Except pyrene all the substrates has close contact with Fe atom with a distance ranges from 1.67 to 2.43 Å and also made a strong π-π stacking interaction with H252 by the distance ranges from 3.40 to 3.90 Å.In conclusion all the substrate except Chrysene is made good interaction with active site residues which in turn suggest that this enzyme is capable of metabolizing these aforementioned compounds.

BIOINFORMATION
[10, 11].The modeling protocol is referred from our previous work [12, 13].Here, the coordinate information of the template FE2+ was also used for the C23O model generation.From the ten models best model is chosen based on their Discrete Optimized Protein Energy (DOPE) and Modeler Objective Function (MOF) scores.The best model was validated using structural analysis and verification server (SAVES), which includes various tools like Ramachandran plot [14], ERRAT [15] and VERIFY_3D [16].

Figure 2 :
Figure 2: The validation results of Catechol 2, 3-dioxygenase using ERRAT and Verify3D a) The overall quality factor predicted by ERRAT shows that 81.311 %. b) The compatibility score for both template and target structure which shows that similar nature of the compatibility.

Figure 3 :
Figure 3: Three dimensional structure of Catechol 2,3dioxygenase with catalytic Fe atom.a) It contains 18 strands and 5 helices named from β1-β17 and α1-α5, in which the active site residues H150, H220 and E271 coordinating with Fe atoms are clearly pictured in the enlarged portion; b) The four repeated βαβββ motifs are shown in each N and C-terminal domain.The figure was prepared using protein molecular viewer PyMol and Chimera.

Figure 4 :
Figure 4: The molecular dynamics simulation results.a) Potential energy graph illustrates the energetically stable nature of Catechol 2,3-dioxygenase; b) RMSD profile explains the equilibration nature of the protein throughout the 10 ns MD simulation; c) RMS fluctuation describes about the three different flexible regions during the production MD run.
), showed it has three flexible regions.The first flexible region is a coil (124-140), which connects the N-terminal and C-terminal domain and fluctuate more with a maximum RMSD deviation of 0.45 nm.The second flexible region is the loop one (160-165), which connects the β9 and α3 with a maximum RMSD ranges from 0.25 to 0.35 nm.Whereas the third one is very important βhairpin motif which connects β10 and β11 (185-193), in combination with C-terminal tail region this motif covers the open side of the C-terminal domain bowl structure which offers a way for the substrate binding near to the Iron atom (Figure 4c).This flexibility furnishes open and close conformation of the C-terminal domain bowl with the RMSD deviation of 0.35 nm.

Figure 6 :Figure 7 :
Figure 6: The molecular interaction of substrate with Catechol-2, 3-dioxygenase.a) Benzene b) O-xylene c) Toluene d) Fluorene e) Napthalene and f) Carbazole interacting with Fe atom in the active site.Fe atom shown in Spear format and the interacting residues are shown in stick format colored by heteroatom, whereas the substrates are shown in yellow color.The hydrogen bond formed between coordinated active site

Table 1 :
Molecular docking of substrate in the active site of Catechol-2,3-dioxygenase.The list of substrate and their HB and VdW interaction residues, free energy of Binding, estimated inhibition constant and substrate distance with Fe atom are given