Residue level quantification of protein stability in living cells

The intracellular milieu differs from the dilute conditions in which most biophysical and biochemical studies are performed. This difference has led both experimentalists and theoreticians to tackle the challenging task of understanding how the intracellular environment affects the properties of biopolymers. Despite a growing number of in-cell studies, there is a lack of quantitative, residue-level information about equilibrium thermodynamic protein stability under nonperturbing conditions. We report the use of NMR-detected hydrogen–deuterium exchange of quenched cell lysates to measure individual opening free energies of the 56-aa B1 domain of protein G (GB1) in living Escherichia coli cells without adding destabilizing cosolutes or heat. Comparisons to dilute solution data (pH 7.6 and 37 °C) show that opening free energies increase by as much as 1.14 ± 0.05 kcal/mol in cells. Importantly, we also show that homogeneous protein crowders destabilize GB1, highlighting the challenge of recreating the cellular interior. We discuss our findings in terms of hard-core excluded volume effects, charge–charge GB1-crowder interactions, and other factors. The quenched lysate method identifies the residues most important for folding GB1 in cells, and should prove useful for quantifying the stability of other globular proteins in cells to gain a more complete understanding of the effects of the intracellular environment on protein chemistry.Proteins function in a heterogeneous and crowded intracellular environment. Macromolecules comprise 20–30% of the volume of an Escherichia coli cell and reach concentrations of 300–400 g/L (1, 2). Theory predicts that the properties of proteins and nucleic acids can be significantly altered in cells compared with buffer alone (3, 4). Nevertheless, most biochemical and biophysical studies are conducted under dilute (<10 g/L macromolecules) conditions. Here, we augment the small but growing list of reports probing the equilibrium thermodynamic stability of proteins in living cells (5–9), and provide, to our knowledge, the first measurement of residue-level stability under nonperturbing conditions.Until recently, the effects of macromolecular crowding on protein stability were thought to be caused solely by hard-core, steric repulsions arising from the impenetrability of matter (4, 10, 11). The expectation was that crowding enhances stability by favoring the compact native state over the ensemble of denatured states. Increased attention to transient, nonspecific protein-protein interactions (12–15) has led both experimentalists (16–19) and theoreticians (20–22) to recognize the effects of chemical interactions between crowder and test protein when assessing the net effect of macromolecular crowding. These weak, nonspecific interactions can reinforce or oppose the effect of hard-core repulsions, resulting in increased or decreased stability depending on the chemical nature of the test protein and crowder (23–26).We chose the B1 domain of streptococcal protein G (GB1) (27) as our test protein because its structure, stability and folding kinetics have been extensively studied in dilute solution (28–38). Its small size (56 aa; 6.2 kDa) and high thermal stability make GB1 well suited for studies by NMR spectroscopy.Quantifying the equilibrium thermodynamic stability of proteins relies on determining the relative populations of native and denatured states. Because the denatured state ensemble of a stable protein is sparsely populated under native conditions, stability is usually probed by adding heat or a cosolute to promote unfolding so that the concentration ratio of the two states can be determined (39). However, stability can be measured without these perturbations by exploiting the phenomenon of backbone amide H/D exchange (40) detected by NMR spectroscopy (41). The observed rate of amide proton (N–H) exchange, k_obs, is related to equilibrium stability by considering a protein in which each N–H exists in an open (exposed, exchange-competent) state, or a closed (protected, exchange-incompetent) state (40, 42):closed(N−H)⇌kclkopopen(N−H)→kintopen(N−D)⇌kopkclclosed(N−D).[1]Each position opens and closes with rate constants, k_op and k_cl (where K_op = k_op/k_cl), and exchange from the open state occurs with intrinsic rate constant, k_int. Values for k_int are based on exchange data from unstructured peptides (43, 44). If the test protein is stable (i.e., k_cl >> k_op), the observed rate becomes:kobs=kopkintkcl+kint.[2]Exchange occurs within two limits (42). At the EX1 limit, closing is rate determining, and k_obs = k_op. This limit is usually observed for less stable proteins and at basic pH (45). Most globular proteins undergo EX2 kinetics, where exchange from the open state is rate limiting (i.e., k_cl >> k_int), and k_obs values can be converted to equilibrium opening free energies, ΔGop°′ (46):kobs=kopkclkint=Kopkint[3]ΔGop°′=−RT⁡lnkobskint,[4]where RT is the molar gas constant multiplied by the absolute temperature.The backbone amides most strongly involved in H-bonded regions of secondary structure exchange only from the fully unfolded state, yielding a maximum value of ΔGop°′ (47–49). For these residues ΔGop°′ approximates the free energy of denaturation, ΔGden°′, providing information on global stability. Lower amplitude fluctuations of the native state can give rise to partially unfolded forms (50), resulting in residues with ΔGop°′ values less than those of the global unfolders.In summary, NMR-detected H/D exchange can measure equilibrium thermodynamic stability of a protein at the level of individual amino acid residues under nonperturbing conditions. Inomata et al. (51) used this technique to measure k_obs values in human cells for four residues in ubiquitin, but experiments confirming the exchange mechanism were not reported and opening free energies were not quantified. Our results fill this void and provide quantitative residue-level protein stability measurements in living cells under nonperturbing conditions.     

Quinary structure modulates protein stability in cells

Protein quinary interactions organize the cellular interior and its metabolism. Although the interactions stabilizing secondary, tertiary, and quaternary protein structure are well defined, details about the protein–matrix contacts that comprise quinary structure remain elusive. This gap exists because proteins function in the crowded cellular environment, but are traditionally studied in simple buffered solutions. We use NMR-detected H/D exchange to quantify quinary interactions between the B1 domain of protein G and the cytosol of Escherichia coli. We demonstrate that a surface mutation in this protein is 10-fold more destabilizing in cells than in buffer, a surprising result that firmly establishes the significance of quinary interactions. Remarkably, the energy involved in these interactions can be as large as the energies that stabilize specific protein complexes. These results will drive the critical task of implementing quinary structure into models for understanding the proteome.The inside of cells is packed with macromolecules whose concentrations reach 300–400 g/L (1). Compared with the ideal (dilute) environments conventionally used to study proteins, crowding inside cells can significantly alter the biophysical landscape of proteins, including their equilibrium thermodynamic stability (2–6). Experimental and computational efforts establish that crowding effects arise from a combination of short-range (steric) repulsions and longer-range (often called soft) interactions between macromolecules (7–13). Despite a growing number of in-cell protein studies (2–6), there is no quantitative information about the energetics of quinary interactions.Amide proton exchange experiments have been used for more than 50 y to measure equilibrium protein stability, defined as the Gibbs free energy required to open the protein and expose individual backbone amide protons to solvent, ΔGop°′ (14). For the B1 domain of protein G (GB1), ΔGop°′ equals −RTln(k_obs/k_uns), where R is the gas constant, T is the absolute temperature, k_obs is the observed rate of exchange, and k_uns is the rate in an unstructured peptide (6). We know that the cytoplasm does not affect k_uns (15). Most importantly, we know that for exchange under these conditions ΔGop°′ approximates the free energy required to denature the protein, ΔGden°′ (6). Therefore, these experiments provide a thermodynamically rigorous measure of equilibrium global protein stability. Using this information, we quantified the stability of GB1 at the residue level in Escherichia coli (6) via NMR-detected backbone amide hydrogen/deuterium exchange (16).Thermodynamic cycles (17) can be used to quantify the energetics of interactions between proteins in specific protein complexes (17, 18) and between side chains in globular proteins (19, 20). Briefly, the individual effects of two single-site amino acid changes are compared with the combined effect of both mutations. If the sites interact, the sum of the contributions from the single-site changes will not equal the contribution from the double mutant. The difference between the two values measures the strength of the interaction.We realized that transferring a variant (“var”) from buffer (“buff”) to cells (“cell”) is analogous to making a second mutation to the protein (Fig. 1 and SI Appendix, Fig. S1). Discrepancies in the horizontal (and vertical) sides of Fig. 1 define the free energy (δΔΔGint°′) associated with quinary interactions. Differences in the residue-level stability change caused by the mutation (ΔΔGop,mut°′) in cells and in buffer are used to calculate δΔΔGop,int°′:δΔΔGint°′=(ΔGcell,var°′−ΔGcell,WT°′)−(ΔGbuff,var°′−ΔGbuff,WT°′)=ΔΔGmut,cell°′−ΔΔGmut,buff°′=ΔΔGcell,var°′−ΔΔGcell,WT°′.A negative value of δΔΔGint°′ indicates the introduction of an attractive interaction (relative to WT) upon transferring the mutant from buffer to cells.Open in a separate windowFig. 1.Thermodynamic cycle used to quantify quinary interactions.     

Label-free probe of HIV-1 TAT peptide binding to mimetic membranes

