Poisson-Boltzmann methods for biomolecular electrostatics

Interactions of Cardiac Proteins with Plasma-Synthesized Polypyrrole (PSPy) to Improve Adult Cardiomyocytes Culture

Plasma-Synthesized Polypyrrole (PSPy) has been reported as a biomaterial suitable for cell growth in vitro and in vivo. An experimental duplicate was carried out that showed the growth of cardiomyocytes with PSPy, following a protocol previously reported by the working group. The cardiomyocytes cultured with the biomaterial retained their native morphological characteristics, a fundamental key to improving cardiac cell therapy procedures. Such observations motivated us to investigate the molecular characteristics of the biomaterial and the type of interactions that could be occurring (mainly electrostatic, hydrogen bonds, and non-polar). Additionally, PSPy has been studied to establish the probable mechanisms of action of the biomaterial, in particular, its action on a group of cell membrane proteins, integrins, which we know participate in the adhesion of cells to the extracellular matrix, in adhesion between cells and as bidirectional signal transducer mechanisms. In this work, we carried out studies of the interactions established between cardiac integrins α2β1 and α5β1 with different PSPy models by molecular docking studies and binding free energies (ΔGb) calculations. The models based on a previously reported PSPy molecule have three variable terminal chemical groups, with the purpose of exploring the differences in the type of interaction that will be established by modifying the position of an amino (-NH2), a hydroxyl (-OH), and a nitrile (C≡N) in (fixed) groups, as well as the length of the terminal chains (a long/short -NH2). A model with short chains for the -OH and -NH2 (lateral) group was the model with the best interactions with cardiac integrins. We experimentally verified the direct interaction of cardiomyocytes with the PSPy biomaterial observed in rat primary cultures, allowing us to validate the favorable interactions predicted by the computational analysis.

Coupling finite and boundary element methods to solve the Poisson-Boltzmann equation for electrostatics in molecular solvation

The Poisson-Boltzmann equation is widely used to model electrostatics in molecular systems. Available software packages solve it using finite difference, finite element, and boundary element methods, where the latter is attractive due to the accurate representation of the molecular surface and partial charges, and exact enforcement of the boundary conditions at infinity. However, the boundary element method is limited to linear equations and piecewise constant variations of the material properties. In this work, we present a scheme that couples finite and boundary elements for the linearised Poisson-Boltzmann equation, where the finite element method is applied in a confined solute region and the boundary element method in the external solvent region. As a proof-of-concept exercise, we use the simplest methods available: Johnson-Nédélec coupling with mass matrix and diagonal preconditioning, implemented using the Bempp-cl and FEniCSx libraries via their Python interfaces. We showcase our implementation by computing the polar component of the solvation free energy of a set of molecules using a constant and a Gaussian-varying permittivity. As validation, we compare against well-established finite difference solvers for an extensive binding energy data set, and with the finite difference code APBS (to 0.5%) for Gaussian permittivities. We also show scaling results from protein G B1 (955 atoms) up to immunoglobulin G (20,148 atoms). For small problems, the coupled method was efficient, outperforming a purely boundary integral approach. For Gaussian-varying permittivities, which are beyond the applicability of boundary elements alone, we were able to run medium to large-sized problems on a single workstation. The development of better preconditioning techniques and the use of distributed memory parallelism for larger systems remains an area for future work. We hope this work will serve as inspiration for future developments that consider space-varying field parameters, and mixed linear-nonlinear schemes for molecular electrostatics with implicit solvent models.

Bridging Eulerian and Lagrangian Poisson-Boltzmann solvers by ESES

The Poisson-Boltzmann (PB) model is a widely used electrostatic model for biomolecular solvation analysis. Formulated as an elliptic interface problem, the PB model can be numerically solved on either Eulerian meshes using finite difference/finite element methods or Lagrangian meshes using boundary element methods. Molecular surface generators, which produce the discretized dielectric interfaces between solutes and solvents, are critical factors in determining the accuracy and efficiency of the PB solvers. In this work, we investigate the utility of the Eulerian Solvent Excluded Surface (ESES) software for rendering conjugated Eulerian and Lagrangian surface representations, which enables us to numerically validate and compare the quality of Eulerian PB solvers, such as the MIBPB solver, and the Lagrangian PB solvers, such as the TABI-PB solver. Furthermore, with the ESES software and its associated PB solvers, we are able to numerically validate an interesting and useful but often neglected source-target symmetric property associated with the linearized PB model.

Variable domain mutational analysis to probe the molecular mechanisms of high viscosity of an IgG1 antibody

Subcutaneous injection is the preferred route of administration for many antibody therapeutics for reasons that include its speed and convenience. However, the small volume limit (typically ≤ 2 mL) for subcutaneous delivery often necessitates antibody formulations at high concentrations (commonly ≥100 mg/mL), which may lead to physicochemical problems. For example, antibodies with large hydrophobic or charged patches can be prone to self-interaction giving rise to high viscosity. Here, we combined X-ray crystallography with computational modeling to predict regions of an anti-glucagon receptor (GCGR) IgG1 antibody prone to self-interaction. An extensive mutational analysis was undertaken of the complementarity-determining region residues residing in hydrophobic surface patches predicted by spatial aggregation propensity, in conjunction with residue-level solvent accessibility, averaged over conformational ensembles from molecular dynamics simulations. Dynamic light scattering (DLS) was used as a medium throughput screen for self-interaction of ~ 200 anti-GCGR IgG1 variants. A negative correlation was found between the viscosity determined at high concentration (180 mg/mL) and the DLS interaction parameter measured at low concentration (2-10 mg/mL). Additionally, anti-GCGR variants were readily identified with reduced viscosity and antigen-binding affinity within a few fold of the parent antibody, with no identified impact on overall developability. The methods described here may be useful in the optimization of other antibodies to facilitate their therapeutic administration at high concentration.

Accurate Boundary Integral Formulations for the Calculation of Electrostatic Forces with an Implicit-Solvent Model

An accurate force calculation with the Poisson-Boltzmann equation is challenging, as it requires the electric field on the molecular surface. Here we present a calculation of the electric field on the solute-solvent interface that is exact for piecewise linear variations of the potential and analyze four different alternatives to compute the force using a boundary element method. We performed a verification exercise for two cases: the isolated and two interacting molecules. Our results suggest that the boundary element method outperforms the finite difference method, as the latter needs a much finer mesh than in solvation energy calculations to get acceptable accuracy in the force, whereas the same surface mesh as in a standard energy calculation is appropriate for the boundary element method. Among the four evaluated alternatives of force calculation, we saw that the most accurate one is based on the Maxwell stress tensor. However, for a realistic application, like the barnase-barstar complex, the approach based on variations of the energy functional, which is less accurate, gives equivalent results. This analysis is useful toward using the Poisson-Boltzmann equation for force calculations in applications where high accuracy is key, for example, to feed molecular dynamics models or to enable the study of the interaction between large molecular structures, like viruses adsorbed onto substrates.

