Electron flow in multiheme bacterial cytochromes is a balancing act between heme electronic interaction and redox potentials

The naturally widespread process of electron transfer from metal reducing bacteria to extracellular solid metal oxides entails unique biomolecular machinery optimized for long-range electron transport. To perform this function efficiently, microorganisms have adapted multiheme c-type cytochromes to arrange heme cofactors into wires that cooperatively span the cellular envelope, transmitting electrons along distances greater than 100 Å. Implications and opportunities for bionanotechnological device design are self-evident. However, at the molecular level, how these proteins shuttle electrons along their heme wires, navigating intraprotein intersections and interprotein interfaces efficiently, remains a mystery thus far inaccessible to experiment. To shed light on this critical topic, we carried out extensive quantum mechanics/molecular mechanics simulations to calculate stepwise heme-to-heme electron transfer rates in the recently crystallized outer membrane deca-heme cytochrome MtrF. By solving a master equation for electron hopping, we estimate an intrinsic, maximum possible electron flux through solvated MtrF of 10⁴–10⁵ s⁻¹, consistent with recently measured rates for the related multiheme protein complex MtrCAB. Intriguingly, our calculations show that the rapid electron transport through MtrF is the result of a clear correlation between heme redox potential and the strength of electronic coupling along the wire: thermodynamically uphill steps occur only between electronically well-connected stacked heme pairs. This observation suggests that the protein evolved to harbor low-potential hemes without slowing down electron flow. These findings are particularly profound in light of the apparently well-conserved staggered cross-heme wire structural motif in functionally related outer membrane proteins.Respiratory electron transfer (ET) is not restricted to the aqueous subunits and membranes inside cells but in specialized cases can also occur across the outer membrane to extracellular space. This possibility is heavily used by dissimilatory metal reducing bacteria (DMRB), which are capable of using extracellular solid metal oxides as terminal respiratory electron sinks, a process that has been suggested to proceed via direct cell-mineral contact (1), extracellular redox shuttles (2), and/or pilus-like appendages (3, 4). Although essential to the survival of the bacterium, extracellular ET also plays an important role in the biogeochemical cycling of transition metals (5–7). It is or could be exploited in a multitude of biotechnological applications ranging from mediator-less biofuel cells (8) to biological waste-to-electricity conversion (1), photocatalytic bioenergy generation (9), and even bioelectronic systems using directional electronic communication between living and nonliving systems (9).The transport of electrons from the inner membrane, where they accumulate as a result of metabolic activity, across the periplasm and outer membrane to the extracellular space relies on an efficient network of ET proteins (10). It has been known for some time that multiheme c-type cytochromes play a central role in this process. Examples of such systems include the MtrCAB and MtrFDE transmembrane complexes of the bacterial strain Shewanella oneidensis MR-1 that form a biological nanowire of 20 c-type hemes, 10 from MtrC(F) and 10 from MtrA(D), wrapped in a β barrel porin MtrB(E) (10) (Fig. 1A). MtrB(E) does not contain any hemes but is supposed to enable contact for ET between the periplasmic MtrA(D) and the outer membrane cytochrome MtrC(F). The latter is assumed to pass electrons on to extracellular substrates either directly or via redox mediators such as flavins (11).Open in a separate windowFig. 1.(A) Model of the protein complex MtrFDE. The outer membrane deca-heme cytochrome MtrF [PDB code 3PMQ (12)] is connected to a model of a periplasmic/membrane cytochrome MtrD via a porin MtrE enabling close contact between the two cytochromes across the outer membrane. (B) Heme cofactor arrangement in MtrF. Arrows denote single ET steps from heme i to heme j with rate constant k_ji (denoted exemplary for pair 1–6), as well as ET steps to/from an external electron acceptor/donor, k_i,out and k_i,in, shown exemplarily for hemes i = 10 and 5. (C) Three different heme pair motifs found in MtrF. From top to bottom: T-shaped, coplanar, and stacked. d_π orbitals involved in the coupling are also depicted.The recently published crystal structure for MtrF (12) (and indeed the first one for any deca-heme cytochrome) reveals hemes arranged side by side in a sequence clearly intended for directional electron flow. However, the arrangement is not simply a linear chain of 10 cofactors; rather, it features a peculiar “staggered-cross” formation of the 10 hemes as shown in Fig. 1B, with a central tetra-heme chain between hemes 2 and 7 and two heme-triples branching off in orthogonal directions to yield an octaheme chain between hemes 5 and 10. The relative orientation of adjacent hemes also varies in three apparent types (see depictions in Fig. 1C), with coplanar pairs within the tetra-heme chain, stacked pairs within the two heme-triples, and T-shaped connections between tetra-heme chain and heme-triples. These motifs as such are not uncommon in biological ET: the stacked heme arrangement bears some similarity with the tightly packed chlorophylls in reaction center proteins (13), and a similar T-shaped connection is found for the heme a-a₃ pair in cytochrome c oxidase (14) (Fig. S1). However, it is thus far an open question why all three heme-heme motifs are present in MtrF, and why the 10 hemes are arranged in a nonlinear cross-wired fashion.Given the central role of multiheme cytochromes like MtrF in the extracellular ET processes of DMRB, elucidating their function on a molecular level is therefore at the heart of coming to understand and possibly adapt these astonishing capabilities. However, elementary aspects of electron transport through these proteins are difficult to assess in experiment. For example, although overall potential windows of operation for the whole protein can be established (12, 15), the cofactors’ almost identical chemical environment impedes deconvolution into redox potentials of single structurally assignable hemes (12). Similarly, although the analysis of tunneling spectroscopy-derived current-voltage curves for MtrC single molecules yielded two individual redox potentials consistent with the whole protein window (16), tunneling experiments do not necessarily involve the same heme-to-heme hopping mechanism of electron transmission as expected in the native protein function.Computational methods fill this accessibility gap to provide molecular-level insight. They are not only able to elucidate properties of individual cofactors in these multiheme cytochromes but also allow for the analysis of structure-function relationships: Our recent previous study of the thermodynamics of electron transfer along hemes in MtrF (17) revealed a roughly symmetric free energy profile featuring two thermodynamic barriers of 0.2 eV, where heme redox potentials were lowered by the close proximity of negatively charged propionate side chains of adjacent hemes. Although this provides the protein with two solvent-exposed heme sites capable of spontaneously reducing the key redox shuttle flavin mononucleotide (FMN), it also raises the puzzling question of how these barriers are incorporated into the chain without affecting through-protein electron transport.To answer this question, we now present extensive simulation work on the kinetics of stepwise and overall ET in MtrF. Specifically, we carry out classical molecular dynamics (MD) and quantum mechanics/molecular mechanics (QM/MM) calculations to obtain electronic couplings for sequential ET between each of the nine adjacent heme pairs in MtrF. Combined with our previous work, this completes the set of quantities necessary to describe ET rates in the framework of nonadiabatic Marcus theory, an approach found to be appropriate to describe ET through MtrF. The simulations reveal that the rapid transport rate through MtrF is a consequence of a subtle balancing act: energetically uphill steps occur only between the tightest packed, electronically best-connected hemes. The potential slowing of ET rates by low potential hemes is thus compensated by high electronic coupling. In fact, the rates for the thermodynamically unfavorable ET steps do not fall below those for the thermodynamically reversible steps. As a result, the maximum intrinsic electron flux through the heme wire is maintained at 10⁴−10⁵ s⁻¹, just slightly higher than recently measured acceptor-limited transport rates through the multiheme protein complex MtrCAB (18). The structural similarity of MtrF with homologs UndA (19) and MtrC (20) allows us to generalize our findings to these cytochromes also, suggesting an important electron transfer strategy in nature well conserved because of its efficiency for long-range electron transport.     