The transacting activator of transduction (TAT) protein plays a key role in the progression of AIDS. Studies have shown that a +8 charged sequence of amino acids in the protein, called the TAT peptide, enables the TAT protein to penetrate cell membranes. To probe mechanisms of binding and translocation of the TAT peptide into the cell, investigators have used phospholipid liposomes as cell membrane mimics. We have used the method of surface potential sensitive second harmonic generation (SHG), which is a label-free and interface-selective method, to study the binding of TAT to anionic 1-palmitoyl-2-oleoyl-sn-glycero-3-phospho-1′-rac-glycerol (POPG) and neutral 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) liposomes. It is the SHG sensitivity to the electrostatic field generated by a charged interface that enabled us to obtain the interfacial electrostatic potential. SHG together with the Poisson–Boltzmann equation yielded the dependence of the surface potential on the density of adsorbed TAT. We obtained the dissociation constants K_d for TAT binding to POPC and POPG liposomes and the maximum number of TATs that can bind to a given liposome surface. For POPC K_d was found to be 7.5 ± 2 μM, and for POPG K_d was 29.0 ± 4.0 μM. As TAT was added to the liposome solution the POPC surface potential changed from 0 mV to +37 mV, and for POPG it changed from −57 mV to −37 mV. A numerical calculation of K_d, which included all terms obtained from application of the Poisson–Boltzmann equation to the TAT liposome SHG data, was shown to be in good agreement with an approximated solution.The HIV type 1 (HIV-1) transacting activator of transduction (TAT) is an important regulatory protein for viral gene expression (1–3). It has been established that the TAT protein has a key role in the progression of AIDS and is a potential target for anti-HIV vaccines (4). For the TAT protein to carry out its biological functions, it needs to be readily imported into the cell. Studies on the cellular internalization of TAT have led to the discovery of the TAT peptide, a highly cationic 11-aa region (protein transduction domain) of the 86-aa full-length protein that is responsible for the TAT protein translocating across phospholipid membranes (5–8). The TAT peptide is a member of a class of peptides called cell-penetrating peptides (CPPs) that have generated great interest for drug delivery applications (ref. 9 and references therein). The exact mechanism by which the TAT peptide enters cells is not fully understood, but it is likely to involve a combination of energy-independent penetration and endocytosis pathways (8, 10). The first step in the process is high-affinity binding of the peptide to phospholipids and other components on the cell surface such as proteins and glycosaminoglycans (1, 9).The binding of the TAT peptide to liposomes has been investigated using a variety of techniques, each of which has its own advantages and limitations. Among the techniques are isothermal titration calorimetry (9, 11), fluorescence spectroscopy (12, 13), FRET (12, 14), single-molecule fluorescence microscopy (15, 16), and solid-state NMR (17). Second harmonic generation (SHG), as an interface-selective technique (18–24), does not require a label, and because SHG is sensitive to the interface potential, it is an attractive method to selectively probe the binding of the highly charged (+8) TAT peptide to liposome surfaces. Although coherent SHG is forbidden in centrosymmetric and isotropic bulk media for reasons of symmetry, it can be generated by a centrosymmetric structure, e.g., a sphere, provided that the object is centrosymmetric over roughly the length scale of the optical coherence, which is a function of the particle size, the wavelength of the incident light, and the refractive indexes at ω and 2ω (25–30). As a second-order nonlinear optical technique SHG has symmetry restrictions such that coherent SHG is not generated by the randomly oriented molecules in the bulk liquid, but can be generated coherently by the much smaller population of oriented interfacial species bound to a particle or planar surfaces. As a consequence the SHG signal from the interface is not overwhelmed by SHG from the much larger populations in the bulk media (25–28).The total second harmonic electric field, E_2ω, originating from a charged interface in contact with water can be expressed as (31–33)E2ω∝∑iχc,i(2)EωEω+∑jχinc,j(2)EωEω+χH2O(3)EωEωΦ,[1]where χc,i(2) represents the second-order susceptibility of the species i present at the interface; χinc,j(2) represents the incoherent contribution of the second-order susceptibility, arising from density and orientational fluctuations of the species j present in solution, often referred to as hyper-Rayleigh scattering; χH2O(3) is the third-order susceptibility originating chiefly from the polarization of the bulk water molecules polarized by the charged interface; Φ is the potential at the interface that is created by the surface charge; and E_ω is the electric field of the incident light at the fundamental frequency ω. The second-order susceptibility, χc,i(2), can be written as the product of the number of molecules, N, at the surface and the orientational ensemble average of the hyperpolarizability α_i⁽²⁾ of surface species i, yielding χc,i(2)=N〈αi(2)〉 (18). The bracket ?? indicates an orientational average over the interfacial molecules. The third term in Eq. 1 depicts a third-order process by which a second harmonic field is generated by a charged interface. This term is the focus of our work. The SHG signal is dependent on the surface potential created by the electrostatic field of the surface charges, often called the χ⁽³⁾ contribution to the SHG signal. The χ⁽³⁾ method has been used to extract the surface charge density of charged planar surfaces and microparticle surfaces, e.g., liposomes, polymer beads, and oil droplets in water (21, 25, 34–39).In this work, the χ⁽³⁾ SHG method is used to explore a biomedically relevant process. The binding of the highly cationic HIV-1 TAT peptide to liposome membranes changes the surface potential, thereby enabling the use of the χ⁽³⁾ method to study the binding process in a label-free manner. Two kinds of liposomes, neutral 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) and anionic 1-palmitoyl-2-oleoyl-sn-glycero-3-phospho-1′-rac-glycerol (POPG), were investigated. The chemical structures of TAT, POPC, and POPG lipids are shown in Scheme 1.Open in a separate windowScheme 1.Chemical structures of HIV-1 TAT (47–57) peptide and the POPC and POPG lipids.     

How amide hydrogens exchange in native proteins

Amide hydrogen exchange (HX) is widely used in protein biophysics even though our ignorance about the HX mechanism makes data interpretation imprecise. Notably, the open exchange-competent conformational state has not been identified. Based on analysis of an ultralong molecular dynamics trajectory of the protein BPTI, we propose that the open (O) states for amides that exchange by subglobal fluctuations are locally distorted conformations with two water molecules directly coordinated to the N–H group. The HX protection factors computed from the relative O-state populations agree well with experiment. The O states of different amides show little or no temporal correlation, even if adjacent residues unfold cooperatively. The mean residence time of the O state is ∼100 ps for all examined amides, so the large variation in measured HX rate must be attributed to the opening frequency. A few amides gain solvent access via tunnels or pores penetrated by water chains including native internal water molecules, but most amides access solvent by more local structural distortions. In either case, we argue that an overcoordinated N–H group is necessary for efficient proton transfer by Grotthuss-type structural diffusion.Before the tightly packed and densely H-bonded structure of globular proteins had been established, Hvidt and Linderstrøm-Lang (1) showed that all backbone amide hydrogens of insulin exchange with water hydrogens, implying that all parts of the polypeptide backbone are, at least transiently, exposed to solvent. In the following 60 y, hydrogen exchange (HX), usually monitored by NMR spectroscopy (2) or mass spectrometry (3), has been widely used to study protein folding and stability (4–10), structure (11, 12), flexibility and dynamics (13–15), and solvent accessibility and binding (16, 17), often with single-residue resolution. However, because the exchange mechanism is unclear, HX data from proteins can, at best, be interpreted qualitatively (18–25).Under most conditions, amide HX is catalyzed by hydroxide ions (26, 27) at a rate that is influenced by inductive and steric effects from adjacent side chains (28). For unstructured peptides, HX is a slow process simply because the hydroxide concentration is low. For example, at 25° C and pH 4, HX occurs on a time scale of minutes. Under similar conditions, amides buried in globular proteins exchange on a wide range of time scales, extending up to centuries. HX can only occur if the amide is exposed to solvent, so conformational fluctuations must be an integral part of the HX mechanism (18).Under sufficiently destabilizing conditions HX occurs from the denatured-state ensemble, but under native conditions few amides exchange by such global unfolding (9, 29–31). For example, in bovine pancreatic trypsin inhibitor (BPTI), 8 amides in the core β-sheet exchange by global unfolding under native conditions (7, 32), whereas the remaining 45 amides require less extensive conformational fluctuations. Much of the debate in the protein HX field over the past half-century has concerned the nature of these subglobal fluctuations and their frequency, duration, amplitude, and cooperativity (18–25).According to the standard HX model (18), each amide can exist in a closed (C) state, where exchange cannot occur, or in an open (O) state, where exchange proceeds at a rate k_int. The kinetic scheme for H exchange into D₂O then reads as(N−H)C⇌kclkop(N−H)O→kint(N−D)Oand the measured steady-state HX rate is k_HX = k_op?k_int/(k_op + k_cl + k_int). To make this phenomenological model practically useful, two auxiliary assumptions are needed to disentangle the conformational and intrinsic parts of the process: (i) The conformational fluctuations (k_op and k_cl) are independent of pH, and (ii) HX from the O state proceeds at the same rate as in model peptides with the same neighboring side chains, so that kint=kHX0.Two HX regimes are distinguished with reference to the pH dependence of k_HX (18). If k_HX is constant in some pH range, it follows that k_int ? k_op + k_cl so that k_HX = k_op. In this so-called EX1 limit, the HX experiment measures the opening rate, or the mean residence time (MRT), of the C state, τ_C = 1/k_op. For BPTI, such pH invariance has only been observed for the eight core amides, and then only in a narrow pH interval (32).More commonly, HX experiments are performed in the EX2 limit, where k_int ? k_op + k_cl. Then k_HX = k_int/(κ + 1), where κ ≡ k_cl/k_op = τ_C/τ_O is the protection factor (PF). At equilibrium, the fractional populations, f_C and f_O, and the rates are linked by detailed balance, k_op?f_C = k_cl?f_O, so the PF may also be expressed as κ = f_C/f_O. Clearly, 1/(κ + 1) is the probability of finding the amide in the O state, 1/κ is the C  ?  O equilibrium constant, and β?ΔG = ln?κ is the free energy difference between the O and C states in units of k_B?T ≡ 1/β. The PF can thus be deduced from the HX rates measured (under EX2 conditions) for the amide in the protein and in a model peptide as κ=kHX0/kHX−1.The vast majority of the available protein HX data pertains to the EX2 regime and thus provides no information about the time scales, τ_C and τ_O, of the conformational fluctuations, except for the EX2 bound: 1/τC+1/τO≫kint≈kHX0. In the typical case where kHX≪kHX0, so that τ_C ? τ_O, we therefore only know that τO≪1/kHX0, which is in the millisecond range at pH 9 (EX2 HX data are usually measured at lower pH, where 1/kHX0 is even longer). Our analysis indicates that τ_O is seven orders of magnitude shorter than this upper bound estimate.The HX experiment is unique in probing sparsely populated conformational states with single-residue resolution. However, the physical significance of the PF is obscured by our ignorance about the structure and dynamics of the O state. Several attempts have been made to correlate experimental PFs with physical attributes of the amides, such as solvent contact (33–37), burial depth (38), intramolecular H-bonds (35, 38–40), packing density (38, 41), or electric field (42). Where significant correlations have been found, they suggest that the chosen attribute can serve as a proxy for the propensity for C → O fluctuations. However, whether based on crystal structures or molecular dynamics (MD) trajectories, these studies examined the time-averaged protein structure, which is dominated by the C state and therefore provides little or no information about the nature of the C → O fluctuations.In principle, the O state can be identified from molecular simulations, but this requires extensive conformational sampling because most C → O transitions are exceedingly rare. To date, this approach has been tried only with coarse-grained and/or empirical protein models without explicit solvent (43–45), or for HX from the denatured-state ensemble (46). The recent availability of ultralong MD simulations with realistic force fields opens up new opportunities in the search for the elusive O state. We have thus analyzed the millisecond MD trajectory of fully solvated native BPTI performed by Shaw et al. (47). Fortunately, BPTI is also among the proteins that have been most thoroughly studied by HX experiments.     

Mechanisms underlying drug-mediated regulation of membrane protein function

The hydrophobic coupling between membrane proteins and their host lipid bilayer provides a mechanism by which bilayer-modifying drugs may alter protein function. Drug regulation of membrane protein function thus may be mediated by both direct interactions with the protein and drug-induced alterations of bilayer properties, in which the latter will alter the energetics of protein conformational changes. To tease apart these mechanisms, we examine how the prototypical, proton-gated bacterial potassium channel KcsA is regulated by bilayer-modifying drugs using a fluorescence-based approach to quantify changes in both KcsA function and lipid bilayer properties (using gramicidin channels as probes). All tested drugs inhibited KcsA activity, and the changes in the different gating steps varied with bilayer thickness, suggesting a coupling to the bilayer. Examining the correlations between changes in KcsA gating steps and bilayer properties reveals that drug-induced regulation of membrane protein function indeed involves bilayer-mediated mechanisms. Both direct, either specific or nonspecific, binding and bilayer-mediated mechanisms therefore are likely to be important whenever there is overlap between the concentration ranges at which a drug alters membrane protein function and bilayer properties. Because changes in bilayer properties will impact many diverse membrane proteins, they may cause indiscriminate changes in protein function.