EPISOL: A software package with expanded functions to perform 3D-RISM calculations for the solvation of chemical and biological molecules

Integral equation theory (IET) provides an effective solvation model for chemical and biological systems that balances computational efficiency and accuracy. We present a new software package, the expanded package for IET-based solvation (EPISOL), that performs 3D-reference interaction site model (3D-RISM) calculations to obtain the solvation structure and free energies of solute molecules in different solvents. In EPISOL, we have implemented 22 different closures, multiple free energy functionals, and new variations of 3D-RISM theory, including the recent hydrophobicity-induced density inhomogeneity (HI) theory for hydrophobic solutes and ion-dipole correction (IDC) theory for negatively charged solutes. To speed up the convergence and enhance the stability of the self-consistent iterations, we have introduced several numerical schemes in EPISOL, including a newly developed dynamic mixing approach. We show that these schemes have significantly reduced the failure rate of 3D-RISM calculations compared to AMBER-RISM software. EPISOL consists of both a user-friendly graphic interface and a kernel library that allows users to call its routines and adapt them to other programs. EPISOL is compatible with the force-field and coordinate files from both AMBER and GROMACS simulation packages. Moreover, EPISOL is equipped with an internal memory control to efficiently manage the use of physical memory, making it suitable for performing calculations on large biomolecules. We demonstrate that EPISOL can efficiently and accurately calculate solvation density distributions around various solute molecules (including a protein chaperone consisting of 120,715 atoms) and obtain solvent free energy for a wide range of organic compounds. We expect that EPISOL can be widely applied as a solvation model for chemical and biological systems. EPISOL is available at https://github.com/EPISOLrelease/EPISOL.

TABI-PB 2.0: An Improved Version of the Treecode-Accelerated Boundary Integral Poisson-Boltzmann Solver

This work describes TABI-PB 2.0, an improved version of the treecode-accelerated boundary integral Poisson-Boltzmann solver. The code computes the electrostatic potential on the molecular surface of a solvated biomolecule, and further processing yields the electrostatic solvation energy. The new implementation utilizes the NanoShaper surface triangulation code, node-patch boundary integral discretization, a block preconditioner, and a fast multipole method based on barycentric Lagrange interpolation and dual tree traversal. Performance-critical portions of the code were implemented on a GPU. Numerical results for protein 1A63 and two viral capsids (Zika, H1N1) demonstrate the code's accuracy and efficiency.

Omicron Spike Protein Has a Positive Electrostatic Surface That Promotes ACE2 Recognition and Antibody Escape

High transmissibility is a hallmark of the Omicron variant of SARS-CoV-2. Understanding the molecular determinants of Omicron’s transmissibility will impact development of intervention strategies. Here we map the electrostatic potential surface of the Spike protein to show that major SARS-CoV-2 variants have accumulated positive charges in solvent-exposed regions of the Spike protein, especially its ACE2-binding interface. Significantly, the Omicron Spike-ACE2 complex has complementary electrostatic surfaces. In contrast, interfaces between Omicron and neutralizing antibodies tend to have similar positively charged surfaces. Structural modeling demonstrates that the electrostatic property of Omicron’s Spike receptor binding domain (S RBD) plays a role in enhancing ACE2 recognition and destabilizing Spike-antibody complexes. Specifically, the Omicron S RBD has favorable electrostatic interaction energy with ACE2 that is 3-5 times greater than the Delta variant over a range of 20 Å, implying efficient recognition of host receptors. Computed binding affinities of six representative S RBD-antibody complexes show that Omicron can escape most antibodies targeting the ACE2-binding region of S RBD. Interestingly, a straightforward assessment of the electrostatic surfaces of 18 neutralizing antibodies correctly predicted the Omicron escape status of 80% of cases. Collectively, our structural analysis implies that Omicron S RBD interaction interfaces have been optimized to simultaneously promote access to human ACE2 receptors and evade antibodies. These findings suggest that electrostatic interactions are a major contributing factor for increased Omicron transmissibility relative to other variants.

Towards optimal boundary integral formulations of the Poisson-Boltzmann equation for molecular electrostatics

The Poisson-Boltzmann equation offers an efficient way to study electrostatics in molecular settings. Its numerical solution with the boundary element method is widely used, as the complicated molecular surface is accurately represented by the mesh, and the point charges are accounted for explicitly. In fact, there are several well-known boundary integral formulations available in the literature. This work presents a generalized expression of the boundary integral representation of the implicit solvent model, giving rise to new forms to compute the electrostatic potential. Moreover, it proposes a strategy to build efficient preconditioners for any of the resulting systems, improving the convergence of the linear solver. We perform systematic benchmarking of a set of formulations and preconditioners, focusing on the time to solution, matrix conditioning, and eigenvalue spectrum. We see that the eigenvalue clustering is a good indicator of the matrix conditioning, and show that they can be easily manipulated by scaling the preconditioner. Our results suggest that the optimal choice is problem-size dependent, where a simpler direct formulation is the fastest for small molecules, but more involved second-kind equations are better for larger problems. We also present a fast Calderón preconditioner for first-kind formulations, which shows promising behavior for future analysis. This work sets the basis towards choosing the most convenient boundary integral formulation of the Poisson-Boltzmann equation for a given problem.

MLL1 is regulated by KSHV LANA and is important for virus latency

Mixed lineage leukemia 1 (MLL1) is a histone methyltransferase. Kaposi's sarcoma-associated herpesvirus (KSHV) is a leading cause of malignancy in AIDS. KSHV latently infects tumor cells and its genome is decorated with epigenetic marks. Here, we show that KSHV latency-associated nuclear antigen (LANA) recruits MLL1 to viral DNA where it establishes H3K4me3 modifications at the extensive KSHV terminal repeat elements during primary infection. LANA interacts with MLL1 complex members, including WDR5, integrates into the MLL1 complex, and regulates MLL1 activity. We describe the 1.5-Å crystal structure of N-terminal LANA peptide complexed with MLL1 complex member WDR5, which reveals a potential regulatory mechanism. Disruption of MLL1 expression rendered KSHV latency establishment highly deficient. This deficiency was rescued by MLL1 but not by catalytically inactive MLL1. Therefore, MLL1 is LANA regulable and exerts a central role in virus infection. These results suggest broad potential for MLL1 regulation, including by non-host factors.