Simulation of the isotropic EXAFS spectra for the S2 and S3 structures of the oxygen evolving complex in photosystem II

Most of the main features of water oxidation in photosystem II are now well understood, including the mechanism for O–O bond formation. For the intermediate S₂ and S₃ structures there is also nearly complete agreement between quantum chemical modeling and experiments. Given the present high degree of consensus for these structures, it is of high interest to go back to previous suggestions concerning what happens in the S₂–S₃ transition. Analyses of extended X-ray adsorption fine structure (EXAFS) experiments have indicated relatively large structural changes in this transition, with changes of distances sometimes larger than 0.3 Å and a change of topology. In contrast, our previous density functional theory (DFT)(B3LYP) calculations on a cluster model showed very small changes, less than 0.1 Å. It is here found that the DFT structures are also consistent with the EXAFS spectra for the S₂ and S₃ states within normal errors of DFT. The analysis suggests that there are severe problems in interpreting EXAFS spectra for these complicated systems.The knowledge of the different steps of water oxidation in photosystem II has increased rapidly the past years. After the first low-resolution X-ray structures appeared ∼10 y ago (1–3), quantum chemical studies using density functional theory (DFT) have played a major role for obtaining a mechanistic understanding. First, an O–O bond formation mechanism was suggested in 2006 (4) in which a terminally bound oxyl radical in the center of the oxygen evolving complex (OEC) was attacked by a manganese-bridging oxo group. Second, an improved structure was suggested in which, most importantly, the outer manganese was placed differently from where it was placed in the previous X-ray structures (5). This position led to an open space in the center of the OEC, which is critical for allowing the low-barrier O–O bond formation suggested earlier (4).In 2011, a major experimental breakthrough occurred when the first high-resolution X-ray structure at 1.9 Å was presented by Umena et al. (6), which essentially confirmed the quantum chemical structure of the OEC. The main difference was that Asp170 was found to bind in a bridging mode between the terminal manganese and calcium instead of only terminally to the manganese as in the quantum chemical structure. The rest of the structure is very similar, including the critical positions of the outer manganese and the oxo groups, and the ligand connections. A minor problem with the X-ray structure is that it is most probably reduced by X-ray radiation (7–9), indicating that it is unlikely to be in the S₁ state as claimed. More recently, spectroscopic studies have played a major role by confirming the most important aspects of the quantum chemical suggestions. On the basis of the new X-ray structure and old DFT(B3LYP) structure (5), using electron paramagnetic resonance (EPR), electron nuclear double resonance (ENDOR), and DFT, a detailed structure of the OEC in the S₂ state was reached (10) that agrees almost perfectly with a structure obtained independently by a DFT(B3LYP) energy minimization (11, 12) (Fig. 1). The positions of the oxo groups and the protonation states, including which ligands are water and which are hydroxides, agree, along with which manganese are Mn(III) and which are Mn(IV) at that stage. Also, the DFT(B3LYP) structure from 2009, before the high-resolution structure, is very similar (13). Two years ago, the substrate oxygen positions were suggested for S₂ using a W-band ¹⁷O ELDOR-detected NMR spectroscopy (14). The position for the slowly exchanging substrate agrees with the one suggested by the DFT studies (4, 11–13), but there is still a minor possible disagreement for the fast-exchanging substrate. Very recently, a combined experimental and theoretical study by Cox et al. (15) used EPR and ⁵⁵Mn–EDNMR spectra to suggest an S₃ structure almost identical to the structure suggested by DFT(B3LYP) 2 y ago (16) (Fig. 2), and again very similar to the one from 2009 (13). It was claimed that only this structural model fits the measured spectra.Open in a separate windowFig. 1.(Left) Previously DFT(B3LYP)-optimized structure for the S₂ state (11, 12). (Right) Structure suggested after a spectroscopic analysis (10).Open in a separate windowFig. 2.(Left) Previously DFT(B3LYP)-optimized structure for the S₃ state (12, 16). (Right) Structure suggested after a spectroscopic analysis (15).Even though the major features of water oxidation can now be claimed to be reasonably well understood, additional studies are required to sort out details of the mechanism. A puzzling observation stems from previous extended X-ray adsorption fine structure (EXAFS) studies of the S₂–S₃ transition. In the EXAFS studies by Yachandra and coworkers (17–19), three short distances of 2.7–2.8 Å were found in S₂. In another EXAFS study by Dau and coworkers (20), only two short Mn–Mn distances of 2.7 Å were suggested. Instead, two of the Mn–Mn distances were proposed to be longer than 3.0 Å. For the S₂–S₃ transition, the discrepancies were even more marked. In the studies by Yachandra and coworkers (17–19), it was concluded that there is a lengthening of one of the 2.7–2.8 Å distances to 3.00 Å. In the study by Dau and coworkers (20) it was instead suggested that there is a shortening of one of the distances, which was >3.0 Å, down to 2.7 Å, indicating a formation of an additional Mn–Mn bis-μ-oxo bridge in S₃. Perhaps the most noteworthy of the differences of the suggested S₃ distances is the one that is 2.80 Å in the Dau and coworkers (20) study, and as long as 3.0 Å in the Yachandra et al. study (17). The suggestions from both these studies give larger deviations to the DFT/spectroscopy structure than are expected from DFT, greater than 0.1 Å on some distances and a topology change. The two different EXAFS interpretations led to different proposals for the water oxidation mechanism. It should in this context be mentioned that for S₂ the raw data from the two groups are the same but not for S₃. In the present work, the EXAFS information from Dau and coworkers (20) has been used. Recently, after the theoretical and spectroscopic consensus structure of S₂ had appeared, a reanalysis of the EXAFS spectra was made by the group of Yachandra and coworkers (21); for the S₂ structure, they now find full agreement between the EXAFS analysis and the DFT/spectroscopic structure in Fig. 1. For the S₃ structure, two alternatives were given, one with essentially four equivalent distances and one where one distance is longer.It is not straightforward to compare DFT and EXAFS data, because the latter are very sensitive to the metal–ligand distances (with an accuracy of 0.01–0.02 Å), whereas DFT calculations often give ∼0.05 Å too-long metal–O bonds and even larger deviations in the metal–metal distances. Therefore, an EXAFS spectrum calculated directly on a DFT structure will be poor, and a direct comparison of DFT and EXAFS distances will also show extensive deviations. Instead, it is better to perform a combined EXAFS/DFT refinement of the EXAFS spectrum, in which the EXAFS raw data (not only the EXAFS distances) are used as a restraint in the DFT geometry optimization (22–24). Thereby, DFT will determine the general structure of the complex, whereas the EXAFS data will determine the detailed distances involving the metals. Such calculations are presented in this study for the S₂ and S₃ states, for a large DFT model used previously (11, 12), and also for a much smaller model. The results are compared with the two different experimental EXAFS analyses. The main purpose of the present study is to investigate whether the discrepancies to experiments for the computational model, concerning the structural changes in S₂ to S₃, really indicate significant differences in the structures or if they are mainly due to minor differences in bond lengths and technical differences in how the spectra are analyzed. The present study agrees with earlier ones (25, 26) in that a major problem of interpreting EXAFS spectra of complicated molecules is that it is possible to fit several different structures to the same spectrum.It should finally be emphasized that a comparison with other theoretical and experimental work is not part of the purpose of the present paper, which is instead focused on EXAFS results. However, other theoretical work on water oxidation in photosystem II has been discussed in detail in recent reviews (12, 27, 28).     

