An Investigation of the BSSE Effect on the Evaluation of Ab Initio Interaction Energies for Cisplatin-Water Complexes

No presente artigo o efeito do pseudopotencial de caroço (ECP) para o átomo de platina no erro de sobreposição de base (BSSE) foi avaliado para as energias de interação entre o complexo cisplatina ([Pt(NH 3 ) 2 Cl 2 ]) e água. Oito orientações espaciais distintas para os monômeros cisplatina-água foram consideradas na determinação do BSSE utilizando o procedimento padrão “counterpoise”. O BSSE foi também avaliado ao longo das curvas de energia potencial, as quais são importantes para a determinação dos parâmetros do potencial de Lennard-Jones. Os cálculos foram feitos em diferentes níveis de teoria, incluindo HF, DFT (funcionais: BLYP, B3LYP, BP86, PBE, PW91) e MP2 com conjunto de funções de base do tipo “split-valence” de Pople e “correlated consistent” de Dunning.


Introduction
6][7][8][9] In a recent study 10 we have emphasized the importance of the intermolecular interactions between cisplatin and water molecule as a pre-step to the hydrolysis reaction which plays an important role in the drug mechanism of action.In that work, the interest was devoted to the trend of theoretical methods used, with the focus on the electronic correlation and electrostatic effects for the interaction of cisplatin with water and their aquated species using ab initio correlated level of theory. 10In the present study the intermolecular interaction of cisplatin with water is investigated using the supermolecule approach and the counterpoise procedure 11 for the evaluation of the basis set superposition error (BSSE) at various levels of calculation, aiming to address its effect on the interaction energy for a molecular complex containing transition metal atom, where the use of effective core potential (ECP) is required.
Generally, finite basis sets are used in quantum mechanical calculations of interaction energy, and the proximity between the monomers involved in the complex formation leads to the so-called BSSE. 12This error is associated with the incompleteness of the monomer's basis set leading to an artificial overestimation of the complex stabilization.One of the most accepted approaches to evaluate the BSSE is the counterpoise (CP) method proposed by Boys and Bernardi. 11In the CP theory the energy of each monomer is calculated with the basis functions of the other subunit (without the nuclei or electrons) using the so-called "ghost-orbital" at the complex equilibrium geometry.A review on the state of the art in counterpoise theory can be found in ref. 13 The counterpoise method was further applied to large interacting systems 14 and complexes greater than dimers. 15In general, the BSSE is more substantial when a small basis set is used and it decreases with the basis set improvement. 14More recently, several authors have been paying attention to the occurrence of intramolecular BSSE involving intramolecular noncovalent interactions. 15,16n our recent work, 10 irregular results concerning the CP BSSE correction was reported for equilibrium structures located on potential energy curves (PECs) for platinum compounds interacting with water, that was tentatively attributed to the ECP used for the platinum atom.Therefore, we found that a more detailed analysis of the BSSE correction for the cisplatin-water complex could be relevant to assess the role played by the ECP and also the BSSE behavior as a function of the level of theory used (HF, DFT and MP2).In this work we report CP BSSE calculations using the standard LANL2DZ pseudopotential and also a modified ECP aiming to improve the platinum atom basis set.These results can indicate the importance of the BSSE correction in calculations of interaction energies of large molecular complexes containing transition metal atoms of interest to biological and catalytic processes.This is an important point since the study of transition metal compounds in solution using classical Molecular Dynamics and/or Monte Carlo methods usually requires the parameterization of pairwise intermolecular parameters.The lack of structural and thermodynamics experimental data to be used in the parameterization procedure forces the use of ab initio pair interaction energies and geometrical parameters.Therefore, in this work we evaluated the effect of the BSSE correction on the calculation of PECs necessary for the adjustment of parameters to be used in Lennard-Jones (12-6) potential which is of great importance to molecular simulation studies.