Cell membranes regulate membrane protein function by providing the necessary lipids and the appropriate bulk bilayer environment for function (1–5). The bilayer-mediated regulation arises because membrane proteins are hydrophobically coupled to their host lipid bilayer. As membrane proteins undergo conformational transitions (change shape), from a state I to a state II, the hydrophobic coupling gives rise to changes in the organization of adjacent lipids, which has an energetic cost (ΔGbilayerI→II=ΔGdef II−ΔGdef I), where ΔGdef I and ΔGdef II are the bilayer deformation energies incurred in each state. The total energetic cost of a transition between two protein conformations (I and II) thus is the sum of contributions from conformational rearrangements within the protein (ΔGprotein I→II) and rearrangements within the bilayer to minimize the exposure of hydrophobic amino acids to the aqueous environment (ΔGbilayer I→II). When the constraints on lipid packing around a protein do not allow for perfect hydrophobic matching, resulting in residual exposure of hydrophobic amino acids, the associated residual exposure energy ΔGres I→II will also contribute to ΔGtotal I→II (6). Thus, the equilibrium constant between I and II is:KI→II=exp{–ΔGtotalI→IIkBT}=exp{–ΔGproteinI→II+ΔGbilayerI→II+ΔGresI→IIkBT},where ΔGbilayer and ΔGres will vary depending on protein shape and conformational transitions.ΔGbilayer varies with changes in bilayer physical properties (3), meaning that changes in bilayer properties will alter the protein’s conformational equilibrium and function. This bilayer-mediated regulation of protein function is nonspecific; any membrane protein is subject to this regulation, which becomes important when ΔGbilayer>kBT, where k_B is Boltzmann’s constant and T temperature in Kelvin. The bilayer’s role in regulating membrane protein function by amphiphiles (including many drugs) is important because amphiphiles reversibly partition into the bilayer/water interface and thereby alter bilayer properties (7). Amphiphilic drugs thus may act through both direct binding to the desired target(s) and bilayer-mediated mechanisms. Distinguishing between direct and indirect bilayer-mediated regulation will help understand the mechanism(s) underlying the pleiotropic effects observed with many biologically active molecules.To develop a framework for understanding the mechanisms of drug-induced effects on membrane proteins, we used KcsA as a prototypical, proton (H⁺)-gated (8) ion channel and examined how changes in bulk bilayer properties and diverse drugs with varying bilayer (7, 9, 10) and clinical off-target effects (11–13) alter its function. We used gramicidin channels, whose bilayer-mediated regulation is well understood (7), to determine the drugs’ bilayer-modifying potency. To calibrate the drugs’ bilayer-modifying effect to its concentration in the bilayer, we estimated the mole fraction of each drug using membrane or octanol partition coefficient. We determined the drugs’ affinity for the purified KcsA to compare to its effective and bilayer-modifying concentrations.We first show that changes in bulk bilayer properties, specifically thickness, regulate KcsA gating and that thin bilayers promote the activated state. A more complex picture emerges when we explore the effects of various (amphiphilic) drugs, in which a drug’s effect(s) varied with bilayer thickness. Yet, the drug mole fraction in the bilayer and changes in bilayer properties predicted the drugs’ effects on KcsA gating, which enabled us to tease apart direct and bilayer-mediated effects. Our results provide a strategy for determining whether a drug’s pleiotropic effects result from direct, nonspecific, and bilayer-mediated mechanisms.     

Modeling the instantaneous normal mode spectra of liquids as that of unstable elastic media

We study the instantaneous normal mode (INM) spectrum of a simulated soft-sphere liquid at different equilibrium temperatures T. We find that the spectrum of eigenvalues ρ(λ) has a sharp maximum near (but not at) λ=0 and decreases monotonically with |λ| on both the stable and unstable sides of the spectrum. The spectral shape strongly depends on temperature. It is rather asymmetric at low temperatures (close to the dynamical critical temperature) and becomes symmetric at high temperatures. To explain these findings we present a mean-field theory for ρ(λ), which is based on a heterogeneous elasticity model, in which the local shear moduli exhibit spatial fluctuations, including negative values. We find good agreement between the simulation data and the model calculations, done with the help of the self-consistent Born approximation (SCBA), when we take the variance of the fluctuations to be proportional to the temperature T. More importantly, we find an empirical correlation of the positions of the maxima of ρ(λ) with the low-frequency exponent of the density of the vibrational modes of the glasses obtained by quenching to T=0 from the temperature T. We discuss the present findings in connection to the liquid to glass transformation and its precursor phenomena.

The investigation of the potential energy surface (PES) V(r1(t)…rN(t)) of a liquid (made up of N particles with positions r1(t)…rN(t) at a time instant t) and the corresponding instantaneous normal modes (INMs) of the (Hessian) matrix of curvatures has been a focus of liquid and glass science since the appearance of Goldstein’s seminal article (1) on the relation between the PES and the liquid dynamics in the viscous regime above the glass transition (2–27).The PES has been shown to form a rather ragged landscape in configuration space (8, 28, 29) characterized by its stationary points. In a glass these points are minima and are called “inherent structures.” The PES is believed to contain important information on the liquid–glass transformation mechanism. For the latter a complete understanding is still missing (28, 30, 31). The existing molecular theory of the liquid–glass transformation is mode-coupling theory (MCT) (32, 33) and its mean-field Potts spin version (28, 34). MCT predicts a sharp transition at a temperature TMCT>Tg, where T_g is the temperature of structural arrest (glass transition temperature). MCT completely misses the heterogeneous activated relaxation processes (dynamical heterogeneities), which are evidently present around and below T_MCT and which are related to the unstable (negative-λ) part of the INM spectrum (28, 30).Near and above T_MCT, apparently, there occurs a fundamental change in the PES. Numerical studies of model liquids have shown that minima present below T_MCT change into saddles, which then explains the absence of activated processes above T_MCT (2–24). Very recently, it was shown that T_MCT is related to a localization–delocalization transition of the unstable INM modes (25, 26).The INM spectrum is obtained in molecular dynamic simulations by diagonalizing the Hessian matrix of the interaction potential, taken at a certain time instant t:Hijαβ(t)=∂2∂xi(α)∂xj(β)V{r1(t)…rN(t)},[1]with ri=(xi(1),xi(2),xi(3)). For large positive values of the eigenvalues λ_j (j=1…N, N being the number of particles in the system) they are related to the square of vibrational frequencies λj=ωj2, and one can consider the Hessian as the counterpart of the dynamical matrix of a solid. In this high-frequency regime one can identify the spectrum with the density of vibrational states (DOS) of the liquid viag(ω)=2 ωρ(λ(ω))=13N∑jδ(ω−ωj) .[2]For small and negative values of λ this identification is not possible. For the unstable part of the spectrum (λ<0) it has become common practice to call the imaginary number λ=iω˜ and define the corresponding DOS asg(ω˜)≡2 ω˜ρ(λ(ω˜)).[3]This function is plotted on the negative ω axis and the stable g(ω), according to [2], on the positive axis. However, the (as we shall see, very interesting) details of the spectrum ρ(λ) near λ = 0 become almost completely hidden by multiplying the spectrum with |ω|. In fact, it has been demonstrated by Sastry et al. (6) and Taraskin and Elliott (7) already 2 decades ago that the INM spectrum of liquids, if plotted as ρ(λ) and not as g(ω) according to [2] and [3], exhibits a characteristic cusp-like maximum at λ = 0. The shape of the spectrum changes strongly with temperature. This is what we find as well in our simulation and what we want to explore further in our present contribution.In the present contribution we demonstrate that the strong change of the spectrum with temperature can be rather well explained in terms of a model, in which the instantaneous harmonic spectrum of the liquid is interpreted to be that of an elastic medium, in which the local shear moduli exhibit strong spatial fluctuations, which includes a large number of negative values. Because these fluctuations are just a snapshot of thermal fluctuations, we assume that they are obeying Gaussian statistics, the variance of which is proportional to the temperature.Evidence for a characteristic change in the liquid configurations in the temperature range above T_g has been obtained in recent simulation studies of the low-frequency vibrational spectrum of glasses, which have been rapidly quenched from a certain parental temperature T*. If T* is decreased from high temperatures toward T_MCT, the low-frequency exponent of the vibrational DOS of the daughter glass (quenched from T* to T = 0) changed from Debye-like g(ω)∝ω2 to g(ω)∝ωs with s > 2. In our numerical investigation of the INM spectra we show a correlation of some details of the low-eigenvalue features of these spectra with the low-frequency properties of the daughter glasses obtained by quenching from the parental temperatures.The stochastic Helmholtz equations (Eq. 7) of an elastic model with spatially fluctuating shear moduli can be readily solved for the averaged Green’s functions by field theoretical techniques (35–37). Via a saddle point approximation with respect to the resulting effective field theory one arrives at a mean-field theory (self-consistent Born approximation [SCBA]) for the self-energy of the averaged Green’s functions. The SCBA predicts a stable spectrum below a threshold value of the variance. Restricted to this stable regime, this theory, called heterogeneous elasticity theory (HET), was rather successful in explaining several low-frequency anomalies in the vibrational spectrum of glasses, including the so-called boson peak, which is an enhancement at finite frequencies over the Debye behavior of the DOS g(ω)∝ω2 (37–41). We now explore the unstable regime of this theory and compare it to the INM spectrum of our simulated soft-sphere liquid.^*We start Results by presenting a comparison of the simulated spectra of the soft-sphere liquid with those obtained by the unstable version of HET-SCBA theory. We then concentrate on some specific features of the INM spectra, namely, the low-eigenvalue slopes and the shift of the spectral maximum from λ = 0. Both features are accounted for by HET-SCBA. In particular, we find an interesting law for the difference between the slopes of the unstable and the stable parts of the spectrum, which behaves as T−2/3, which, again, is accounted for by HET-SCBA.In the end we compare the shift of the spectral maximum with the low-frequency exponent of the DOS of the corresponding daughter glasses and find an empirical correlation. We discuss these results in connection with the saddle to minimum transformation near T_MCT.     

Transient conformational fluctuation of TePixD during a reaction