Geometric and electronic structure differences between the type 3 copper sites of the multicopper oxidases and hemocyanin/tyrosinase

The coupled binuclear “type 3” Cu sites are found in hemocyanin (Hc), tyrosinase (Tyr), and the multicopper oxidases (MCOs), such as laccase (Lc), and play vital roles in O₂ respiration. Although all type 3 Cu sites share the same ground state features, those of Hc/Tyr have very different ligand-binding properties relative to those of the MCOs. In particular, the type 3 Cu site in the MCOs (Lc^T3) is a part of the trinuclear Cu cluster, and if the third (i.e., type 2) Cu is removed, the Lc^T3 site does not react with O₂. Density functional theory calculations indicate that O₂ binding in Hc is ≈9 kcal mol⁻¹ more favorable than for Lc^T3. The difference is mostly found in the total energy difference of the deoxy states (≈7 kcal mol⁻¹), where the stabilization of deoxy Lc^T3 derives from its long equilibrium Cu–Cu distance of ≈5.5–6.5 Å, relative to ≈4.2 Å in deoxy Hc/Tyr. The O₂ binding in Hc is driven by the electrostatic destabilization of the deoxy Hc site, in which the two Cu(I) centers are kept close together by the protein for facile 2-electron reduction of O₂. Alternatively, the lack of O₂ reactivity in Lc^T3 reflects the flexibility of the active site, capable of minimizing the electrostatic repulsion of the 2 Cu(I)s. Thus, the O₂ reactivity of the MCOs is intrinsic to the trinuclear Cu cluster, leading to different O₂ intermediates as required by its function of irreversible reduction of O₂ to H₂O.     

DFT Insight into Conductive and Magnetic Properties of Heterostructures with BaTiO3 Overlayer

The ab initio calculations of a heterostructure based on the ferroelectric phase of barium titanate and dielectrics lanthanum manganese (LaMnO₃) or silicon (Si) are presented. We analyze structures of BaTiO₃/LaMnO₃ and BaTiO₃/Si interfaces, investigate magnetic properties and the impact of ferroelectric polarization. The use of ferroelectrics in the heterostructure plays a crucial role; in particular, ferroelectric polarization leads to the appearance of the conducting state at the interface and in the layers close to it. We show that defects (here, oxygen vacancies) incorporated into the system may change the electronic and magnetic properties of a system. Experimental results of magnetic susceptibility measurements for the Ba_0.8Sr_0.2TiO₃/LaMnO₃ heterostructure are also presented. It is shown that a correlation between the behavior of the ferromagnetic ordering and the resistance takes place. In addition, the ferromagnetic ordering at the interface of the heterostructure can be associated with the exchange interaction through current carriers that appear in high carrier concentration regions.     

Stilbestrols and stilbestrol derivatives: estrogenic potency and temporal relationships between estrogen receptor binding and uterine growth

