comparison example lysozyme

[Marr575]

Hello again

I’m trying to do a comparison with Pygbe and APBS. I took the example of lysozyme I took the file lys1_charges.pqr and I submit it directly to APBS and I put the next inputs to the APBS webser

amber single Debye-Huckel

ion1:+1 ion2:-1

concentration ions 0.15M Radious ions Na: 1.16 Cl: 1.67

biomolecular dielectric constant: 4 dielectric constant solvent: 80

Radius solvent 1.4 Temperature 298.15

I think are the same inputs used in Pygbe but I get to complete different results

with pygbe I got

Calculate Solvation Energy (E_solv) Calculating solvation energy for region 2, stored in E_solv[0] 4998.481481481482 of 14846 analytical integrals for phi_reac calculation in region 2 Time E_solv: 0.24305462837219238s Region 2: E_solv = -581.2948612575798 kcal/mol = -2432.137699501714 kJ/mol

Calculate Surface Energy (E_surf) Time E_surf: 0.0133514404296875s

Calculate Coulomb Energy (E_coul) Calculate Coulomb energy for region 2 Region 2: E_coul = -7858.158720994905 kcal/mol = -32878.536088642686 kJ/mol Time E_coul: 0.008977651596069336s

---

Totals: E_solv = -581.2948612575798 kcal/mol E_surf = 0 kcal/mol E_coul = -7858.158720994905 kcal/mol

Time = 27.12166738510132 s (pygbe) bash-4.2$

but with APBS I got

CALCULATION #2: MULTIGRID Setting up problem... Vpbe_ctor: Using max ion radius (1.67 A) for exclusion function Debye length: 7.92928 A Current memory usage: 385.531 MB total, 418.972 MB high water Using cubic spline charge discretization. Grid dimensions: 129 x 97 x 129 Grid spacings: 0.475 x 0.584 x 0.514 Grid lengths: 60.840 x 56.103 x 65.775 Grid center: (0.649, 21.997, 19.351) Multigrid levels: 4 Molecule ID: 1 Linearized traditional PBE Boundary conditions from focusing 2 ion species (0.150 M ionic strength): 1.160 A-radius, 1.000 e-charge, 0.150 M concentration 1.670 A-radius, -1.000 e-charge, 0.150 M concentration Solute dielectric: 4.000 Solvent dielectric: 80.000 Using "molecular" surface definition;harmonic average smoothing Solvent probe radius: 1.400 A Temperature: 298.150 K Electrostatic energies will be calculated Potential to be written to 15550611003-pot.dx Solving PDE (see io.mc for details)... Calculating energy (see io.mc for details)... Total electrostatic energy = 4.459905281715E+04 kJ/mol Calculating forces... Writing potential to 15550611003-pot-PE0.dx

CLEANING UP AND SHUTTING DOWN... Destroying force arrays. No energy arrays to destroy. Destroying multigrid structures. Destroying 1 molecules Final memory usage: 0.001 MB total, 418.972 MB high water

I would like to ask why, or what I’m doing wrong?

Regards

Manuel Ruiz

[cdcooper84]

Hi Manuel,

Are you comparing the results from pygbe with the 4.4599e4 kJ/mol of APBS? If so, those are not the right numbers to compare. We are calculating the solvation energy, which is the energetic difference of having the protein in vacuum and inside the solvent. With BEM, this difference can be computed in one calculation, however, in a volumetric approach (like APBS) you need to do two simulations: one in the solvent and one in vacuum, and subtract those two. From the text, I see you only performed the solvent calculation. Take a look at the input files of the examples in APBS under apbs/examples/solv, and you'll see they call APBS twice and subtract the results in the end.

Hope this helps!

Chris

[Marr575]

thanks for the patience of replay to my questions

mmmmm better to be sure, I indicate 1(vacuum) for the relative permittivity in APBS and I got this

CALCULATION #2: MULTIGRID Setting up problem... Vpbe_ctor: Using max ion radius (0 A) for exclusion function Debye length: 0 A Current memory usage: 501.142 MB total, 538.905 MB high water Using cubic spline charge discretization. Grid dimensions: 129 x 129 x 129 Grid spacings: 0.478 x 0.443 x 0.516 Grid lengths: 61.163 x 56.751 x 65.997 Grid center: (0.649, 21.997, 19.351) Multigrid levels: 6 Molecule ID: 1 Linearized traditional PBE Boundary conditions from focusing 0 ion species (0.000 M ionic strength): Solute dielectric: 4.000 Solvent dielectric: 1.000 Using "molecular" surface definition;harmonic average smoothing Solvent probe radius: 1.400 A Temperature: 298.150 K Electrostatic energies will be calculated Potential to be written to 15550759479-pot.dx Solving PDE (see io.mc for details)... Calculating energy (see io.mc for details)... Total electrostatic energy = 7.063923431265E+04 kJ/mol Calculating forces... Writing potential to 15550759479-pot-PE0.dx

CLEANING UP AND SHUTTING DOWN... Destroying force arrays. No energy arrays to destroy. Destroying multigrid structures. Destroying 1 molecules Final memory usage: 0.001 MB total, 538.905 MB high water

so 4.459905281715E+04 kJ/mol - 7.063923431265E+04 kJ/mol = -26040.1814955 kJ/mol

is ok?

[cdcooper84]

That sounds like right. One thing: the idea of the 'vacuum' calculation is to substract off the 'self-energy' (ie. the Coulomb component). To get the Coulomb component only, the protein and solvent dielectric constants should be the same, and you are using different dielectrics there. Try with solvent dielectric = 4.

[Marr575]

ok, need to read more, some references that you could recommend me about salvation energy and in general to understand better the topic?

Other question in your article you put some nice graphs of the protein coloured with the charges, (I think are VMD images) but how you generate the file with the charges when the protein is under the influence of the electrode? The pygbe output doesn't give me something that I know I could use in VMD