Knowledge of the dynamical behavior of proteins, and in particular their conformational fluctuations, is essential to understanding the mechanisms underlying their reactions. Here, transient enhancement of the isothermal partial molar compressibility, which is directly related to the conformational fluctuation, during a chemical reaction of a blue light sensor protein from the thermophilic cyanobacterium Thermosynechococcus elongatus BP-1 (TePixD, Tll0078) was investigated in a time-resolved manner. The UV-Vis absorption spectrum of TePixD did not change with the application of high pressure. Conversely, the transient grating signal intensities representing the volume change depended significantly on the pressure. This result implies that the compressibility changes during the reaction. From the pressure dependence of the amplitude, the compressibility change of two short-lived intermediate (I₁ and I₂) states were determined to be +(5.6 ± 0.6) × 10⁻² cm³⋅mol⁻¹⋅MPa⁻¹ for I₁ and +(6.6 ± 0.7)×10⁻² cm³⋅mol⁻¹⋅MPa⁻¹ for I₂. This result showed that the structural fluctuation of intermediates was enhanced during the reaction. To clarify the relationship between the fluctuation and the reaction, the compressibility of multiply excited TePixD was investigated. The isothermal compressibility of I₁ and I₂ intermediates of TePixD showed a monotonic decrease with increasing excitation laser power, and this tendency correlated with the reactivity of the protein. This result indicates that the TePixD decamer cannot react when its structural fluctuation is small. We concluded that the enhanced compressibility is an important factor for triggering the reaction of TePixD. To our knowledge, this is the first report showing enhanced fluctuations of intermediate species during a protein reaction, supporting the importance of fluctuations.Proteins often transfer information through changes in domain–domain (or intermolecular) interactions. Photosensor proteins are an important example. They have light-sensing domains and function by using the light-driven changes in domain–domain interactions (1). The sensor of blue light using FAD (BLUF) domain is a light-sensing module found widely among the bacterial kingdom (2). The BLUF domain initiates its photoreaction by the light excitation of the flavin moiety inside the protein, which changes the domain–domain interaction, causing a quaternary structural change and finally transmitting biological signals (3, 4). It has been an important research topic to elucidate how the initial photochemistry occurring in the vicinity of the chromophore leads to the subsequent large conformation change in other domains, which are generally apart from the chromophore.It may be reasonable to consider that the conformation change in the BLUF domain is the driving force in its subsequent reaction; that is, the change in domain–domain interaction. However, sometimes, clear conformational changes have not been observed for the BLUF domain; its conformation is very similar before and after photo-excitation (5–13). The circular dichroism (CD) spectra of BLUF proteins AppA and PixD from thermophilic cyanobacterium Thermosynechococcus elongatus BP-1 (TePixD) did not change on illumination (5, 13). Similarly, solution NMR studies of AppA and BlrB showed only small chemical shifts on excitation (9, 10). The solution NMR structure of BlrP1 showed a clear change, but this was limited in its C-terminal extension region and not core BLUF (11). Furthermore, the diffusion coefficient (D) of the BLUF domain of YcgF was not changed by photo-excitation (12), although D is sensitive to global conformational changes. These results imply that a minor structural change occurs in the BLUF domain. In such cases, how does the BLUF domain control its interdomain interaction? Recently, a molecular dynamics (MD) simulation on another light-sensing domain, the light-oxygen-voltage (LOV) sensing domain, suggested that fluctuation of the LOV core structure could be a key to understanding the mechanism of information transfer (14–16).Because proteins work at room temperature, they are exposed to thermal fluctuations. The importance of such structural fluctuations for biomolecular reactions has been also pointed out: for example, enzymatic activity (17–20). Experimental detections of such conformation fluctuations using single molecular detection (21) or NMR techniques such as the hydrogen-deuterium (H-D) exchange, relaxation dispersion method, and high-pressure NMR (22–24) have succeeded. However, these techniques could not detect the fluctuation of short-lived transient species. Indeed, single molecule spectroscopy can trace the fluctuation in real time, but it is still rather difficult to detect rapid fluctuations for a short-lived intermediate during a reaction. Therefore, information about the fluctuation of intermediates is thus far limited.A thermodynamic measurement is another way to characterize the fluctuation of proteins. In particular, the partial molar isothermal compressibility [K¯T=−(∂V¯/∂P)T] is essential, because this property is directly linked to the mean-square fluctuations of the protein partial molar volume by 〈(V¯−〈V¯〉)2〉≡〈δV¯2〉=kBTK¯T (25). (Here, <X> means the averaged value of a quantity of X.) Therefore, isothermal compressibility is thought to reflect the structural fluctuation of molecules (26). However, experimental measurement of this parameter of proteins in a dilute solution is quite difficult. Indeed, this quantity has been determined indirectly from the theoretical equation using the adiabatic compressibility of a protein solution, which was determined by the sound velocity in the solution (26–31). Although the relation between volume fluctuations and isothermal compressibility is rigorously correct only with respect to the intrinsic part of the volume compressibility, and not the partial molar volume compressibility (32), we considered that this partial molar volume compressibility is still useful for characterizing the fluctuation of the protein structure including its interacting water molecules. In fact, the relationship between β¯T and the volume fluctuation has been often used to discuss the fluctuation of proteins (17, 26–28), and the strong correlation of β¯T of reactants with the functioning for some enzymes (17, 33, 34) has been reported. These studies show the functional importance of the structural fluctuation represented by β¯T. However, thermodynamic techniques lack time resolution, and it has been impossible to measure the fluctuations of short-lived intermediate species.Recently, we developed a time-resolving method for assessing thermodynamic properties using the pulsed laser induced transient grating (TG) method. Using this method, we thus far succeeded in measuring the enthalpy change (ΔH) (35–38), partial molar volume change (ΔV¯) (12, 35, 37), thermal expansion change (Δα¯th) (12, 37), and heat capacity change (ΔC_p) (36–38) for short-lived species. Therefore, in principle, the partial molar isothermal compressibility change (ΔK¯T) of a short-lived intermediate become observable if we conduct the TG experiment under the high-pressure condition and detect ΔV¯ with varying external pressure.There are several difficulties in applying the traditional high-pressure cell to the TG method to measure thermodynamic parameters quantitatively. The most serious problem is ensuring the quantitative performance of the intensity of TG signals measured under the high-pressure condition. On this point, our group has developed a new high-pressure cell specially designed for TG spectroscopy (39) and overcome this problem. In this paper, by applying this high-pressure TG system to the BLUF protein TePixD, we report the first measurement, to our knowledge, of ΔK¯T of short-lived intermediates to investigate the mechanism underlying signal transmission by BLUF proteins, from the view point of the transient fluctuation.TePixD is a homolog of the BLUF protein PixD, which regulates the phototaxis of cyanobacterium (40) and exists in a thermophilic cyanobacterium Thermocynechococcus elongates BP-1 (Tll0078). TePixD is a relatively small (17 kDa) protein that consists only of the BLUF domain with two extended helices in the C-terminal region. In crystals and solutions, it forms a decamer that consists of two pentameric rings (41). The photochemistry of TePixD is typical among BLUF proteins (42–45); on blue light illumination, the absorption spectrum shifts toward red by about 10 nm within a nanosecond. The absorption spectrum does not change further, and the dark state is recovered with a time constant of ∼5 s at room temperature (40, 43). The spectral red shift was explained by the rearrangement of the hydrogen bond network around the chromophore (6, 46–48). The TG method has revealed the dynamic photoreaction mechanism, which cannot be detected by conventional spectroscopic methods. The TG signal of TePixD (Fig. S1) showed that there are two spectrally silent reaction phases: a partial molar volume expansion with the time constant of ∼40 μs and the diffusion coefficient (D) change with a time constant of ∼4 ms. Furthermore, it was reported that the pentamer and decamer states of TePixD are in equilibrium and that the final photoproduct of the decamer is pentamers generated by its dissociation (13, 49). On the basis of these studies, the reaction scheme has been identified as shown in Fig. 1. Here, I₁ is the intermediate of the spectrally red-shifted species (generated within a nanosecond) and I₂ is the one created on the subsequent volume expansion process of +4 cm³⋅mol⁻¹ (∼40 μs). Furthermore, an experiment of the excitation laser power dependence of its TG signal revealed that the TePixD decamer undergoes the original dissociation reaction when only one monomer in the decamer is excited (50). In this study, we investigated the transient compressibility of the intermediates I₁ and I₂ of the photoreaction of TePixD and found a direct link between their fluctuation and reactivity.Open in a separate windowFig. 1.Schematic illustration of the photoreaction of TePixD. Yellow circles represent the TePixD monomer in the ground state, which constructs the decamer and pentamer states. In the dark state, these two forms are in equilibrium. The excited, spectral red-shifted state of the TePixD monomer is indicated by a red circle. The square represents the I₂ state of the monomer, which is created by the volume expansion process.     

Nature of dynamic gradients,glass formation,and collective effects in ultrathin freestanding films

Molecular, polymeric, colloidal, and other classes of liquids can exhibit very large, spatially heterogeneous alterations of their dynamics and glass transition temperature when confined to nanoscale domains. Considerable progress has been made in understanding the related problem of near-interface relaxation and diffusion in thick films. However, the origin of “nanoconfinement effects” on the glassy dynamics of thin films, where gradients from different interfaces interact and genuine collective finite size effects may emerge, remains a longstanding open question. Here, we combine molecular dynamics simulations, probing 5 decades of relaxation, and the Elastically Cooperative Nonlinear Langevin Equation (ECNLE) theory, addressing 14 decades in timescale, to establish a microscopic and mechanistic understanding of the key features of altered dynamics in freestanding films spanning the full range from ultrathin to thick films. Simulations and theory are in qualitative and near-quantitative agreement without use of any adjustable parameters. For films of intermediate thickness, the dynamical behavior is well predicted to leading order using a simple linear superposition of thick-film exponential barrier gradients, including a remarkable suppression and flattening of various dynamical gradients in thin films. However, in sufficiently thin films the superposition approximation breaks down due to the emergence of genuine finite size confinement effects. ECNLE theory extended to treat thin films captures the phenomenology found in simulation, without invocation of any critical-like phenomena, on the basis of interface-nucleated gradients of local caging constraints, combined with interfacial and finite size-induced alterations of the collective elastic component of the structural relaxation process.

Spatially heterogeneous dynamics in glass-forming liquids confined to nanoscale domains (1–7) play a major role in determining the properties of molecular, polymeric, colloidal, and other glass-forming materials (8), including thin films of polymers (9, 10) and small molecules (11–15), small-molecule liquids in porous media (2, 4, 16, 17), semicrystalline polymers (18, 19), polymer nanocomposites (20–22), ionomers (23–25), self-assembled block and layered (26–33) copolymers, and vapor-deposited ultrastable molecular glasses (34–36). Intense interest in this problem over the last 30 y has also been motivated by the expectation that its understanding could reveal key insights concerning the mechanism of the bulk glass transition.Considerable progress has been made for near-interface altered dynamics in thick films, as recently critically reviewed (1). Large amplitude gradients of the structural relaxation time, τ(z,T), converge to the bulk value, τ_bulk(T), in an intriguing double-exponential manner with distance, z, from a solid or vapor interface (13, 37–42). This implies that the corresponding effective activation barrier, F_total(z,T,H) (where H is film thickness), varies exponentially with z, as does the glass transition temperature, T_g (37). Thus the fractional reduction in activation barrier, ε(z,H), obeys the equation ε(z,H)≡1−Ftotal(z,T,H)/Ftotal,bulk(T)=ε0⁡exp(−z/ξF), where F_total,bulk(T) is the bulk temperature-dependent barrier and ξ_F a length scale of modest magnitude. Although the gradient of reduction in absolute activation barriers becomes stronger with cooling, the amplitude of the fractional reduction of the barrier gradient, quantified by ε₀, and the range ξ_F of this gradient, exhibit a weak or absent temperature dependence at the lowest temperatures accessed by simulations (typically with the strength of temperature dependence of ξ_F decreasing rather than increasing on cooling), which extend to relaxation timescales of order 10⁵ ps. This finding raises questions regarding the relevance of critical-phenomena–like ideas for nanoconfinement effects (1). Partially due to this temperature invariance, coarse-grained and all-atom simulations (1, 37, 42, 43) have found a striking empirical fractional power law decoupling relation between τ(z,T) and τ_bulk(T):τ(T,z)τbulk(T)∝(τbulk(T))−ε(z).[1]Recent theoretical analysis suggests (44) that this behavior is consistent with a number of experimental data sets as well (45, 46). Eq. 1 also corresponds to a remarkable factorization of the temperature and spatial location dependences of the barrier:Ftotal(z,T)=[1−ε(z)]Ftotal,bulk(T).[2]This finding indicates that the activation barrier for near-interface relaxation can be factored into two contributions: a z-dependent, but T-independent, “decoupling exponent,” ε(z), and a temperature-dependent, but position-insensitive, bulk activation barrier, F_total,bulk(T). Eq. 2 further emphasizes that ε(z) is equivalent to an effective fractional barrier reduction factor (for a vapor interface), 1−Ftotal(z,T,H)/Ftotal,bulk(T), that can be extracted from relaxation data.In contrast, the origin of “nanoconfinement effects” in thin films, and how much of the rich thick-film physics survives when dynamic gradients from two interfaces overlap, is not well understood. The distinct theoretical efforts for aspects of the thick-film phenomenology (44, 47–50) mostly assume an additive summation of one-interface effects in thin films, thereby ignoring possibly crucial cooperative and whole film finite size confinement effects. If the latter involve phase-transition–like physics as per recent speculations (14, 51), one can ask the following: do new length scales emerge that might be truncated by finite film size? Alternatively, does ultrathin film phenomenology arise from a combination of two-interface superposition of the thick-film gradient physics and noncritical cooperative effects, perhaps in a property-, temperature-, and/or thickness-dependent manner?Here, we answer these questions and establish a mechanistic understanding of thin-film dynamics for the simplest and most universal case: a symmetric freestanding film with two vapor interfaces. We focus on small molecules (modeled theoretically as spheres) and low to medium molecular weight unentangled polymers, which empirically exhibit quite similar alterations in dynamics under “nanoconfinement.” We do not address anomalous phenomena [e.g., much longer gradient ranges (29), sporadic observation of two distinct glass transition temperatures (52, 53)] that are sometimes reported in experiments with very high molecular weight polymers and which may be associated with poorly understood chain connectivity effects that are distinct from general glass formation physics (54–56).We employ a combination of molecular dynamics simulations with a zero-parameter extension to thin films of the Elastically Cooperative Nonlinear Langevin Equation (ECNLE) theory (57, 58). This theory has previously been shown to predict well both bulk activated relaxation over up to 14 decades (44–46) and the full single-gradient phenomenology in thick films (1). Here, we extend this theory to treat films of finite thickness, accounting for coupled interface and geometric confinement effects. We compare predictions of ECNLE theory to our previously reported (37, 43) and new simulations, which focus on translational dynamics of films comprised of a standard Kremer–Grest-like bead-spring polymer model (see SI Appendix). These simulations cover a wide range of film thicknesses (H, from 4 to over 90 segment diameters σ) and extend to low temperatures where the bulk alpha time is ∼0.1 μs (10⁵ Lennard Jones time units τ_LJ).The generalized ECNLE theory is found to be in agreement with simulation for all levels of nanoconfinement. We emphasize that this theory does not a priori assume any of the empirically established behaviors discovered using simulation (e.g., fractional power law decoupling, double-exponential barrier gradient, gradient flattening) but rather predicts these phenomena based upon interfacial modifications of the two coupled contributions to the underlying activation barrier– local caging constraints and a long-ranged collective elastic field. It is notable that this strong agreement is found despite the fact the dynamical ideas are approximate, and a simple hard sphere fluid model is employed in contrast to the bead-spring polymers employed in simulation. The basic unit of length in simulation (bead size σ) and theory (hard sphere diameter d) are expected to be proportional to within a prefactor of order unity, which we neglect in making comparisons.As an empirical matter, we find from simulation that many features of thin-film behavior can be described to leading order by a linear superposition of the thick-film gradients in activation barrier, that is:ε(z,H)=1−Ftotal(z,T,H)/Ftotal,bulk(T)≈ε0[exp(−z/ξF)+exp(−(H−z)/ξF)],[3]where the intrinsic decay length ξ_F is unaltered from its thick-film value and where ε₀ is a constant that, in the hypothesis of literal gradient additivity, is invariant to temperature and film thickness. We employ this functional form [originally suggested by Binder and coworkers (59)], which is based on a simple superposition of the two single-interface gradients, as a null hypothesis throughout this study: this form is what one expects if no new finite-size physics enters the thin-film problem relative to the thick film.However, we find that the superposition approximation progressively breaks down, and eventually entirely fails, in ultrathin films as a consequence of the emergence of a finite size confinement effect. The ECNLE theory predicts that this failure is not tied to a phase-transition–like mechanism but rather is a consequence of two key coupled physical effects: 1) transfer of surface-induced reduction of local caging constraints into the film, and 2) interfacial truncation and nonadditive modifications of the collective elastic contribution to the activation barrier.     