Genetic, cellular, and structural characterization of the membrane potential-dependent cell-penetrating peptide translocation pore

Cell-penetrating peptides (CPPs) allow intracellular delivery of bioactive cargo molecules. The mechanisms allowing CPPs to enter cells are ill-defined. Using a CRISPR/Cas9-based screening, we discovered that KCNQ5, KCNN4, and KCNK5 potassium channels positively modulate cationic CPP direct translocation into cells by decreasing the transmembrane potential (V_m). These findings provide the first unbiased genetic validation of the role of V_m in CPP translocation in cells. In silico modeling and live cell experiments indicate that CPPs, by bringing positive charges on the outer surface of the plasma membrane, decrease the V_m to very low values (–150 mV or less), a situation we have coined megapolarization that then triggers formation of water pores used by CPPs to enter cells. Megapolarization lowers the free energy barrier associated with CPP membrane translocation. Using dyes of varying dimensions in CPP co-entry experiments, the diameter of the water pores in living cells was estimated to be 2 (–5) nm, in accordance with the structural characteristics of the pores predicted by in silico modeling. Pharmacological manipulation to lower transmembrane potential boosted CPP cellular internalization in zebrafish and mouse models. Besides identifying the first proteins that regulate CPP translocation, this work characterized key mechanistic steps used by CPPs to cross cellular membranes. This opens the ground for strategies aimed at improving the ability of cells to capture CPP-linked cargos in vitro and in vivo.

Before a drug can have its desired effect, it must reach its target tissue or organ, and enter its cells. This is not easy because cells are surrounded by the plasma membrane, a fat-based barrier that separates the cell from its external environment. The plasma membrane contains proteins that act as channels, shuttling specific molecules in and out of the cell, and it also holds charge, with its inside surface being more negatively charged than its outside surface.

Cell-penetrating peptides are short sequences of amino acids (the building blocks that form proteins) that carry positive charges. These positive charges allow them to cross the membrane easily, but it is not well understood how.

To find out how cell-penetrating peptides cross the membrane, Trofimenko et al. attached them to dyes of different sizes. This revealed that the cell-penetrating peptides enter the cell through temporary holes called water pores, which measure about two nanometres across. The water pores form when the membrane becomes ‘megapolarized’, this is, when the difference in charge between the inside and the outside of the membrane becomes greater than normal. This can happen when the negative charge on the inside surface or the positive charge on the outer surface of the membrane increase. Megapolarization depends on potassium channels, which transport positive potassium ions outside the cell, making the outside of the membrane positive. When cell-penetrating peptides arrive at the outer surface of the cell near potassium channels, they make it even more positive. This increases the charge difference between the inside and the outside of the cell, allowing water pores to form. Once the peptides pass through the pores, the charge difference between the inside and the outside of the cell membrane dissipates, and the pores collapse.

Drug developers are experimenting with attaching cell-penetrating peptides to drugs to help them get inside their target cells. Currently there are several experimental medications of this kind in clinical trials. Understanding how these peptides gain entry, and what size of molecule they could carry with them, provides solid ground for further drug development.

Quantifying the Separation of Positive and Negative Areas in Electrostatic Potential for Predicting Feasibility of Ammonium Sulfate for Protein Crystallization

Ammonium sulfate (AS) and poly(ethylene glycol) (PEG) are the most popular precipitants in protein crystallization. Some proteins are preferably crystallized by AS, while some are by PEG. The electrostatic potential is related to the preference of the precipitant agents. The iso-surfaces of the electrostatic potentials for the AS-crystallized proteins display a common shape and a distinct separation between the positive and negative areas. In contrast, the PEG-crystallized proteins show unclear positive and negative separation. In this work, we propose schemes to quantitatively evaluate the separation for predicting which precipitant is favorable for crystal growth between AS or PEG. Three methods were attempted to quantify the amplitude of the separation, separation distance, dipole moment, and shape regularity. The positive and negative areas are approximated to the spherical potentials caused by point charges. The first method is a measurement of the distance between the positive and negative point charges. The second one is an assessment including the quantity of electric charge into the distance. The last one is an approach monitoring the clarity of the positive and negative separation. The average value for 25 kinds of AS-preferring proteins was higher than that for the PEG-preferring ones in all three methods. Therefore, every method can distinguish the proteins preferring AS for crystal growth from those preferring PEG. These methods require an iso-surface of the electrostatic potential depicted at a certain contouring value. The shape of the iso-surface depends on the contouring value. The dependency on contour was examined by depicting the iso-surfaces of electrostatic potential with three values at ±0.8, ±0.5, and ±0.2 kT/e. While reducing the contouring value leads to the increase in separation distance and the decrease in shape regularity, dipole moment is independent of the alteration of contouring value. While the AS-preferring proteins are distinguishable from the PEG-preferring ones in any contouring values, the iso-surface at ±0.5 kT/e seems adequate for regular use. The dipole moment assessment is feasible for the choice of potent precipitants for crystal growth in experiments.

Comparison of the MSMS and NanoShaper molecular surface triangulation codes in the TABI Poisson-Boltzmann solver

The Poisson-Boltzmann (PB) implicit solvent model is a popular framework for studying the electrostatics of solvated biomolecules. In this model the dielectric interface between the biomolecule and solvent is often taken to be the molecular surface or solvent-excluded surface (SES), and the quality of the SES triangulation is critical in boundary element simulations of the model. This work compares the performance of the MSMS and NanoShaper surface triangulation codes for a set of 38 biomolecules. While MSMS produces triangles of exceedingly small area and large aspect ratio, the two codes yield comparable values for the SES surface area and electrostatic solvation energy, where the latter calculations were performed using the treecode-accelerated boundary integral (TABI) PB solver. However we found that NanoShaper is computationally more efficient and reliable than MSMS, especially when parameters are set to produce highly resolved triangulations.

Efficient mesh refinement for the Poisson-Boltzmann equation with boundary elements

The Poisson-Boltzmann equation is a widely used model to study electrostatics in molecular solvation. Its numerical solution using a boundary integral formulation requires a mesh on the molecular surface only, yielding accurate representations of the solute, which is usually a complicated geometry. Here, we utilize adjoint-based analyses to form two goal-oriented error estimates that allow us to determine the contribution of each discretization element (panel) to the numerical error in the solvation free energy. This information is useful to identify high-error panels to then refine them adaptively to find optimal surface meshes. We present results for spheres and real molecular geometries, and see that elements with large error tend to be in regions where there is a high electrostatic potential. We also find that even though both estimates predict different total errors, they have similar performance as part of an adaptive mesh refinement scheme. Our test cases suggest that the adaptive mesh refinement scheme is very effective, as we are able to reduce the error one order of magnitude by increasing the mesh size less than 20% and come out to be more efficient than uniform refinement when computing error estimations. This result sets the basis toward efficient automatic mesh refinement schemes that produce optimal meshes for solvation energy calculations.