This study analyzes some of the differences between the estrogenic potencies of two homo logous stilbestrols, the potent estiogen diethylstilbestrol (DES) and the weak estrogen dimethylstilbestrol (DMS). The action of these compounds and their corresponding dimethyl ethers is compared in terms of the duration of their interaction with the estrogen receptor in immature rat uterus and the time-course of responses elicited in this tissue.Dose-response curves of 3-day uterotrophic assays indicate that etherification of DMS, which is only weakly uterotrophic, converts it into a compound, DMS-(OMe)₂, that has enhanced uterotrophic activity, while etherification of the more active estrogen, DES, diminishes its potency. Only at high doses is DES-(OMe)₂ as effective a uterotrophic agent as DES. After a single injection, DMS (20 μg) and DES (10 μg) both rapidly translocate estrogen receptor from the cytoplasmic to the nuclear compartment, but while uterine weight (by 24 h) and nuclear receptor levels (by 6 h) rapidly return to control levels after DMS, they remain elevated for a more prolonged period after DES. Likewise, DMS stimulates only an early (2 h) wave of uterine deoxyglucose phosphorylation. In contrast to DMS, DMS-(OMe)₂ (20 μg) shows a gradual movement of receptor to the nucleus after 1 h, with moderate but above control levels of receptor being maintained for at least 36–48 h. This retention of nuclear receptor correlates with prolonged elevation of uterine weight (beyond 60 h) and stimulation of deoxyglucose metabolism (beyond 24 h). Likewise a high (10 μg) dose of DES-(OMe)₂ evokes a slower but more protracted elevation of nuclear receptor levels, and a more prolonged elevation of uterine weight than does DES (10 μg).The weak activity of DMS appears to be due to its short duration of interaction with receptor, and conversion to the methyl ether prolongs nuclear receptor occupancy and increases its biological potency. DES is potent because it is itself long-acting. Methylation of DES further extends its period of nuclear receptor occupancy; this increases its duration of action at high doses, but reduces its potency at low doses.     

Free energy landscape for the binding process of Huperzine A to acetylcholinesterase

Drug-target residence time (t = 1/k_off, where k_off is the dissociation rate constant) has become an important index in discovering better- or best-in-class drugs. However, little effort has been dedicated to developing computational methods that can accurately predict this kinetic parameter or related parameters, k_off and activation free energy of dissociation (). In this paper, energy landscape theory that has been developed to understand protein folding and function is extended to develop a generally applicable computational framework that is able to construct a complete ligand-target binding free energy landscape. This enables both the binding affinity and the binding kinetics to be accurately estimated. We applied this method to simulate the binding event of the anti-Alzheimer’s disease drug (−)−Huperzine A to its target acetylcholinesterase (AChE). The computational results are in excellent agreement with our concurrent experimental measurements. All of the predicted values of binding free energy and activation free energies of association and dissociation deviate from the experimental data only by less than 1 kcal/mol. The method also provides atomic resolution information for the (−)−Huperzine A binding pathway, which may be useful in designing more potent AChE inhibitors. We expect this methodology to be widely applicable to drug discovery and development.     

INAUGURAL ARTICLE by a Recently Elected Academy Member:Long-ranged contributions to solvation free energies from theory and short-ranged models

Long-standing problems associated with long-ranged electrostatic interactions have plagued theory and simulation alike. Traditional lattice sum (Ewald-like) treatments of Coulomb interactions add significant overhead to computer simulations and can produce artifacts from spurious interactions between simulation cell images. These subtle issues become particularly apparent when estimating thermodynamic quantities, such as free energies of solvation in charged and polar systems, to which long-ranged Coulomb interactions typically make a large contribution. In this paper, we develop a framework for determining very accurate solvation free energies of systems with long-ranged interactions from models that interact with purely short-ranged potentials. Our approach is generally applicable and can be combined with existing computational and theoretical techniques for estimating solvation thermodynamics. We demonstrate the utility of our approach by examining the hydration thermodynamics of hydrophobic and ionic solutes and the solvation of a large, highly charged colloid that exhibits overcharging, a complex nonlinear electrostatic phenomenon whereby counterions from the solvent effectively overscreen and locally invert the integrated charge of the solvated object.Solvation thermodynamics underlies a vast array of important processes, ranging from protein folding (1, 2) and ligand binding (3) to self-assembly at interfaces (4). Thus, understanding solvation, and driving forces rooted in solvation, has been a focus of chemistry and physics for over a century (5, 6).Quantitatively successful theories of self-solvation and solvophobic solvation in simple fluids have been developed (7–16). However, a generally useful analytic approach for solvation in complex charged and polar environments is lacking, and solvation is typically studied with computer simulations. Contributions from the long-ranged components of Coulomb interactions in periodic images of the simulation cell are typically evaluated using computationally intense Ewald and related lattice summation techniques (17). These methods generate distorted, system size-dependent interaction potentials (18) and do not scale well in massively parallel simulations (19), adding considerable computational overhead. Moreover, artifacts can arise from spurious interactions between the periodic images of solutes, as observed for proteins in water (20).The local molecular field (LMF) theory of nonuniform fluids is a promising avenue for substantially improving free energy calculations by removing many of the computational and conceptual burdens associated with long-ranged interactions (14, 21). LMF theory prescribes a way to accurately determine the structure of a full system with long-ranged intermolecular interactions in a general single particle field by studying a simpler mimic system wherein particles interact with short-ranged intermolecular interactions only. An effective field in the mimic system accounts for the averaged effects of the long-ranged “far-field” interactions in the full system.This approach is especially powerful for studying solvation in charged and polar solvents, where in the simplest case the effective field can represent the interactions between a fixed solute and the solvent. In this paper, we show that when the effective field and induced density around the solute are accurately determined by LMF theory it is very easy to integrate over the solvent structure and accurately compute the far-field contributions to the solvation free energy as well, using quantities determined solely in the short-ranged mimic system, where simulations scale linearly with system size.LMF theory presents a general conceptual framework that gives qualitative as well as quantitative insight into many other problems. Its treatment of long- and short-ranged forces makes suggestive connections to other well-established theoretical methods, such as perturbation theory for uniform simple fluids (6, 7), classical density functional theory (DFT) of nonuniform fluids (11), and the successful quasichemical approach for solvation (16). Although our focus in this paper is on the quantitative determination of the solvation free energy, many of these connections will be touched upon in our discussion here and in Supporting Information. The treatment of solvation free energies we present here can be readily generalized to determine more complex free energies, including alchemical transformations and potentials of mean force (3), and extended to more general charged and polar mixtures (21) with mobile solutes.The conceptual development of the LMF approach to solvation thermodynamics is introduced in the next section, with derivations and other technical details given in Materials and Methods, Derivation of the Far-Field Solvation Free Energy and Supporting Information. We first focus on the solvophobic solvation of a repulsive, spherical solute in a Lennard-Jones (LJ) fluid, where most of the ideas can be understood in their simplest form and the basic physics is well understood. We then turn to more challenging and experimentally relevant problems involving the length-scale transition in hydrophobic solvation of an apolar solute in water and its effect on the solvation free energies; the hydration of single ions is discussed in Supporting Information. Finally we discuss the solvation and “overcharging” of a large, highly charged colloid in an ionic fluid, a highly nontrivial process involving ion correlations (22) that is completely missed in classic mean field treatments of ionic solutions (23).     