Fracture of model end-linked networks

Advances in polymer chemistry over the last decade have enabled the synthesis of molecularly precise polymer networks that exhibit homogeneous structure. These precise polymer gels create the opportunity to establish true multiscale, molecular to macroscopic, relationships that define their elastic and failure properties. In this work, a theory of network fracture that accounts for loop defects is developed by drawing on recent advances in network elasticity. This loop-modified Lake–Thomas theory is tested against both molecular dynamics (MD) simulations and experimental fracture measurements on model gels, and good agreement between theory, which does not use an enhancement factor, and measurement is observed. Insight into the local and global contributions to energy dissipated during network failure and their relation to the bond dissociation energy is also provided. These findings enable a priori estimates of fracture energy in swollen gels where chain scission becomes an important failure mechanism.

Models that link materials structure to macroscopic behavior can account for multiple levels of molecular structure. For example, the statistical, affine deformation model connects the elastic modulus E to the molecular structure of a polymer chain,Eaff=3νkbT(ϕo13Roϕ13R)2,[1]where ν is density of chains, ϕ is polymer volume fraction, R is end-to-end distance, ϕo and R_o represent the parameters taken in the reference state that is assumed to be the reaction concentration in this work, and k_bT is the available thermal energy where k_b is Boltzmann’s constant and T is temperature (1–6). Refinements to this model that account for network-level structure, such as the presence of trapped entanglements or number of connections per junction, have been developed (7–11). Further refinements to the theory of network elasticity have been developed to account for dynamic processes such as chain relaxation and solvent transport (12–17). Together these refinements link network elasticity to chain-level molecular structure, network-level structure, and the dynamic processes that occur at both size scales.While elasticity has been connected to multiple levels of molecular structure, models for network fracture have not developed to a similar extent. The fracture energy G_c typically relies upon the large strain deformation behavior of polymer networks, making it experimentally difficult to separate the elastic energy released upon fracture from that dissipated through dynamic processes (18–26). In fact, most fracture theories have been developed at the continuum scale and have focused on modeling dynamic dissipation processes (27). An exception to this is the theory of Lake and Thomas that connects the elastic energy released during chain scission to chain-level structure,Gc,LT=ChainsArea×Energy DissipatedChain=νRoNU,[2]where NU is the total energy released when a chain ruptures in which N represents the number of monomer segments in the chain and U the energy released per monomer (26).While this model was first introduced in 1967, experimental attempts to verify Lake–Thomas theory as an explicit model, as summarized in SI Appendix, have been unsuccessful. Ahagon and Gent (28) and Gent and Tobias (29) attempted to do this on highly swollen networks at elevated temperature but found that, while the scalings from Eq. 2 work well, an enhancement factor was necessary to observe agreement between theory and experiment. This led many researchers to conclude that Lake–Thomas theory worked only as a scaling argument. In 2008, Sakai et al. (30) introduced a series of end-linked tetrafunctional, star-like poly(ethylene glycol) (PEG) gels. Scattering measurements indicated a lack of nanoscale heterogeneities that are characteristic of most polymer networks (30–32). Fracture measurements on these well-defined networks were performed and it was again observed that an enhancement factor was necessary to realize explicit agreement between experiment and theory (33). Arora et al. (34) recently attempted to address this discrepancy by accounting for loop defects; however, different assumptions were used when inputting U to calculate Lake–Thomas theory values that again required the use of an enhancement factor to achieve quantitative agreement. In this work we demonstrate that refining the Lake–Thomas theory to account for loop defects while using the full bond dissociation energy to represent U yields excellent agreement between the theory and both simulation and experimental data without the use of any adjustable parameters.PEG gels synthesized via telechelic end-linking reactions create the opportunity to build upon previous theory to establish true multiscale, molecular to macroscopic relationships that define the fracture response of polymer networks. This paper combines pure shear notch tests, molecular dynamics (MD) simulations, and theory to quantitatively extend the concept of network fracture without the use of an enhancement factor. First, the control of molecular-level structure in end-linked gel systems is discussed. Then, the choice of molecular parameters used to estimate chain- and network-level properties is discussed. Experimental and MD simulation methods used when fracturing model end-linked networks are then presented. A theory of network fracture that accounts for loop defects is developed, in the context of other such models that have emerged recently, and tested against data from experiments and MD simulations. Finally, a discussion of the local and global energy dissipated during failure of the network is presented.     

Discontinuous shear thickening in Brownian suspensions by dynamic simulation

Dynamic particle-scale numerical simulations are used to show that the shear thickening observed in dense colloidal, or Brownian, suspensions is of a similar nature to that observed in noncolloidal suspensions, i.e., a stress-induced transition from a flow of lubricated near-contacting particles to a flow of a frictionally contacting network of particles. Abrupt (or discontinuous) shear thickening is found to be a geometric rather than hydrodynamic phenomenon; it stems from the strong sensitivity of the jamming volume fraction to the nature of contact forces between suspended particles. The thickening obtained in a colloidal suspension of purely hard frictional spheres is qualitatively similar to experimental observations. However, the agreement cannot be made quantitative with only hydrodynamics, frictional contacts, and Brownian forces. Therefore, the role of a short-range repulsive potential mimicking the stabilization of actual suspensions on the thickening is studied. The effects of Brownian and repulsive forces on the onset stress can be combined in an additive manner. The simulations including Brownian and stabilizing forces show excellent agreement with experimental data for the viscosity η and the second normal stress difference N₂.The rheology of dense suspensions is of considerable theoretical and technological importance, yet the shear rheology of even the simplest case of a suspension of hard spheres in a Newtonian suspending fluid is incompletely understood (1). Many of the features observed in these suspensions, including shear thinning (2) or thickening (3, 4) and the magnitudes and even the algebraic signs of normal stress differences (5), are at best understood at a qualitative level, and a general theoretical framework is lacking. Furthermore, there has been a tendency to treat the rheology of Brownian (colloidal) suspensions and non-Brownian suspensions as distinct.Recently, a picture has emerged in which central aspects of the rheology of non-Brownian dense suspensions are interpreted as manifestations of proximity to jamming transitions in the parameter space. These transitions are singularities whose locations in the volume fraction ϕ and shear stress σ-plane depend on the details of the microscopic interactions (shape of the particles, friction, interparticle forces). In turn, the locations of these singularities shape the large ϕ-portion of the rheological landscape, i.e., the effective viscosity and the normal stresses as functions of ϕ and σ. In particular, in the “stress-induced friction” scenario (4, 6–13), shear thickening is a transition from a rheological response dominated by frictionless jamming to one controlled by frictional jamming upon increase of the shear stress. This transition is argued to be due to the creation of frictional contacts between particles at high stresses; the contacts are prevented at low stresses by a short-range stabilizing repulsive force, as would be expected to be present to stabilize a colloidal dispersion against aggregation by, for example, attractive van der Waals forces (14, 15). This picture contrasts with previous models for Brownian suspensions (16–18), in which frictional contacts are neglected based on idealized lubrication hydrodynamics.Suspensions are said to be colloidal, or Brownian, when the immersed particles are small enough: a commonly accepted upper bound for Brownian motion to be significant is a diameter of 1 µm (2). For these systems, the Brownian forces have been seen to be an essential factor in non-Newtonian behavior (19). Physically, for a system of strictly hard colloidal spheres under shear in the Stokes regime, there are only two independent time scales: the inverse shear rate γ˙−1 and the diffusion time a²/D₀; here, a is the sphere radius and D₀ is the single-particle diffusion coefficient, which is related to the suspending fluid viscosity η₀ and thermal energy k_BT through the Stokes–Einstein relation D₀ = k_BT/6πη₀a. The shear rate dependence of the rheology can be stated in terms of a competition between advection and diffusion described by the Péclet number Pe≡6πη0a3γ˙/kT. Smooth spheres with ideal lubrication resulting from hydrodynamic interactions would be noncontacting, and hence exhibit the following rheology: a shear-rate–independent regime close to thermal equilibrium [that is for γ˙≲τα−1, where τ_α is the “caging” time, or the typical time for which the thermal motion leads to a structural reorganization (20)] where most forces are Brownian, followed by a shear thinning regime at intermediate values of Pe, over which the Brownian forces become progressively less important relative to the hydrodynamic ones, and finally a shear thickening regime at large Pe that is dominated by the hydrodynamic lubrication forces due to increasingly smaller interparticle gaps (16–18). τ_α diverges at the glass transition ? = ?_G, the system develops a yield stress σ_y above ?_G, and the low Pe viscosity plateau yields to an asymptotic shear thinning regime η∼σy/γ˙ for γ˙→0 (21–23). Although numerical simulations including only hydrodynamic interactions agree well with experimental data in the shear thinning regime, the simulated shear thickening regime is much weaker than is often experimentally observed, with the disagreement increasing as the volume fraction increases (18).Most of the experimentally studied thickening suspensions are in the colloidal size range (4). (Notable exceptions include cornstarch suspensions.) It is thus essential to address the validity of the stress-induced friction scenario for these systems. In this scenario, the number of frictional contacts directly depends on the ratio of the shear stress to the Brownian stress scale σa³/k_BT. At small stresses, i.e., when σa³/k_BT ? 1, the thermal motion keeps particles separated and makes contacts unlikely. [In the equilibrium limit Pe → 0, the average contact number per particle is zero at volume fractions below the jamming transition, otherwise the pressure would diverge as required by the virial equation for hard particles (24).] When σa³/k_BT ? 1, however, the Brownian forces are not strong enough to overcome the forces bringing particles together due to the shear, and contacts are created. For dense suspensions, this shear-activated friction mechanism can be related to a jamming transition framework: at the largest shear stress for which shear thinning occurs, the rheology is dominated by the proximity to the frictionless jamming transition point at ϕJ0, in the sense that rheological properties are (roughly) diverging functions of the form (ϕ−ϕJ0)−λ0, whereas at shear stresses above shear thickening, the frictional jamming transition at ϕJμ dominates and rheological properties scale with (ϕ−ϕJμ)−λμ, with λ₀ and λ_μ two positive exponents (whose exact values are still debated). Because ϕJμ<ϕJ0, this leads to shear thickening, which can be continuous or discontinuous depending on the proximity to ϕJμ (11, 25).In this work, we show that simulations of frictional colloidal suspensions can reproduce both continuous and discontinuous shear thickening, hence demonstrating that the stress-induced friction scenario extends to the Brownian case. Quantitative agreement with experiments cannot, however, be achieved with only frictional contacts and Brownian motion combined with hydrodynamic lubrication interactions. Instead, one must also consider the repulsive force that is induced between immersed colloids by the stabilization process. We study the qualitative influence of the range and amplitude of the repulsive force on the shear thickening. In particular, we show that the effects of the Brownian and stabilizing forces on the onset stress are additive. For a suitable choice of amplitude and range of the repulsive force, the simulations of the relative viscosity and the second normal stress difference (which is large relative to the first) are in good agreement with recent experiments by Cwalina and Wagner (26).     