Coupling Monte Carlo, Variational Implicit Solvation, and Binary Level-Set for Simulations of Biomolecular Binding

We develop a hybrid approach that combines the Monte Carlo (MC) method, a variational implicit-solvent model (VISM), and a binary level-set method for the simulation of biomolecular binding in an aqueous solvent. The solvation free energy for the biomolecular complex is estimated by minimizing the VISM free-energy functional of all possible solute-solvent interfaces that are used as dielectric boundaries. This functional consists of the solute volumetric, solute-solvent interfacial, solute-solvent van der Waals interaction, and electrostatic free energy. A technique of shifting the dielectric boundary is used to accurately predict the electrostatic part of the solvation free energy. Minimizing such a functional in each MC move is made possible by our new and fast binary level-set method. This method is based on the approximation of surface area by the convolution of an indicator function with a compactly supported kernel and is implemented by simple flips of numerical grid cells locally around the solute-solvent interface. We apply our approach to the p53-MDM2 system for which the two molecules are approximated by rigid bodies. Our efficient approach captures some of the poses before the final bound state. All-atom molecular dynamics simulations with most of such poses quickly reach the final bound state. Our work is a new step toward realistic simulations of biomolecular interactions. With further improvement of coarse graining and MC sampling, and combined with other models, our hybrid approach can be used to study the free-energy landscape and kinetic pathways of ligand binding to proteins.

Computing Protein pKas Using the TABI Poisson–Boltzmann Solver

A common approach to computing protein pKas uses a continuum dielectric model in which the protein is a low dielectric medium with embedded atomic point charges, the solvent is a high dielectric medium with a Boltzmann distribution of ionic charges, and the pKa is related to the electrostatic free energy which is obtained by solving the Poisson–Boltzmann equation. Starting from the model pKa for a titrating residue, the method obtains the intrinsic pKa and then computes the protonation probability for a given pH including site–site interactions. This approach assumes that acid dissociation does not affect protein conformation aside from adding or deleting charges at titratable sites. In this work, we demonstrate our treecode-accelerated boundary integral (TABI) solver for the relevant electrostatic calculations. The pKa computing procedure is enclosed in a convenient Python wrapper which is publicly available at the corresponding author’s website. Predicted results are compared with experimental pKas for several proteins. Among ongoing efforts to improve protein pKa calculations, the advantage of TABI is that it reduces the numerical errors in the electrostatic calculations so that attention can be focused on modeling assumptions.

Constant-pH Brownian Dynamics Simulations of a Protein near a Charged Surface

We have developed a rigid-body Brownian dynamics algorithm that allows for simulations of a globular protein suspended in an ionic solution confined by a charged planar boundary, with an explicit treatment of pH-dependent protein protonation equilibria and their couplings to the electrostatic potential of the plane. Electrostatic interactions are described within a framework of the continuum Poisson-Boltzmann model, whereas protein-plane hydrodynamic interactions are evaluated based on analytical expressions for the position- and orientation-dependent near-wall friction tensor of a spheroid. The algorithm was applied to simulate near-surface diffusion of lysozyme in solutions having pH in the range 4-10 and ionic strengths of 10 and 150 mM. As a reference, we performed Brownian dynamics simulations in which the protein is assigned a fixed, most probable protonation state, appropriate for given solution conditions and unaffected by the presence of the charged plane, and Brownian dynamics simulations in which the protein probes possible protonation states with the pH-dependent probability, but these variations are not coupled to the electric field generated by the boundary. We show that electrostatic interactions with the negatively charged plane substantially modify probabilities of different protonation states of lysozyme and shift protonation equilibria of both acidic and basic amino acid side chains toward higher pH values. Consequently, equilibrium energy distributions, equilibrium position-orientation distributions, and functions that characterize rotational dynamics, which for a protein with multiple ionization sites, such as lysozyme, in the presence of a charged obstacle are pH-dependent, are significantly affected by the approach taken to incorporate the solution pH into simulations.

The Structure of an AAV5-AAVR Complex at 2.5 Å Resolution: Implications for Cellular Entry and Immune Neutralization of AAV Gene Therapy Vectors

Adeno-Associated Virus is the leading vector for gene therapy. Although it is the vector for all in vivo gene therapies approved for clinical use by the US Food and Drug Administration, its biology is still not yet fully understood. It has been shown that different serotypes of AAV bind to their cellular receptor, AAVR, in different ways. Previously we have reported a 2.4Å structure of AAV2 bound to AAVR that shows ordered structure for only one of the two AAVR domains with which AAV2 interacts. In this study we present a 2.5Å resolution structure of AAV5 bound to AAVR. AAV5 binds to the first polycystic kidney disease (PKD) domain of AAVR that was not ordered in the AAV2 structure. Interactions of AAV5 with AAVR are analyzed in detail, and the implications for AAV2 binding are explored through molecular modeling. Moreover, we find that binding sites for the antibodies ADK5a, ADK5b, and 3C5 on AAV5 overlap with the binding site of AAVR. These insights provide a structural foundation for development of gene therapy agents to better evade immune neutralization without disrupting cellular entry.

pKa Calculations with the Polarizable Drude Force Field and Poisson-Boltzmann Solvation Model

Electronic polarization effects have been suggested to play an important role in proton binding to titratable residues in proteins. In this work, we describe a new computational method for pK_a calculations, using Monte Carlo (MC) simulations to sample protein protonation states with the Drude polarizable force field and Poisson-Boltzmann (PB) continuum electrostatic solvent model. While the most populated protonation states at the selected pH, corresponding to residues that are half-protonated at that pH, are sampled using the exact relative free energies computed with Drude particles optimized in the field of the PB implicit solvation model, we introduce an approximation for the protein polarization of low-populated protonation states to reduce the computational cost. The highly populated protonation states used to compute the polarization and pK_a's are then iteratively improved until convergence. It is shown that for lysozyme, when considering 9 of the 18 titratable residues, the new method converged within two iterations with computed pK_a's differing only by 0.02 pH units from pK_a's estimated with the exact approach. Application of the method to predict pK_a's of 94 titratable side chains in 8 proteins shows the Drude-PB model to produce physically more correct results as compared to the additive CHARMM36 (C36) force field (FF). With a dielectric constant of two assigned to the protein interior the Root Mean Square (RMS) deviation between computed and experimental pK_a's is 2.07 and 3.19 pH units with the Drude and C36 models, respectively, and the RMS deviation using the Drude-PB model is relatively insensitive to the choice of the internal dielectric constant in contrast to the additive C36 model. At the higher internal dielectric constant of 20, pK_a's computed with the additive C36 model converge to the results obtained with the Drude polarizable force field, indicating the need to artificially overestimate electrostatic screening in a nonphysical way with the additive FF. In addition, inclusion of both syn and anti orientations of the proton in the neutral state of acidic groups is shown to yield improved agreement with experiment. The present work, which is the first example of the use of a polarizable model for the prediction of pK_a's in proteins, shows that the use of a polarizable model represents a more physically correct model for the treatment of electrostatic contributions to pK_a shifts in proteins.