Leveraging nonstructural data to predict structures and affinities of protein–ligand complexes

Over the past five decades, tremendous effort has been devoted to computational methods for predicting properties of ligands—i.e., molecules that bind macromolecular targets. Such methods, which are critical to rational drug design, fall into two categories: physics-based methods, which directly model ligand interactions with the target given the target’s three-dimensional (3D) structure, and ligand-based methods, which predict ligand properties given experimental measurements for similar ligands. Here, we present a rigorous statistical framework to combine these two sources of information. We develop a method to predict a ligand’s pose—the 3D structure of the ligand bound to its target—that leverages a widely available source of information: a list of other ligands that are known to bind the same target but for which no 3D structure is available. This combination of physics-based and ligand-based modeling improves pose prediction accuracy across all major families of drug targets. Using the same framework, we develop a method for virtual screening of drug candidates, which outperforms standard physics-based and ligand-based virtual screening methods. Our results suggest broad opportunities to improve prediction of various ligand properties by combining diverse sources of information through customized machine-learning approaches.

Binding of small-molecule ligands to proteins is one of the most fundamental processes in biology, and the great majority of drugs exert their effects by binding to a target protein. Predicting properties of protein–ligand interactions—including three-dimensional (3D) structures, binding affinities, binding kinetics, selectivity, and functional effects—is critical both for the rational design of effective medicines and for addressing important questions in molecular biology. A great deal of work has thus focused on the development of computational methods to predict these properties (1, 2).Such computational methods generally fall into two categories. “Physics-based” approaches use a 3D structure of the target protein and exploit an understanding of the physics of protein–ligand interactions (3). “Ligand-based” approaches use experimental measurements of a given property (e.g., affinity at a particular target) for many ligands and employ pattern-matching to predict the corresponding property for other ligands (4, 5).Can one combine these two paradigms and the orthogonal sources of information they leverage in a systematic, principled manner? This has proven challenging, particularly when making predictions for ligands substantially different from those for which experimental data are available. It is especially difficult when one wishes to predict properties different from those measured experimentally (e.g., to predict ligand properties that are difficult to determine experimentally by exploiting experimental data that is easy to collect).Here, we present a rigorous statistical framework to combine the distinct sources of information exploited by physics-based and ligand-based approaches. Using this framework, we develop ComBind, a method to improve prediction of a ligand’s binding pose at a target protein by exploiting readily available nonstructural data. We use the same framework to develop ComBindVS, a virtual-screening method that leverages both structural and nonstructural data to predict ligand binding affinities.Determining a ligand’s binding pose—the 3D coordinates of the ligand’s atoms when bound to the target—is critical for structure-based optimization of the ligand’s pharmacological properties as well as for understanding how the ligand influences its target. Medicinal chemists have long used the binding pose of a lead compound—when available—as an intuitive guide in choosing which analogs to synthesize and assay (6–9). Binding poses also serve as starting points for computational methods that predict ligand properties such as affinity and selectivity (10–15). Indeed, knowledge of a ligand’s binding pose is so advantageous that researchers in industry and academia often spend months or years to solve an experimental structure of a particular ligand in complex with a target protein.Because experimental structure determination is time consuming, expensive, and sometimes intractable, tremendous effort has been invested in the development of in silico “docking” methods for predicting ligand binding poses (16–25). These methods are physics based: given a structure of the target protein, they sample many candidate poses of a ligand and rank these poses using scoring functions that approximate the energetic favorability of each pose, typically by capturing interatomic interactions such as hydrogen bonds and van der Waals forces (Fig. 1A). Despite the development of dozens of docking software packages over the past five decades, binding pose predictions are typically correct less than half the time for ligands substantially different from those in the experimental structures used for docking (SI Appendix, Table S1).Open in a separate windowFig. 1.ComBind leverages nonstructural data to improve ligand binding pose predictions. (A) Standard docking methods take as input the chemical structure of the query ligand and the 3D structure of the target protein and predict a binding pose using a per-ligand scoring function. (B) ComBind additionally considers other ligands known to bind the target protein (whose binding poses are not known), resulting in more-accurate predictions. For clarity, hydrogen and fluorine atoms are omitted from the 3D renderings.ComBind improves binding pose prediction by exploiting a widely available type of nonstructural data: the identities of other ligands known to bind the same target (Fig. 1B). Collecting such data is typically far easier than structure determination. Indeed, such data are routinely collected in drug development campaigns and are already available in public databases such as ChEMBL for most recognized drug targets (26).How can a list of other ligands that bind to the target protein—but whose binding poses are unknown—be used to improve pose prediction? Medicinal chemists have long recognized that distinct ligands tend to bind a given protein in similar poses. Even ligands sharing no common substructure often form similar interactions with the target protein (Fig. 2A). This intuition has a sound basis in physics. For example, the energetic favorability of a protein–ligand hydrogen bond depends on the mobility of the protein atoms involved and their ability to form hydrogen bonds with water in the absence of the ligand (9). These factors contribute similarly to binding of different ligands but are difficult to predict from a static protein structure alone.Open in a separate windowFig. 2.Distinct ligands that bind to a given target protein often adopt similar binding poses and do so more frequently than predicted by a state-of-the-art per-ligand docking method. (A) Chemically distinct ligands share key interactions with the mineralocorticoid receptor (Protein Data Bank IDs: 2AA2, 5L7E, and 5MWP). (B) Across a set of 3,115 ligand pairs, interaction similarities are generally higher in pairs of correct poses than in pairs of poses ranked highly by a per-ligand scoring function. Shading depicts the per-target SEM. A.U.: arbitrary units. (C) Across a set of 690 ligand pairs with a shared substructure, the substructure tends to be placed more similarly in correct poses than in other poses ranked highly by a per-ligand scoring function (SI Appendix, Supplementary Text).We use a large set of experimentally determined structures to quantify the medicinal chemist’s intuition—in particular, to determine the probability that binding poses for different ligands will share various features. We use the results to define the ComBind scoring function, which predicts the favorability of a set of binding poses comprising one pose for each ligand known to bind the target protein. By contrast, scoring functions typically utilized by docking software assign a score to the pose of a single ligand at a time; we thus refer to them as per-ligand scoring functions. The ComBind scoring function takes into account the similarities and differences between the poses of different ligands as well as the energetic favorability of each ligand’s pose, as evaluated by a per-ligand scoring function. By using this scoring function to predict the binding poses of a set of ligands simultaneously, we can predict the pose of each ligand more accurately, even when the ligands share no common scaffold and none of their binding poses are known in advance.ComBindVS uses the same sources of information—a structure of the target protein and a set of ligands known to bind the target—for virtual screening. Here, we use the ComBind scoring function not only to predict binding poses of known binders but also to predict binding affinities of unrelated molecules.We benchmark ComBind pose prediction by comparing its results to 248 experimentally determined ligand binding poses across 30 proteins representing all major families of drug targets. ComBind improves the pose prediction accuracy of state-of-the-art docking software for all major drug target families.We benchmark ComBindVS for virtual screening using the Directory of Useful Decoys, Enhanced (DUD-E) benchmark set (27), which includes diverse protein targets. ComBindVS outperforms state-of-the-art structure-based docking and ligand-based virtual-screening methods, as well as approaches that combine the results of docking and ligand-based methods. ComBindVS yields performance improvements even when the candidate molecules are very different from the known binders, making it suitable for discovery of novel chemotypes.We also illustrate the use of ComBind to predict the previously unknown binding poses of several antipsychotics at the D₂ dopamine receptor (D₂R), an important drug target for which experimental structure determination has proven difficult. We validate ComBind’s predictions—which differ from those of state-of-the-art docking software—using mutagenesis experiments. These results reveal a structural motif that influences the subtype selectivity of D₂R-targeted drugs and may thus prove useful in optimization of these ligands. ComBindVS also enables improved prediction of the effects of ligand modifications on binding affinity.Our approach provides a principled manner to integrate physics-based structural modeling with inference based on experimental data for other ligands, including ligands that share no common scaffold or substructure. Similar methods may prove useful in combining physics-based modeling with ligand-based approaches to improve prediction of various ligand properties by exploiting diverse sources of data.     