Thermodynamic controls on rates of iron oxide reduction by extracellular electron shuttles

Anaerobic microbial respiration in suboxic and anoxic environments often involves particulate ferric iron (oxyhydr-)oxides as terminal electron acceptors. To ensure efficient respiration, a widespread strategy among iron-reducing microorganisms is the use of extracellular electron shuttles (EES) that transfer two electrons from the microbial cell to the iron oxide surface. Yet, a fundamental understanding of how EES–oxide redox thermodynamics affect rates of iron oxide reduction remains elusive. Attempts to rationalize these rates for different EES, solution pH, and iron oxides on the basis of the underlying reaction free energy of the two-electron transfer were unsuccessful. Here, we demonstrate that broadly varying reduction rates determined in this work for different iron oxides and EES at varying solution chemistry as well as previously published data can be reconciled when these rates are instead related to the free energy of the less exergonic (or even endergonic) first of the two electron transfers from the fully, two-electron reduced EES to ferric iron oxide. We show how free energy relationships aid in identifying controls on microbial iron oxide reduction by EES, thereby advancing a more fundamental understanding of anaerobic respiration using iron oxides.

The use of iron oxides as terminal electron acceptors in anaerobic microbial respiration is central to biogeochemical element cycling and pollutant transformations in many suboxic and anoxic environments (1–6). To ensure efficient electron transfer to solid-phase ferric iron, Fe(III), at circumneutral pH, metal-reducing microorganisms from diverse phylae use dissolved extracellular electron shuttle (EES), including quinones (7–9), flavins (10–16), and phenazines (17–19), to transfer two electrons per EES molecule from the respiratory chain proteins in the outer membrane of the microbial cell to the iron oxide (17, 20, 21). The oxidized EES can diffuse back to the cell surface for rereduction, thereby completing the catalytic redox cycle involving the EES.The electron transfer from the reduced EES to Fe(III) is considered a key step in overall microbial Fe(III) respiration. Several lines of evidence suggest that the free energy of the electron transfer reaction, ΔrG, controls Fe(III) reduction rates (15, 17, 22, 23). For instance, microbial Fe(III) oxide reduction by dissolved model quinones as EES was accelerated only for quinones with standard two-electron reduction potentials, EH,1,20, that fell into a relatively narrow range of −180±80 mV at pH 7 (24). Furthermore, in abiotic experiments, Fe(III) reduction rates by EES decreased with increasing ΔrG that resulted from increasing either EH,1,20 of the EES (25, 26), the concentration of Fe(II) in the system (27), or solution pH (25, 26, 28). However, substantial efforts to relate Fe(III) reduction rates for different EES species, iron oxides, and pH to the EH,1,20 averaged over both electrons transferred from the EES to the iron oxides were only partially successful (25, 28). Reaction free energies of complex redox processes involving the transfer of multiple electrons can readily be calculated using differences in the reduction potentials averaged over all electrons transferred, and this approach is well established in biogeochemistry and microbial ecology. For kinetic considerations, however, the use of averaged reduction potentials is inappropriate.Herein, we posit that rates of Fe(III) reduction by EES instead relate to the ΔrG of the less exergonic first one-electron transfer from the two-electron reduced EES species to the iron oxide, following the general notion that reaction rates scale with reaction free energies (29). Our hypothesis is based on the fact that, at circumneutral to acidic pH and for many EES, the reduction potential of the first electron transferred to the fully oxidized EES to form the one-electron reduced intermediate semiquinone species, EH,1, is lower than the reduction potential of the second electron transferred to the semiquinone to form the fully two-electron reduced EES species, EH,2 [i.e., EH,1<EH,2 (30–33)]. This difference in one-electron reduction potentials implies that the two-electron reduced EES (i.e., the hydroquinone) is the weaker one-electron reductant for Fe(III) as compared to the semiquinone species. We therefore expect that rates of iron oxide reduction relate to the ΔrG of the first electron transferred from the hydroquinone to Fe(III). The ΔrG of this first electron transfer may even be endergonic provided that the two-electron transfer is exergonic.We verified our hypothesis in abiotic model systems by demonstrating that reduction rates of two geochemically important crystalline iron oxides, goethite and hematite, by two-electron reduced quinone- and flavin-based EES over a wide pH range, and therefore thermodynamic driving force for Fe(III) reduction, correlate with the ΔrG of the first electron transferred from the fully reduced EES to Fe(III). We further show that rates of goethite and hematite reduction by EES reported in the literature are in excellent agreement with our rate data when comparing rates on the basis of the thermodynamics of the less exergonic first of the two electron transfers.     

Submesoscale dispersion in the vicinity of the Deepwater Horizon spill

Reliable forecasts for the dispersion of oceanic contamination are important for coastal ecosystems, society, and the economy as evidenced by the Deepwater Horizon oil spill in the Gulf of Mexico in 2010 and the Fukushima nuclear plant incident in the Pacific Ocean in 2011. Accurate prediction of pollutant pathways and concentrations at the ocean surface requires understanding ocean dynamics over a broad range of spatial scales. Fundamental questions concerning the structure of the velocity field at the submesoscales (100 m to tens of kilometers, hours to days) remain unresolved due to a lack of synoptic measurements at these scales. Using high-frequency position data provided by the near-simultaneous release of hundreds of accurately tracked surface drifters, we study the structure of submesoscale surface velocity fluctuations in the Northern Gulf of Mexico. Observed two-point statistics confirm the accuracy of classic turbulence scaling laws at 200-m to 50-km scales and clearly indicate that dispersion at the submesoscales is local, driven predominantly by energetic submesoscale fluctuations. The results demonstrate the feasibility and utility of deploying large clusters of drifting instruments to provide synoptic observations of spatial variability of the ocean surface velocity field. Our findings allow quantification of the submesoscale-driven dispersion missing in current operational circulation models and satellite altimeter-derived velocity fields.The Deepwater Horizon (DwH) incident was the largest accidental oil spill into marine waters in history with some 4.4 million barrels released into the DeSoto Canyon of the northern Gulf of Mexico (GoM) from a subsurface pipe over ∼84 d in the spring and summer of 2010 (1). Primary scientific questions, with immediate practical implications, arising from such catastrophic pollutant injection events are the path, speed, and spreading rate of the pollutant patch. Accurate prediction requires knowledge of the ocean flow field at all relevant temporal and spatial scales. Whereas ocean general circulation models were widely used during and after the DwH incident (2–6), such models only capture the main mesoscale processes (spatial scale larger than 10 km) in the GoM. The main factors controlling surface dispersion in the DeSoto Canyon region remain unclear. The region lies between the mesoscale eddy-driven deep water GoM (7) and the wind-driven shelf (8) while also being subject to the buoyancy input of the Mississippi River plume during the spring and summer months (9). Images provided by the large amounts of surface oil produced in the DwH incident revealed a rich array of flow patterns (10) showing organization of surface oil not only by mesoscale straining into the loop current “Eddy Franklin,” but also by submesoscale processes. Such processes operate at spatial scales and involve physics not currently captured in operational circulation models. Submesoscale motions, where they exist, can directly influence the local transport of biogeochemical tracers (11, 12) and provide pathways for energy transfer from the wind-forced mesoscales to the dissipative microscales (13–15). Dynamics at the submesoscales have been the subject of recent research (16–20). However, the investigation of their effect on ocean transport has been predominantly modeling based (13, 21–23) and synoptic observations, at adequate spatial and temporal resolutions, are rare (24, 25). The mechanisms responsible for the establishment, maintenance, and energetics of such features in the Gulf of Mexico remain unclear.Instantaneous measurement of all representative spatiotemporal scales of the ocean state is notoriously difficult (26). As previously reviewed (27), traditional observing systems are not ideal for synoptic sampling of near-surface flows at the submesoscale. Owing to the large spacing between ground tracks (28) and along-track signal contamination from high-frequency motions (29), gridded altimeter-derived sea level anomalies only resolve the largest submesoscale motions. Long time-series ship-track current measurements attain similar, larger than 2 km, spatial resolutions, and require averaging the observations over evolving ocean states (30). Simultaneous, two-point accoustic Doppler current profiler measurements from pairs of ships (25) provide sufficient resolution to show the existence of energetic submesoscale fluctuations in the mixed layer, but do not explicitly quantify the scale-dependent transport induced by such motions at the surface. Lagrangian experiments, centered on tracking large numbers of water-following instruments, provide the most feasible means of obtaining spatially distributed, simultaneous measurements of the structure of the ocean’s surface velocity field on 100-m to 10-km length scales.Denoting a trajectory by x(a, t), where x(a, t₀) = a, the relative separation of a particle pair is given by D(t,D0)=x(a1,t)−x(a2,t)=D0+∫t0tΔv(t′,D0)dt′, where the Lagrangian velocity difference is defined by Δv(t, D₀) = v(a₁, t) − v(a₂, t). The statistical quantities of interest, both practically and theoretically, are the scale-dependent relative dispersion D²(t) = 〈D ⋅ D〉 (averaged over particle pairs) and the average longitudinal or separation velocity, Δv(r), at a given separation, r. The velocity scale is defined by the second order structure function Δv(r)=〈δv2〉, where δv(r) = (v(x + r) − v(x)) ⋅ r/∥r∥ (31, 32) where the averaging is now conditioned on the pair separation r.The applicability of classical dispersion theories (32–34) developed in the context of homogeneous, isotropic turbulence with localized spectral forcing, to ocean flows subject to the effects of rotation, stratification, and complex forcing at disparate length and time scales remains unresolved. Turbulence theories broadly predict two distinct dispersion regimes depending upon the shape of the spatial kinetic energy spectrum, E(k) ∼ k^−β, of the velocity field (35). For sufficiently steep spectra (β ≥ 3) the dispersion is expected to grow exponentially, D ∼ e^λt with a scale-independent rate. At the submesoscales (∼ 100 m–10 km), this nonlocal growth rate will then be determined by the mesoscale motions currently resolved by predictive models. For shallower spectra (1 < β < 3), however, the dispersion is local, D ∼ t^2/(3−β), and the growth rate of a pollutant patch is dominated by advective processes at the scale of the patch. Accurate prediction of dispersion in this regime requires resolution of the advecting field at smaller scales than the mesoscale.Whereas compilations of data from dye measurements broadly support local dispersion in natural flows (36), the range of scales in any particular dye experiment is limited. A number of Lagrangian observational studies have attempted to fill this gap. LaCasce and Ohlmann (37) considered 140 pairs of surface drifters on the GoM shelf over a 5-y period and found evidence of a nonlocal regime for temporally smoothed data at 1-km scales. Koszalka et al. (38) using ??(100) drifter pairs with D₀ < 2 km launched over 18 mo in the Norwegian Sea, found an exponential fit for D²(t) for a limited time (t = 0.5 − 2 d), although the observed longitudinal velocity structure function is less clearly fit by a corresponding quadratic. They concluded that a nonlocal dispersion regime could not be identified. In contrast, Lumpkin and Elipot (39) found evidence of local dispersion at 1-km scales using 15-m drogued drifters launched in the winter-time North Atlantic. It is not clear how the accuracy of the Argos positioning system (150–1,000 m) used in these studies affects the submesoscale dispersion estimates. Schroeder et al. (40), specifically targeting a coastal front using a multiscale sampling pattern, obtained results consistent with local dispersion, but the statistical significance (maximum 64 pairs) remained too low to be definitive.     

