DOI:
10.1039/C3MB70300B
(Paper)
Mol. BioSyst., 2014,
10, 117-127
Hydration of porphyrin and Mg–porphyrin: ab initio quantum mechanical charge field molecular dynamics simulations
Received
22nd July 2013
, Accepted 14th October 2013
First published on 31st October 2013
Abstract
Ab initio QMCF-MD simulations were performed for porphyrin (POR) and magnesium–porphyrin (Mg–POR) immersed in water to study their structural and dynamical properties. The observed hydration behaviour of these solutes representing biomimetic models is in fair agreement with structural and dynamical features of their biological analogues, protoporphyrin IX (PPIX) and chlorophyll (CHl). Structural data obtained from the radial, angular and spatial distribution functions as well as the angular–radial distributions have a consensus on possessing a contrasting hydration behaviour of POR and Mg–POR. Flexibility of the ring in both solutes described by the improper torsional distribution and root mean square fluctuation showed an influence on H-bond interactions between the nitrogen atoms and water molecules that are also reflected in the respective dynamics. An axial water molecule coordinated to the Mg(II) ion indicates the penta-coordinated Mg–POR to be stable along the simulation. It was also shown that complexation of the Mg(II) ion to the porphyrin influences the hydration patterns significantly compared to the porphyrin itself, which is further supported by the vibrational power spectra evaluated for both solutes. Free energy of binding and solvent accessible surface area calculations also confirmed that these two solutes have distinct hydration behaviour. Detailed knowledge of the individual hydration patterns is expected to be of particular benefit.
1 Introduction
Porphyrins are biological molecules which are abundantly found in nature,1 for instance protoporphyrin IX (PPIX), one of the porphyrins that widely acts as a carrier molecule for divalent cations such as Fe(II) in heme,2 myoglobin3 and many other heme-containing enzymes like cytochrome c4 and catalase,5 as well as Mg(II) in chlorophylls,6etc. PPIX is generated within the mitochondria of cells and acts as a marker in photodynamic therapy against different forms of cancer.7 In particular, the medicinally relevant dysfunctions related to PPIX are known as porphyrias, which are a group of rare inherited or acquired disorders of certain enzymes participating in the production of porphyrins.8 For instance, the biosynthesis of all porphyrins is generally carried out by their first common precursor, δ-aminolevulinic acid catalysed by aminolevulinic acid synthase that is thought to be impaired in the porphyrial condition.9 In plants, the enzyme magnesium chelatase converts protoporphyrin IX to Mg–protoporphyrin IX during chlorophyll biosynthesis. This is the first unique step in the synthesis of bacteriochlorophyll,10 but the name chlorophyll (CHl) is generally thought to be associated with all green plants, cynobacteria and in the chloroplasts of algae. CHl is a chlorin pigment which has structural similarity to other porphyrin pigments such as heme and all porphyrins are produced via the same metabolic pathways. One of the CHl variants, CHl A, plays a major role in the light-harvesting complex (LHC) involved in the photosynthetic processes,11,12 which is influenced by the structure, dynamics and electronic properties of the CHl. The structure of PPIX is composed of tetrapyrroles interconnected via methine (–CH
) bridges forming a porphyrin ring and of functional groups commonly referred to as the tail attached to the pyrrole rings, which are considered to have hydrophobic and hydrophilic profiles, respectively.13 On the other hand, CHl A has a chlorin ring (a modified porphyrin) and a phytol tail with hydrophobic character, while magnesium(II) is coordinatively unsaturated enabling binding of nucleophilic polar molecules.14 The C
O group in the phytol tail saturates the Mg leading to dimerization and aggregation of CHl whereas in an aqueous environment large aggregates are expected due to intermolecular hydrogen bonding.11 Oligomerizations or aggregations of CHl A can be prevented by other polar solvents such as methanol.15 It was also found that in a water–acetonitrile mixture, the hydrophobic effects of the phytol tail restrict the CHl A molecules to forming smaller aggregates.16
As described briefly about the key biological roles played by PPIX and CHl A in the animal and plant kingdoms, the structural and dynamical properties of porphyrin rings of these molecules in aqueous solutions will help to understand the behaviour of these molecules in biological systems. They have been investigated for many decades and have attracted an immense number of experimental and theoretical approaches from physical, chemical and biological sciences, to establish the structure–function relationship of porphyrins and their derivatives. Experimental techniques applied to porphyrins experience severe limitations to provide insights into their structures and dynamics at an atomic level that can be complemented by the use of theoretical approaches. Molecular mechanical dynamics simulations have been extensively applied to characterise porphyrins, which served to unravel their veritable properties related to biological functions.17 Structure and dynamics of porphyrins strongly depend on their electronic properties, which need to be explored by engaging more sophisticated theoretical methods based on ab initio quantum mechanical theories18 to obtain reliable interpretations of porphyrin–water interactions and their respective dynamics. Since polarisation, charge transfer and many-body effects are known to be of particular importance in the case of metal-containing systems, the combination of force field parameters developed for the treatment of the isolated species in an aqueous environment can be expected to impose limits on the achievable accuracy. On the other hand such effects have not been focus of previous, purely classical simulation studies. Therefore, an ab initio quantum mechanical charge field molecular dynamics (QMCF-MD)19,20 study of dianionic porphyrin (POR) and magnesium–porphyrin (Mg–POR) in aqueous solution was performed which is expected to serve as a biomimetic model for PPIX and CHl molecules in aqueous environments whereas the computational effort required to treat the entire PPIX and CHl molecules in water15,21 at a quantum mechanical level is beyond the capabilities of present computational facilities. Fig. 1 depicts a sketch of the POR and Mg–POR molecules showing numbering and labelling of atoms employed throughout this manuscript as well as a snapshot illustrating the partitioning of the simulation box into the QM and MM regions.
![image file: c3mb70300b-f1.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f1.gif) |
| Fig. 1 Sketch showing numbering and labelling of atoms employed in this study for (a) POR and (b) Mg–POR. (c) Snapshot obtained from the POR-simulation illustrating the partitioning of the simulation box into the QM and MM regions. | |
2 Methods
QMCF-MD formalism
The ab initio quantum mechanical charge field molecular dynamics (QMCF-MD)19,20 formalism follows the same partition scheme as implemented for most of the quantum mechanical molecular mechanical (QM/MM) frameworks. The simulation of complex chemical systems requires a multi-level description which has been successfully fulfilled by the QMCF-MD framework. In the QMCF procedure, the system is typically partitioned into a region of special interest treated via a computationally demanding quantum mechanical (QM) approach, while the remaining environment is described by a simpler (and thus faster) molecular mechanical (MM) technique. The further division of the QM part into an inner core zone and an outer layer region makes the formalism more efficient in the case of investigations of solvated systems. If the size of the QM region is sufficiently large, the only required interaction between particles located in the inner core region and atoms located in the MM region are Coulombic contributions. The QMCF-MD code was interfaced to the parallel version of the program “TURBOMOLE”22 for the quantum chemical calculations. The Hartree–Fock (HF) method at the QM level of theory along with 6-31G** basis sets23–25 were employed for all atoms present in the QM region. Further information related to the QMCF-MD formalism has been reported here.19,20
Simulation protocols
Ab initio QMCF-MD simulations were performed for porphyrin (POR) and magnesium porphyrin (Mg–POR) inserted in cubic boxes (side lengths of 39.28 Å) with ∼2000 explicit water molecules using periodic boundary conditions. The porphyrin was modelled with the generalised Amber force field (GAFF),26 restrained electrostatic potential (RESP) charges27–29 were obtained by the ab initio calculation as described in the literature. The charges of the water molecules were adopted from the SPC/E water model (extended simple point charge model) which constitutes the MM region.30,31 The system was run in the NPT (isothermal–isobaric) ensemble using the Velocity–Verlet algorithm with a time step of 2 fs to integrate the Newtonian equations of motion. The Berendsen thermostat and manostat32 were employed to keep the temperature at 298.15 K using a relaxation time of 0.1 ps, the pressure was controlled using a reference of 1 atm and a relaxation time of 0.5 ps. For the treatment of long-range electrostatic interactions, the reaction field method was employed using a cutoff distance of 18.0 Å. The M-RATTLE algorithm,33 similar to the M-SHAKE34 technique was used to constrain bond lengths involving hydrogen atoms of the systems. In order to equilibrate the system, classical MD simulations were performed for 150 picoseconds followed by the ab initio simulations. In the case of POR, the geometric mean of the four nitrogen atoms served as the centroid for the QM zone, while the Mg marked the centre for the Mg–POR simulation. The radii were set to 0.5 and 8.5 Å for the core and layer zone, respectively. In the case of the POR simulation, the core zone was empty whereas the Mg atom is assigned to the core zone in the case of Mg–POR. On average the QM region contains in addition to the solute 57 and 62 water molecules in the case of POR and Mg–POR, respectively. A sketch illustrating the partitioning of the simulation box into the QM and MM regions for the POR system is shown in Fig. 1c. A smoothing function between 8.3 and 8.5 Å was set to ensure smooth transitions of particles between the QM and MM region.19,20 Finally, the equilibration of both systems was performed at 298.15 K prior to sampling of simulation trajectories for 10 ps, recording data every 2 fs.
Analysed properties
Structural and dynamical properties were analysed for POR and Mg–POR in aqueous solutions; the structural behaviour was evaluated via distributions of water oxygen (Ow) and hydrogen (Hw) around the POR and Mg–POR molecule. Radial distribution functions (RDF), g(r) were computed for oxygen and hydrogen of water molecules with reference to the Mg and nitrogen atoms as well as the centroid (CEN) in the case of POR.
Being a macrocycle, the porphyrin ring was found to possess conformational freedom in various environments including biological systems where porphyrins are involved in a number of processes linked to their conformational flexibility. Therefore, improper torsional distributions were evaluated for ring atoms arranged in a plane with the centroid and the Mg(II) ion, respectively. To further analyse the conformational flexibility of the porphyrin systems, averaged atomic root mean square fluctuations (RMSF) were calculated for each non-hydrogen atom of the porphyrin ring. A comparison of the RMSF plots enables a evaluation of the flexibility of the porphyrin ring in the Mg–POR and POR systems, respectively.
Angular distribution function (ADF) and coordination number distribution (CND) for hydrated Mg–POR enabled the investigation of the geometry and the number of coordinating atoms around the Mg(II) ion. The number of hydrogen atoms forming H-bonds with nitrogens of the porphyrin ring were also computed via CNDs.
The distribution of solvent molecules around both solutes was investigated via angular–radial distribution functions (ARD) combining the radial and angular information into a 3D density plot. The spatial distribution function (SDF) provided further information on the spatial arrangements of water molecules with reference to a predefined set of atoms arranged in a planar structure. The SDF spans both the radial and the angular spherical coordinates and characterises the local 3D packing of molecules.
Besides structural characteristics, the dynamics of water molecules in the vicinity of the solutes were evaluated. The properties of the water molecules involved in H-bonding to POR as well as in direct coordination to Mg(II) were quantified via ligand mean residence times (MRTs) using the direct method.35 The MRT, (τ) was calculated for t* values of 0.0 and 0.5 ps, represented as τ0.0 and τ0.5 respectively, where t* stands for the minimum time span that is required to register an exchange of a ligand between its hydration shell. A value of 0.5 ps is the time interval chosen as the experimental average lifetime of hydrogen bonds in pure water.36,37 The number of attempts, Rex, corresponds to the number of required ligand migrations that are on average required, to achieve one lasting exchange event (t > t*).
| ![image file: c3mb70300b-t1.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-t1.gif) | (1) |
where
Nex0.0 and
Nex0.5 are the number of exchanges registered for
t* 0.0 and 0.5 ps, respectively.
The underlying nature of the dynamics of POR and Mg–POR was also evaluated by calculating velocity autocorrelation function (VACF), C(t),38 expressed as
| ![image file: c3mb70300b-t2.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-t2.gif) | (2) |
where
N is the number of particles,
Nt corresponds to the number of time origins
ti, and
vj denotes a given velocity component of the particle
j. The power spectrum of the VACF was calculated by Fourier transformation,
39 using a correlation length of 2.0 ps with 2000 averaged time origins.
The potential of mean force (PMF),40W(r) of the interaction of water molecules with the centroid of POR and Mg of Mg–POR, as a function of radial separation, were calculated from the RDFs, using the equation;
| ![image file: c3mb70300b-t3.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-t3.gif) | (3) |
where
R is the universal gas constant and
T is the temperature.
The time-averaged solvent accessible surface area (SASA) was computed for both solutes using VMD41 to distinguish the hydrophilic and hydrophobic contribution of atoms in the porphyrin ring, enabling the determination of the relative magnitude of the H-bonding potential.
3 Results and discussion
Structure
Radial distribution function.
Strongly contrasting hydration properties of POR and Mg–POR was observed that can be attributed to the metallic centre occupied by Mg(II) in the latter case. The CEN–Ow and Mg–Ow RDF are illustrated in Fig. 2 and depict well-structured first peaks centred at ≃2.3 and ≃2.07 Å, respectively, indicating axially arranged water molecules at the centroid and the magnesium ion. The CEN–Hw RDF depicts a broad first peak between 1.09 Å and 2.88 Å corresponding to a mean CEN–Hw distance of 1.71 Å. Both CEN–Ow and CEN–Hw RDFs provide insight into the orientation of water molecules. This also signifies the contribution of the four nitrogen atoms, each attracting hydrogen atoms of water molecules due to opposite polarities. In the g(r) profile of POR, the non-zero first peak minimum reveals axial waters interacting with other water molecules above and below the plane of the POR ring that is visible via a broad second peak in CEN–Ow RDF centred at ≃4.8 Å. A very prominent peak from 6.40 Å to 9.68 Å centred at 8.36 Å is associated with a full hydration layer that covers the entire POR molecule including axial and nitrogen bound water molecules. Beyond the third peak, the distribution of water molecules is homogeneous. Integration of the CEN–Ow RDF up to the first and second peak minima indicates the presence of ∼2 and ∼8 water molecules, respectively which are equally arranged at each side of the plane. Fig. 2b displays the RDFs for atoms of Mg–POR with Ow and Hw: the Mg–Ow RDF exhibits a very sharp first peak at a distance of ≃2.07 Å which is in fair agreement with the crystal structure of ethyl-chlorophyllide A dihydrate,42 whereas the first peak of the Mg–Hw RDF was observed at a larger distance of 2.74 Å. In contrast to the case of POR, the first peak in the Mg–Ow RDF is well-separated from the second peak found between 3.37 Å and 5.12 Å, suggesting a strong coordination of water molecules to the Mg(II) ion. As mentioned earlier, Mg–POR is assumed to have one or more open coordination sites where ligands can be accommodated.14 This phenomenon, termed coordinative unsaturation, is satisfied by coordination of a water molecule in the case of hydrated Mg–POR.11 A broad second peak is located at 4.5 Å separated from the first peak indicating strong localisation of water molecules at either side of the plane at a larger distance. Mg–POR shows contrasting hydration properties compared to the POR molecule: Hw atoms around Mg–POR have a more diffuse localisation in contrast to the water oxygens as reflected by the non-zero first peak minimum in the Mg–Hw RDF. This also indicates interactions between water molecules arranged at different positions with respect to the Mg(II) atom. Integration of the first peak in the Mg–Ow RDF corresponds to one water molecule coordinated axially to the Mg(II) ion. This suggests that the Mg(II) ion has a predominantly penta-coordinated structure formed by the four nitrogen atoms of the porphyrin system and an axial water molecule, which correlates with other studies. An MP2- and DFT-based theoretical study reported that the binding of a second axial water ligand to Mg–POR is not exergonic, thus exclusively exhibiting five coordination of the Mg(II) atom.43 Our results are also in good agreement with crystallographic data44 revealing that Chl A tends to bind the fifth ligand in particular water which was shown to be an important stabilising factor in the special Chl pairs. A further study related to the diastereotopic ligation of magnesium atoms of chlorophylls in Photosystem I also reflects the penta-coordination of the Mg(II) atom.45
![image file: c3mb70300b-f2.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f2.gif) |
| Fig. 2 Radial distribution function of oxygen (Ow) (full lines) and hydrogen (Hw) (dashed lines) atoms of water, with respect to (a) CEN of the nitrogen atoms in the POR and with respect to (b) Mg in Mg–POR. | |
To investigate the arrangement and orientation of the water molecules toward POR, N–Ow and N–Hw RDFs were plotted to ascertain possible interactions between nitrogen atoms and water molecules via H-bonds as reflected by the N–Hw RDF peaks centred at ∼2.5 Å shown in Fig. 3a. Integration of each N–Hw RDF plot corresponds to approx. one Hw atom bound to each nitrogen of POR. The low occurrence of Hw at the nitrogen atoms could be inferred as a high degree of veritable flexibility and collective motion of the porphyrin ring.46,47 The conformational freedom includes twisting, stretching and even bending of the whole porphyrin ring as observed in the improper torsional distributions and enhanced orthogonal rotation with respect to the plane that was further confirmed from visual inspection. This also led to an assumption of dynamical H-bond formation between water molecules and nitrogens at one side of the plane which is the consequence of non-planarity of the POR molecule. In the case of Mg–POR, N–Ow and N–Hw RDFs were plotted to further investigate the hydration at each side of the ring plane (cf.Fig. 3b). Similar patterns were anticipated for the nitrogen–water g(r) profiles as in the POR–water system, but the N–Hw distribution functions reveal very low probability of H-bond formation at each nitrogen, thus indicating that H-bonding only occurs on the opposite of POR with respect to the Mg(II)-coordinating water molecule. In addition to a disruption of H-bonding resulting from rotation and/or bending of the porphyrin, a competing influence of the Mg(II)–water interaction can also be envisaged for the low occurrence of Hw atoms at the nitrogen atoms. A well-structured third peak between 6.39 Å to 9.00 Å separated by a small flat shoulder peak in the Mg–Ow RDF reflects a full hydration layer around the whole Mg–POR molecule, thus giving a clear indication of water oxygens oriented towards Mg(II).
![image file: c3mb70300b-f3.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f3.gif) |
| Fig. 3 RDFs of oxygen (Ow) (dashed lines) and hydrogen (Hw) (full lines) atoms of water, with respect to the nitrogen atoms of the porphyrin ring in (a) POR and (b) Mg–POR. | |
Angular distribution function and coordination number distributions.
The coordination geometry of hydrated Mg–POR was assessed via ADFs for Mg(II) coordinating with the nitrogen atoms of the porphyrin ring and water molecules (cf.Fig. 4a). The distribution function of the X–Mg(II)–X angle (X = N, Ow) exhibits two main peaks centred at ∼88.7° and ∼156.5°, which agree with the angular criteria for a square pyramidal arrangement. In the ADF profile, a side-peak results from the angular distortion which can be attributed to the conformational flexibility of the porphyrin ring (see Fig. 8–10). The CND for Mg shown in Fig. 4b further confirms the Mg(II) as being penta-coordinated in Mg–POR in excellent agreement with the interpretation of interactions between chlorophyll and polar solvents such as water and methanol.49 The probability of H-bond formation between nitrogen atoms and water molecules was also assessed by the CNDs for the N–Hw bond in hydrated POR and Mg–POR. From the CNDs shown in Fig. 5a, it is clear that on average at least one H-bond is formed with each nitrogen in hydrated POR, resulting in a total of 4 hydrogen bonds, whereas in the case of Mg–POR, the probability of finding H-bonds in the vicinity of nitrogen atoms was very low as revealed by the highest peaks at zero value (cf.Fig. 5b).
![image file: c3mb70300b-f4.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f4.gif) |
| Fig. 4 (a) Angular distribution function and (b) coordination number distribution for Mg coordinating with porphyrin nitrogens and the Ow. | |
![image file: c3mb70300b-f5.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f5.gif) |
| Fig. 5 Coordination number distribution of Hw involved in possible H-bonding with porphyrin nitrogens in (a) POR and (b) Mg–POR. | |
Spatial distribution functions.
SDF is a useful tool to depict the spatial arrangement of ligands, since both radial and angular contributions are visible. The probability of observing water oxygens in different regions of space around the centroid (CEN) and Mg(II) in the case of the POR and Mg–POR molecules is shown in Fig. 6. The resulting CEN–Ow SDF profiles demonstrate two axial water molecules above and below the centroid (cf.Fig. 6a; lateral view) that correlates well with the CEN–Ow RDF. The hydration layer represented as a wireframe on each side of the plane corresponds to those water ligands that are H-bonded to the nitrogen atoms of POR. Moreover, one can also observe that the entire POR molecule is surrounded by a spherical hydration layer as observed in the CEN–Ow RDF (see Fig. 2a). This cage results from the H-bond network formed between water molecules, causing them to arrange in a specific pattern due to residual hydrophobicity of the porphyrin nucleus.13 For Mg–POR a volume with high density near Mg(II) is visible in the SDF depicted in Fig. 6b, which represents the spatial density of an axial water molecule on one side of the ring (lateral view), thus highlighting the penta-coordination of the Mg(II) ion42,48 and the uneven hydration pattern on each side of the ring (cf.Fig. 6b). The presence of the axial water on one side of the ring perturbs the H-bonding between the nitrogen atoms and water molecules, whereas H-bond formation occurring on the opposite side of the plane is not influenced.
![image file: c3mb70300b-f6.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f6.gif) |
| Fig. 6 Isosurface of the spatial distribution function of Ow with respect to (a) CEN and (b) Mg of POR and Mg–POR molecule, respectively. | |
Angular–radial distribution.
Fig. 7a depicts the ARD plots of Ow surrounding POR and Mg–POR. A strong correlation with the structural data presented above can be observed. Contour plots of the hydrated solutes provide information on distinct features of the ligand density distribution of solvent molecules in different regions. Two intense peaks near the centre of the plot verify two water ligands arranged axially at the centroid of the POR molecule whereas a lower particle density parallel to these two intense peaks indicates a probability of water molecules on average showing a symmetric distribution in each hemisphere. The density distribution of water ligands at ∼8.40 Å corresponds to the hydration layer surrounding the entire solute including axial water molecules and the hydration layer formed at the nitrogen atoms above and below the plane of the porphyrin ring. Overall, the hydration pattern on both sides of the plane of the POR molecule is symmetrical with some fluctuations that are attributed to the inherent flexibility of the porphyrin system. In the case of Mg–POR, different features for the ARD functions were expected corresponding to the distinct structural features evaluated so far (see Fig. 7b and c). The preferential localisation of one water molecule just above the Mg(II) ion is evident from a very intense peak, thus further suggesting a very strong coordination between Mg and this water molecule. Since this high density suppresses the density distributions for all other surrounding water molecules, a separate ARD for water molecules excluding the axial one was generated, enabling the detailed investigation of the hydration pattern. The contour plots in Fig. 7c clearly indicate the formation of different layers of water molecules surrounding the Mg–POR system. The presence of the axial water only at one side of the ring plane influences the water density, being unevenly distributed above and below the plane. This is illustrated by the varying density distributions in the contour plots: the peaks in the lower half correspond to water molecules H-bonded to nitrogen atoms while peaks in the upper half represent water molecules interacting with the Mg(II) coordinated water molecule via H-bonds.
![image file: c3mb70300b-f7.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f7.gif) |
| Fig. 7 Angular radial distribution contour plots of Ow with respect to (a) CEN of POR, (b) Mg of Mg–POR showing the highest probability of axial water and (c) without axial water. | |
Non-planarity of porphyrins.
Fig. 8 illustrates improper torsional distributions of the porphyrin rings in POR and Mg–POR, respectively, that demonstrate the conformational flexibility of the macrocyclic ring in aqueous solutions. Hydration also seems to alter the ring flexibility leading to a lowered probability of H-bond formation (see Fig. 3a and b). The distribution functions for the improper torsions in the POR ring exhibit larger variations reflected by the peaks at 3.0 ± 40° and 180 ± 25°, whereas in the case of Mg–POR, the improper torsional distributions have narrow peaks at 6.0 ± 20° and 180 ± 20°. In both cases the peak maxima are displaced from zero, highlighting the bent structure of the porphyrin ring in solution (cf.Fig. 8a and c). In the case of Mg–POR, the coordination between Mg(II) and nitrogen atoms causes them to remain planar as observed in the improper torsional distribution shown in Fig. 8d. On the other hand, the water molecule strongly bound to the Mg(II) ion leads to a distortion of the porphyrin ring, resulting in displacement of Mg from the ring plane similar to that observed in the crystal structure of ethyl-chlorophyllide A dihydrate.42 A large variation of the improper torsional distribution of the nitrogen atoms in POR was observed, demonstrating a large degree of conformational flexibility. The relative degree of the flexibility of the porphyrin ring in POR and Mg–POR was further elaborated by computing RMSFs of all heavy atoms, which are depicted in Fig. 9. The RMSF plots show that atoms at the exterior of the pyrrole rings have larger fluctuations compared to those atoms interconnecting these rings. The nitrogen atoms of the porphyrin ring have less flexibility, as can be seen from the lowered RMSF values. Complexation of Mg(II) by the porphyrin ring resulted in a similar sequence in the RMSF profile, however all atomic fluctuations were found to be lower compared to un-complexed porphyrin which again indicates the lower flexibility of Mg–POR. This behaviour could be implicated in the structure–function relationship of metal-coordinated porphyrins.48 The observed non-planarity of porphyrins is associated with numerous modes of conformational flexibility that are represented as snapshots of the porphyrin ring taken from the simulation trajectory shown in Fig. 10.
![image file: c3mb70300b-f8.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f8.gif) |
| Fig. 8 Distribution functions of improper torsion formed between (a) ring atoms including CEN and (b) nitrogen atoms in POR, and between (c) ring atoms including Mg and (d) nitrogen atoms in Mg–POR (atoms forming improper torsion are highlighted; see inserted figures). | |
![image file: c3mb70300b-f9.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f9.gif) |
| Fig. 9 Atomic root mean square fluctuations of heavy atoms of the porphyrin ring in POR and Mg–POR (see Fig. 1 for atom numbering). | |
![image file: c3mb70300b-f10.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f10.gif) |
| Fig. 10 Different snapshots of the POR showing conformational flexibility possessed by porphyrins in aqueous solution. | |
Dynamics
Fluctuation of charges.
The fluctuations of the partial charges of atoms in the porphyrin ring for both POR and Mg–POR molecules were monitored throughout the simulation and are depicted in Fig. 11. The nitrogen atoms of POR carry varying negative charges with an average value of −0.63 ± 0.2 a.u. This value indicates that the nitrogen atoms serve as hydrogen-bond acceptor between POR and water molecules (cf.Fig. 11a). In comparison to the POR molecule, contrasting characteristics demonstrated by Mg–POR in water give a clear indication of the reduced polarity of HW towards the nitrogen atoms resulting from the presence of the Mg(II) ion. Fig. 11b depicts the variation of the atomic charge of Mg fluctuating about a mean value of 1.18 ± 0.1 a.u., which is significantly lower than its formal charge, thereby indicating a large degree of charge delocalization over the porphyrin ring. In particular, the charges of the nitrogen and the CA carbon atoms were influenced due to complexation of Mg(II) by the ring system as compared to the POR molecule having no metallic center. Table 1 shows the average atomic charges of different atom types of the porphyrin ring of POR and Mg–POR. In the case of Mg–POR, an average value of −0.87 ± 0.15 a.u. was observed for the nitrogen atoms whereas the average atomic charge of the carbon atoms (CA) increases from 0.23 to 0.32. The other carbon atoms (CB and CC) of both solutes have similar atomic charges indicating that the charges are not significantly perturbed by the presence of the Mg(II) ion. The combined effect of charge delocalisation and the presence of the Mg(II) ion in the porphyrin ring restricts Hw atoms to be oriented toward the nitrogen atoms on one side of the ring plane whereas the increased negative charges of the nitrogen atoms attract Hw on the opposite side of the plane, where no axial water molecules are bound to Mg(II).
![image file: c3mb70300b-f11.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f11.gif) |
| Fig. 11 Fluctuation of partial atomic charges of (a) nitrogens of POR, and of (b) Mg and (c) nitrogens of the Mg–POR molecule. | |
Table 1 Average partial atomic charges on different atom types of POR and Mg–POR
Types |
POR |
Mg–POR |
Mg |
— |
1.18 |
N |
−0.63 |
−0.87 |
CA |
0.23 |
0.32 |
CB |
−0.21 |
−0.21 |
CC |
−0.25 |
−0.25 |
Mean residence time and hydrogen bonding.
MRTs for the H-bonds between nitrogen atoms and water molecules determined for hydrated POR and Mg–POR are listed in Table 2 along with their relevant dynamical data. The average MRT values for the N–Hw bond in hydrated POR and Mg–POR were computed as 4.7 and 0.6 ps, respectively, using a cutoff distance of 2.5 Å, which generally serves as a basic geometric criterion for the hydrogen bonding. The MRT value of H-bonds in hydrated POR was much higher than that of pure water (1.3 ps),50 whereas the hydrogen bonds in Mg–POR seem to be perturbed by the presence of Mg which attracts Ow more strongly than Hw. Besides MRTs, other relevant dynamical data also demonstrate the contrasting character of H-bond formation in POR and Mg–POR. The residence time of nucleophilic ligands such as methanol coordinated to Mg in Chl A was estimated as spanning from picoseconds up to a few nanoseconds assuming exponential decay approximating the distribution for the residences times obtained from the force field based MD simulation.15 A similar approach was adopted for the estimation of residence times of the Mg(II) coordinating water molecule in the chlorin ring and, among three values, an intermediate time scale (hundreds of ps) of residence time was also reported that was attributed to the interaction of other water molecules hydrogen bonded to the Mg-coordinated molecule.15 On the other hand, dynamical data representing the Mg–Ow coordination obtained from the ab initio QMCF-MD exhibit no exchange, which demonstrates enhanced stability of the water molecule coordinated to Mg(II). Hydrogen-bond formation between the Mg(II) coordinated water molecule and other water molecules was further assessed by the Oaxial1–Owg(r) profile (data not shown here). This also correlates with a recent simulation study15 and the proposed crystal structure of ethyl-chlorophyllide A dihydrate.42 The presence of H-bonds in the crystal structure is chained through the C
O bond of the phytol tail that leads to aggregation of chlorophylls in aqueous medium.14,51
Table 2 Dynamics data for the water molecules making H-bonds with nitrogens of the porphyrin ring in POR and Mg–POR, and of pure water
H-bond (system) |
CN |
N
ex
0.5
|
N
ex
0.0
|
τ
0.5 (ps) |
R
ex
|
N–Hw (POR) |
∼1 |
∼1 |
264 |
4.7 |
202 |
N–Hw (Mg-POR) |
0.1 |
— |
39 |
0.6 |
7.7 |
H2O 50 |
4.1 |
20 |
131 |
1.3 |
6.6 |
Power spectra via Fourier transformed velocity autocorrelation functions.
Fig. 12 displays the power spectra representing all vibrational modes associated to POR and Mg–POR molecules immersed in aqueous solutions observed in the QMCF MD simulations (since-bonds containing hydrogen atoms were treated as rigid via the application of the M-Rattle algorithm,33 the respective modes cannot be determined in this simulation study). The highest peak observed for the POR molecule was at approx. 940 cm−1 along with a number of other peaks ranging between 36 cm−1 and 2800 cm−1. An interesting feature of great importance was the spectral shift observed for the Mg–POR molecule which correlates with the infrared spectra for the divalent metalloporphyrins demonstrating the frequency shifts based on the metal–ligand coordinate bond strength.52 The frequency shifts to higher values confirm the strong coordination of nitrogens with Mg and the most intense peak was observed at approx. 1115 cm−1 which is larger than that observed for the POR molecule, thus demonstrating a crucial property for a reliable interpretation of the flexibility and stability of porphyrins in the presence and absence of a metal such as Mg.
![image file: c3mb70300b-f12.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f12.gif) |
| Fig. 12 Power spectra in cm−1 of the POR and Mg–POR molecule. | |
Potential of mean force and free energy
Energetics of the interaction between the solute and water in different hydration shells can be quantified via the PMF W(r). Calculation of W(r) associated with the CEN–Owg(r) profile of a POR molecule yields the free energy of water molecules with respect to the centroid of the ring (cf.Fig. 13a). The calculated free energy was obtained as approximately −0.76 kJ mol−1, which indicates a less attractive character of Ow from 1.82 to 3.28 Å (see Fig. 2a), thus revealing preferential orientation of Hw involved in H-bond interactions with ring nitrogens. On the other hand, the W(r) profile of the interaction between Mg and water molecules in hydrated Mg–POR illustrated in Fig. 13b yielded a free energy of binding of ∼−2.27 kJ mol−1. The free energy of solvation of the Mg atom from force field based MD simulation of Chl A in water was calculated as −1.92 kJ mol−1,15 which was considered reasonable on the basis of an estimation from QM calculations that the binding energy decreases with increasing permittivity of the medium.53 The free energy of binding based on the MM-MD study is similar to the value obtained from this QMCF-MD study.
![image file: c3mb70300b-f13.tif](/image/article/2014/MB/c3mb70300b/c3mb70300b-f13.gif) |
| Fig. 13 Potential of mean force of the solvent molecules derived from the CEN–Ow and Mg–Ow RDFs in (a) POR and (b) Mg–POR, respectively. | |
Solvent accessible surface area
Porphyrins are known to have residual hydrophobicity and hydrophilicity due to varying polarities attained by the atoms in the ring system. Nitrogens of porphyrins have the potential to form several hydrogen bonds with water molecules, acting as H-bond acceptors, whereas carbons are unable to act as either hydrogen donor or acceptor. Structural and dynamical data for hydrated POR and Mg–POR molecules showed that H-bonds are perturbed by the complexation of porphyrins with metal, as also evaluated by the SASA (cf.Table 3). The POR molecule is shown to have a higher hydrophilicity than Mg–POR, showing a greater capacity to establish H-bonds with water molecules, which is decreased in the case of Mg–POR due to complexation.
Table 3 Time-averaged solvent accessible surface area in nm2. Hydrophobic and hydrophilic surfaces are decomposed to show their individual contributions
|
POR |
Mg–POR |
Hydrophobic |
4.58 |
4.62 |
Hydrophilic |
0.29 |
0.21 |
Total |
4.87 |
4.83 |
4 Conclusion
Structural and dynamical properties of the hydrated apo-porphyrin and the Mg(II)–porphyrin complex obtained via QMCF-MD simulations demonstrate dramatically contrasting hydration behaviour. Complexation of porphyrin to Mg(II) results in Mg–POR, which is characterised by a notably different coordination of water molecules compared to the uncomplexed POR, leading to different H-bonding patterns and corresponding dynamics. The estimation of the free energy resulting from the potential of mean force revealed preferential orientation of water hydrogens involved in H-bond interactions with ring nitrogens in the POR molecule and a strong coordination between the Mg(II) ion and water in the case of Mg–POR.
The QMCF-MD simulations of POR and Mg–POR provide more insights for a reliable interpretation of the hydration properties of porphyrins and the respective change upon complexation with a metal ion, enabling the characterisation of the contrasting hydration behaviours of the apo-form. In particular since in addition to the QM treatment of the solute no empirical parameters for the description of the solute–solvent interaction are required in the QMCF procedure, all related intra- and intermolecular forces are obtained on the first principle basis. Other features of the QMCF-MD approach such as the accurate description of charge transfer leading to varying partial charges of the involved atoms underline the capabilities of quantum chemical simulation methods in complex chemistry and related disciplines. The findings from the QMCF MD simulations of both solutes nicely mimic the hydration properties of protoporphyrin IX and chlorophyll molecules to a large extent. For this reason, all results obtained from the simulations will provide a useful reference for studies of porphyrin in large systems via theoretical as well as through experimental methods.
As a consequence of the successful application of the QMCF-MD method to the simple porphyrins presented in this study, investigations of the theoretically more challenging species containing other ions such as zinc(II), copper(II) and iron(II/III) are envisaged for the future.
Acknowledgements
Financial support for this work by an Ernst-Mach-Stipendium from BMWF ASEA-UNINET Austria for Syed Tarique Moin is gratefully acknowledged. This work was supported by the Austrian Ministry of Science BMWF UniInfrastrukturprogramm as part of the Research Focal Point Scientific Computing at the University of Innsbruck.
References
- A. R. Battersby, C. J. Fookes, G. W. Matcham and E. McDonald, Nature, 1980, 285, 17–21 CrossRef CAS.
- R. Hardison, Am. Sci., 1999, 87, 126–137 Search PubMed.
- T. Takano, J. Mol. Biol., 1977, 110(3), 569–584 CrossRef CAS.
- J. P. Collman, R. Boulatov, C. J. Sunderland and L. Fu, Chem. Rev., 2004, 104, 561–588 CrossRef CAS PubMed.
- P. Chelikani, I. Fita and P. C. Loewen, Cell. Mol. Life Sci., 2004, 61, 192–208 CrossRef CAS PubMed.
- J. J. Katz, W. Oettmeier, J. R. Norris and S. I. Beale, Philos. Trans. R. Soc. London, 1976, 273, 227–253 CrossRef CAS.
- R. R. Allison, G. H. Downie, R. Cuenca, X. Hu, C. J. H. Childs and C. H. Sibata, Photodiagn. Photodyn., 2004, 1, 27–42 CrossRef CAS.
-
C. M. Lourenço, C. Lee and K. E. Anderson, Disorders of Haem Biosynthesis: Inborn Metabolic Diseases, 2012, pp. 519–539 Search PubMed.
- M. Doss, F. Sixel-Dietrich and F. Verspohl, Clin. Chem. Lab. Med., 1985, 23, 505–514 CAS.
- D. W. Bollivar, J. Y. Suzuki and C. E. Bauer, Annu. Rev. Genet., 1997, 31, 61–89 CrossRef PubMed.
-
J. J. Katz, Naturwissenschaften, 1973, pp. 32–39 Search PubMed.
- A. Agostiano, P. Cosma, M. Trotta, L. Monsù-Scolaro and N. Micali, J. Phys. Chem. B, 2002, 106, 12820–12829 CrossRef CAS.
- A. Wiehe, E. J. Simonenko, M. O. Senge and B. Röder, J. Porphyrins Phthalocyanines, 2001, 5, 758–761 CrossRef CAS.
- J. J. Katz, J. R. Norris, L. L. Shipman, M. C. Thurnauer and M. R. Wasielewski, Annu. Rev. Biomed. Eng., 1978, 7, 393–434 CAS.
- K. Karki and D. Roccatano, J. Chem. Theor. Comput., 2011, 7, 1131–1140 CrossRef CAS.
- A. Agostiano, M. D. Monica, G. Palazzo and M. Trotta, Biophys. Chem., 1993, 47, 193–202 CrossRef CAS.
- H. M. Marques and K. L. Brown, Coord. Chem. Rev., 2002, 225, 123–158 CrossRef CAS.
-
A. Szabo and N. S. Ostlund, Modern quantum chemistry: introduction to advanced electronic structure theory, Courier Dover Publications, 1989 Search PubMed.
- T. S. Hofer, A. B. Pribil, B. Randolf and B. Rode, Adv. Quantum Chem., 2010, 213–246 CrossRef CAS.
- B. M. Rode, T. S. Hofer, B. R. Randolf, C. F. Schwenk, D. Xenides and V. Vchirawongkwin, Theor. Chem. Acc., 2006, 115(2), 77–85 CrossRef CAS PubMed.
- T. R. Guizado, S. R. W. Louro and C. Anteneodo, J. Chem. Phys., 2011, 134, 055103 CrossRef PubMed.
- R. Ahlrichs, M. Bär, M. Häser, H. Horn and C. Kölmel, Chem. Phys. Lett., 1989, 162, 165–169 CrossRef CAS.
- W. J. Hehre, R. Ditchfield and J. A. Pople, J. Chem. Phys., 1972, 56, 2257–2261 CrossRef CAS.
- M. M. Francl, W. J. Pietro, W. J. Hehre, J. S. Binkley, M. S. Gordon, D. J. Defrees and J. A. Pople, J. Chem. Phys., 1982, 77, 3654–3665 CrossRef CAS.
- P. C. Hariharan and J. A. Pople, Theor. Chem. Acc., 1973, 3, 213–222 Search PubMed.
- J. Wang, R. M. Wolf, J. W. Caldwell, P. A. Kollman and D. A. Case, J. Comput. Chem., 2004, 25, 1157–1174 CrossRef CAS PubMed.
- C. I. Bayly, P. Cieplak, W. Cornell and P. A. Kollman, J. Phys. Chem., 1993, 97, 10269–10280 CrossRef CAS.
- W. D. Cornell, P. Cieplak, C. I. Bayly and P. A. Kollmann, J. Am. Chem. Soc., 1993, 115, 9620–9631 CrossRef CAS.
- J. Wang, P. Cieplak and P. A. Kollman, J. Comput. Chem., 2000, 21(12), 1049–1074 CrossRef CAS.
- H. J. C. Berendsen, J. R. Grigera and T. P. Straatsma, J. Phys. Chem., 1987, 91, 6269–6271 CrossRef CAS.
- S. Chatterjee, P. G. Debenedetti, F. H. Stillinger and R. M. Lynden-Bell, J. Chem. Phys., 2008, 128, 124511 CrossRef PubMed.
- H. J. C. Berendsen, J. P. M. Postma, W. F. V. Gunsteren, A. DiNola and J. R. Haak, J. Chem. Phys, 1984, 81(8), 3684–3690 CrossRef CAS.
- H. W. Horn, W. C. Swope, J. W. Pitera, J. D. Madura, T. J. Dick, G. L. Hura and T. Head-Gordon, J. Chem. Phys., 2004, 120, 9665–9678 CrossRef CAS PubMed.
- V. Kräutler, W. F. van Gunsteren and P. H. Hünenberger, J. Comput. Chem., 2001, 22, 501–508 CrossRef.
- T. S. Hofer, H. T. Tran, C. F. Schwenk and B. M. Rode, J. Comput. Chem., 2004, 25(2), 211 CrossRef CAS PubMed.
- A. J. Lock, S. Woutersen and H. J. Bakker, J. Phys. Chem. A, 2001, 105(8), 1238 CrossRef CAS.
-
A. J. Lock, S. Woutersen and H. J. Bakker, Femtochemistry and Femtobiology, Word Scientific, Singapore, 2001 Search PubMed.
- U. Balucani, J. P. Brodholt and R. Vallauri, J. Phys.: Condens. Matter, 1996, 8, 6139–6144 CrossRef CAS.
- R. N. Bracewell, Sci. Am., 1989, 260, 86–89 CrossRef CAS PubMed.
-
D. Chandler, Introduction to modern statistical mechanics, Oxford University Press, 1987 Search PubMed.
- W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33–38 CrossRef CAS.
- H. Chow, R. Serlin and C. E. Strouse, J. Am. Chem. Soc., 1975, 97, 7230–7237 CrossRef CAS.
- A. B. Fredj, Z. B. Lakhdar and M. F. Ruiz-López, Chem. Phys. Lett., 2009, 472, 243–247 CrossRef PubMed.
- T. Oba and H. Tamiaki, Photosynth. Res., 2002, 74, 1–10 CrossRef CAS.
- T. S. Balaban, P. Fromme, A. R. Holzwarth, N. Krauß and V. I. Prokhorenko, Biochim. Biophys. Acta, 2002, 1556, 197–207 CrossRef CAS.
-
W. R. Scheidt and Y. J. Lee, Metal Complexes with Tetrapyrrole Ligands I, 1987, pp. 1–70 Search PubMed.
- M. O. Senge, Chem. Commun., 2006, 243–256 RSC.
- J. A. Shelnutt, X. Song, J. Ma, S. Jia, W. Jentzen and C. J. Medforth, Chem. Soc. Rev, 1998, 27, 31–42 RSC.
- J. J. Katz, H. H. Strain, D. L. Leussing and R. C. Dougherty, J. Am. Chem. Soc., 1968, 90, 784–791 CrossRef.
- T. S. Hofer, B. M. Rode, A. B. Pribil and B. R. Randolf, Adv. Inorg. Chem., 2010, 62, 143–175 CrossRef CAS.
- K. Ballschmiter and J. J. Katz, J. Am. Chem. Soc., 1969, 91, 2661–2677 CrossRef CAS.
- L. J. Boucher and J. J. Katz, J. Am. Chem. Soc., 1967, 89, 1340–1345 CrossRef CAS.
- A. B. Fredj and M. F. Ruiz-Lopez, J. Phys. Chem. B, 2009, 114, 681–687 CrossRef PubMed.
|
This journal is © The Royal Society of Chemistry 2014 |