High-resolution structures of the M2 channel from influenza A virus reveal dynamic pathways for proton stabilization and transduction

The matrix 2 (M2) protein from influenza A virus is a proton channel that uses His37 as a selectivity filter. Here we report high-resolution (1.10 Å) cryogenic crystallographic structures of the transmembrane domain of M2 at low and high pH. These structures reveal that waters within the pore form hydrogen-bonded networks or “water wires” spanning 17 Å from the channel entrance to His37. Pore-lining carbonyl groups are well situated to stabilize hydronium via second-shell interactions involving bridging water molecules. In addition, room temperature crystallographic structures indicate that water becomes increasingly fluid with increasing temperature and decreasing pH, despite the higher electrostatic field. Complementary molecular dynamics simulations reveal a collective switch of hydrogen bond orientations that can contribute to the directionality of proton flux as His37 is dynamically protonated and deprotonated in the conduction cycle.Proton transport and conduction is essential to life. Proteins conduct protons over long distances through membranes to facilitate proton-coupled electron transfer and the formation and utilization of proton gradients. The M2 proton channel from the influenza A virus (1) is not only a medically important protein but also a simple, well-defined system for studying proton transport through confined spaces (2–4). This channel is the target of the anti-flu drug amantadine. M2 is activated at low pH by protonation of His37, which also participates in proton conduction by shuttling protons into the interior of the virus (5–7). His37 lies near the center of the bilayer, where it is connected to the viral exterior by a water-filled pore through which protons must pass to gain access to the viral interior (8–13).Visualizing the flow of protons within protein channels such as M2 is one of the long-standing challenges in molecular biophysics. Based on computational studies (9, 14–19) it has been suggested that protons reach His37 through “water wires” via the Grotthuss mechanism, but there is little high-resolution information concerning the path by which protons are conducted. A previously solved 1.65-Å crystal structure (9) showed six ordered waters immediately above the His37 tetrad, but ordered waters spanning the entire aqueous pore of M2 have not been observed until now. Previous MD simulations suggested a pore with mobile waters (12, 15), whereas the results of NMR and IR experiments are more consistent with an environment that is more similar to bulk water at low pH (13, 19, 20). However, it is difficult to deconvolute the changes in the water structure and dynamics when the protonation of His37 is raised from those induced indirectly via the conformation of the protein’s main chain.The M2 channel is known to have at least two conformational states that are populated to differing extents at low versus high pH (1, 10, 12). One, seen primarily at high pH, has been characterized extensively by solution NMR (21, 22), solid-state NMR (SSNMR) (10, 12), and X-ray crystallography (9). A second form is observed in dynamic equilibrium at lower pH (21–23), as evidenced by a broadening of magnetic resonances that thus far has made it impractical to determine a high-resolution structure of the protein in this state by SSNMR or solution NMR. X-ray crystallographic studies, however, have provided structures of both states (8, 9), which differ primarily in the conformation of the C terminus where protons exit the channel. Here we have obtained crystals that diffract to high resolution (1.10 Å) at both low and high pH, allowing visualization of water wires leading to His37 as a function of pH. The conformations of the backbone at the two pH values are essentially identical, permitting us to isolate changes in the organization of the water without any confounding factors.     

A Theoretical Study of the Binding of [Re6Se8(OH)2(H2O)4] Rhenium Clusters to DNA Purine Base Guanine