Glasses denser than the supercooled liquid

When aged below the glass transition temperature, Tg, the density of a glass cannot exceed that of the metastable supercooled liquid (SCL) state, unless crystals are nucleated. The only exception is when another polyamorphic SCL state exists, with a density higher than that of the ordinary SCL. Experimentally, such polyamorphic states and their corresponding liquid–liquid phase transitions have only been observed in network-forming systems or those with polymorphic crystalline states. In otherwise simple liquids, such phase transitions have not been observed, either in aged or vapor-deposited stable glasses, even near the Kauzmann temperature. Here, we report that the density of thin vapor-deposited films of N,N′-bis(3-methylphenyl)-N,N′-diphenylbenzidine (TPD) can exceed their corresponding SCL density by as much as 3.5% and can even exceed the crystal density under certain deposition conditions. We identify a previously unidentified high-density supercooled liquid (HD-SCL) phase with a liquid–liquid phase transition temperature (TLL) ∼35 K below the nominal glass transition temperature of the ordinary SCL. The HD-SCL state is observed in glasses deposited in the thickness range of 25 to 55 nm, where thin films of the ordinary SCL have exceptionally enhanced surface mobility with large mobility gradients. The enhanced mobility enables vapor-deposited thin films to overcome kinetic barriers for relaxation and access the HD-SCL state. The HD-SCL state is only thermodynamically favored in thin films and transforms rapidly to the ordinary SCL when the vapor deposition is continued to form films with thicknesses more than 60 nm.

Glasses are formed when the structural relaxations in supercooled liquids (SCLs) become too slow, causing the system to fall out of equilibrium at the glass transition temperature (Tg). The resulting out-of-equilibrium glass state has a thermodynamic driving force to evolve toward the SCL state through physical aging (1). At temperatures just below Tg, the extent of equilibration is limited by the corresponding SCL state, while at much lower temperatures, equilibration is limited by the kinetic barriers for relaxation. As such, the degree of thermodynamic stability achieved through physical aging is limited (2).Physical vapor deposition (PVD) is an effective technique to overcome kinetic barriers for relaxation to produce thermodynamically stable glasses (3–10). The accelerated equilibration in these systems is due to their enhanced surface mobility (11–14). During PVD, when the substrate temperature is held below Tg, molecules or atoms can undergo rearrangements and adopt more stable configurations at the free surface and proximate layers underneath (13). After the molecules are buried deeper into the film, their relaxation dynamics significantly slow down, which prevents further equilibration. Through this surface-mediated equilibration process, stable glasses can achieve low-energy states on the potential energy landscape that would otherwise require thousands or millions of years of physical aging (2, 3, 15, 16).As such, the degree of enhanced surface mobility and mobility gradients are critical factors in the formation of stable glasses (3, 11, 17, 18). While the effect of film thickness on the surface mobility and gradients of liquid-quenched (LQ) glasses has been studied in the past (19, 20), there are limited data on the role of film thickness in the stability of vapor-deposited glasses. In vapor-deposited toluene, it has been shown that decreasing the film thickness from 70 to 5 nm can increase the thermodynamic stability but decrease the apparent kinetic stability (5, 6). In contrast, thin films covered with a top layer of another material do not show a significant evidence of reduced kinetic stability (21), indicating the nontrivial role of mobility gradients in thermal and kinetic stability.Stable glasses of most organic molecules, with short-range intramolecular interactions, have properties that are indicative of the same corresponding metastable SCL state as LQ and aged glasses, without any evidence of the existence of generic liquid–liquid phase transitions that can potentially provide a resolution for the Kauzmann entropy crisis (22). The Kauzmann crisis occurs at the Kauzmann temperature (TK), where the extrapolated SCL has the same structural entropy as the crystal, producing thermodynamically impossible states just below this temperature. Recently, Beasley et al. (16) showed that near-equilibrium states of ethylbenzene can be produced using PVD down to 2 K above TK and hypothesized that any phase transition to an “ideal glass” state to avoid the Kauzmann crisis must occur at TK.In some glasses of elemental substances (23, 24) and hydrogen-bonding compounds (25, 26), liquid–liquid phase transitions can occur between polyamorphic states with distinct local packing structures that correspond to polymorphic crystalline phases. For example, at high pressures, high- and low-density supercooled water phases are interconvertible through a first-order phase transition (27, 28). Recent studies have demonstrated that such polyamorphic states can also be accessed through PVD in hydrogen-bonding systems with polymorphic crystal states at depositions above the nominal Tg (29, 30). However, these structure-specific transitions do not provide a general resolution for the Kauzmann crisis.Here, we report the observation of a liquid–liquid phase transition in vapor-deposited thin films of N,N′-bis(3-methylphenyl)-N,N′-diphenylbenzidine (TPD). TPD is a molecular glass former with only short-range intermolecular interactions. When thin films of TPD are vapor deposited onto substrates held at deposition temperatures (Tdep) below the nominal glass transition temperature of bulk TPD, Tg (bulk), films in the thickness range of 25 nm<h<55 nm achieve a high-density supercooled liquid (HD-SCL) state, which has not been previously observed. The liquid–liquid phase transition temperature (TLL) between the ordinary SCL and HD-SCL states is measured to be TLL≃Tg(bulk)−35 K. The density of thin films deposited below TLL tangentially follows the HD-SCL line, which has a stronger temperature dependence than the ordinary SCL. When vapor deposition is continued to produce thicker films (h>60 nm), the HD-SCL state transforms into the ordinary SCL state, indicating that the HD-SCL is only thermodynamically favored in the thin-film geometry. This transition is qualitatively different from the previously reported liquid–liquid phase transitions, as it is not related to a specific structural motif in TPD crystals, and it can only be observed in thin films, indicating that the energy landscape of thin films is favoring this high-density state.We observe an apparent correlation between enhanced mobility gradients in LQ thin films of TPD and the thickness range where HD-SCL states are produced during PVD. We hypothesize that enhanced mobility gradients are essential in providing access to regions of the energy landscape corresponding to the HD-SCL state, which are otherwise kinetically inaccessible. This hypothesis should be further investigated to better understand the origin of this phenomenon.     

Understanding the kinetic mechanism of RNA single base pair formation

RNA functions are intrinsically tied to folding kinetics. The most elementary step in RNA folding is the closing and opening of a base pair. Understanding this elementary rate process is the basis for RNA folding kinetics studies. Previous studies mostly focused on the unfolding of base pairs. Here, based on a hybrid approach, we investigate the folding process at level of single base pairing/stacking. The study, which integrates molecular dynamics simulation, kinetic Monte Carlo simulation, and master equation methods, uncovers two alternative dominant pathways: Starting from the unfolded state, the nucleotide backbone first folds to the native conformation, followed by subsequent adjustment of the base conformation. During the base conformational rearrangement, the backbone either retains the native conformation or switches to nonnative conformations in order to lower the kinetic barrier for base rearrangement. The method enables quantification of kinetic partitioning among the different pathways. Moreover, the simulation reveals several intriguing ion binding/dissociation signatures for the conformational changes. Our approach may be useful for developing a base pair opening/closing rate model.RNAs perform critical cellular functions at the level of gene expression and regulation (1–4). RNA functions are determined not only by RNA structure or structure motifs [e.g., tetraloop hairpins (5, 6)] but also by conformational distributions and dynamics and kinetics of conformational changes. For example, riboswitches can adopt different conformations in response to specific conditions of the cellular environment (7, 8). Understanding the kinetics, such as the rate and pathways for the conformational changes, is critical for deciphering the mechanism of RNA function (9–19). Extensive experimental and theoretical studies on RNA folding kinetics have provided significant insights into the kinetic mechanism of RNA functions (19–36). However, due to the complexity of the RNA folding energy landscape (37–46) and the limitations of experimental tools (47–55), many fundamental problems, including single base flipping and base pair formation and fraying, remain unresolved. These unsolved fundamental problems have hampered our ability to resolve other important issues, such as RNA hairpin and larger structure folding kinetics. Several key questions remain unanswered, such as whether the hairpin folding is rate-limited by the conformational search of the native base pairs, whose formation leads to fast downhill folding of the whole structure, or by the breaking of misfolded base pairs before refolding to the native structure (18, 19, 54–73).Motivated by the need to understand the basic steps of nucleic acids folding, Hagan et al. (74) performed forty-three 200-ps unfolding trajectories at 400 K and identified both on- and off-pathway intermediates and two dominant unfolding pathways for a terminal C-G base pair in a DNA duplex. In one of the pathways, base pairing and stacking interactions are broken concomitantly, whereas in the other pathway, base stacking is broken after base pairing is disrupted. Furthermore, the unfolding requires that the Cyt diffuse away from the pairing Gua to a distance such that the C-G hydrogen bond cannot reform easily. More recently, Colizzi and Bussi (75) performed molecular dynamics (MD) pulling simulations for an RNA duplex and construct free energy landscape from the pulling simulation. The simulation showed that the base pair opening reaction starts with the unbinding of the 5′-base, followed by the unbinding of the 3′-base (i.e., the 5′-base is less stable than the 3′-base). These previous unfolding simulations offered significant insights into the pathways and transition states. However, as shown below, several important issues remain.One intriguing problem is the rate model for base pairing. There are currently three main types of models. In the first type of model, the barrier ΔG+‡ for closing a base pair is dominated by the entropic cost ΔS for positioning the nucleotides to the base-paired configuration and the barrier ΔG−‡ for opening a base pair is the enthalpic cost ΔH for disrupting the hydrogen bonds and base stacking interactions (18, 59, 60). In the second type of model, ΔG+‡ is the net free energy change for base pairing ΔG = ΔH ? TΔS and ΔG−‡ is zero (76, 77). In the third type of model, ΔG±‡=±ΔG/2 is used (78). In addition to the above three main types, other models, such as more sophisticated hybrid rate models, have been proposed (29).In this paper, we report a hybrid method (see Fig. 1) to investigate the single base pairing process. In contrast to the previous simulations for temperature- or force-induced unfolding reactions, we directly model the folding process here (i.e., the base pair closing process). Specifically, we use MD simulations to identify the conformational clusters. Based on the network of the conformational clusters as a reduced conformational ensemble, we apply kinetic Monte Carlo (KMC) and master equation (ME) methods to elucidate the detailed roles of base pairing and stacking interactions, as well as the roles of water and ions (79–82). The study reveals previously unidentified kinetics pathways, misfolded states, and rate-limiting steps. A clear understanding of the microscopic details of the elementary kinetic move is a prerequisite for further rigorous study of large-scale RNA kinetic studies. The method described here may provide a feasible way to develop a rate model for the base pair/stack-based kinetic move set. Furthermore, the mechanism of RNA single base folding may provide useful insights into many biologically significant processes, such as nucleotide flipping (83) in helicases and base pair fraying (84) (as the possible first step for nucleic duplex melting in nucleic acid enzymatic processes).Open in a separate windowFig. 1.(A) Folding of a single nucleotide (G1, red) from the unfolded (Left) to the native folded (Right) state. (B) Exhaustive sampling for the (discrete) conformations of the G1 nucleotide (Right) through enumeration of the torsion angles (formed by the blue bonds). (C) Schematic plot shows the trajectories on the energy landscape (depicted with two reaction coordinates for clarity) explored by the MD simulations. The lines, open circles, and hexagons denote the trajectories; the initial states; and the (centroid structures of the) clusters, respectively. (D) Conformational network based on six clusters. (E) The rmsds to the different clusters provide information about the structural changes in a MD trajectory.     