Distinct 3-disulfide-bonded isomers of tridegin differentially inhibit coagulation factor XIIIa: The influence of structural stability on bioactivity

Tridegin is a 66mer cysteine-rich coagulation factor XIIIa (FXI-IIa) inhibitor from the giant amazon leech Haementeria ghilianii of yet unknown disulfide connectivity. This study covers the structural and functional characterization of five different 3-disulfide-bonded tridegin isomers. In addition to three previously identified isomers, one isomer containing the inhibitory cystine knot (ICK, knottin) motif, and one isomer with the leech antihemostatic protein (LAP) motif were synthesized in a regioselective manner. A fluorogenic enzyme activity assay revealed a positive correlation between the constriction of conformational flexibility in the N-terminal part of the peptide and the inhibitory potential towards FXI-IIa with clear differences between the isomers. This observation was supported by molecular dynamics (MD) simulations and subsequent molecular docking studies. The presented results provide detailed structure-activity relationship studies of different tridegin disulfide isomers towards FXI-IIa and reveal insights into the possibly existing native linkage compared to non-native disulfide tridegin species.

A size-modified poisson-boltzmann ion channel model in a solvent of multiple ionic species: Application to voltage-dependent anion channel

We present a new size-modified Poisson-Boltzmann ion channel (SMPBIC) model and use it to calculate the electrostatic potential, ionic concentrations, and electrostatic solvation free energy for a voltage-dependent anion channel (VDAC) on a biological membrane in a solution mixture of multiple ionic species. In particular, the new SMPBIC model adopts a membrane surface charge density and a natural Neumann boundary condition to reflect the charge effect of the membrane on the electrostatics of VDAC. To avoid the singularity difficulties caused by the atomic charges of VDAC, the new SMPBIC model is split into three submodels such that the solution of one of the submodels is obtained analytically and contains all the singularity points of the SMPBIC model. The other two submodels are then solved numerically much more efficiently than the original SMPBIC model. As an application of this SMPBIC submodel partitioning scheme, we derive a new formula for computing the electrostatic solvation free energy. Numerical results for a human VDAC isoform 1 (hVDAC1) in three different salt solutions, each with up to five different ionic species, confirm the significant effects of membrane surface charges on both the electrostatics and ionic concentrations. The results also show that the new SMPBIC model can describe well the anion selectivity property of hVDAC1, and that the new electrostatic solvation free energy formula can significantly improve the accuracy of the currently used formula. © 2019 Wiley Periodicals, Inc.

Short solvent model for ion correlations and hydrophobic association

Coulomb interactions play a major role in determining the thermodynamics, structure, and dynamics of condensed-phase systems, but often present significant challenges. Computer simulations usually use periodic boundary conditions to minimize corrections from finite cell boundaries but the long range of the Coulomb interactions generates significant contributions from distant periodic images of the simulation cell, usually calculated by Ewald sum techniques. This can add significant overhead to computer simulations and hampers the development of intuitive local pictures and simple analytic theory. In this paper, we present a general framework based on local molecular field theory to accurately determine the contributions from long-ranged Coulomb interactions to the potential of mean force between ionic or apolar hydrophobic solutes in dilute aqueous solutions described by standard classical point charge water models. The simplest approximation leads to a short solvent (SS) model, with truncated solvent-solvent and solute-solvent Coulomb interactions and long-ranged but screened Coulomb interactions only between charged solutes. The SS model accurately describes the interplay between strong short-ranged solute core interactions, local hydrogen-bond configurations, and long-ranged dielectric screening of distant charges, competing effects that are difficult to capture in standard implicit solvent models.

DelPhi Suite: New Developments and Review of Functionalities

Electrostatic potential, energies, and forces affect virtually any process in molecular biology, however, computing these quantities is a difficult task due to irregularly shaped macromolecules and the presence of water. Here, we report a new edition of the popular software package DelPhi along with describing its functionalities. The new DelPhi is a C++ object-oriented package supporting various levels of multiprocessing and memory distribution. It is demonstrated that multiprocessing results in significant improvement of computational time. Furthermore, for computations requiring large grid size (large macromolecular assemblages), the approach of memory distribution is shown to reduce the requirement of RAM and thus permitting large-scale modeling to be done on Linux clusters with moderate architecture. The new release comes with new features, whose functionalities and applications are described as well. © 2019 The Authors. Journal of Computational Chemistry published by Wiley Periodicals, Inc.

Sensor I Regulated ATPase Activity of FleQ Is Essential for Motility to Biofilm Transition in Pseudomonas aeruginosa

Members of the AAA+ (ATPase associated with various cellular activities) family of ATPases couple chemical energy derived from ATP hydrolysis for generation of mechanical force, resulting in conformational changes. The hydrolysis is brought about by highly conserved domains and motifs. The sensor I motif is critical for sensing and hydrolysis of the nucleotide. Pseudomonas aeruginosa FleQ is an ATPase that is a positive regulator of flagellar gene expression. We have determined the crystal structures of the ATPase domain of wild-type FleQ and sensor I mutants H287N and H287A in complex with ATPγS and Mg²⁺ to 2.4, 1.95, and 2.25 Å resolution, respectively. The structural data highlight the role of sensor I in regulating the ATPase activity. The in vitro and in vivo data demonstrate that the moderate ATPase activity of FleQ due to the presence of histidine in sensor I is essential for maintaining the monotrichous phenotype and for the rapid motility to biofilm transition.

A Boundary-Integral Approach for the Poisson-Boltzmann Equation with Polarizable Force Fields

