Predicting the binding properties of cibacron blue F3GA in affinity separation systems

International Journal of Biological Macromolecules 41 (2007) 430–438

Predicting the binding properties of cibacron blue F3GA in affinity separation systems Cenk A. Andac a,b , Muge Andac c , Adil Denizli c,∗ a

Dicle University, Department of Chemistry, Diyarbakir, Turkey b Gazi University, School of Pharmacy, Ankara, Turkey c Hacettepe University, Department of Chemistry, Ankara, Turkey Received 19 March 2007; received in revised form 4 June 2007; accepted 4 June 2007 Available online 26 June 2007

Abstract The binding properties of cibacron blue F3GA (CB-F3GA) bound to a model NAD(P)H/FAD(H2 )-dependent protein system, namely cytosolic quinone reductase (QR), was characterized by AMBER in an attempt to address the binding properties of immobilized CB-F3GA used in the separation of serum albumin. A favorable binding free energy of −4.52 kcal/mol (KD = 5.09 × 10−4 kcal/mol) was determined for CB-F3GA binding by MM-PBSA method, which was found to be a ballpark estimate of empirical values reported in literature (G ≈ −6 kcal/mol). We propose that CB-F3GA primarily follows a class III binding motif in presence of FAD in the binding site of QR in solution, while a class II binding motif is observed in the crystal form. It was found that favorable van der Waals/hydrophobic interactions take place in the binding site making a major contribution to a favorably dominating enthalpy of binding (Htot = −25.87 kcal/mol) as compared to a disfavorable binding entropy term (TStot = −21.35 kcal/mol). Additional MM-PBSA experiments in the absence of FAD gave rise to a disfavorable binding free energy for CB in complex with QR, suggesting that FAD is an essential determinant of CB-F3GA binding. This is in contrast to an earlier observation of Denizli et al. on separation of human serum albumin (HSA) by immobilized CB-F3GA in the absence of FAD. Therefore, a class I binding model for CB-F3GA is proposed here to account for the efficient separation of HSA in affinity chromatography systems. © 2007 Elsevier B.V. All rights reserved. Keywords: Cibacron blue F3GA; Quinone reductase; HSA; Molecular dynamics; Affinity; AMBER; MM-PBSA; Electrostatic surface map

1. Introduction A number of textile dyes, known as reactive dyes, have been used for protein purification in dye–ligand affinity systems, as they bind a variety of proteins in a selective and reversible manner. Dyes are commercially available, inexpensive and can be easily immobilized, especially on matrices bearing hydroxyl groups [1]. Immobilized textile triazine dyes, particularly cibacron blue F3GA (CB-F3GA), have been used as affinity chromatography tools for protein and enzyme purification for about 30 years [2]. CB-F3GA is a derivative of monochlorotriazine dye (Fig. 1), containing three titratable acidic sulfonate groups and four basic primary and secondary aromatic amine groups, which binds with considerable specificity and affinity to nucleotide dependent enzymes and to a series of other proteins

∗

Corresponding author. E-mail address: [email protected] (A. Denizli).

0141-8130/$ – see front matter © 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.ijbiomac.2007.06.008

including serum albumin [3]. For instance, Odabası and Denizli [4] observed an increase in the level of human serum albumin (HSA) adsorption onto a CB-F3GA-immobilized PHEMA adsorbent. However, the exact binding mechanism of immobilized CB-F3GA and the origin of its selectivity towards serum albumin currently remain mysterious. Thus, seeking ways to address a binding mechanism for HSA + CB-F3GA mutual relationship has been of primary importance in the development of more efficient affinity separation systems. There has been only two bound structures of CB-F3GA ever determined thus far by X-ray crystallographic studies, whose coordinates were deposited at the Protein Data Bank for CBF3GA in complex with a glutathione S-transferase (GST) (PDB registry code: 2OGS) [5] and CB-F3GA in complex with cytosolic quinone reductase (QR) (PDB registry code: 1QRD) [6]. The GSTs are a family of dimeric enzymes which play a crucial rule in the detoxification of carcinogenic, mutagenic, toxic and pharmacologically active compounds by conjugating them to the thiol group of the cellular nucleophile glutathione (GSH, ␥-

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

431

NAD(P)H/FAD(H2 )-dependent enzymes, which provides further insight into the likelihood of what class of binding motif CB-F3GA immobilized polymers [8–11] may exhibit in the separation of NAD(P)H/FAD(H2 )-independent plasma proteins such as HSA. 2. Materials and methods Fig. 1. The chemical structure of CB-F3GA (CB-F3GA). The aromatic functional groups are numbered as I for the 1,4-diamino-2-sulfonyl-9,10diketoanthracenyl (anthraquinonyl) group, II and IV for the benzenesulfonyl groups and III for the 4-chloro-1,3,5-triazabenzenyl group, which are all linked by secondary amino groups. All sulfonate groups are deprotonated at physiological pH values, yielding a total of 3 negative charges overall.

Glu–Cys–Gly). In the X-ray structure of the GST + CB-F3GA complex, only the anthraquinone domain of CB-F3GA (ring I in Fig. 1) is bound to each subunit of GST as no density was observed for the rest of the ligand [5], indicating that the anthraquinone moiety of CB-F3GA is buried in the binding site while other rings (rings II–IV in Fig. 1) of the ligand stick out of the binding site. Thus the bound CB-F3GA exhibits, what we refer to as, a class I binding motif in complex with GST type proteins. QR on the other hand is also a dimeric enzyme but, unlike GSTs, its activity is strongly dependent on cofactors FAD(H2 ) or NAD(P)H to catalyze two-electron reductions to protect cells against the toxic and neoplastic effects of quinones, quinoneimines, azo dyes and nitro groups [6]. Moreover, it is known that the QR level is elevated in tumor cells and also QR plays a role in the reductive activation of chemotherapeutic quinones such as mitomycins and aziridylbenzoquinones [7]. Each subunit in the X-ray structure of the QR + CB-F3GA complex (PDB registry code: 1QRD) [6] is associated with flavin adenine dinucleotide (FAD) and duroquinone (DQ). In this XRay structure, CB-F3GA follows, again we classify another interaction mode, a class II binding motif in which, unlike the class I binding mode, the anthraquinone moiety this time sticks out of the binding site and all other rings of CB-F3GA (rings II–IV in Fig. 1) interacts with the binding site of protein. The class I and II binding modes of bound CB-F3GA are roughly generalized here based only on the two aforementioned X-ray structures available at the Protein Data Bank. In addition to class I and II binding modes, a third type of, namely class III, binding mode in general would be useful to define in the context of bound CB-F3GA studies, which is utilized in Section 3, in that all rings of CB-F3GA interact with the binding site of protein. Although no precedents of class III binding motif is available yet in a crystal form of bound CB-F3GA, we believe this binding mode represents some solution structures of bound CB-F3GA. The primary aim of this paper is to computationally characterize molecular dynamics (MD) and binding properties of bound CB-F3GA by Assisted Model Building with Energy Refinement (AMBER) suite of programs using the X-ray structure of cytosolic quinone reductase in complex with CB-F3GA (PDB registry number: 1QRD) [6] as a model. In addition, the function of FAD in the binding process of cibacron blue-F3GA is thermodynamically investigated here. Structural and thermodynamic information gained in this AMBER study is utilized to determine the binding mode of CB-F3GA in complex with