Methodology
The geometry of the cisplatin and water monomers were optimized, without geometrical constrains, always yielding the expected C 2v symmetries, and interaction energies calculated considering eight distinct monomer orientations used for the evaluation of the PECs which are shown in Figure 1 (equilibrium structures).It is important to state that the cisplatin and water geometries were kept rigid during the PEC calculation; only the intermolecular distance, dPt•••O w is varied, what makes the computational procedure of practical use.Various levels of theory were used: Hartree-Fock (HF), Density Functional Theory (DFT) with the PBE1PBE, 17 PW91, 18,19 BLYP, [20][21][22] B3LYP [21][22][23] and BP86 22,24 functional and Møller-Plesset second-order perturbation theory (MP2) 25 , employing the 6-31G(d,p), 6-311++G(2df,p) and 6-311++G(2df,2p) 26 basis sets for the water and ligand atoms.The 60 inner shell electrons of platinum were replaced by the LANL2DZ effective core potential of Hay and Wadt. 27The 18 valence electrons were explicitly included in the calculation and were described by two different approaches: (i) using the associated LANL2DZ double-ζ basis set, (ii) LANL2DZ double-ζ basis set augmented by a set of f polarization function (exponent: α = 0.9930). 6An alternative ECP proposed by Burda, 28 named here ECP*, was also used.The improved ECP* core potential includes three sets of f and two sets of g polarization functions on the Pt valence shell.With the changes made by Burda, 28 the ECP* provides 29 basis set functions for the valence electrons.The core potential still evaluates the 60 platinum core electrons.More details of the procedure can be found in reference 28.
The quantum mechanical interaction energy (∆E) was evaluated by means of equation ( 1), (1)   where E cpx (r) stands for the electronic plus nuclear repulsion energy of the cisplatin(A)-water(B) complex, which is a function of the distance between the interacting monomers, and E A and E B correspond to the total energies of the free monomers A and B, calculated at the specified level of theory.The BSSE was accounted for using the standard CP approach.The BSSE corrected interaction energy (∆E BSSEc ) is evaluated by the following equation: The terms E A (B) -cpx and E B (A) -cpx are evaluated through a calculation with ghost orbital on monomers A and B indicated in the superscript parenthesis, with zero nuclear charge and no electrons, but containing the basis functions of the respective monomer (A or B), with both monomers being placed at the corresponding equilibrium geometry of the complex (i.e., using the basis set of the complex and deformed monomer geometries (E A-cpx , E B-cpx ) upon complex formation).The BSSE value can be obtained by the difference between the BSSE corrected (∆E BSSEc ) and uncorrected (∆E) energies: (3) Using the equation ( 3) with ( 1) and (2), one can write an equation in terms of the individual contributions to BSSE according to equation ( 4): (4)   As already mentioned, no geometry deformation due to complex formation was considered for the BSSE calculation, so E A-cpx = E A and E B-cpx = E B .The cisplatinwater complex potential energy curves were calculated for various monomer orientations, depicted in the Figure 1.The quantum mechanical calculations were carried out using the Gaussian 2003 suite program, 29 with the default "frozen core" mode being used in MP2 calculations (inner-shells are excluded from the correlation calculation).