Implicit-solvent models are widely used to study the electrostatics in dissolved biomolecules, which are parameterized using force fields. Standard force fields treat the charge distribution with point charges; however, other force fields have emerged which offer a more realistic description by considering polarizability. In this work, we present the implementation of the polarizable and multipolar force field atomic multipole optimized energetics for biomolecular applications (AMOEBA), in the boundary integral Poisson-Boltzmann solver PyGBe. Previous work from other researchers coupled AMOEBA with the finite-difference solver APBS, and found difficulties to effectively transfer the multipolar charge description to the mesh. A boundary integral formulation treats the charge distribution analytically, overlooking such limitations. This becomes particularly important in simulations that need high accuracy, for example, when the quantity of interest is the difference between solvation energies obtained from separate calculations, like happens for binding energy. We present verification and validation results of our software, compare it with the implementation on APBS, and assess the efficiency of AMOEBA and classical point-charge force fields in a Poisson-Boltzmann solver. We found that a boundary integral approach performs similarly to a volumetric method on CPU. Also, we present a GPU implementation of our solver. Moreover, with a boundary element method, the mesh density to correctly resolve the electrostatic potential is the same for standard point-charge and multipolar force fields. Finally, we saw that for binding energy calculations, a boundary integral approach presents more consistent results than a finite difference approximation for multipolar force fields. © 2019 Wiley Periodicals, Inc.

Effect of Conformational Diversity on the Bioactivity of µ-Conotoxin PIIIA Disulfide Isomers

Cyclic µ-conotoxin PIIIA, a potent blocker of skeletal muscle voltage-gated sodium channel NaV1.4, is a 22mer peptide stabilized by three disulfide bonds. Combining electrophysiological measurements with molecular docking and dynamic simulations based on NMR solution structures, we investigated the 15 possible 3-disulfide-bonded isomers of µ-PIIIA to relate their blocking activity at NaV1.4 to their disulfide connectivity. In addition, three µ-PIIIA mutants derived from the native disulfide isomer, in which one of the disulfide bonds was omitted (C4-16, C5-C21, C11-C22), were generated using a targeted protecting group strategy and tested using the aforementioned methods. The 3-disulfide-bonded isomers had a range of different conformational stabilities, with highly unstructured, flexible conformations with low or no channel-blocking activity, while more constrained molecules preserved 30% to 50% of the native isomer's activity. This emphasizes the importance and direct link between correct fold and function. The elimination of one disulfide bond resulted in a significant loss of blocking activity at NaV1.4, highlighting the importance of the 3-disulfide-bonded architecture for µ-PIIIA. µ-PIIIA bioactivity is governed by a subtle interplay between an optimally folded structure resulting from a specific disulfide connectivity and the electrostatic potential of the conformational ensemble.

Fludarabine-Specific Molecular Interactions with Maltose-Modified Poly(propyleneimine) Dendrimer Enable Effective Cell Entry of the Active Drug Form: Comparison with Clofarabine

Fludarabine is an anticancer antimetabolite essential for modern chemotherapy, but its efficacy is limited due to the complex pharmacokinetics. We demonstrated the potential use of maltose-modified poly(propyleneimine) dendrimer as drug delivery agent to improve the efficiency of therapy with fludarabine. In this study, we elaborated a novel synthesis technique for radioactively labeled fludarabine triphosphate to prove for the first time the direct ability of nucleotide-glycodendrimer complex to enter and kill leukemic cells, without the involvement of membrane nucleoside transporters and intracellular kinases. This will potentially allow to bypass the most common drug resistance mechanisms observed in the clinical setting. Further, we applied surface plasmon resonance and molecular modeling to elucidate the properties of the drug-dendrimer complexes. We showed that clofarabine, a more toxic nucleoside analogue drug, is characterized by significantly different molecular interactions with poly(propyleneimine) dendrimers than fludarabine, leading to different cellular outcomes (decreased rather than increased treatment efficiency). The most probable mechanistic explanation of uniquely dendrimer-enhanced fludarabine toxicity points to a crucial role of both an alternative cellular uptake pathway and the avoidance of intracellular phosphorylation of nucleoside drug form.

The Role of Tertiary Structure in MicroRNA Target Recognition

Translational repression and degradation of transcripts by microRNAs (miRNAs) is mediated by a ribonucleoprotein complex called the miRNA-induced silencing complex (miRISC, or RISC). Advances in experimental determination of RISC structures have enabled detailed analysis and modeling of known miRNA targets, yet a full appreciation of the structural factors influencing target recognition remains a challenge, primarily because target recognition involves a combination of RNA-RNA and RNA-protein interactions that can vary greatly among different miRNA-target pairs. In this chapter, we review progress toward understanding the role of tertiary structure in miRNA target recognition using computational approaches to assemble RISC complexes at known targets and physics-based methods for computing target interactions. Using this framework to examine RISC structures and dynamics, we describe how the conformational flexibility of Argonautes plays an important role in accommodating the diversity of miRNA-target duplexes formed at canonical and noncanonical target sites. We then discuss applications of tertiary structure-based approaches to emerging topics, including the structural effects of SNPs in miRNA targets and cooperative interactions involving Argonaute-Argonaute complexes. We conclude by assessing the prospects for genome-scale modeling of RISC structures and modeling of higher-order Argonaute complexes associated with miRNA biogenesis, mRNA regulation, and other functions.

Dynamics of metal ions around an RNA molecule

Noncoding RNA molecules take part in many biological processes, while metal ions play crucial roles in helping RNAs to perform their functions. However, the statics and dynamics of these metal ions around RNA molecules are still not well understood. In this work, we report a detailed molecular dynamics study of the type-I preQ_{1}-bound riboswitch aptamer domain (PRAD) at different ionic conditions (K^{+}, Na^{+}, and Mg^{2+}). The results show that the structural properties and flexibility of the PRAD molecule greatly influence the distributions and dynamics of metal ions around it. Simultaneously, Na^{+} ions show a stronger competitiveness with Mg^{2+} ions than K^{+} ions, and the three types of metal ions have different modes of interaction with the RNA molecule. Furthermore, we have also investigated specific binding sites of metal ions on the PRAD molecule and found that the dynamics and hydration structures of metal ions located at the ion-binding sites were obviously affected by the RNA structure near these ion-binding sites. These results may be useful to understand the role of the metal ions in noncoding RNA functions.

A simple model for electrical charge in globular macromolecules and linear polyelectrolytes in solution