2.1. System set-up and initial structures Molecular dynamics (MM) and thermodynamic computations were implemented by AMBER v9 (2006) suite of programs [12] running mainly under a 64 bit parallel Linux-PC system. The initial structure of CB-F3GA, shown in Fig. 1, in complex with cytosolic quinone reductase associated with flavin adenine dinucleotide and duroquinone was obtained from crystal structure PDB ID 1QRD [11]. The published X-ray structure reflects two identical complex units non-covalently attached to each other. Thus, only the first complex unit listed in the original PDB file was selected for molecular dynamics studies. Water molecules were manually removed form the selected complex unit, yielding a cytosolic QR protein (with 273 amino acid residues giving a total molecular weight (M.W.) of 30838 Da) in complex with FAD (M.W. 786.58 Da.), DQ (M.W. 164.2 Da.) and CB-F3GA (M.W. 753.76 Da.). Protons of the crude complex structure were added by the LeaP [13] module of AMBER v9 using AMBER1999SB force field [14], general AMBER force field (GAFF) [15], contributed parameters for FAD [16] and GAFF parameters for CB-F3GA and DQ (see GAFF parameters section below). 2.2. Partial atomic charges and GAFF parameters Cornell’s [17] partial atomic charges (determined at the HF/6-31G* basis set [18,19] via the RESP [20] method) in AMBER1999SB force field were assigned to the amino acid residues on QR. Atomic partial charges for CB-F3GA and DQ were computed at the AM1-BCC [21] semiempirical hemiltonian level by Divcon 2005 quantum mechanics software [22] incorporated into the antechamber module of AMBER v9 [23]. A total of three negative charges was assigned to CB-F3GA and a zero charge was given to DQ, which were then optimized at 300 K using a conjugate gradient method to converge the structures for AM1-BCC partial atomic charge computations. The computed AM1-BCC charges and the optimized coordinates for CB-F3GA and DQ were saved in GAFF mol2 format, for which GAFF parameters were prepared by the parmchk module of AMBER v9 [12] to be used as inputs in Leap. Partial atomic charges (determined at the HF/6-31G* basis set) and AMBER parameters for FAD were obtained from Antony et al. [16] as posted at http://amber.scripps.edu. 2.3. Molecular mechanics Addition of counterions, solvation, and preparation of parameter/topology and coordinate files were all implemented in LeaP. The molecular total charge of the complex molecule (QR + CB-

432

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