Results and Discussion
The optimized intramolecular geometrical parameters for cisplatin and water, calculated with the 6-31G(d,p)/ LANL2DZ and 6-311++G(2df,p)/LANL2DZ(f) Pople's basis set, show a good agreement for all levels of calculation.The use of an improved effective potential for the platinum atom (ECP*) 28 combined with the 6-311++G(2df,p) and 6-311++G(2df,2p) basis sets does not lead to significant changes in the MP2 optimized structures.The HF and DFT results for the complex equilibrium distances, interaction energies and BSSE corrections are given in Table 1, and the MP2 results given in Table 2.The other DFT and MP2 results with the smaller basis set are given as supplementary material available upon request.
The overall trend in the DFT energies is similar to the HF results, with a small variation in the energy values among the exchange-correlation functional used.For the relevant low energy structures the stability order is IV > VI > II > V at the HF level and II > IV > VI > V for all DFT methods.The DFT BSSE corrections for these four spatial arrangements are, on average, 16, 13, 54 and 12%, respectively, being smaller at the PW91 level (15, 9, 28 and 10%).The average DFT result is very close to the HF prediction (14, 11, 45 and 10%).The relative contributions to BSSE from cisplatin and water found at the DFT level are also similar to the HF values with the same basis set.It can be seen from Table 1 (structures II, IV, V and VI) that that the BSSE contribution is mainly due to the cisplatin basis set (> 50%) for the arrangements where water O-H bond is pointed towards the cisplatin plane (II and VI).For the other structures (IV and V), more than 50% of the BSSE error comes from the water basis set.The difference between HF and DFT results is more pronounced for those arrangements where water interacts with cisplatin directly through the oxygen atom (IV and V).
A similar analysis of the BSSE contribution was done at the MP2 level using the LANL2DZ(f) and an improved ECP* for the Pt atom.The calculated values for equilibrium distances and energies are shown in Table 2.By analyzing the data reported in Table 2 it can be noted that the equilibrium distances do not change significantly when the ECP* is used.However, the potential energy minimum (∆E) is deeper with the standard LANL2DZ(f) basis set.As found in DFT calculations, structure II was predicted as the lowest energy complex at the MP2 level with both basis sets, followed by IV, VI and V.For these equilibrium structures the BSSE were considerable larger than the values calculated at the HF and DFT levels, being equal to 52, 30, 65 and 27% of the interaction energy for structures II, IV, V and VI, respectively, when the LANL2DZ(f) core potential is used for Pt (see Table 2).The improvement of basis set to ECP* changes these BSSE values to 27, 22, 54 and 20%.It is noticeable the decreasing of BSSE when the basis set for Pt is enhanced, especially for structure II where it drops from 52 to 27%.The behavior of individual contributions is also interesting at the MP2 level.It can be seen that a better description of the heavy atom, seems to equilibrate the BSSE contributions from the cisplatin and water.This is observed even for structures that shows a major contribution of water molecule as Ddp-IV which presents 31% for cisplatin contribution (69% for water molecule).
The use of improved ECP* reduces substantially the BSSE error, but at an increased computational cost.However, as can be seen from Table 2, the BSSE corrected energies (∆E BSSEc ) using the modified ECP* do not differ significantly from the corresponding values obtained using standard ECP (LANL2DZ) with the same triple-zeta basis set.We have also found that the same energetic order is observed for the eight complexes investigated, with or without BSSE correction, with a discrepancy of less than 2 kcal mol -1 between the MP2/6-311++G(2df,2p)/ECP* and PBE/6-311++G(2df,p)/LANL2DZ(f) results, with similar behavior being also observed for the Pt•••O equilibrium intermolecular distance calculated at the PBE and MP2 levels (supplementary material available upon request).
][32] The results are reported in Figure 2 at the DFT (PBE functional) and MP2 levels for structures I, II, IV, V, VI and VIII.It can be promptly seen that distinct complex spatial orientations have a dissimilar sensitivity to the basis set size, the behavior being more pronounced at the MP2 level, both for ∆E and BSSE correction (Figures 2a and 2c respectively).It can be seen that at the MP2 level the BSSE increases substantially with the basis set size for structure II, an odd behavior dictated by the use of the ECP, as already pointed out.A regular behavior of the PBE calculated BSSE value with the improvement of the Dunning's basis set is quite evident, including structure II.Nevertheless, the BSSE corrected interaction energies (∆E BSSEc ) shows an almost linear profile both at the DFT and MP2 levels (Figure 2b).From the profile shown in Figure 2b it can be said that the interaction energy can be very satisfactorily described using a small basis set (6-31G(d,p) or cc-pVDZ) including the BSSE correction (∆E BSSEc ).This result is particularly useful for calculation of interaction energy involving big molecular systems where the use of large basis sets is precluded.For reason of comparison only results for water dimer are given in Figure 3, where it can be seen that when there is no ECP a smooth expected behavior with the basis set improvement is observed.The non-corrected (∆E) and BSSE corrected (∆E BSSEc ) interaction energies are virtually the same when the cc-pV5Z basis set is used.These water dimer results are also in agreement with experimental 33 and theoretical 34,35 data available.The results reported in Figure 2 also provide support for the use of DFT based methods, as a lower cost computational approach, for further studies involving the interaction of platinum compounds (or other transition metal species) with biological targets, i.e., intermolecular interactions of biological relevance involving large size complexes, and also molecular systems relevant to catalytic processes (organometallic chemistry) where large size ligands are usually present.Our results indicate that the PBE functional may be recommended as a computational affordable quantum mechanical method to treat large interacting systems, at a tolerable uncertainty.Nevertheless, we still need to be careful about the quality of the DFT interaction energies, in what electron correlation effect is concerned, regardless the BSSE correction, particularly when dispersion forces are predominant.
A last point that deserves our attention concerns the calculation of the entire potential energy curve, far from equilibrium distances, as for example in the determination of empirical force field parameters required in classical computational simulation studies, where BSSE plays a role on the calculation of the energy curve in the short and long distance region.A good example is the use of BSSE corrected potential energy curves for the generation of empirical parameters that feed the Lennard-Jones (12-6) potential function used in classical Monte Carlo simulation of liquids and also solvent effects.We have reported ε and σ Lennard-Jones parameters obtained for cisplatin using BSSE corrected 9 and uncorrected 10 potential energy curves in the fitting procedure.The results are presented in Table 3 (the TIP3P potential was used for the water molecule), 36 where it can be promptly seen that the difference is quite substantial.
Figure 4 shows the Lennard-Jones (12-6) plus Coulomb parameterized potential energy curve (solid line) with and without BSSE correction and the corresponding MP2/6-31G(d,p)/LANL2DZ ab initio curves for the cisplatin-water complex (DDP-I), where a relative better agreement for the BSSE uncorrected curves can be inferred.In Figure 5 the quantum mechanical (LJ') and true Lennard-Jones (12-6) potential energy curves (also including BSSE correction) for the monoaquo (Madp-I) and diaquo (Dap-I) structures 10 are presented, showing the transferability of the LJ parameters obtained for the cisplatin species.Here the LJ is the true Lennard-Jones classical curve and LJ' obtained as a difference between the calculated ab initio curve (∆E) and the electrostatic contribution (U EL ) evaluated through the classical Coulomb potential using ab initio CHELPG atomic charges. 10It can be promptly seen that the BSSE uncorrected curves exhibit a better accordance along the distance variation.What remains to be seen is the effect of using the BSSE corrected PECs for the generation of empirical parameters on the results of a Monte Carlo simulation in order to validate or not the use of BSSE corrected Lennard-Jones parameters.
Preliminary results 37 shows that the Radial Distribution Function (RDF) for cis-DDP and water, obtained using BSSE corrected LJ parameters, yield Pt•••H and Pt•••O distances larger than the corresponding BSSE uncorrected values.A direct comparison between calculated bulk properties, using computer simulation methods, and experimentally measured values for the cisplatin-water system would provide a definite answer to the relevance of BSSE correction for the determination of Lennard-Jones parameters to be used in classical Monte Carlo (MC) and Molecular Dynamics (MD) simulations.Finally, it should be mentioned that so far we have investigated cisplatin-water complexes for eight chosen spatial orientations, keeping the equilibrium monomer geometries frozen at their isolated fully optimized values (and also the intermolecular angles unchanged), which were relevant for the construction of potential energy curves aiming at the determination of Lennard-Jones parameters for further studies of cisplatin compounds in aqueous solution.Therefore, it is also important to analyze the true minima on the potential energy surface for the cisplatin-water complex, through a complete geometry optimization of the intra and intermolecular geometrical parameters for the most relevant complexes, and so the effect of monomer geometry relaxation on the BSSE correction.In order to assess the effect of the monomer geometry deformation (∆E Mon-Def ) on complex formation the geometry of the six bound cisplatin-water complex structures (I, II, IV, V, VI and VIII) were fully optimized without any symmetry or geometrical constraint at the MP2 and PBE levels using the 6-311++G(2df,p)/LANL2DZ(f) basis set.Three equilibrium structures survived after full geometry optimization, with only structure VI resembling the original spatial orientation on the frozen-monomer potential energy curves (structure I and IV optimize to II, and VIII optimizes to V).The fully optimized structures are depicted in Figure 6, with MP2 optimized intermolecular geometrical parameters and interaction energies given in Table 4.The monomer deformation energy is defined as ( 5) with (6)   The ab initio MP2 fully optimized results reported in Table 4 shows that the monomer geometry deformation effect for the cisplatin-water complexes is small corresponding to less than 5% of the interaction energy (less than 20% of the BSSE correction).The MP2 BSSE correction stays virtually unchanged for complexes DDP-II and DDP-VI when full geometry optimization is performed, and increasing by 1.5 kcal mol -1 (almost doubled) for DDP-V, due to a substantial decrease of the equilibrium intermolecular distance.
The PBE BSSE values have a corresponding smaller variation.The energy trend predicted by the PBE functional, when compared to the MP2 results, is    quite similar to the previously reported relative energy results for the potential energy curves with no angular dependence (see values in parenthesis).It can also be seen from Figure 6 that there is a significant structural re-organization when full geometry optimization is performed.The water molecule tends to stay above the cisplatin plane, except for structure VI where the water molecule stays in line with the cisplatin with both hydrogen atoms pointed towards the chlorine atoms.
The ∆E values for structures II, V and VI evaluated at the CCSD(T) level of theory (coupled cluster with single, double and perturbative triple excitations) with the MP2 fully optimized geometries are respectively −13.92, −13.32 and −7.15 kcal mol -1 .The corresponding MP4(SDTQ) (Mfller-Plesset fourth order perturbation theory with single, double, triple and quadruple excitations) values are respectively −14.41, −13.49 and −7.16 kcal mol -1 .The CCSD(T) BSSE corrections for structures II, V and VI are respectively 4.67, 2.97 and 1.55 kcal mol -1 (the corresponding MP4(SDTQ) values are respectively 4.77, 3.03 and 1.58 kcal mol -1 ) and the CCSD(T) BSSE contribution due to the cisplatin monomer respectively 3.54, 1.84 and 1.10 kcal mol -1 (the corresponding MP4(SDTQ) values are respectively 3.62, 1.88 and 1.12 kcal mol -1 ).It can be seen that the MP2 level of theory is quite satisfactory for the calculation of interaction energies, and most important, it confirms that the BSSE behavior shown in previous Tables and Figures is not a particular characteristic of the MP2 level of theory, since the MP4(SDTQ) and CCSD(T) BSSE results follow closely the same trend as the MP2 data.
Lastly, the MP2 fully optimized equilibrium geometries, including the BSSE correction during the geometry optimization procedure as implemented in the Gaussian 2003 package, 29 are given in Table 4, along with the BSSE correction for equilibrium structures located on a CPcorrected PES, i.e., using an a priori CP correction scheme.It can be seen that there is a small geometry deformation and increase in the Pt-O distance, with a consequent decrease of the BSSE value for the CP-corrected structures.However, the interaction energies corrected for BSSE (∆E B M S P S 2 E-c ), calculated for the BSSE corrected equilibrium geometries (a priori CP correction), do not differ significantly from the corresponding values evaluated in a usual way as single point calculations using the BSSE uncorrected optimized geometries (a posteriori scheme).Similar results are found for the water dimer where the BSSE value using the MP2 CP-corrected optimized geometry is 0.12 kcal mol -1 smaller than the a posteriori scheme value, and the ∆E B M S P S 2 E-c value calculated using a priori scheme is −4.72 kcal mol -1 comparing to −4.66 kcal mol -1 evaluated using a posteriori BSSE correction.In addition, a small geometrical distortion is also observed with the MP2 CP-corrected optimized R O-O distance being 0.07 Å longer and the intermolecular angles having a maximum deviation of 1.5 degrees.

Conclusions
The cisplatin-water complex may be viewed as a model system that can be treated by accurate ab initio methods, where a comparative study can provide a guide for theoretical investigations involving larger molecular complexes containing transition metal atoms, where more modest theoretical approaches are to be used.The BSSE results reported here illustrate how the counterpoise correction works for molecular complexes involving transition metal species, where the use of effective core potential is required.The analysis of the individual contributions of cisplatin and water to the BSSE revealed that the water monomer shows a regular behavior as the basis set is improved, with the cisplatin subunit carrying most of the BSSE correction as the basis set is enlarged, showing even an increase of the BSSE for some spatial orientations at the MP2 level.Improving the ECP for the platinum atom (ECP*), 28 lead to a significant decrease in the MP2 BSSE results at a substantial computational cost, nevertheless, being still approximately twice the corresponding DFT values.To further improve the ECP* many additional polarization functions should be added to the metal valence shell what would make the calculations very expensive and perhaps of no practical use to supramolecular and organometallic chemistry.As already known in the literature for hydrogen bonded and weakly bound dimers, 35 the BSSE is relatively quite small at the HF and DFT levels with triple-zeta quality basis sets augmented by polarization and diffuse functions compared to the corresponding post-HF values, what corroborates the results presented here for cisplatin-water complexes.
A last point concerns the relevance of BSSE corrected PECs for the determination of Lennard-Jones (12-6)  parameters used in classical Monte Carlo simulation studies.It can be seen from Figures 4 and 5 that the evaluation of interaction energies at short, middle and long intermolecular distances is affected by the BSSE correction in a distinct extension (also depending strongly on a specific monomers spatial orientation), with the DFT correction expected to be substantially smaller.Based on theoretical grounds, we should calculate BSSE corrected PECs for the determination of LJ parameters.However, the results presented here for cisplatin-water, and monoaquo and diaquo species (structure I), show that the use of BSSE uncorrected PEC for the generation of LJ parameters leads to a better accordance between the classical and quantum mechanical calculated curves.However, to make a more sound comparison it would be necessary to calculate bulk properties using simulation methods that use the LJ potential function and make comparisons with the corresponding experimental values for cisplatin-water system in order to decide which set of LJ parameters really lead to the best agreement with experiment.

Figure 3 .
Figure 3. MP2 and PBE interaction energy calculation results (kcal mol -1 ) for the C 2v equilibrium structure of the water dimer.The results reported with the Dunning's basis sets were obtained as single point energy calculation at the equilibrium geometry determined with the cc-pVDZ basis set.

Figure 4 .
Figure 4.The Lennard-Jones (12-6) plus Coulomb parameterized potential energy curve (solid line) with and without BSSE correction and the corresponding MP2/6-31G(d,p)/LANL2DZ ab initio curves for the cisplatin-water complex (DDP-I) showing the relative better agreement for the BSSE uncorrected curves.

Table 4 .a
MP2 and PBE/6-311++G(2df,p)/LANL2DZ(f) fully optimized geometrical parameters a and ∆E/ kcal mol -1 values for the relevant structures located on the PES for the cisplatin-water complex bThe intermolecular geometrical parameters are specified in Figure6.The frozen monomer equilibrium values for geometry and energy are given in parenthesis.The values calculated with MP2 geometry including a priori BSSE correction during full geometry optimization (CP-corrected PES) is shown in brackets.b During the geometry optimization procedure it was found that DDP−I and DDP-IV optimize to DDP-II, and DDP-VIII optimizes to DDP-V.

Figure 6 .
Figure 6.(a) Geometric arrangement of the two interacting monomers.R AB is the intermolecular distance, θ A and θ B the intermolecular angles, and f the dihedral angle between the two C 2v axes.(b) The MP2/6-311++G(2df,p)/LANL2DZ(f) fully optimized structures (true minima) located on the PES for the cisDDP-water complex (side view) and (c) Top view of the MP2 optimized structures.

Table 2 .
Equilibrium distances (d(Pt … O w )/Å), interaction energy (∆E / kcal mol -1 ) and BSSE / kcal mol -1 calculated at MP2 level for the cisplatin-water system.The individual BSSE contribution for cisplatin (BSSE cis-DDP / % ) is also included as well as the percentage contributions in parenthesis