Hexanuclear rhenium complexes are promising candidates for use as antitumor drugs. However, to date, there has been no investigation into the nature of their binding to DNA. In this study, density functional theory (DFT) was used to examine the binding of [Re₆Se₈(OH)₂(H₂O)₄] to the DNA purine base guanine. The geometrical structures of cluster-guanine adducts in water were modeled at the zero order regular approximation (ZORA)-PW91 level. Calculating the bond energies allowed us to compare the cis and trans forms of the cluster, and a possible manners of interaction between [Re₆Se₈(OH)₂(H₂O)₃] clusters and DNA was obtained and explained.     

Effect of Impurity Adsorption on the Electronic and Transport Properties of Graphene Nanogaps

Graphene stands out as a versatile material with several uses in fields that range from electronics to biology. In particular, graphene has been proposed as an electrode in molecular electronics devices that are expected to be more stable and reproducible than typical ones based on metallic electrodes. In this work, we study by means of first principles, simulations and a tight-binding model the electronic and transport properties of graphene nanogaps with straight edges and different passivating atoms: Hydrogen or elements of the second row of the periodic table (boron, carbon, nitrogen, oxygen, and fluoride). We use the tight-binding model to reproduce the main ab-initio results and elucidate the physics behind the transport properties. We observe clear patterns that emerge in the conductance and the current as one moves from boron to fluoride. In particular, we find that the conductance decreases and the tunneling decaying factor increases from the former to the latter. We explain these trends in terms of the size of the atom and its onsite energy. We also find a similar pattern for the current, which is ohmic and smooth in general. However, when the size of the simulation cell is the smallest one along the direction perpendicular to the transport direction, we obtain highly non-linear behavior with negative differential resistance. This interesting and surprising behavior can be explained by taking into account the presence of Fano resonances and other interference effects, which emerge due to couplings to side atoms at the edges and other couplings across the gap. Such features enter the bias window as the bias increases and strongly affect the current, giving rise to the non-linear evolution. As a whole, these results can be used as a template to understand the transport properties of straight graphene nanogaps and similar systems and distinguish the presence of different elements in the junction.     

Improved Hydrogen Generation of Al-H2O Reaction by BiOX (X = Halogen) and Influence Rule

In this work, three additives BiOX (BiOI, BiOBr, and BiOF) for Al-H₂O reaction have been synthesized using chemical methods. SEM analysis shows that the structure of BiOF is nanoparticles, while BiOBr and BiOI have flower-like structures composed of nanosheets. Then, Al-BiOI, Al-BiOBr, and Al-BiOF composites have been prepared using the ball milling method. The effect of halogen ions on the performance of hydrogen generation from Al hydrolysis has been explored. The results indicate that the conversion yields of Al-BiOBr, Al-BiOI, and Al-BiOF for hydrogen generation are 96.3%, 95.3%, and 8.9%, respectively. In particular, the maximum hydrogen generation rate (MHGR) of Al-BiOI is as high as 3451.8 mL g⁻¹ min⁻¹, eight times higher than that of Al-BiOBr. Furthermore, the influence rule of BiOX (X = F, Cl, Br, and I) on Al-H₂O reaction has been studied using density functional theory. The results illustrate that HI can be more easily adsorbed on the Al surface as compared with HF, HCl, and HBr. Meanwhile, the bond length between halogen ions and the Al atom increased in the order of F⁻, Cl⁻, Br⁻, and I⁻. Therefore, the dissociation of I⁻ from the Al surface becomes easier and will expose more active sites to enhance the reaction activity of Al. In summary, the BiOI has the most favorable performance to Al-H₂O reaction.     

First-Principles Calculations of Thermoelectric Transport Properties of Quaternary and Ternary Bulk Chalcogenide Crystals

Chalcogenide crystals have a wide range of applications, especially as thermoelectric materials for energy conversion. Thermoelectric materials can be used to generate an electric current from a temperature gradient based on the Seebeck effect and based on the Peltier effect, and they can be used in cooling applications. Using first-principles calculations and semiclassical Boltzmann theory, we have computed the Seebeck coefficient, electrical conductivity, electronic thermal conductivity, power factor, and figure of merit of 30 chalcogenide crystals. A Quantum Espresso package is used to calculate the electronic properties and locate the Fermi level. The transport properties are then calculated using the BoltzTraP code. The 30 crystals are divided into two groups. The first group has four crystals with quaternary composition (A₂BCQ₄) (A = Tl; B = Cd, Hg; C = Si, Ge, Sn; Q = S, Se, Te). The second group contains 26 crystals with the ternary composition (A’B’Q₂) (A’ = Ag, Cu, Au, Na; B’ = B, Al, Ga, In; Q = S, Se, Te). Among these 30 chalcogenide crystals, the results for 11 crystals: Tl₂CdGeSe₄, Tl₂CdSnSe₄, Tl₂HgSiSe₄, Tl₂HgSnS₄, AuBSe₂, AuBTe₂, AuAlTe₂, AuGaTe₂, AuInTe₂, AgAlSe₂, and AgAlTe₂ are revealed for the first time. In addition, temperature-dependent transport properties of pure and doped AgSbSe₂ and AgSbTe₂ crystals with dopant compositions of AgSb_0.94Cd_0.06Te₂ and AgSbTe_1.85Se_0.15 were explored. These results provide an excellent database for bulk chalcogenides crucial for a wide range of potential applications in renewable energy fields.     

Electronic,Structural, and Optical Properties of Mono-Doped and Co-Doped (210) TiO2 Brookite Surfaces for Application in Dye-Sensitized Solar Cells—A First Principles Study

Titanium dioxide (TiO₂) polymorphs have recently gained a lot of attention in dye-sensitized solar cells (DSSCs). The brookite polymorph, among other TiO₂ polymorphs, is now becoming the focus of research in DSSC applications, despite the difficulties in obtaining it as a pure phase experimentally. The current theoretical study used different nonmetals (C, S and N) and (C-S, C-N and S-N) as dopants and co-dopants, respectively, to investigate the effects of mono-doping and co-doping on the electronic, structural, and optical structure properties of (210) TiO₂ brookite surfaces, which is the most exposed surface of brookite. The results show that due to the narrowing of the band gap and the presence of impurity levels in the band gap, all mono-doped and co-doped TiO₂ brookite (210) surfaces exhibit some redshift. In particular, the C-doped, and C-N co-doped TiO₂ brookite (210) surfaces exhibit better absorption in the visible region of the electromagnetic spectrum in comparison to the pure, S-doped, N-doped, C-S co-doped and N-S co-doped TiO₂ brookite (210) surfaces.     