F3GA + FAD + DQ) is 2 (−). Therefore, 2 Na+ cations were added to neutralize the complex system. The complex molecule was solvated by putting together small units of TIP3 water boxes ´˚ × 18 A ´˚ × 18 A ´˚ to constitute an octahedral in dimensions of 18 A ´˚ × 91 A ´˚ × 91 A. ´˚ 11025 water box in dimensions of about 91 A TIP3 water molecules were added to solvate the complex yielding a MD system with a total of 37617 atoms (including water molecules, ions and the complex). The distance between the outer boundary of the octahedral solvent box and the solute sur´˚ A space of 0.4 A ´˚ was used to set water face was set to 8 A. molecules at the solute/solvent boundary. Parameter/topology and coordinate files for the neutralized and solvated complex system was prepared in Leap. Temperature equilibration and MD routines were conducted for the complex molecule by the sander module of AMBER v9 [13]. The coordinates of water molecules, and the Na + ions in the solvated starting structure were initially relaxed to remove bad close contacts through a short restrained minimization routine for 0.2 ps using a steepest descent method with 1 fs time steps over 200 iterations, during which the atom coordinates for the solute species were restrained to their initial coordinates ´˚ The temperature of the with a force constant of 100 kcal/mol A. relaxed system was then equilibrated at 300 K through 2 ps of MD using 2 fs time steps over 1000 iterations. A constant volume periodic boundary was set to equilibrate the temperature of the system by the Langevin dynamics [24] using a collision frequency of 10/ps and a velocity limit of five temperature units. The complex structure in the solvated system was restrained to its ´˚ initial coordinates with a weak force constant of 20 kcal/mol A during the temperature equilibration routine. The final coordinates of the temperature equilibration routine (after 2 ps) was then used to complete a 1100 ps (1.1 ns) molecular dynamics run using 2 fs time steps through 550,000 iterations, during which the temperature was kept at 300 K using the same Langevin dynamics parameters as applied in the restrained–relaxation routines, and the pressure of the solvated system was equilibrated at 1 bar at a certain density in a constant pressure periodic boundary by an isotropic pressure scaling method employing a pressure relaxation time of 2 ps. During the temperature equili´˚ bration and MD routines, a non-bonded cutoff distance of 9 A was applied to handle electrostatic interactions in periodic boxes by the Particle Mesh Ewald method [25] and SHAKE method [26] was applied to keep the bond lengths of protons attached to heteroatoms constant. Coordinates and energy outputs for the relaxation and the molecular dynamics routines were saved every 200 iterations. 2.4. Thermodynamic computations Thermodynamic binding computations were conducted by the mm-pbsa [27] module of AMBER v9. For the formation of the complex, as generically shown in Eq. (1), MM-PBSA [28] computations including absolute enthalpy (H), where absolute means single-state, and absolute entropy (S) were separately implemented for an isolated complex including

QR + FAD + CB-F3GA, its receptor (QR + FAD) and the ligand (CB-F3GA). Receptor + Ligand ↔ Complex

(1)

In order to evaluate the binding energy in the absence of FAD, a second set of MM-PBSA computations was carried out for the components of the complex, in which case FAD was not treated as parts of the receptor and the complex (QR + CB-F3GA) species. DQ was excluded in all thermodynamic computations because it becomes fully detached from the protein during the course of the MD experiment. The absolute free energy (G) of the complex systems, their receptors, and the ligand was computed in a classical manner as in Eq. (2), in which T is the temperature of the system at 300 K. G = H − TS

(2)

The binding free energies (G) of the complex systems were computed as in Eq. (3) where Gcomp is the absolute free energy of the complex, Grec is the absolute free G = Gcomp − [Grec + Glig ]

(3)

energy of the receptor, and Glig is the absolute free energy of the ligand. Fifty snapshots were extracted for the coordinates of the solute species (complex, receptor and ligand) at 2 ps time intervals between 1.0 ns and 1.1 ns of the trajectory. H subenergy terms were averaged over 100 snapshots to constitute mean absolute H terms. The enthalpy term in Eq. (2) is dissected into subenergy terms as in Eq. (4) Htot = Hgas + Htans/rot + Gsolv

(4)

where Hgas is the potential energy of the solute which is determined as a sum of van der Waals (Evdw ), electrostatic (Eel ) and internal energies (Eint ) in gas phase by sander using the molecular mechanics algorithm in Cornell et al. [17] force field, Htans/rot is due to six translational and rotational degrees of freedom (6*1/2RT) and is equal to 1.79 kcal/mol at 300 K [27], and Gsolv is the solvation free energy for transferring the solute from vacuum into solvent and is a sum of electrostatic (Gel ) and nonelectrostatic (hydrophobic) contributions (Gnonel ) as seen in Eq. (5). Gsolv = Gel + Gnonel

(5)

Gel in Eq. (5) was computed at 0.15 M salt concentration by the pbsa module of AMBER v9 by dividing implicitly solvated ´˚ cubic grid points and summing up the solute species into 0.4 A electrostatic potentials computed at each grid point. Electrostatic potential ∅(r) at a grid point r that is not at the solvent–solute boundary was computed by a linear Poisson Boltzmann (PB) equation [28], which is a three dimensional vector differential equation as in Eq. (6)  ∇ε(r) ∇∅(r) = −4πρ(r) − 4π zi exp(−zi ∅(r) /kB T ) (6) in which ε(r) is the dielectric constant (ε = 1 for the solute interior and ε = 80 for implicit PB water), ρ(r) is the charge density, zi is

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

the charge of ion type i, kB is the Boltzman constant and T is the temperature; the summation is over all ion types. An atomic electrostatic potential ∅i at the solute–solvent boundary was computed by pbsa as implemented in Delphi v.2 [29] and takes the simple form as in Eq. (7)   exp(−kr) ∅i = qi (7) εw r in which qi is the partial charge of atom i at the boundary, r is the distance of the atom i from the boundary, k is the Debye salt screening parameter [30], and εw is the solvent dielectric constant (εw = 80 for water). A sum of the atomic electrostatic potentials for individual atoms (∅i ) as a factor of their partial atomic charges qi at the solute–solvent boundary and outside this boundary constituted Gel as seen in Eq. (8) Gel = 1/2Σqi ∅i

(8)

The PB implicit solvent molecules at the solute-solvent boundary were allowed to energetically converge over 1000 iterations before the single-point PB computations were implemented by ´˚ pbsa for each snapshot. A spherical solvent probe (radii) of 1.6 A and atomic radii provided by AMBER force field were used for the implicit solvent molecules and solute atoms, respectively, during the pbsa computations. Gnonel in Eq. (5) is proportional to the solvent accessible surface area (SASA) of the solute and was computed by the molsurf module of AMBER v9 [29] according to Eq. (9) Gnonel = γSASA + b

(9)