Specific heat and sound velocity at the relevant competing phase of high-temperature superconductors

Recent highly accurate sound velocity measurements reveal a phase transition to a competing phase in YBa₂Cu₃O_6+δ that is not identified in available specific heat measurements. We show that this signature is consistent with the universality class of the loop current-ordered state when the free-energy reduction is similar to the superconducting condensation energy, due to the anomalous fluctuation region of such a transition. We also compare the measured specific heat with some usual types of transitions, which are observed at lower temperatures in some cuprates, and find that the upper limit of the energy reduction due to them is about 1/40th the superconducting condensation energy.A significant step forward toward understanding high-temperature superconductivity is the variety of experimental results that have led to the widespread acceptance of the idea (1) that a phase with a broken symmetry competes with superconductivity in the underdoped region, often called the pseudogap region. There are a plethora of suggested phases (2–15). However, experimental results (16) consistent with transition to only one of them are observed at the pseudogap temperature T^?, where major changes in transport and thermodynamics occur in all cuprates. Charge density waves (CDWs) are observed in some compounds at lower temperatures. The idea that a broken symmetry phase competes with superconductivity makes thermodynamic sense only if the energy gained due to it is comparable to that gained through the superconducting transition in their coexistence region. The energy gained is related to the specific heat associated with the transition. Extraordinarily, however, no specific heat signature of a phase transition has been identified in the available measurements at the pseudogap temperature T^?. The much more accurately measured sound velocity singularity near the transition temperature is proportional to the heat capacity and can be used to find the symmetry class of the phase transition. In this paper, we use the recent highly accurate sound velocity measurements (17) and the best available specific heat measurements in YBa₂Cu₃O_6+δ (18–20) to show that phase transitions to the universality class of the loop current-ordered state with free-energy reduction similar to the measured superconducting condensation are consistent with the sound velocity and with lack of identifiable observation in the specific heat.Sound velocity changes near a phase transition, as shown below, are proportional to ?γ(T) = ?C_v/T, where C_v is the specific heat, if the transition temperature depends linearly on the strain. Since they are measured with a factor of O(10^?2) greater accuracy than the deductions from the best available specific heat measurements, they can be used to decipher the universality class to which the transition at T^? belongs and therefore their specific heat. The free-energy reduction can be calculated from the specific heat expected for the class of broken symmetry and its magnitude.The free energy due to a phase transition at temperature T_x is a homogeneous function of (T ? T_x). The elastic constants of a solid are given by the second derivative of the free energy with respect to the relevant strain. Therefore, the isothermal sound velocity variation δc_λ in a polarization λ associated with the phase transition, normalized to the background smoothly varying sound velocity c_0λ for δc_λ ? c_0λ, is given byδcλ(T−Tx)c0λ=12ρc0λ2[−Cv(T−Tx)T(dTxduλ)2+S(T−Tx)d2Txduλ2].[1]C_v(T ? T_x) and S(T ? T_x) are the specific heat at constant volume and entropy associated with the part of the free energy associated with the transition at T_x, i.e., the part that is a homogeneous function of (T ? T_x). Here ρ is the density. In mean field phase transitions, such as the superconducting transition, this reduces to the relation commonly used. Noting that the second contribution above is much smoother than the first and typically (1/Tx)(dTxduλ)2 is similar or larger than d2Txduλ2, we need to consider only the first term. Comparing the sound velocity variations at two different transitions, one at T_c and the other at T^?,δcλ(T−T∗)δcλ(T−Tc)=Cv(T−T∗)Cv(T−Tc)(dT∗/duλdTc/duλ)2.[2]     

Complex pathways and memory in compressed corrugated sheets

The nonlinear response of driven complex materials—disordered magnets, amorphous media, and crumpled sheets—features intricate transition pathways where the system repeatedly hops between metastable states. Such pathways encode memory effects and may allow information processing, yet tools are lacking to experimentally observe and control these pathways, and their full breadth has not been explored. Here we introduce compression of corrugated elastic sheets to precisely observe and manipulate their full, multistep pathways, which are reproducible, robust, and controlled by geometry. We show how manipulation of the boundaries allows us to elicit multiple targeted pathways from a single sample. In all cases, each state in the pathway can be encoded by the binary state of material bits called hysterons, and the strength of their interactions plays a crucial role. In particular, as function of increasing interaction strength, we observe Preisach pathways, expected in systems of independently switching hysterons; scrambled pathways that evidence hitherto unexplored interactions between these material bits; and accumulator pathways which leverage these interactions to perform an elementary computation. Our work opens a route to probe, manipulate, and understand complex pathways, impacting future applications in soft robotics and information processing in materials.

The response of complex media to external driving is intermittent, featuring smooth reversible episodes, associated with a single (meta)stable state of the system, punctuated by sharp irreversible steps between states that together form a multistep pathway (1–7). These steps are typically hysteretic and for several systems, such as amorphous media, can be associated with local rearrangements that act as two-state degrees of freedom. The ensuing complex pathways are often modeled by collections of hysteretic, two-state elements called hysterons (8). These two-state hysteretic elements switch up and down between internal states s = 0 and s = 1 when a driving field U passes through the upper and lower switching fields U⁺ or U^–(with U+>U−); the state of the hysteron for U−≤U≤U+ depends on its driving history. One can think of these as material bits (9–12) that collectively label the (meso)state of the physical system. Properties such as memory are then determined by the sequences of bit switches as function of a global driving U, which can be encoded in so-called transition graphs (t graphs), whose nodes represent the mesostates and edges represent their transitions (13, 14).Collections of n uncoupled hysterons form the Preisach model (8), which has been studied extensively in the context of complex hysteresis and memory effects. The absence of coupling implies that hysteron i changes state at switching fields Ui+ and Ui−, independent of the state of the other hysterons. As a result, the sequence of bit switches in response to sweeping U is given by the ordering of the 2n switching fields. This restricts the type of pathways that are possible, with the t graphs featuring a hierarchical structure of loops within loops and exhibiting return point memory (RPM), the widespread ability of complex systems to remember their extremal driving, i.e., to return to a previous state when the driving revisits an extremum (15–19).However, interactions between hysterons can break the no-passing (NP) property that underlies RPM (5, 18). Recent simulations of models of interacting hysterons, as well as amorphous media, have presented examples for complex pathways and transition graphs featuring, e.g., avalanches, transient memories, and multiperiodic orbits, which cannot be captured by models of noninteracting hysterons (5, 20–22). Unfortunately, distinguishing, observing, and manipulating individual hysterons and their interactions is experimentally challenging for most complex systems. Moreover, we lack a conceptual framework that organizes the distinct impacts that hysteron interactions have on the phenomenology. Hence, both the connection between hysteron models and experimentally observable pathways and the relevance of hysteron interactions for driven complex media remain unclear.Here we introduce mechanical compression of curved, corrugated elastic sheets to directly observe mechanical hysterons, their interactions, and their concomitant nontrivial pathways (Fig. 1). We experimentally observe that the driving value where a given hysteron switches is modified by the states of the other hysterons, thus evidencing interactions between hysterons. To organize the resulting phenomenology, we distinguish between two characteristics of the pathways that are impacted by hysteron interactions. Most strikingly, interactions can modify the topology of the transition graph, and we identify the first three steps in a hierarchy of increasingly complex t graphs and give concrete examples of each. In addition, we show that even for a given t-graph topology, interactions can have a more subtle effect depending on the precise ordering of the switching fields. The strict hierarchy of t-graph topologies and the more subtle effects of the relative ordering are experimentally observable and testable and provide a conceptual tool to organize the plethora of pathways observed in driven frustrated matter. Together, our work shows how hysteron interactions bring sequences of bit flips that encode forms of information processing within reach, creating opportunities for soft robotics (23–26) and information processing in materials (9–11, 27).Open in a separate windowFig. 1.Robust pathways in a cyclically compressed corrugated sheet. (A) Our samples are corrugated elastic cylindrical shells of height H, thickness t, radius of curvature R, and N sinusoidal corrugations of pitch p and amplitude A₀. For sample A, shown here, N = 3 and {H,t,R,p,A0}={35,0.2,1.0,8,3} mm. (B) Compression U, force F, and bottom plate tilt angle α. (C) Upon compression at α = 0 mrad, sample A reaches four different mesostates, associated with sudden snapping of distinct regions Ω_i (colored strips; Movie S1). (D) The force F exhibits three sharp jumps during compression (red) and decompression (blue). (E) Each force jump is associated with a sudden deformation, evidenced by spikes in Δi2, the sum of the squared differences between subsequent digital images restricted to Ω_i (Materials and Methods). (F) Schematic representation of a hysteron, its two states (here gray corresponds to state 1, white corresponds to state 0, and dashed corresponds to the hysteretic range where the state is either 0 or 1 depending on the history), and the switching fields U±, where U−<U<U+ is the hysteretic range. (G) Evolution of hysteron state during a compression cycle. (H) Our samples behave as collections of parallel hysterons with distinct thresholds. (I) Force–displacement curve corresponding to increasing compression cycles (Inset). The mechanical response of the system features connected hysteresis loops and multiple pathways. (J) The transition graph of sample A at α = 0 mrad contains four states (nodes) labeled by the state of each hysteron. Red (blue) arrows correspond to up (down) transitions at (de)compression U_c as indicated in mm.