We present a model for calculating the net and effective electrical charge of globular macromolecules and linear polyelectrolytes such as proteins and DNA, given the concentration of monovalent salt and pH in solution. The calculation is based on a numerical solution of the non-linear Poisson-Boltzmann equation using a finite element discretized continuum approach. The model simultaneously addresses the phenomena of charge regulation and renormalization, both of which underpin the electrostatics of biomolecules in solution. We show that while charge regulation addresses the true electrical charge of a molecule arising from the acid-base equilibria of its ionizable groups, charge renormalization finds relevance in the context of a molecule’s interaction with another charged entity. Writing this electrostatic interaction free energy in terms of a local electrical potential, we obtain an “interaction charge” for the molecule which we demonstrate agrees closely with the “effective charge” discussed in charge renormalization and counterion-condensation theories. The predictions of this model agree well with direct high-precision measurements of effective electrical charge of polyelectrolytes such as nucleic acids and disordered proteins in solution, without tunable parameters. Including the effective interior dielectric constant for compactly folded molecules as a tunable parameter, the model captures measurements of effective charge as well as published trends of pKa shifts in globular proteins. Our results suggest a straightforward general framework to model electrostatics in biomolecules in solution. In offering a platform that directly links theory and experiment, these calculations could foster a systematic understanding of the interrelationship between molecular 3D structure and conformation, electrical charge and electrostatic interactions in solution. The model could find particular relevance in situations where molecular crystal structures are not available or rapid, reliable predictions are desired.

A small molecule activator of SIRT3 promotes deacetylation and activation of manganese superoxide dismutase