in which γ is the surface tension proportionality constant ´˚ 2 ), b is the free energy of non-polar solva(0.0072 kcal/mol A tion for a point solute (0 kcal/mol). SASA was computed by molsurf via a fast algorithm called the Linear Combinations of Pairwise Overlaps (LCPO) which simply computes a surface area of two overlapping atoms by substracting their overlapping region from the total area of their van der Waals surfaces [31]. Van der Waals radii in AMBER and GAFF force fields were used to define solute atoms which were then augmented by a ˚ for the molsurf computations. probe sphere of 1.4 A The absolute entropy was computed for each solute species by normal-mode analysis [32] integrated in the nmode module of AMBER v9 [12]. Because entropy computations are very costly in terms of time and entropy fluctuations between snapshots are usually negligible, the entropy of each solute species was computed for one snapshot only. The total entropy (Stot ) computed, as formulated in Eq. (10), arose from changes in the degree Stot = Strans + Srot + Svib

which can be further simplified as Eq. (12). Gbinding = Htot − TStot

(12)

Parameter/topology files used in mm-pbsa computations were prepared for the complex, the receptor and the ligand species by implicit treatment of these species in Leap. Snapshots extracted from trajectories were pre-minimized in the gas phase by sander using a conjugate gradient method with a distance dependent dielectric of 4r (with r being the distance between two atoms) until the root-mean-square-deviation (RMSD) of the elements ˚ frequenof the gradient vector was less than 10−4 kcal/mol A, cies of the vibrational modes were computed at 300 K for these minimized structures including all snapshot atoms and using a harmonic approximation of the energies [27]. 2.5. Trajectory analysis RMSD, inter-residual H-bond analysis and average PDB structure computations throughout the trajectory were implemented by the ptraj module of AMBER v9 [12]. 3. Results and discussion 3.1. Cibacron Blue F3GA (CB-F3GA) structure in complex with protein CB-F3GA possesses four aromatic structures linked by secondary aromatic amino groups (Fig. 1). CB-F3GA in the bound form rather tends to adopt a rigid dynamic structure throughout the molecular dynamics experiment over 1100 ps (1.1 ns). The root-mean-square-deviation plot in Fig. 2 (lower line) shows that CB-F3GA in the bound form reaches an equilibrium state beyond 400 ps of molecular dynamics and this steady-state is ˚ RMSD maintained through the rest of the experiment around 2 A referenced to the crystal structure coordinates, suggesting that the binding mode of CB-F3GA is quite specific to the binding site on QR. Fig. 3 is a superposed image of the bound CB-F3GA in the X-ray structure (red) and after 1.1 ns of molecular dynamics (blue), indicating that the observed RMSD with CB-F3GA mainly arises from a slightly greater, and yet limited, degree of flexibility about the secondary aromatic amino linkages (Fig. 3).

(10)

of freedom [translational (Strans ), rotational (Srot ), and vibrational (Svib )] of each species [27]. When all absolute energy terms are put together in Eq. (3), the binding free energy G takes the following form, Gbinding = [Hgas + Hsolv − 1.79] − TStot

433

(11)

˚ of cibacron blue (lower line) and Fig. 2. Root-mean-square deviations (A) quinone reductase (upper line) for 1100 ps (1.1 ns) simulation referenced to the X-ray structure coordinates.

434

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

Fig. 3. Superposed structures of bound CB-F3GA in the X-ray structure (red) and after 1.1 ns of molecular dynamics (blue). Nomenclature for the numbered aromatic functional groups is given in Fig. 1. The aromatic amino linkages are shown with asterisks. RMSD of the AMBER structure (blue) referenced to the ˚ (For interpretation of the references to colour X-Ray structure (red) is 1.99 A. in this figure legend, the reader is referred to the web version of the article.)

In comparison to the bound CB-F3GA, the equilibration of cytosolic quinone reductase in the complex takes less time during molecular dynamics simulation and seems to reach a steady-state beyond 800 ps as observed on the RMSD plot in Fig. 2 (upper line). A superposed view of the QR backbone traces for the X-ray (red) and 1100 ps simulation (last snapshot) coordinates in Fig. 4 suggests that QR does not sample a major conformational change at the steady state in solution, which is ˚ RMSD during the course of also supported by less than 3 A simulation (Fig. 2, upper line). The binding process of CB-F3GA to QR follows an interesting trend in the course of MD in aqueous media. As proposed in the introduction section, CB-F3GA exhibits a class II binding motif in the binding site of QR in the crystal form. That is, the anthraquinone moiety (ring I in Fig. 5A) sticks out of the binding site while rings II–IV interacts with QR binding site in the X-ray structure (Fig. 5A). It is clearly seen in the X-ray structure that DQ (orange surface in Fig. 5A) stacks over the pteridine ring of the flavin moiety of FAD (green surface in Fig. 5A), which is the reaction site of FAD, and also interacts with ring III of CB-F3GA (cyan in Fig. 5A), while a part of DQ is bound to QR. Moreover, the sulphonate group on ring IV of CB-F3GA electrostatically interacts with the hydroxyl groups on the penthanetetraol linker of FAD in the X-ray structure. However, DQ dislocates during the course of MD and strongly interacts with the anthraquinone moiety (ring I) of CB-F3GA

Fig. 4. A backbone-trace view of the superimposed structures for the X-ray structure of QR in the complex form (red, PDB ID 1QRD) and the last snapshot of a 1100 ps molecular dynamics simulation for the QR part of the complex (blue). RMSD of the AMBER structure (blue) referenced to the X-Ray structure ˚ (red) is 2.82 A.

(Fig. 5B), becoming fully detached from the protein. Rings III and IV of CB-F3GA then shift towards FAD, finally stacking over the flavin moiety of bound FAD in the MD structure. This is a strong indication of a competitive inhibition of QR by CBF3GA against the DQ ligand. In that case, CB-F3GA plays two functional roles: one is to inhibit the active site (pteridine ring) of FAD and the other is to capture the ligand. The MD solution structure of bound CB-F3GA exhibits a more comprehensive binding motif, different than classes I and II binding modes as proposed in the introduction section, in that the anthraquinone moiety of CB-F3GA, which is found to be free in the X-ray structure, now tends to interact with the amino acid residues on the surface of the protein in the MD structure. This binding motif of CB-F3GA fully interacts with the components of the receptor and thus we refer to this as a third type of binding motif, class III binding mode, which we believe it is only observed in the binding site of NAD(P)H/FAD(H2 )-dependent proteins in solution. Electrostatic surface map of QR, shown in Fig. 6, clearly depicts that the anthraquinone moiety (I in Fig. 6) and the

Fig. 5. Surface views for (A) the X-ray (PDB registry number: 1QRD) [6] and (B) MD solution structures of CB-F3GA (cyan) in complex with QR (grey), FAD (green) and DQ (orange). Hydrogen is not included in the surface views of the bound components for simplicity.

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

Fig. 6. Delphi electrostatic surface map for the final coordinates of 1100 ps molecular dynamics simulation of cytosolic QR in complex with FAD (green stick view), DQ (orange stick view) and CB-F3GA (atom-color stick view). Aromatic groups of CB-F3GA are numbered according to Fig. 1. Redish spots on the protein surface indicate electronegative surface areas, blueish spots point out electropositive surface areas and whitish spots stress mainly on van der Waals/hydrophobic surface areas. The Delphi electrostatic surface map was determined by Delphi v.4 (2002) [33,34] and was visualized by Chimera (UCSF) [35].

terminal benzenesulfonyl group (IV in Fig. 6) of CB-F3GA mainly undergo van der Waals (vdW)/hydrophobic interactions within the binding site of QR (the whitish area in the CBF3GA binding site in Fig. 6), while the second benzenesulfonyl group (II in Fig. 6) and the 1,3,5-triaza-monochlorobenzene group (III in Fig. 6) seem to sample both electrostatic (the corresponding bluish and redish areas on the protein surface) and vdW/hydrophobic interactions, suggesting that CB-F3GA primarily makes specific vdw/hydrophobic interactions in the binding site of proteins while some degrees of inter-molecular electrostatic interactions also contribute to the binding affinity and specificity. In addition, ptraj analyses for equilibrated snapshots of CB-F3GA in complex with QR yields a moderate H-bond formation between the 9-oxygen of the anthraquinone moiety of CB-F3GA and the ε-nitrogen of the His161 residue of the protein. This indicates that the free anthraquinone moiety of CB-F3GA in the X-ray structure has drifted in the course of MD to make a specific interaction with the protein. Interestingly, despite the three negatively charged aromatic sulfonate groups on CB-F3GA, none of these charged groups is associated with electrostatic interactions with the protein in the MD solution structure, but they rather stick out of the binding site to perhaps increase solubility in the outer water shell. In affinity chromatography systems, CB-F3GA is usually immobilized on to a solid support (matrix) via a nucleophilic substitution reaction between the Cl atom at position 4 of the 1,3,5-triazabenzenyl group (ring III in Figs. 1 and 6) and an aliphatic linker bound to the matrix. This in turn suggests that the position 4 of the 1,3,5-triazabenzenyl group should essentially face the outer water shell in the bound form. This notion is, indeed, in contrast with the binding motif of CB-F3GA in complex with QR, in which position 4 of the 1,3,5-triazabenzenyl group (ring III) is deeply immersed in the binding site while rings III and IV of CB-F3GA are able to interact with FAD by stacking over the flavin moiety of FAD. Therefore, it would

435

be plausible to consider that separation of QR by such CBF3GA-immobilized affinity chromatography systems may not be as applicable because the immobilization site on CB-F3GA is required to properly contact in the binding site of QR. We believe this rule should apply to class II and class III binding motifs of CB-F3GA observed in the binding sites of NAD(P)H/FAD(H2 )dependent enzymes/proteins. In terms of providing a rational explanation for the separation of serum albumin by CB-F3GA-immobilized affinity chromatography systems, we propose that the immobilized CBF3GA interacts with the binding site of serum albumin through a class I binding motif, similar in spirit to a bound motif observed in the binding site of NAD(P)H/FAD(H2 )-independent GST proteins [5]. We reach this conclusion based on an earlier experimental work of Odabası and Odabasi and Denizli [4] on an efficient separation of serum albumin in the absence of NAD(P)H2 or FAD(H2 ). 3.2. Thermodynamic properties of bound CB-F3GA In regards to finding more evidence to illuminate the binding properties of CB-F3GA in complex with NAD(P)H/FAD(H2 )dependent proteins, in which case our protein is the cytosolic quinone reductase, thermodynamic computations by the MMPBSA method was conducted on the complex according to experimental details given in Section 2.4, in which QR + FAD constitutes the receptor while CB-F3GA is the ligand. A second set of MM-PBSA computations were also carried out in the absence of FAD to determine the influence of FAD on CB-F3GA binding, in which case QR is the receptor and CB-F3GA is the ligand. Thermodynamic binding energy results are listed in 1 of Table 1 for the complex with FAD and in 2 of Table 1 for the complex without FAD. The internal energy changes for the binding processes in gas phase (Eint ) amounts to zero in Table 1 because conformations of the reacting solute species originate from the same complex structure, and thus any structural changes upon binding is neglected here. In the presence of FAD (1 of Table 1), the electrostatic interaction energy in gas phase (Eel ) between the receptor and CB-F3GA is favorable by −103.29 ± 13.20 kcal/mol. On the other hand, electrostatic contributions to the solvation energy (Gel ) is disfavorable by 120.57 ± 13.43 kcal/mol, which altogether sum to a disfavorable electrostatic contribution (Eel + Gel ) to the enthalpy of binding by 17.28 kcal/mol, suggesting that the binding of CB-F3GA displaces all water molecules at the binding interface and energy is consumed to align the ligand and the receptor during the binding process. Compared to the disfavorable electrostatic interactions in the bound state in solution (Eel + Gel ), the van der Waals interaction energy in gas phase (Evdw = −36.40 ± 3.07 kcal/mol) associated with non-polar contributions to the solvation energy (Gnonel = −4.96 ± 0.26 kcal/mol) seems to favorably manage the binding process, yielding a total vdW/hydrophobic energy (Evdw + Gnonel ) of −41.36 kcal/mol. The favorable vdW/hydrophobic energy compensates the disfavorable effects of electrostatic interactions at the binding interface, turning the case into a favorable binding process predominantly driven by

436

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

Table 1 MM-PBSA binding energies. All energy values are given in kcal/mol Hgas

1 2

1 2

Gsolv

H

Eel

Evdw

Eint

Gel

Gnonel

(Hgas + Gsolv )

−103.29 ± 13.20 −110.66 ± 9.40

−36.40 ± 3.07 −22.34 ± 3.07

0.00 ± 0.00 0.00 ± 0.00

120.57 ± 13.43 123.09 ± 10.39

−4.96 ± 0.26 −3.66 ± 0.29

−24.08 ± 3.33 −13.57 ± 3.45

TStrans

TSrot

TSvib

TStot

Htot

Gbinding

−13.71 −13.71

−12.18 −12.18

4.55 −6.40

−21.35 −32.29

−25.87 −15.35

−4.52 +16.93

Mean energy values were averaged over 100 snapshots. 1 represents the complex with FAD [QR + FAD + CB-F3GA] and 2 refers to the complex without FAD [QR + CB-F3GA]. Eel , Evdw and Eint are electrostatic, van der Waals and internal, respectively, energies of binding in gas phase and a sum of them yield the energy of binding in the gas phase (Hgas in Eq. (4)). Gel (Eq. (8)) and Gnonel (Eq. (9)) represent electrostatic and non-electrostatic, respectively, contributions to the solvation free energy (Gsolv , Eq. (5)). Htot is the enthalpy of binding and is equal to [Hgas + Gsolv −1.79], where 1.79 is due to Htrans/rot (Section 2.4). Stot is the entropy of binding and is a product of translational (Strans ), rotational (Srot ) and vibrational (Svib ) entropy changes (Eq. (10)). Gbinding is the binding free energy and is equal to Eq. (11).

vdW/hydrophobic forces in solution. This in a way makes sense because, as indicated in Section 3.1, the MD solution structure of bound CB-F3GA in complex with QR suggests that none of the charged groups on CB-F3GA electrostatically interacts with the components of the receptor. Not surprisingly though, such dominant vdW/hydrophobic binding forces have some precedents for charged ligands, i.e. netropsin, a naturally produced antibiotic which binds AT rich regions of DNA, has + 2 charges and its association with the binding site of an AATT DNA duplex reveals a binding motif which is primarily driven by van der Waals forces despite the polyionic natures of netropsin and DNA [36]. The overall enthalpy of binding (Htot , Eq. (12)) was computed by Eq. (11) to be −25.87 kcal/mol, which makes a favorable contribution to the energy of binding. The total entropy term, TStot , presented in Table 1, corresponds to a binding entropy term and is a product of temperature (T) at 300K and the total entropy change (Stot ) which arises from changes in translational (Strans ), rotational (Srot ) and vibrational (Svib ) degrees of freedom upon binding. The binding entropy term in 1 of Table 1 is disfavorable (−21.35 kcal/mol) due primarily to disfavorable entropy contributions from translational and rotational entropy terms, TStrans and TSrot , respectively, while the vibrational entropy term TSvib , is favorable, but too weak to compensate the disfavorable effects of the translational and rotational entropy changes. The binding free energy for the complex (Gbinding ) computed according to Eq. (11), which gave rise to a favorable −4.52 kcal/mol (1 of Table 1). The dissociation constant (KD ) for CB-F3GA in complex with QR + FAD was computed in a classical manner [G◦ = −R.T.ln 1/KD , where G◦ is the binding free energy in cal/mol, T is the temperature at 300 K and R is the ideal gas constant, 1.987 cal/mol K] to be 5.09 × 10−4 kcal/mol, which suggests a moderate binding affinity towards QR. It appears that the binding process of CB-F3GA is favorably enthalpy driven as the disfavorable entropy term remains insufficient to counterbalance the favorable enthalpy effect. As a conclusion, the binding of CB-F3GA to QR in the presence of FAD is enthalpy driven, in which van der Waals/hydrophobic forces play a dominant role in the binding process.

MM-PBSA data in 2 of Table 1 for the complex formation in the absence of FAD resulted in a highly disfavorable Gbinding (+16.93 kcal/mol) emphasizing the importance of FAD on binding. That is, the best complex between CB-F3GA and QR forms only when FAD is in the vicinity of the binding site of rings III and IV of CB-F3GA. This strongly suggests that the highest affinity of CB-F3GA towards NAD(P)H/FAD(H2 )-dependent enzymes/proteins is reached only in association with a reducing cofactor such as NAD(P)H or FAD(H2 ). However, compared to CB-F3GA, we believe that the thermodynamics of immobilized CB-F3GA would be very similar to that of CB-F3GA binding observed in the absence of FAD (2 of Table 1), which would be disfavorable because, as pointed out in Section 3.1, ring III of CB-F3GA is the immobilization site and thus it should have sticked out of the binding site whether or not FAD is present. Even if an FAD is present, immobilized CB-F3GA would disrupt the required stacking interaction in between the flavin moiety of FAD and rings III and IV of CB-F3GA, thus, significantly lowering the binding affinity by increasing a disfavorable entropy of binding. This notion is also supported with an increased disfavorable entropy term (by −10.94 kcal/mol) in 2 of Table 1 for the loss of mutual interactions between FAD and CB-F3GA in the complex. Usually empirical binding free energies for CB-F3GA in complex with several proteins are reported in literature to be around −6 kcal/mol [37]. Therefore, the computational binding free energy result (−4.52 kcal/mol) presented in 1 of Table 1 is a ballpark estimate for the binding affinity of CB-F3GA towards QR. 4. Conclusion Based on our Delphi electrostatic surface map (Fig. 5) and MM-PBSA data (Table 1), it is concluded that bound CB-F3GA primarily faces favorable van der Waals/hydrophobic binding interactions in the binding site of QR in the presence of FAD, in which we propose that CB-F3GA adopts a class III binding motif in solution, while a class II binding motif is observed in crystal form. The binding of CB-F3GA in the presence of FAD

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

is favorably driven by dominating enthalpy contributions to the energy of binding. We believe that the class III binding motif inferred here from the dynamics and MM-PBSA data in the presence of FAD represents, in a general sense, a binding motif of CB-F3GA in complex with NAD(P)H/FAD(H2 )-dependent enzymes/proteins in solution. We determined that the affinity of CB-F3GA towards QR becomes highly disfavorable when FAD is not involved in binding. This is associated with the absence of a stacking interaction between rings III and IV of CB-F3GA and the flavin moiety of FAD in the binding site. Thus, we think that immobilized CB-F3GA would have a reduced affinity towards NAD(P)H/FAD(H2 )-dependent enzymes/proteins because ring III of CB-F3GA is conjugated to a matrix in affinity separation systems which would presumably disrupt the required interactions with FAD in the binding site. Interestingly, this finding opens a new research discussion as to whether or not CB-F3GA can be immobilized from a different site that would provide the same binding interactions with NAD(P)H/FAD(H2 )dependent enzymes/proteins as observed with non-immobilized CB-F3GA. For instance; a suitable linkage from ring I (the anthraquinonyl group) of CB-F3GA to a matrix system would be interesting to try to selectively separate NAD(P)H/FAD(H2 )dependent enzymes/proteins in solution. It has been known for quite sometime that CB-F3GAimmobilized matrixes can efficiently separate HSA out of plasma in the absence of FAD [8–11]. In similar spirit to CBF3GA binding to GST proteins [5], which is referred to as a class I binding model in the introduction section, we believe that the mechanism of separating HSA via immobilized CB-F3GA follows the principals of class I binding model in which only the anthraquinonyl moiety of CB-F3GA is involved in binding while the immobilized site (ring III) remains free in solution. This suggests that rings II–IV of CB-F3GA may be replaced with a suitable linker without affecting the affinity and selectivity of binding in the affinity separation of NAD(P)H/FAD(H2 )independent proteins. The proposed class I binding mode for the anthraquinone moiety of immobilized CB-F3GA in well-identified binding sites of HSA has been currently an ongoing research in our laboratory in order to develop novel CB-F3GA-immobilized matrixes which would presumably separate HSA out of plasma in a more efficient way. However, we do not know yet if the anthraquinonyl moiety of CB-F3GA can be further modified to increase its selectivity and affinity towards serum albumin. Perhaps the best way to test this would be to initially determine dominant forces constituting the binding interactions in a class I binding mode of CB-F3GA in complex with NAD(P)H/FAD(H2 )-independent proteins/enzymes. Acknowledgements We would like to express our best gratitudes to Prof. David Case at the Scripps Research Institute, La Jolla, CA, USA for providing us with a license agreement for AMBER v.9 (2006), and to Prof. Barry Honig at Columbia University, New York, NY, USA for making Delphi v.4 (2002) available to our laboratory.

437

We also would like to thank the TR-Grid e-Infrastructure of Turkey for allowing us to use their cluster computing facilities. This work was initially presented at the ICS-UNIDO meeting at Trieste, Italy on 9.5.2006. We also would like to thank all folks at the ICS-UNIDO meeting for their valuable criticisms and suggestions on improving the CB-F3GA binding study. References [1] A. Denizli, E. Pis¸kin, J. Biochem. Biophys. Methods 49 (2001) 391. [2] Y.D. Clonis, N.E. Labrou, V.Ph. Kotsira, C. Mazitsos, S. Melissis, G. Gogolas, J. Chromatogr. A 891 (2000) 33. [3] G. Kopperschlager, H.J. Bohme, E. Hofmann, Adv. Biochem. Eng. 25 (1982) 101. [4] M. Odabasi, A. Denizli, Polym. Int. 53 (2004) 332. [5] A.J. Oakley, M. Lo Bello, M. Nuccetelli, A.P. Mazzetti, M.W. Parker, J. Mol. Biol. 291 (1999) 913. [6] R. Li, M.A. Bianchet, P. Talalay, L.M. Amzel, Proc. Natl. Acad. Sci. U.S.A. 92 (1995) 8846. [7] D. Ross, D. Siegel, H. Beall, A.S. Prakash, R.T. Mulcahy, N.W. Gibson, Cancer Metastasis Rev. 12 (1993) 83. [8] E. Ruckenstein, X. Zeng, J. Membr. Sci. 142 (1998) 13. [9] A. Denizli, G. Kokturk, H. Yavuz, E. Pis¸kin, J. Appl. Polym. Sci. 74 (1999) 2803. [10] E. Altıntas, A. Denizli, J. Chromatogr. B 832 (2006) 216. [11] Z. Ma, Y. Guan, H. Liu, H. React, Funct. Polym. 66 (2006) 618. [12] D.A. Case, T.A. Darden, T.E. Cheatham III, C.L. Simmerling, J. Wang, R.E. Duke, R. Luo, K.M. Merz, D.A. Pearlman, M. Crowley, R.C. Walker, W. Zhang, B. Wang, S. Hayik, A. Roitberg, G. Seabra, K.F. Wong, F. Paesani, X. Wu, S. Brozell, V. Tsui, H. Gohlke, L. Yang, C. Tan, J. Mongan, V. Hornak, G Cui, P. Beroza, D.H Mathews, C. Schafmeister, W.S. Ross, P.A. Kollman, AMBER 9, University of California, San Francisco, 2006, http://amber.scripps.edu. [13] D.A. Case, T.E. Cheatham III, T. Darden, H. Gohlke, R. Luo, K.M. Merz Jr., A. Onufriev, C. Simmerling, B. Wang, R. Woods, J. Comput. Chem. 26 (2005) 1668. [14] V. Hornak, R. Abel, A. Okur, B. Strockbine, A. Roitberg, C. Simmerling, Proteins 65 (2006) 712. [15] Y. Duan, C. Wu, S. Chowdhury, M.C. Lee, G. Xiong, W. Zhang, R. Yang, P. Cieplak, R. Luo, T. Lee, J. Comput. Chem. 24 (2003) 1999. [16] J. Antony, D.M. Medvedev, A.A. Stuchebrukhov, J. Am. Chem. Soc. 122 (2000) 1057. [17] W.D. Cornell, P. Cieplak, C.I. Bayly, I.R. Gould, K.M. Merz Jr., D.M. Ferguson, D.C. Spellmeyer, T. Fox, J.W. Caldwell, P.A. Kollman, J. Am. Chem. Soc. 117 (1995) 5179. [18] Y. Osamura, Y. Yamaguchi, D.J. Fox, M.A. Vincent, H.F. Schaefer, J. Mol. Struct. 103 (1983) 183. [19] M. Duran, Y. Yamaguchi, H.F. Schaefer, J. Phys. Chem. 92 (1988) 3070. [20] C.I. Bayly, P. Cieplak, W.D. Cornell, P.A. Kollman, J. Phys. Chem. 97 (1993) 10269. [21] A. Jakalian, B.L. Bush, B.D. Jack, C.I. Bayly, J. Comp. Chem. 21 (2000) 132. [22] S.L. Dixon, A. van der Vaart, B. Wang, V. Gogonea, J.J. Vincent, E.N. Brothers, D. Su´arez, L.M. Westerhoff, K.M. Merz Jr., DivCon Lite, The Pennsylvania State University, University Park, PA, 2004, 16802. [23] J. Wang, W. Wang, P.A. Kollman, D.A. Case, J. Mol. Graph. Model. 25 (2006) 247. [24] R.W. Pastor, B.R. Brooks, A. Szabo, Mol. Phys. 65 (1988) 1409. [25] T. Darden, D. York, L. Pedersen, J. Chem. Phys. 103 (1993) 8577. [26] J.P. Ryckaert, G. Ciccotti, H.J.C. Berendsen, J. Comput. Phys. 23 (1977) 327. [27] H. Gohlke, D.A. Case, J. Comput. Chem. 25 (2004) 238. [28] J. Srinivasan, T.E. Cheatham III, P. Cieplak, P.A. Kollman, D.A. Case, J. Am. Chem. Soc. 120 (1998) 9401. [29] D. Sitkoff, K.A. Sharp, B. Honig, J. Phys. Chem. 98 (1994) 1978.

438

C.A. Andac et al. / International Journal of Biological Macromolecules 41 (2007) 430–438

[30] M.K. Gilson, K.A. Sharp, B. Honig, J. Comput. Chem. 9 (1987) 327. [31] J. Weiser, P.S. Shenkin, W.C. Still, J. Comput. Chem. 20 (1999) 217. [32] D.A. Case, in: M.F. Thorpe (Ed.), Rigidity Theory and Applications, Plenum, 1999, pp. 329–344. [33] W. Rocchia, E. Alexov, B. Honig, J. Phys. Chem. B 105 (2001) 6507.

[34] W. Rocchia, S. Sridharan, A. Nicholls, E. Alexov, A. Chiabrera, B. Honig, J. Comp. Chem. 23 (2002) 128. [35] E.F. Pettersen, T.D. Goddard, C.C. Huang, G.S. Couch, D.M. Greenblatt, E.C. Meng, T.E. Ferrin, J. Comput. Chem. 25 (2004) 1605. [36] S.B. Singh, P.A. Kollman, J. Am. Chem. Soc. 121 (1999) 3267. [37] G. Finette, Q. Mao, M. Hearn, J. Chromatogr. A 763 (1997) 71.