Dynamics, transition states, and timing of bond formation in Diels-Alder reactions

The time-resolved mechanisms for eight Diels–Alder reactions have been studied by quasiclassical trajectories at 298 K, with energies and derivatives computed by UB3LYP/6-31G(d). Three of these reactions were also simulated at high temperature to compare with experimental results. The reaction trajectories require 50–150 fs on average to transverse the region near the saddle point where bonding changes occur. Even with symmetrical reactants, the trajectories invariably involve unequal bond formation in the transition state. Nevertheless, the time gap between formation of the two new bonds is shorter than a C─C vibrational period. At 298 K, most Diels–Alder reactions are concerted and stereospecific, but at high temperatures (approximately 1,000 K) a small fraction of trajectories lead to diradicals. The simulations illustrate and affirm the bottleneck property of the transition state and the close connection between dynamics and the conventional analysis based on saddle point structure.     

First-Principles Study on Mechanical and Optical Behavior of Plutonium Oxide under Typical Structural Phases and Vacancy Defects

The chemical corrosion aging of plutonium is a very important topic. It is easy to be corroded and produces oxidation products of various valence states because of its 5f electron orbit between local and non-local. On the one hand, the phase diagram of plutonium and oxygen is complex, so there is still not enough research on typical structural phases. On the other hand, most of the studies on plutonium oxide focus on PuO₂ and Pu₂O₃ with stoichiometric ratio, while the understanding of non-stoichiometric ratio, especially for Pu₂O_3-x, is not deep enough. Based on this, using the DFT + U theoretical scheme of density functional theory, we have systematically studied the structural stability, lattice parameters, electronic structure, mechanical and optical properties of six typical high temperature phases of β-Pu₂O₃, α-Pu₂O_3, γ-Pu₂O_3, PuO, α-PuO_2, γ-PuO₂. Further, the mechanical properties and optical behavior of Pu₂O_3-x under different oxygen vacancy concentrations are analyzed and discussed in detail. The result shows that the elasticity modulus of single crystal in mechanical properties is directly related to the oxygen/plutonium ratio and crystal system. As the number of oxygen vacancies increases, the mechanical constants continue to increase. In terms of optical properties, PuO has the best optical properties, and the light absorption rate decreases with the increase of oxygen vacancy concentration.     

INAUGURAL ARTICLE by a Recently Elected Academy Member:From the Cover: Atomic-scale insight and design principles for turbine engine thermal barrier coatings from theory

To maximize energy efficiency, gas turbine engines used in airplanes and for power generation operate at very high temperatures, even above the melting point of the metal alloys from which they are comprised. This feat is accomplished in part via the deposition of a multilayer, multicomponent thermal barrier coating (TBC), which lasts up to approximately 40,000 h before failing. Understanding failure mechanisms can aid in designing circumvention strategies. We review results of quantum mechanics calculations used to test hypotheses about impurities that harm TBCs and transition metal (TM) additives that render TBCs more robust. In particular, we discovered a number of roles that Pt and early TMs such as Hf and Y additives play in extending the lifetime of TBCs. Fundamental insight into the nature of the bonding created by such additives and its effect on high-temperature evolution of the TBCs led to design principles that can be used to create materials for even more efficient engines.     

Cigarette smoking, alcohol intake, and oral contraceptives: relationships to lipids and lipoproteins in adolescent school-children.

The effects of cigarette smoking, alcohol intake, and oral contraceptives on plasma cholesterol, triglyceride, high density lipoprotein cholesterol (C-HDL), and low density lipoprotein cholesterol (C-LDL) were assessed in 965 12--19-year-old school-children in the Cincinnati Lipid Research Clinic's Princeton school survey. After pair matching for age, sex, race, and total plasma cholesterol, adolescent children who smoked had mean C-HDL 6.1 mg/dl lower, and mean C-LDL 4.1 mg/dl higher, than nonsmokers (p less than 0.01). These findings for C-HDL were replicated by covariance analysis, adjusting for age, race, sex, alcohol intake, and triglyceride levels. Adolescents who drank alcohol had higher C-HDL and triglyceride levels and lower C-LDL than nondrinkers, but the differences were not significant. Adolescent females taking oral contraceptives had higher triglyceride, C-HDL, and C-LDL levels than matched controls, but the differences were not significant. If a portion of smoking's contribution to coronary heart disease risk is mediated through its inverse association with C-HDL, and if smoking habits initiated in adolescence continue into adulthood, this report provides additional physiologic data relevant to programs designed to prevent, reduce, or stop cigarette smoking in the adolescent years.     

Half Metallic Ferromagnetism and Transport Properties of Zinc Chalcogenides ZnX2Se4 (X = Ti,V, Cr) for Spintronic Applications

In ferromagnetic semiconductors, the coupling of magnetic ordering with semiconductor character accelerates the quantum computing. The structural stability, Curie temperature (T_c), spin polarization, half magnetic ferromagnetism and transport properties of ZnX₂Se₄ (X = Ti, V, Cr) chalcogenides for spintronic and thermoelectric applications are studied here by density functional theory (DFT). The highest value of T_c is perceived for ZnCr₂Se₄. The band structures in both spin channels confirmed half metallic ferromagnetic behavior, which is approved by integer magnetic moments (2, 3, 4) μ_B of Ti, V and Cr based spinels. The HM behavior is further measured by computing crystal field energy ΔE _crystal, exchange energies Δ _x(d), Δ _x (pd) and exchange constants (N_o α and N_o β ). The thermoelectric properties are addressed in terms of electrical conductivity, thermal conductivity, Seebeck coefficient and power factor in within a temperature range 0–400 K. The positive Seebeck coefficient shows p-type character and the PF is highest for ZnTi2Se4 (1.2 × 10¹¹ W/mK²) among studied compounds.