The modulation of protein acetylation network is a promising strategy for life span extension and disease treatment (Sabari et al., 2016; Giblin et al., 2014) [1,2]. A variety of small molecules have been developed to target deacetylases, but extremely few of these molecules are capable of activating the mitochondrial NAD-dependent deacetylase sirtuin-3 (SIRT3) (Gertz and Steegborn, 2016; Scholz et al., 2015) [3,4]. Manganese superoxide dismutase (MnSOD) is the major superoxide scavenger in mitochondria, whose activity is regulated by SIRT3-mediated deacetylation, particularly at the Lys68 site (Chen et al., 2011) [5]. To investigate the influence of Lys68 acetylation on MnSOD activity, we produced a mutant MnSOD protein-bearing N-acetyllysine (AcK) at its Lys68 position through the genetic code expansion approach. We solved the crystal structure of this acetylated MnSOD (MnSODK68AcK), thus revealing the structural and electrostatic basis for the significant activity decrease upon Lys68 acetylation. On the basis of an assay we developed for the SIRT3-mediated deacetylation of MnSODK68AcK, we identified a novel SIRT3 activator, 7-hydroxy-3-(4'-methoxyphenyl) coumarin (C12), which binds to SIRT3 with high affinity and can promote the deacetylation and activation of MnSOD. C12 adds to the current repertoire of extremely few SIRT3 activators, which are potentially valuable for treating a wide array of diseases via modulating the cellular acetylome.

Free energy landscape of siRNA-polycation complexation: Elucidating the effect of molecular geometry, polymer flexibility, and charge neutralization

The success of medical threatments with DNA and silencing interference RNA is strongly related to the design of efficient delivery technologies. Cationic polymers represent an attractive strategy to serve as nucleic-acid carriers with the envisioned advantages of efficient complexation, low cost, ease of production, well-defined size, and low polydispersity index. However, the balance between efficacy and toxicity (safety) of these polymers is a challenge and in need of improvement. With the aim of designing more effective polycationic-based gene carriers, many parameters such as carrier morphology, size, molecular weight, surface chemistry, and flexibility/rigidity ratio need to be taken into consideration. In the present work, the binding mechanism of three cationic polymers (polyarginine, polylysine and polyethyleneimine) to a model siRNA target is computationally investigated at the atomistic level. In order to better understand the polycationic carrier-siRNA interactions, replica exchange molecular dynamic simulations were carried out to provide an exhaustive exploration of all the possible binding sites, taking fully into account the siRNA flexibility together with the presence of explicit solvent and ions. Moreover, well-tempered metadynamics simulations were employed to elucidate how molecular geometry, polycation flexibility, and charge neutralization affect the siRNA-polycations free energy landscape in term of low-energy binding modes and unbinding free energy barriers. Significant differences among polymer binding modes have been detected, revealing the advantageous binding properties of polyarginine and polylysine compared to polyethyleneimine.

Divalent Ion-Mediated DNA-DNA Interactions: A Comparative Study of Triplex and Duplex

Ion-mediated interaction between DNAs is essential for DNA condensation, and it is generally believed that monovalent and nonspecifically binding divalent cations cannot induce the aggregation of double-stranded (ds) DNAs. Interestingly, recent experiments found that alkaline earth metal ions such as Mg²⁺ can induce the aggregation of triple-stranded (ts) DNAs, although there is still a lack of deep understanding of the surprising findings at the microscopic level. In this work, we employed all-atom dynamic simulations to directly calculate the potentials of mean force (PMFs) between tsDNAs, between dsDNAs, and between tsDNA and dsDNA in Mg²⁺ solutions. Our calculations show that the PMF between tsDNAs is apparently attractive and becomes more strongly attractive at higher [Mg²⁺], although the PMF between dsDNAs cannot become apparently attractive even at high [Mg²⁺]. Our analyses show that Mg²⁺ internally binds into grooves and externally binds to phosphate groups for both tsDNA and dsDNA, whereas the external binding of Mg²⁺ is much stronger for tsDNA. Such stronger external binding of Mg²⁺ for tsDNA favors more apparent ion-bridging between helices than for dsDNA. Furthermore, our analyses illustrate that bridging ions, as a special part of external binding ions, are tightly and positively coupled to ion-mediated attraction between DNAs.

Unravelling the Molecular Basis of High Affinity Nanobodies against HIV p24: In Vitro Functional, Structural, and in Silico Insights

Preventing the spread of infectious diseases remains an urgent priority worldwide, and this is driving the development of advanced nanotechnology to diagnose infections at the point of care. Herein, we report the creation of a library of novel nanobody capture ligands to detect p24, one of the earliest markers of HIV infection. We demonstrate that these nanobodies, one tenth the size of conventional antibodies, exhibit high sensitivity and broad specificity to global HIV-1 subtypes. Biophysical characterization indicates strong 690 pM binding constants and fast kinetic on-rates, 1 to 2 orders of magnitude better than monoclonal antibody comparators. A crystal structure of the lead nanobody and p24 was obtained and used alongside molecular dynamics simulations to elucidate the molecular basis of these enhanced performance characteristics. They indicate that binding occurs at C-terminal helices 10 and 11 of p24, a negatively charged region of p24 complemented by the positive surface of the nanobody binding interface involving CDR1, CDR2, and CDR3 loops. Our findings have broad implications on the design of novel antibodies and a wide range of advanced biomedical applications.

Heterogeneous dielectric generalized Born model with a van der Waals term provides improved association energetics of membrane-embedded transmembrane helices

The heterogeneous dielectric generalized Born (HDGB) implicit membrane formalism is extended by the addition of a van der Waals dispersion term to better describe the nonpolar components of the free energy of solvation. The new model, termed HDGBvdW, improves the energy estimates in the hydrophobic interior of the membrane, where polar and charged species are rarely found and nonpolar interactions become significant. The implicit van der Waals term for the membrane environment extends the model from Gallicchio et al. (J. Comput. Chem. 2004, 25, 479) by combining separate contributions from each of the membrane components. The HDGBvdW model is validated with a series of test cases ranging from membrane insertion and pair association profiles of amino acid side chain analogs and transmembrane helices. Overall, the HDGBvdW model leads to increased agreement with explicit membrane simulation results and experimental data. © 2016 Wiley Periodicals, Inc.

Accurate estimation of electrostatic binding energy with Poisson-Boltzmann equation solver DelPhi program

Poisson–Boltzmann (PB) model is a widely used implicit solvent approximation in biophysical modeling because of its ability to provide accurate and reliable PB electrostatic salvation free energies ([Formula: see text] as well as electrostatic binding free energy ([Formula: see text] estimations. However, a recent study has warned that the 0.5[Formula: see text]Å grid spacing which is normally adopted can produce unacceptable errors in [Formula: see text] estimation with the solvent excluded surface (SES) (Harris RC, Boschitsch AH and Fenley MO, Influence of grid spacing in Poisson–Boltzmann equation binding energy estimation, J Chem Theory Comput 19: 3677–3685, 2013). In this work, we investigate the grid dependence of the widely used PB solver DelPhi v6.2 with molecular surface (MS) for estimating both electrostatic solvation free energies and electrostatic binding free energies. Our results indicate that, for the molecular complex and components the absolute errors of [Formula: see text] are smaller than that of [Formula: see text], and grid spacing of 0.8[Formula: see text]Å with DelPhi program ensures the accuracy and reliability of [Formula: see text]; however, the accuracy of [Formula: see text] largely relies on the order of magnitude of [Formula: see text] itself rather than that of [Formula: see text] or [Formula: see text]. Our findings suggest that grid spacing of 0.5[Formula: see text]Å is enough to produce accurate [Formula: see text] for molecules whose [Formula: see text] are large, but finer grids are needed when [Formula: see text] is very small.

The Role of Correlation and Solvation in Ion Interactions with B-DNA

The ionic atmospheres around nucleic acids play important roles in biological function. Large-scale explicit solvent simulations coupled to experimental assays such as anomalous small-angle x-ray scattering can provide important insights into the structure and energetics of such atmospheres but are time- and resource intensive. In this article, we use classical density functional theory to explore the balance among ion-DNA, ion-water, and ion-ion interactions in ionic atmospheres of RbCl, SrCl2, and CoHexCl3 (cobalt hexamine chloride) around a B-form DNA molecule. The accuracy of the classical density functional theory calculations was assessed by comparison between simulated and experimental anomalous small-angle x-ray scattering curves, demonstrating that an accurate model should take into account ion-ion correlation and ion hydration forces, DNA topology, and the discrete distribution of charges on the DNA backbone. As expected, these calculations revealed significant differences among monovalent, divalent, and trivalent cation distributions around DNA. Approximately half of the DNA-bound Rb(+) ions penetrate into the minor groove of the DNA and half adsorb on the DNA backbone. The fraction of cations in the minor groove decreases for the larger Sr(2+) ions and becomes zero for CoHex(3+) ions, which all adsorb on the DNA backbone. The distribution of CoHex(3+) ions is mainly determined by Coulomb and steric interactions, while ion-correlation forces play a central role in the monovalent Rb(+) distribution and a combination of ion-correlation and hydration forces affect the Sr(2+) distribution around DNA. This does not imply that correlations in CoHex solutions are weaker or stronger than for other ions. Steric inaccessibility of the grooves to large CoHex ions leads to their binding at the DNA surface. In this binding mode, first-order electrostatic interactions (Coulomb) dominate the overall binding energy as evidenced by low sensitivity of ionic distribution to the presence or absence of second-order electrostatic correlation interactions.

Charged Substrate and Product Together Contribute Like a Nonreactive Species to the Overall Electrostatic Steering in Diffusion-Reaction Processes

The Debye-Hückel limiting law is used to study the binding kinetics of substrate-enzyme system as well as to estimate the reaction rate of a electrostatically steered diffusion-controlled reaction process. It is based on a linearized Poisson-Boltzmann model and known for its accurate predictions in dilute solutions. However, the substrate and product particles are in nonequilibrium states and are possibly charged, and their contributions to the total electrostatic field cannot be explicitly studied in the Poisson-Boltzmann model. Hence the influences of substrate and product on reaction rate coefficient were not known. In this work, we consider all the charged species, including the charged substrate, product, and mobile salt ions in a Poisson-Nernst-Planck model, and then compare the results with previous work. The results indicate that both the charged substrate and product can significantly influence the reaction rate coefficient with different behaviors under different setups of computational conditions. It is interesting to find that when substrate and product are both considered, under an overall neutral boundary condition for all the bulk charged species, the computed reaction rate kinetics recovers a similar Debye-Hückel limiting law again. This phenomenon implies that the charged product counteracts the influence of charged substrate on reaction rate coefficient. Our analysis discloses the fact that the total charge concentration of substrate and product, though in a nonequilibrium state individually, obeys an equilibrium Boltzmann distribution, and therefore contributes as a normal charged ion species to ionic strength. This explains why the Debye-Hückel limiting law still works in a considerable range of conditions even though the effects of charged substrate and product particles are not specifically and explicitly considered in the theory.

Nonlocal Poisson-Fermi model for ionic solvent

We propose a nonlocal Poisson-Fermi model for ionic solvent that includes ion size effects and polarization correlations among water molecules in the calculation of electrostatic potential. It includes the previous Poisson-Fermi models as special cases, and its solution is the convolution of a solution of the corresponding nonlocal Poisson dielectric model with a Yukawa-like kernel function. The Fermi distribution is shown to be a set of optimal ionic concentration functions in the sense of minimizing an electrostatic potential free energy. Numerical results are reported to show the difference between a Poisson-Fermi solution and a corresponding Poisson solution.

Continuum Electrostatics Approaches to Calculating pKas and Ems in Proteins

Proteins change their charge state through protonation and redox reactions as well as through binding charged ligands. The free energy of these reactions is dominated by solvation and electrostatic energies and modulated by protein conformational relaxation in response to the ionization state changes. Although computational methods for calculating these interactions can provide very powerful tools for predicting protein charge states, they include several critical approximations of which users should be aware. This chapter discusses the strengths, weaknesses, and approximations of popular computational methods for predicting charge states and understanding the underlying electrostatic interactions. The goal of this chapter is to inform users about applications and potential caveats of these methods as well as outline directions for future theoretical and computational research.