Stagnation points control chaotic fluctuations in viscoelastic porous media flow

Viscoelastic flows through porous media become unstable and chaotic beyond critical flow conditions, impacting widespread industrial and biological processes such as enhanced oil recovery and drug delivery. Understanding the influence of the pore structure or geometry on the onset of flow instability can lead to fundamental insights into these processes and, potentially, to their optimization. Recently, for viscoelastic flows through porous media modeled by arrays of microscopic posts, Walkama et al. [D. M. Walkama, N. Waisbord, J. S. Guasto, Phys. Rev. Lett. 124, 164501 (2020)] demonstrated that geometric disorder greatly suppressed the strength of the chaotic fluctuations that arose as the flow rate was increased. However, in that work, disorder was only applied to one originally ordered configuration of posts. Here, we demonstrate experimentally that, given a slightly modified ordered array of posts, introducing disorder can also promote chaotic fluctuations. We provide a unifying explanation for these contrasting results by considering the effect of disorder on the occurrence of stagnation points exposed to the flow field, which depends on the nature of the originally ordered post array. This work provides a general understanding of how pore geometry affects the stability of viscoelastic porous media flows.

Unlike viscous Newtonian liquids (e.g., water), many fluids exhibit an elastic response to an applied strain. Such “viscoelastic” fluids are widespread in biology (blood, mucus, synovial fluid) and industry (paints, coatings, fracking fluids). The elasticity is imparted by the presence of a microstructure (formed by, e.g., polymers, proteins, or self-assemblies of lipids or surfactants) that relaxes after deformation (1). The strength of the elastic response of the fluid to an imposed deformation (or flow) is quantified by the Weissenberg number Wi=τγ˙, with τ the fluid relaxation time and γ˙ the rate of strain. While flows of Newtonian fluids become unstable and turbulent due to the onset of inertial effects at high Reynolds number, Re≫1, viscoelastic flows can become unstable and exhibit so-called “elastic turbulence” even for Re≪1, purely due to elastic effects that arise at high Wi (2–7).Viscoelastic porous media flow occurs in diverse processes ranging from enhanced oil recovery (EOR) and filtration to drug delivery (8, 9). Porous media flow subjects a fluid to a complex cycle of deformation with high shear rates through the pore throats or between obstacles and high elongational rates at points of constriction or at stagnation points, leading to stretching of the fluid microstructure if Wi≳1 (10, 11). Stagnation points (which occur at the front and rear poles of obstacles in a flow) are particularly effective at causing high stretching and large tensile stresses due to the combination of zero flow velocity and finite velocity gradient that exists in such regions (10–15). Elastic tensile stresses due to stretching on curvilinear streamlines (as through porous media) are conditions well established to lead to linear instabilities in viscoelastic fluids (16–19), which can be precursors to elastic turbulence as Wi is further increased (7, 20–23). The chaotic fluctuations that result are expected to greatly enhance the pressure loss and the dispersion in porous media, with positive impacts on, for example, removing oil ganglia from the pore space in EOR or improving the distribution of drugs throughout a tumor (24–27).There have been various recent advances in modeling viscoelastic porous media flows both experimentally and numerically (14, 27–36). However, the complexity of the problem has limited numerical simulations to extremely simplified regular geometries (32, 36) and/or small computational domains and/or regimes of low Wi (31). Experimentally, Browne and coworkers (27, 33) have achieved the detailed characterization of the pore-scale dynamics in model porous media formed by three-dimensional (3D) random packings of spherical glass particles, correlating a global increase in the pressure drop across the media with the onset of elastic turbulence in the pores. Importantly, due to the complexity of the random sphere packings of Browne and coworkers (27, 33), fluid arriving at each pore experiences a unique flow history, and the flow through different pores becomes unstable at different values of the nominal Wi (computed based on macroscopic flow conditions). Fundamental questions remain over how the details of the pore-space geometry affect the onset and strength of the chaotic fluctuations that arise.While randomly packed beds of polydisperse spheres provide a good model for the complex pore geometries that arise in real media such as sandstone or carbonate rock (37, 38), ordered and regular geometries enable investigation of the role of different packing structures and hence, pore shape (29, 39, 40). This is most conveniently achieved by arrangements of posts forming either linear models of the interconnecting capillary network through the pore space (e.g., refs. 36, 41, and 42) or two-dimensional (2D) arrays that represent the tortuous flow paths around closely spaced grains (e.g., refs. 14, 24, 25, 28, 30, 35, 43, and 44).An outstanding open question concerns how the chaotic dynamics of viscoelastic flows are affected when geometric disorder (inherent in real heterogeneous systems) is introduced to a regular model porous medium. In a recent attempt to address this issue, Walkama et al. (28) performed experiments in a series of 2D microfluidic post arrays using shear thinning viscoelastic polymeric test solutions. They examined how the introduction of increasing random disorder to a hexagonal post array (arranged as shown in Fig. 1A) affected the onset and strength of the chaotic fluctuations observed for Wi≳1. Their results led to the broad general conclusion that “disorder suppresses chaos in viscoelastic flows.” However, other works have shown that instabilities and fluctuations in viscoelastic flows through 2D ordered post arrays strongly depend on the orientation of the array relative to the flow direction (30). Thus, different behavior might be anticipated from an ordered array of posts that are staggered along the flow direction (Fig. 1A) (as employed in ref. 28) than from an identical array rotated by 30° such that the posts become aligned (Fig. 1B). Indeed, as shown in Fig. 1 C and ​andD,D, even the low-Re flow of a simple Newtonian fluid shows qualitatively different flow patterns in the two contrasting post arrangements. Notably, in Fig. 1C, it is clear that each post presents both an upstream point and a downstream stagnation point that are accessible to the flow field. However, in the rotated arrangement in Fig. 1D, the flow is concentrated between the aligned rows of posts, largely bypassing the stagnation points. Given the known role of stagnation point regions in driving the onset of instabilities and fluctuations in viscoelastic flows (e.g., refs. 15, 17, and 45–49), we question the generality of the conclusions drawn by Walkama et al. (28), based on modifications made to a single-ordered geometry like in Fig. 1 A and ​andCC.Open in a separate windowFig. 1.(A and B) Unit cell representations of two contrasting ordered hexagonal arrays of posts used in the flow experiments. In A, the posts are staggered along the x direction in which the flow is imposed. The post radius is R, and lattice spacing is S. Rotating the array by 30° aligns the posts in the flow direction (B). Disordered aligned arrays are generated by the random displacement of each post within a hexagon of circumradius βS, as described in ref. 28. (C and D) Streamlines determined by flow velocimetry (Materials and Methods) with a Newtonian fluid in the staggered (C) and aligned (D) arrays at Re≈10−3. The red crossed circles in C and D indicate the locations of the leading and trailing-edge stagnation points on one of the circular posts.Here, we show by microfluidic experiments with a viscoelastic wormlike micelle (WLM) solution that a rotation of the hexagonal post array in Fig. 1 A and ​andCC in order to align the posts with the flow direction (Fig. 1 B and ​andD)D) strongly suppresses the chaotic fluctuations for a range of Wi≳1, consistent with our expectation based on the removal of stagnation points. Subsequently, following the methods of Walkama et al. (28), we introduce random disorder to the aligned array of posts (Fig. 1B). In this case, contrary to Walkama et al. (28), disorder does not further suppress but rather, promotes chaotic fluctuations over a wide range of Wi. Although our results appear to contradict those recently reported in Walkama et al. (28), both are simply explained by considering how disorder affects the prominence of stagnation points in the flow field (which is opposite, depending on the originally ordered geometric arrangement). Furthermore, we significantly extend the range of imposed Wi beyond that studied by Walkama et al. (28), showing that at sufficiently high Wi, the nature of the flow fluctuations becomes essentially geometry independent. Our work reaches an intuitive and general understanding of the role of geometry (specifically the importance of stagnation points) in controlling the onset and strength of chaotic fluctuations in viscoelastic porous media flows.     

Nonlinear rotational spectroscopy reveals many-body interactions in water molecules

Because of their central importance in chemistry and biology, water molecules have been the subject of decades of intense spectroscopic investigations. Rotational spectroscopy of water vapor has yielded detailed information about the structure and dynamics of isolated water molecules, as well as water dimers and clusters. Nonlinear rotational spectroscopy in the terahertz regime has been developed recently to investigate the rotational dynamics of linear and symmetric-top molecules whose rotational energy levels are regularly spaced. However, it has not been applied to water or other lower-symmetry molecules with irregularly spaced levels. We report the use of recently developed two-dimensional (2D) terahertz rotational spectroscopy to observe high-order rotational coherences and correlations between rotational transitions that were previously unobservable. The results include two-quantum (2Q) peaks at frequencies that are shifted slightly from the sums of distinct rotational transitions on two different molecules. These results directly reveal the presence of previously unseen metastable water complexes with lifetimes of 100 ps or longer. Several such peaks observed at distinct 2Q frequencies indicate that the complexes have multiple preferred bimolecular geometries. Our results demonstrate the sensitivity of rotational correlations measured in 2D terahertz spectroscopy to molecular interactions and complexation in the gas phase.

Water has attracted extensive spectroscopic interest because of its critical implications for theoretical and applied sciences (1, 2). Water shows anomalous properties because of complicated fluxional hydrogen-bonded networks and has been investigated by Raman and infrared spectroscopy (3, 4), sum frequency spectroscopy (5), optical Kerr effect spectroscopy (6), vibration-rotation-tunneling spectroscopy (7), and recently by two-dimensional (2D) infrared spectroscopy (8, 9) and 2D Raman-terahertz (THz) spectroscopy (10). In the gas phase, water is of utmost importance for atmospheric science, astrophysics, combustion research, and fundamental chemistry and physics (1, 11, 12). Although the pure rotation spectrum of water vapor has been well known for decades (13, 14), nonlinear THz spectroscopy of water rotational dynamics has not been previously reported. Nonlinear rotational spectroscopy in the microwave spectral range is well established (15), but because of the small moments of inertia of water, most of its rotational transition frequencies lie in the THz frequency range (Fig. 1). Nonlinear THz rotational spectroscopy was reported only recently (16, 17), and 2D THz rotational spectroscopy (10, 18, 19) more recently still. As in 2D spectroscopy of vibrational, electronic, and other degrees of freedom (9, 20–25), 2D rotational spectroscopy can reveal correlations between rotational states, many-body effects, and distinct multiple-field interactions that cannot be observed by linear spectroscopy (26). The large dipole moment of water, manifest in strong atmospheric absorption in the THz window (27–31), and the existence of water dimers and larger clusters with complex structures and dynamics (7, 32, 33), suggest that 2D spectroscopy of water could generate previously elusive insights (34–37).Open in a separate windowFig. 1.Overview of the experiment. (A) Water molecule in the laboratory frame, showing the dipole moment μ at an angle θ from the Z axis (THz polarization direction). (B) Water molecule in the molecule-fixed frame with the three moments of inertia I_a, I_b, and I_c along the corresponding axes. (C) Relative population distribution as a function of the J rotational quantum number. All relevant K_a and K_c components are included in the population distribution. (Inset) The relative population distribution within the state J = 2. (D) Rotational energy levels of para-H₂O and ortho-H₂O molecules. Red arrows illustrate rotational transitions and transition frequencies involved in this work. (E) Measured THz FID (Top) and Fourier transform (Bottom) showing rotational transitions (marked by dashed vertical lines) of water vapor in ambient air. (F) Schematic illustration of the 2D THz experimental setup. Linear THz spectra (example in E) are measured with only one THz pump pulse.Two-dimensional THz rotational spectroscopy has not been extended previously to water or any asymmetric-top molecules, although such molecules, whose rotational spectra are complicated because all three of their moments of inertia are unequal, are the majority of naturally occurring molecular species. Unlike a linear or symmetric-top molecule, in which the spectroscopic transitions between successive total rotational angular momentum levels denoted by the quantum number J are spaced by even-integer multiples of a common factor, the rotational constant B, the asymmetric-top nature of water molecules leads to irregularly spaced rotational energy levels. These levels are described approximately by quantum numbers J, K_a, and K_c that indicate the total angular momentum and its symmetry-axis projections (1, 2). The spectrum of water vapor consists of many transitions, with ΔJ=−1,0,+1 (P, Q, and R branches) all allowed and, for each, changes ΔKa=±1,ΔKc=±1 (1, 2). The large centrifugal distortion of water and the distinct sets of rotational states occupied by its nuclear spin isomers (even symmetry for para, odd symmetry for ortho) further complicate its rotational spectrum (1, 2, 38). A typical THz time-domain free-induction decay (FID) signal from water vapor at ambient conditions, induced by a weak single-cycle THz pulse, is shown in Fig. 1E. The irregular oscillations arise from more than 15 transitions (Fig. 1D) that contribute significantly to the water vapor absorption spectrum in the 0.1- to 2-THz region.     

Electrophilic aromatic substitution reactions of compounds with Craig-Möbius aromaticity

Electrophilic aromatic substitution (EAS) reactions are widely regarded as characteristic reactions of aromatic species, but no comparable reaction has been reported for molecules with Craig-Möbius aromaticity. Here, we demonstrate successful EAS reactions of Craig-Möbius aromatics, osmapentalenes, and fused osmapentalenes. The highly reactive nature of osmapentalene makes it susceptible to electrophilic attack by halogens, thus osmapentalene, osmafuran-fused osmapentalene, and osmabenzene-fused osmapentalene can undergo typical EAS reactions. In addition, the selective formation of a series of halogen substituted metalla-aromatics via EAS reactions has revealed an unprecedented approach to otherwise elusive compounds such as the unsaturated cyclic chlorirenium ions. Density functional theory calculations were conducted to study the electronic effect on the regioselectivity of the EAS reactions.

Aromaticity, a core concept in chemistry, was initially introduced to account for the bonding, stability, reactivity, and other properties of many unsaturated organic compounds. There have been many elaborations and extensions of the concept of aromaticity (1, 2). The concepts of Hückel aromaticity and Möbius aromaticity are widely accepted (Fig. 1A). A π-aromatic molecule of the Hückel type is planar and has 4n + 2 conjugated π-electrons (n = 0 or an integer), whereas a Möbius aromatic molecule has one twist of the π-system, similar to that in a Möbius strip, and 4n π-electrons (3, 4). Since the discovery of naphthalene in 1821, aromatic chemistry has developed into a rich field and with a variety of subdisciplines over the course of its 200-y history, and the concept of aromaticity has been extended to other nontraditional structures with “cyclic delocalization of mobile electrons” (5). For example, benzene-like metallacycles—predicted by Hoffmann et al. as metallabenzenes—in which a metal replaces a C–H group in the benzene ring (6), have garnered extensive research interest from both experimentalists and theoreticians (7–12). As paradigms of the metalla-aromatic family, most complexes involving metallabenzene exhibit thermodynamic stability, kinetic persistence, and chemical reactivity associated with the classical aromaticity concept (13–15). Typically, like benzene, metallabenzene can undergo characteristic reactions of aromatics such as electrophilic aromatic substitution (EAS) reactions (16–18) (Fig. 1B, I) and nucleophilic aromatic substitution reactions (19–21).Open in a separate windowFig. 1.Schematic representations of aromaticity classification (A) and EAS reactions (B) of benzene, metallabenzene, and polycyclic metallacycles with Craig-Möbius aromaticity.The incorporation of transition metals has also led to an increase in the variety of the aromatic families (22–25). We have reported that stable and highly unusual bicyclic systems, metallapentalenes (osmapentalenes), benefit from Craig-Möbius aromaticity (26–30). In contrast to other reported Möbius aromatic compounds with twisted topologies, which are known as Heilbronner-Möbius aromatics (31–34), the involvement of transition metal d orbitals in π-conjugation switches the Hückel anti-aromaticity of pentalene into the planar Craig-Möbius aromaticity of metallapentalene (35–38) (Fig. 1A, III). Both the twisted topology and the planar Craig-Möbius aromaticity are well established and have been accepted as reasonable extensions of aromaticity (39–43). There has been no experimental evidence, however, as to whether these Möbius aromatic molecules can undergo classical aromatic substitution reactions, such as EAS reactions, instead of addition reactions. Given the key role of EAS in aromatic chemistry to obtain various derivatives, we sought to extend the understanding of the reactivity paradigm in the metalla-aromatic family.Our recent synthetic efforts associated with the metallapentalene system prompted us to investigate whether typical EAS reactions could proceed in these Craig-Möbius aromatics. If so, how could substitution be achieved in the same way that it is with traditional Hückel aromatics such as benzenes? In this paper, we present EAS reactions, mainly the halogenation of osmapentalene, osmafuran-fused osmapentalene, and osmabenzene-fused osmapentalene, which follow the classic EAS mechanistic scheme (Fig. 1B). With the aid of density functional theory (DFT) calculations, we characterized the effects on EAS reactivity and regioselectivity.     

Conjugated polymers based on metalla-aromatic building blocks

Conjugated polymers usually require strategies to expand the range of wavelengths absorbed and increase solubility. Developing effective strategies to enhance both properties remains challenging. Herein, we report syntheses of conjugated polymers based on a family of metalla-aromatic building blocks via a polymerization method involving consecutive carbyne shuttling processes. The involvement of metal d orbitals in aromatic systems efficiently reduces band gaps and enriches the electron transition pathways of the chromogenic repeat unit. These enable metalla-aromatic conjugated polymers to exhibit broad and strong ultraviolet–visible (UV–Vis) absorption bands. Bulky ligands on the metal suppress π–π stacking of polymer chains and thus increase solubility. These conjugated polymers show robust stability toward light, heat, water, and air. Kinetic studies using NMR experiments and UV–Vis spectroscopy, coupled with the isolation of well-defined model oligomers, revealed the polymerization mechanism.

Conjugated polymers are macromolecules usually featuring a backbone chain with alternating double and single bonds (1–3). These characteristics allow the overlapping p-orbitals to form a system with highly delocalized π-electrons, thereby giving rise to intriguing chemical and physical properties (4–6). They have exhibited many applications in organic light-emitting diodes, organic thin film transistors, organic photovoltaic cells, chemical sensors, bioimaging and therapies, photocatalysis, and other technologies (7–10). To facilitate the use of solar energy, tremendous efforts have been devoted in recent decades to developing previously unidentified conjugated polymers exhibiting broad and strong absorption bands (11–13). The common strategies for increasing absorption involve extending π-conjugation by incorporating conjugated cyclic moieties, especially fused rings; modulating the strength of intramolecular charge transfer between donor and acceptor units (D–A effect); increasing the coplanarity of π conjugation through weak intramolecular interactions (e.g., hydrogen bonds); and introducing heteroatoms or heavy atoms into the repeat units of conjugated polymers (11–16). Additionally, appropriate solubility is a prerequisite for processing and using polymers and is usually achieved with the aid of long alkyl or alkoxy side chains (12, 17).Aromatic rings are among the most important building blocks for conjugated polymers. In addition to aromatic hydrocarbons, a variety of aromatic heterocycles composed of main-group elements have been used as fundamental components. These heteroatom-containing conjugated polymers show unique optical and electronic properties (4–10). However, while metalla-aromatic systems bearing a transition metal have been known since 1979 due to the pioneering work by Thorn and Hoffmann (18), none of them have been used as building blocks for conjugated polymers. The HOMO–LUMO gaps (E_g) of metalla-aromatics are generally narrower (Fig. 1) than those of their organic counterparts (19–22). We reasoned that this feature should broaden the absorption window if polymers stemming from metalla-aromatics are achievable.Open in a separate windowFig. 1.Comparison of traditional organic skeletons with metalla-aromatic building blocks (the computed energies are in eV). (A) HOMO–LUMO gaps of classic aromatic skeletons. (B) Carbolong frameworks as potential building blocks for novel conjugated polymers with broad absorption bands and improved solubility.In recent years, we have reported a series of readily accessible metal-bridged bicyclic/polycyclic aromatics, namely carbolong complexes, which are stable in air and moisture (23–25). The addition of osmium carbynes (in carbolong complexes) and alkynes gave rise to an intriguing family of d_π–p_π conjugated systems, which function as excellent electron transport layer materials in organic solar cells (26, 27). These observations raised the following question: Can this efficient addition reaction be used to access metalla-aromatic conjugated polymers? It is noteworthy that incorporation of metalla-aromatic units into conjugated polymers is hitherto unknown. In this contribution, we disclose a polymerization reaction involving M≡C analogs of C≡C bonds, which involves a unique carbyne shuttling strategy (Fig. 2A). This led to examples of metalla-aromatic conjugated polymers (polycarbolongs) featuring metal carbyne units in the main chain. On the other hand, the development of polymerization reactions plays a crucial role in involving certain building blocks in conjugated polymers (28–32). These efficient, specific, and feasible polymerizations could open an avenue for the synthesis of conjugated polymers.Open in a separate windowFig. 2.Design of polymers and synthesis of monomers. (A) Schematic illustration of the polymerization strategy. (B) Preparation of carbolong monomers. Insert: X-ray molecular structure for the cations of complex 3. Ellipsoids are shown at the 50% probability level; phenyl groups in PPh₃ are omitted for clarity.     

Electrochemical borylation of carboxylic acids

A simple electrochemically mediated method for the conversion of alkyl carboxylic acids to their borylated congeners is presented. This protocol features an undivided cell setup with inexpensive carbon-based electrodes and exhibits a broad substrate scope and scalability in both flow and batch reactors. The use of this method in challenging contexts is exemplified with a modular formal synthesis of jawsamycin, a natural product harboring five cyclopropane rings.

Boronic acids are among the most malleable functional groups in organic chemistry as they can be converted into almost any other functionality (1–3). Aside from these versatile interconversions, their use in the pharmaceutical industry is gaining traction, resulting in approved drugs such as Velcade, Ninlaro, and Vabomere (4). It has been shown that boronic acids can be rapidly installed from simple alkyl halides (5–19) or alkyl carboxylic acids through the intermediacy of redox-active esters (RAEs) (Fig. 1A) (20–24). Our laboratory has shown that both Ni (20) and Cu (21) can facilitate this reaction. Conversely, Aggarwal and coworkers (22) and Li and coworkers (23) demonstrated photochemical variations of the same transformation. While these state-of-the-art approaches provide complementary access to alkyl boronic acids, each one poses certain challenges. For example, the requirement of excess boron source and pyrophoric MeLi under Ni catalysis is not ideal. Although more cost-effective and operationally simple, Cu-catalyzed borylation conditions can be challenging on scale due to the heterogeneity resulting from the large excess of LiOH•H₂O required. In addition to its limited scope, Li and coworkers’ protocol requires 4 equivalence of B₂pin₂ and an expensive Ir photocatalyst. The simplicity of Aggarwal and coworkers’ approach is appealing in this regard and represents an important precedent for the current study.Open in a separate windowFig. 1.(A) Prior approaches to access alkyl boronic esters from activated acids. (B) Inspiration for initiating SET events electrochemically to achieve borylation. (C) Summary of this work.At the heart of each method described above, the underlying mechanism relies on a single electron transfer (SET) event to promote decarboxylation and form an alkyl radical species. In parallel, the related borylation of aryl halides via a highly reactive aryl radical can also be promoted by SET. While numerous methods have demonstrated that light can trigger this mechanism (Fig. 1B) (16, 25–31), simple electrochemical cathodic reduction can elicit the same outcome (32–35). It was postulated that similar electrochemically driven reactivity could be translated to alkyl RAEs. The development of such a transformation would be highly enabling, as synthetic organic electrochemistry allows the direct addition or removal of electrons to a reaction, representing an incredibly efficient way to forge new bonds (36–40). This disclosure reports a mild, scalable, and operationally simple electrochemical decarboxylative borylation (Fig. 1C) not reliant on transition metals or stoichiometric reductants. In addition to mechanistic studies of this interesting transformation, applications to a variety of alkyl RAEs, comparison to known decarboxylative borylation methods, and a formal synthesis of the polycyclopropane natural product jawsamycin [(–)-FR-900848] are presented.     

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.     

Chlorination of arenes via the degradation of toxic chlorophenols

Aryl chlorides are among the most versatile synthetic precursors, and yet inexpensive and benign chlorination techniques to produce them are underdeveloped. We propose a process to generate aryl chlorides by chloro-group transfer from chlorophenol pollutants to arenes during their mineralization, catalyzed by Cu(NO₃)₂/NaNO₃ under aerobic conditions. A wide range of arene substrates have been chlorinated using this process. Mechanistic studies show that the Cu catalyst acts in cooperation with NO_x species generated from the decomposition of NaNO₃ to regulate the formation of chlorine radicals that mediate the chlorination of arenes together with the mineralization of chlorophenol. The selective formation of aryl chlorides with the concomitant degradation of toxic chlorophenol pollutants represents a new approach in environmental pollutant detoxication. A reduction in the use of traditional chlorination reagents provides another (indirect) benefit of this procedure.

Chlorophenols are widely encountered moieties present in herbicides, drugs, and pesticides (1). Owing to the high dissociation energies of carbon‒chloride bonds, chlorophenols biodegrade very slowly, and their prolonged exposure leads to severe ecological and environmental problems (Fig. 1A) (2–4). Several well-established technologies have been developed for the treating of chlorophenols, including catalytic oxidation (5–11), biodegradation (12–15), solvent extraction (16, 17), and adsorption (18–20) Among these methods, adsorption is the most versatile and widely used method due to its high removal efficiency and simple operation, but the resulting products are of no value, and consequently, these processes are not viable.Open in a separate windowFig. 1.Background and reaction design. (A) Examples of chlorophenol pollutants. (B) Examples of aryl chlorides. (C) The chlorination process reported herein was based on chloro-group transfer from chlorophenol pollutants.With the extensive application of substitution reactions (21, 22), transfunctionalizations (23, 24), and cross-coupling reactions (25, 26), aryl chlorides are regarded as one of the most important building blocks widely used in the manufacture of polymers, pharmaceuticals, and other types of chemicals and materials (Fig. 1B) (27–31). Chlorination of arenes is usually carried out with toxic and corrosive reagents (32–34). Less toxic and more selective chlorination reagents tend to be expensive [e.g., chloroamides (35, 36)] and are not atom economic (37–39). Consequently, from the perspective of sustainability, the ability to transfer a chloro group from unwanted chlorophenols to other substrates would be advantageous.Catalytic isofunctional reactions, including transfer hydrogenation and alkene metathesis, have been widely exploited in organic synthesis. We hypothesized that chlorination of arenes also could be achieved by chloro-group transfer, and since stockpiles of chlorophenols tend to be destroyed by mineralization and chlorophenol pollutants may be concentrated by adsorption (18–20), they could be valorized as chlorination reagents via transfer of the chloro group to arene substrates during their mineralization, thereby adding value to the destruction process (Fig. 1C). Although chlorophenol pollutants are not benign, their application as chlorination reagents, with their concomitant destruction to harmless compounds, may be considered as not only meeting the criteria of green chemistry but also potentially surpassing it. Herein, we describe a robust strategy to realize chloro-group transfer from chlorophenol pollutants to arenes and afford a wide range of value-added aryl chlorides.     

Energy penalties enhance flexible receptor docking in a model cavity

Protein flexibility remains a major challenge in library docking because of difficulties in sampling conformational ensembles with accurate probabilities. Here, we use the model cavity site of T4 lysozyme L99A to test flexible receptor docking with energy penalties from molecular dynamics (MD) simulations. Crystallography with larger and smaller ligands indicates that this cavity can adopt three major conformations: open, intermediate, and closed. Since smaller ligands typically bind better to the cavity site, we anticipate an energy penalty for the cavity opening. To estimate its magnitude, we calculate conformational preferences from MD simulations. We find that including a penalty term is essential for retrospective ligand enrichment; otherwise, high-energy states dominate the docking. We then prospectively docked a library of over 900,000 compounds for new molecules binding to each conformational state. Absent a penalty term, the open conformation dominated the docking results; inclusion of this term led to a balanced sampling of ligands against each state. High ranked molecules were experimentally tested by T_m upshift and X-ray crystallography. From 33 selected molecules, we identified 18 ligands and determined 13 crystal structures. Most interesting were those bound to the open cavity, where the buried site opens to bulk solvent. Here, highly unusual ligands for this cavity had been predicted, including large ligands with polar tails; these were confirmed both by binding and by crystallography. In docking, incorporating protein flexibility with thermodynamic weightings may thus access new ligand chemotypes. The MD approach to accessing and, crucially, weighting such alternative states may find general applicability.

Proteins interchange between conformational states of varying probabilities (1). These rearrangements, naturally, also alter its physicochemical properties (2, 3). Exploiting these varying features can benefit ligand discovery (4–7) but also presents several challenges. Key among them is weighting the different states by their energies, which has been shown to be crucial for docking success (4, 8); without such weights, high-energy protein conformations, often better suited to ligand complementarity but harder to access, can dominate docking results, acting effectively as decoy conformations.Structural models of proteins in distinct conformational states can be obtained from experiments like X-ray crystallography, NMR, or cryoelectron microscopy. The choice of the single structure used for a docking campaign contributes to its likelihood of success and choosing any single conformation inevitably leads to false negatives, even in successful campaigns. A solution to this problem is to consider multiple protein conformations, often referred to as ensemble docking or flexible receptor docking (4, 9–15). Yet incorporating multiple protein conformations only increases accuracy in ligand discovery when they are weighted according to their ensemble probabilities (10, 16, 17). When such energies have been incorporated in docking campaigns, they have enabled the discovery of ligands that are inaccessible to single-state docking, often with high fidelity to the subsequent structure determination of ligand–protein complexes. However, incorporating these weights has relied on experimental observables, such as occupancies from high-resolution structures. This has both limited the range of states that may be used—since states higher in energy than a few kilocalories per mole above the ground state will not be observed experimentally—and cannot be generalized to the vast number of targets for which such information is unavailable. Even when alternative conformational states can be observed in complex with different ligands (4, 11, 18, 19), their thermodynamic weights in the apo ensemble are typically unknown. It would be useful to have a general method of sampling and energy-weighting conformational states that would enable their exploitation in ligand discovery, in general, and for molecular docking in particular.In principle, computationally modeled conformations, such as those derived from molecular dynamics (MD) simulations, can sample such states (9, 20, 21) and can estimate their thermodynamic weighting (22–24). Encouraging studies on how MD simulation can be leveraged to explore the flexibility of ligand binding sites include work from the Bowman group, in which exhaustive MD simulations aided the discovery of allosteric binders (6, 10). More recently, exascale simulations of proteins central to SARS-CoV-2 immune evasion were able to explain and predict cryptic binding sites (25, 26). In practice, however, challenges with many MD simulations include insufficient sampling and the difficulty in weighting states by relative energies. The free energy minima, representing conformational states of a protein, are often separated by high-energy barriers, which are rarely overcome on time scales covered by conventional MD (cMD) simulations (1). Enhanced sampling algorithms, such as accelerated MD (aMD) (27), introduce a bias potential to lower the barriers between individual conformational states. This makes the sampling of a diverse, conformational ensemble, including higher-energy conformational states, more efficient by increasing the sampling by up to three orders of magnitude (28–30). A core question is whether the assumptions and approximations made in aMD affect its ability to usefully weight the conformations sampled.Here, we test energetic weights from MD simulation for ligand discovery in the engineered cavity site of T4 lysozyme L99A (L99A). This hydrophobic cavity was first introduced by Eriksson, Morton, Baase, and Matthews (31–34), as a model system to explore ligand binding and thermodynamics. While binding to this site is not thought to affect the enzymes function (it is over 20 Å from the catalytic aspartate and does not overlap with the muramyl peptide binding site; SI Appendix, Fig. S7), it has important advantages for exploring terms in ligand binding and docking (here, protein flexibility). The cavity site is relatively small, only 150 ų in its apo state, and is completely enclosed from solvent in that conformation (Fig. 1A). Combined with its dominance by apolar interactions, this simplifies the determinants of ligand binding. Despite its small size, there are still many hundreds of likely ligands that are readily available and testable from within docking libraries, enabling prospective predictions to test new docking terms and methods (11, 32, 33, 35–37). Previous studies have revealed at least 68 ligands for this cavity, many of which have protein-bound crystal structures determined (31, 38), enabling detailed retrospective studies. Despite its simplicity, L99A has complexities that make it interesting and relevant as a model site, and its thermodynamics (34, 39–41), dynamics (28, 42–50), and ligand (un)binding (33, 36, 51–59) have been intensely studied.Open in a separate windowFig. 1.Three conformations of the L99A cavity binding site. (A) Crystal structures of T4 lysozyme L99A in its apo state show a small, buried, and entirely apolar cavity. (B) Structures of the protein in complex with ligands of increasing size show three major conformational states of the binding site: closed (purple), intermediate (blue), and open (green). (C) Workflow.Particularly germane to this study, the cavity undergoes a conformational change as larger and larger ligands bind to it, adopting three principal conformations termed closed (150 ų), intermediate (∼200 ų), and open (<300 ų) (35) (Fig. 1B). As larger ligands bind, the cavity opens owing to the unwinding of helix F from an α- to a 3 to 10-helix. In the most voluminous state of the cavity (twice that of the closed state), it opens to form a channel between bulk solvent and the hydrophobic cavity. Thus, despite being a simplified model system, L99A exhibits substantial structural rearrangements, making it a useful site to test flexible receptor docking (11, 60).For this study, we derive conformational state definitions from apo and holo crystal structures of L99A in its closed, intermediate, and open state (Fig. 1C). Removing the ligands, we perform aMD and cMD simulations for exhaustive and efficient sampling. We then construct a Markov state model (MSM) (61–64) to estimate the relative probability of each crystallographic conformational state in the apo ensemble. Converted into a conformational energy penalty Ep (Eq. 1), we incorporate the state probabilities into our flexible receptor docking scoring function (4).Ep=−m⋅kB⋅T⋅ln(p),[1]where k_B is the Boltzmann constant, T is the temperature in K, P is the population, and m is the weighting multiplier.The multiplier m weights the conformational penalty energy to bring it into balance with the other terms in the DOCK3.7 scoring function, which are typically higher in magnitude than true ligand-binding energies. As in earlier studies that used crystallographic occupancies to measure populations, this m value may need to be optimized for each system studied, at least for DOCK3.7 (4); for scoring functions whose energies are already aligned with experimental binding energies, this may not be necessary. While this is admittedly a weakness, we test the reliability of the applied penalties in retrospective screens based on the known ligands and their property-matched decoys, as we do with the normal scoring function (see Results). Thus, while this weighting may change from system to system, doing so fits with the retrospective control calculations that are already typical in docking.Crucially, we evaluate the ability of the approach to predict ligands with new chemotypes selective for each of the three relevant conformations of the cavity in a prospective docking screen. We consider the usefulness of this approach to the general problem of predicting weighted, conformational ensembles of proteins for docking and ligand discovery.     

Empirical social triad statistics can be explained with dyadic homophylic interactions

The remarkable robustness of many social systems has been associated with a peculiar triangular structure in the underlying social networks. Triples of people that have three positive relations (e.g., friendship) between each other are strongly overrepresented. Triples with two negative relations (e.g., enmity) and one positive relation are also overrepresented, and triples with one or three negative relations are drastically suppressed. For almost a century, the mechanism behind these very specific (“balanced”) triad statistics remained elusive. Here, we propose a simple realistic adaptive network model, where agents tend to minimize social tension that arises from dyadic interactions. Both opinions of agents and their signed links (positive or negative relations) are updated in the dynamics. The key aspect of the model resides in the fact that agents only need information about their local neighbors in the network and do not require (often unrealistic) higher-order network information for their relation and opinion updates. We demonstrate the quality of the model on detailed temporal relation data of a society of thousands of players of a massive multiplayer online game where we can observe triangle formation directly. It not only successfully predicts the distribution of triangle types but also explains empirical group size distributions, which are essential for social cohesion. We discuss the details of the phase diagrams behind the model and their parameter dependence, and we comment on to what extent the results might apply universally in societies.

Recognizing the fundamental role of triadic interactions in shaping social structures, Heider (1) introduced the notion of balanced and unbalanced triads. A triad (triangle) of individuals is balanced if it includes zero or two negative links; otherwise, it is unbalanced. Heider (1) hypothesized that social networks have a tendency to reduce the number of unbalanced triangles over time such that balanced triads would dominate in a stationary situation. This theory of “social balance” has been confirmed empirically in many different contexts, such as schools (2), monasteries (3), social media (4), or computer games (5). Social balance theory and its generalizations (6–8) have been studied extensively for more than a half century for their importance in understanding polarization of societies (9), global organization of social networks (10), evolution of the network of international relations (11), opinion formation (12, 13), epidemic spreading (14, 15), government formation (16), and decision-making processes (17).Following Heider’s intuition (18–41), current approaches toward social balance often account for the effect of triangles on social network formation in one way or another. For example, the models in refs. 22 and 23 consider a reduction of the number of unbalanced triads either in the neighborhood of a node or in the whole network. The latter process sometimes leads to imbalance due to the existence of so-called jammed states (42). In order to reach social balance, individuals can also update their links according to their relations to common neighbors (18–21) or adjust link weights via opinion updates (24, 25) or via a minimization of social stress based on triadic interactions (37–44). These works not only ignore the difficulty of individuals to know the social interactions beyond their direct neighbors in reality, so far, they also have not considered the detailed statistical properties of the over- or underrepresentation of the different types of triads, such as those reported in refs. 4 and 5, with the exception of refs. 43 and 44.It is generally believed that the similarity of individuals plays a crucial role in the formation of social ties in social networks, something that has been called homophily (45–48). This means that to form a positive or negative tie with another person, people compare only pairwise overlaps in their individual opinions (dyadic interaction). It has also been argued that social link formation takes into account a tendency in people to balance their local interaction networks in the sense that they introduce friends to each other, that they do give up friendships if two mutual friends have negative attitudes toward each other, and that they tend to avoid situations where everyone feels negatively about the others. This is the essence of social balance theory (1). Obviously, link formation following social balance is cognitively much more challenging than homophily-based link formation since in the former, one has to keep in mind the many mutual relations between all your neighbors in a social network. While social balance–driven link formation certainly occurs in the context of close friendships, it is less realistic to assume that this mechanism is at work in social link formation in general. In Fig. 1, we schematically show the situation in a portion of a social network. It is generally hard for node i to know all the relations between his neighbors j, k, and l.Open in a separate windowFig. 1.Schematic view of opinion and link updates in a society. Every individual has an opinion vector whose components represent (binary) opinions on G=5 different subjects. Red (blue) links denote positive (negative) relationships. The question marks denote unknown relationships between i’s neighbors. As an agent i flips one of its opinions (red circle), si1, from 1 to –1, i can either decrease or increase its individual stress, H(i), depending on the value of the parameter α (Eq. 1). For instance, H(i) would increase if α=1 but would decrease for α=0. For high “rationality” values of individuals w.r.t. social stress, as quantified by β, the latter is more likely to be accepted, resulting in a reduction of the number of unbalanced triads in i’s neighborhood.Here, assuming that it is generally unrealistic for individuals to know their social networks at the triadic level, we aim to understand the emergence and the concrete statistics of balanced triads on the basis of dyadic or one-to-one interactions. Therefore, we use a classic homophily rule (45, 46) to define a “stress level” between any pair of individuals based on the similarity (or overlap) of their individual opinions. Here, the opinions of an individual i are represented by a vector with G components, si, that we show in Fig. 1. Homophily implies that i and j tend to become friends if the overlap (e.g., scalar product of their opinion vectors) is positive, and they become enemies if the overlap is negative. Such a specification of homophily is often referred to as an attraction–repulsion or assimilation–differentiation rule (49, 50). Assuming that, generally, social relations rearrange such as to minimize individual social stress on average, we will show that balanced triads naturally emerge from purely dyadic homophilic interactions without any explicit selection mechanisms for specific triads. We formulate the opinion link dynamics leading to social balance within a transparent physics-inspired framework. In particular, we observe a dynamic transition between two different types of balanced steady states that correspond to different compositions of balanced triads.Explaining the empirical statistics of triangles in social systems is a challenge. Early works considered groups of a few monks in a monastery (3) or a few students in classrooms (51). The studies suffered from limited data and small network sizes. Large-scale studies were first performed in online platforms (4) and in the society of players of the massive multiplayer online game (MMOG) Pardus. Players in Pardus engage in a form of economic life, such as trade and mining, and in social activities, such as communication on a number of channels, forming friendships and enmities (details are in refs. 5, 52, and 53). In the social networks of this game, balanced triads were once more confirmed to be overrepresented compared with what is expected by chance. Similar patterns of triad statistics were also observed in Epinion, Slashdot, and Wikipedia (4). More details on the Pardus society are in Materials and Methods. This dataset gives us the unique possibility to validate the model and compare the predictions with actual triangle statistics and formation of positively connected groups that are foundational to social cohesion.     

SOD1 aggregation in ALS mice shows simplistic test tube behavior

A longstanding challenge in studies of neurodegenerative disease has been that the pathologic protein aggregates in live tissue are not amenable to structural and kinetic analysis by conventional methods. The situation is put in focus by the current progress in demarcating protein aggregation in vitro, exposing new mechanistic details that are now calling for quantitative in vivo comparison. In this study, we bridge this gap by presenting a direct comparison of the aggregation kinetics of the ALS-associated protein superoxide dismutase 1 (SOD1) in vitro and in transgenic mice. The results based on tissue sampling by quantitative antibody assays show that the SOD1 fibrillation kinetics in vitro mirror with remarkable accuracy the spinal cord aggregate buildup and disease progression in transgenic mice. This similarity between in vitro and in vivo data suggests that, despite the complexity of live tissue, SOD1 aggregation follows robust and simplistic rules, providing new mechanistic insights into the ALS pathology and organism-level manifestation of protein aggregation phenomena in general.So far, the difficulty to experimentally measure protein aggregation in live tissue has focused many researchers to infer mechanistic details of neurodegenerative disease from molecular studies in vitro. An important outcome of this in vitro development is the establishment of rational protocols for quantitative assessment of protein aggregation data (1–4), which now start to consolidate our view of what is happening (5). Protein aggregation follows general and simplistic rules dictated by the amino acid sequence. However, the sheer number of competing aggregation sites within a typical protein chain (6) makes the process intrinsically malleable and dependent on experimental conditions (7). The nagging concern is then to what extent these already complex in vitro data are transferable to the even more complex situation in vivo? Here, we shed light on this question by comparing directly in vitro aggregation kinetics with corresponding data from transgenic mice using a recently developed in vivo quantification strategy based on antibodies (8). Our model system is the aggregation of superoxide dismutase 1 (SOD1) associated with the motor neuron disease ALS (8) (Fig. 1). A key feature of this system is that the immature apoSOD1 monomer, which is also implicated as a precursor in human pathology (9–12), needs to be globally unfolded to fibrillate in vitro (7) (Fig. 1). This simplistic behavior presents the experimental advantage that the fibrillation kinetics of apoSOD1 show simple dependence on structural stability (13, 14):ΔGD−N=−RTlnKD−N=−RTln[N][D],[1]where N is the soluble native structure, and D is the aggregation-competent unfolded state. Accordingly, it has been shown that the in vitro fibrillation of apoSOD1 displays the characteristic fingerprint of fragmentation-assisted growth (15) with a square root dependence on [D] (7), consistent with the requirement of sample agitation to expedite the reaction (1–4, 10). Analogous fibrillation behavior is found for β2-microglobulin (2), yeast prions Sup35 (16) and Ure2p (17), insulin (18), WW domain (19), TI 127 (20), and α-synuclein (21). The main difference between these proteins seems to be that some are intrinsically disordered and constantly aggregation-competent by lacking the ability to hide sticky sequence material by folding. In this study, we see that this simplistic in vitro behavior also translates to the more complex conditions in live tissue: the survival times of ALS mice expressing SOD1 variants of different stabilities are directly correlated with the in vivo levels of globally unfolded protein. Also, spinal cords of mice expressing the human SOD1 mutation G93A show exponential buildup of SOD1 aggregates with a square root dependence on log[D] indistinguishable from the fibrillation kinetics observed in agitated test tubes. The data raise fundamental questions about not only the striking resemblance between mouse and test tube aggregation but also, the apparent differences with human ALS pathology, which seems to have less ordered progression. Clues to the latter, however, are hinted in data from homozygous D90A mice showing two strains of structurally distinct SOD1 aggregates.Open in a separate windowFig. 1.SOD1 aggregation in vitro and in ALS mice. (A) Aggregation of SOD in test tubes yields fibrillar structures similar to those of other proteins (7). (B) Immunohistochemistry of the ventral horn in the terminal hSOD1^G93A mouse showing characteristics of aggresomes (44). (C) Competition between SOD1 folding and fibrillation in vitro, where elongation occurs by unfolded monomers through an encounter complex (7). The question that we ask is how do the in vitro and in vivo aggregations compare mechanistically. (D) Agitation-induced fibrillation in vitro with representative data from an SOD1 mutant in 0 (blue) and 5 M (red) urea with the associated statistics of τ_1/2 for repeated measures. To account for this statistical variation, we use the distribution average (Table S1). (E) Log plot of ν_max vs. τ_1/2 for all individual measures in this study showing uniform behavior of the various SOD1 mutants and a slope of one characteristic for exponential growth (1–4). ALS-associated SOD1 mutations examined in ALS mice (red) (Table S1), other ALS-associated mutations (blue) (Table S1), and SOD1 control mutations (black) (Table S1).     

LinearTurboFold: Linear-time global prediction of conserved structures for RNA homologs with applications to SARS-CoV-2

The constant emergence of COVID-19 variants reduces the effectiveness of existing vaccines and test kits. Therefore, it is critical to identify conserved structures in severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) genomes as potential targets for variant-proof diagnostics and therapeutics. However, the algorithms to predict these conserved structures, which simultaneously fold and align multiple RNA homologs, scale at best cubically with sequence length and are thus infeasible for coronaviruses, which possess the longest genomes (∼30,000 nt) among RNA viruses. As a result, existing efforts on modeling SARS-CoV-2 structures resort to single-sequence folding as well as local folding methods with short window sizes, which inevitably neglect long-range interactions that are crucial in RNA functions. Here we present LinearTurboFold, an efficient algorithm for folding RNA homologs that scales linearly with sequence length, enabling unprecedented global structural analysis on SARS-CoV-2. Surprisingly, on a group of SARS-CoV-2 and SARS-related genomes, LinearTurboFold’s purely in silico prediction not only is close to experimentally guided models for local structures, but also goes far beyond them by capturing the end-to-end pairs between 5′ and 3′ untranslated regions (UTRs) (∼29,800 nt apart) that match perfectly with a purely experimental work. Furthermore, LinearTurboFold identifies undiscovered conserved structures and conserved accessible regions as potential targets for designing efficient and mutation-insensitive small-molecule drugs, antisense oligonucleotides, small interfering RNAs (siRNAs), CRISPR-Cas13 guide RNAs, and RT-PCR primers. LinearTurboFold is a general technique that can also be applied to other RNA viruses and full-length genome studies and will be a useful tool in fighting the current and future pandemics.

RNA plays important roles in many cellular processes (1, 2). To maintain their functions, secondary structures of RNA homologs are conserved across evolution (3–5). These conserved structures provide critical targets for diagnostics and treatments. Thus, there is a need for developing fast and accurate computational methods to identify structurally conserved regions.Commonly, conserved structures involve compensatory base pair changes, where two positions in primary sequences mutate across evolution and still conserve a base pair; for instance, an AU or a CG pair replaces a GC pair in homologous sequences. These compensatory changes provide strong evidence for evolutionarily conserved structures (6–10). Meanwhile, they make it harder to align sequences when structures are unknown. Initially, the process of determining a conserved structure, termed comparative sequence analysis, was manual and required substantial insight to identify the conserved structure. A notable early achievement was the determination of the conserved transfer RNA (tRNA) secondary structure (11). Comparative analysis was also demonstrated to be 97% accurate compared to crystal structures for ribosomal RNAs, where the models were refined carefully over time (12).To automate comparative analysis, three distinct algorithmic approaches were developed (13, 14). The first, “joint fold-and-align” method, seeks to simultaneously predict structures and a structural alignment for two or more sequences. This was first proposed by Sankoff (15) using a dynamic programming algorithm. The major limitation of this approach is that the algorithm runs in O(n3k) against k sequences with the average sequence length n. Several software packages provide implementations of the Sankoff algorithm (16–21) that use simplifications to reduce runtime. The second, “align-then-fold” approach, is to input a sequence alignment and predict the conserved structure that can be identified across sequences in the alignment. This was described by Waterman (22) and was subsequently refined and popularized by RNAalifold (23). The third, “fold-then-align” approach, is to predict plausible structures for the sequences and then align the structures to determine the sequence alignment and the optimal conserved structures. This was described by Waterman (24) and implemented in RNAforester (25) and MARNA (26) (SI Appendix, Fig. S1).As an alternative, TurboFold II (27), an extension of TurboFold (28), provides a more computationally efficient method to align and fold sequences. Taking multiple unaligned sequences as input, TurboFold II iteratively refines alignments and structure predictions so that they conform more closely to each other and converge on conserved structures. TurboFold II is significantly more accurate than other methods (16, 18, 23, 29, 30) when tested on RNA families with known structures and alignments.However, the cubic runtime and quadratic memory usage of TurboFold II prevent it from scaling to longer sequences such as full-length severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) genomes, which contain ∼30,000 nucleotides; in fact, no joint-align-and-fold methods can scale to these genomes, which are the longest among RNA viruses. As a (not very principled) workaround, most existing efforts for modeling SARS-CoV-2 structures (31–36) resort to local folding methods (37, 38) with sliding windows plus a limited pairing distance, abandoning all long-range interactions, and only consider one SARS-CoV-2 genome (Fig. 1 B and C), ignoring signals available in multiple homologous sequences. To address this challenge, we designed a linearized version of TurboFold II, LinearTurboFold (Fig. 1A), which is a global homologous folding algorithm that scales linearly with sequence length. This linear runtime makes it, to our knowledge, the first joint-fold-and-align algorithm scale to full-length coronavirus genomes without any constraints on window size or pairing distance, taking about 13 h to analyze a group of 25 SARS-CoV homologs. It also leads to significant improvement on secondary structure prediction accuracy as well as an alignment accuracy comparable to or higher than all benchmarks.Open in a separate windowFig. 1.(A) The LinearTurboFold framework. Like TurboFold II, LinearTurboFold takes multiple unaligned homologous sequences as input and outputs a secondary structure for each sequence and a multiple-sequence alignment (MSA). But unlike TurboFold II, LinearTurboFold employs two linearizations to ensure linear runtime: a linearized alignment computation (module 1) to predict posterior coincidence probabilities (red squares) for all pairs of sequences (first four sections in Methods) and a linearized partition function computation (module 2) to estimate base-pairing probabilities (yellow triangles) for all the sequences (Methods, Extrinsic Information Calculation and Methods, LinearPartition for Base Pairing Probabilities Estimation with Extrinsic Information). These two modules take advantage of information from each other and iteratively refine predictions (SI Appendix, Fig. S2). After several iterations, module 3 generates the final multiple-sequence alignments (Methods, MSA Generation and Secondary Structure Prediction), and module 4 predicts secondary structures. Module 5 can stochastically sample structures. (B and C) Prior studies (31–36) [except for the purely experimental work by Ziv et al. (39)] used local folding methods with limited window size and maximum pairing distance. B shows the local folding of the SARS-CoV-2 genome by Huston et al. (32), which used a window of 3,000 nt that was advanced 300 nt. It also limited the distance between nucleotides that can form base pair at 500. Some studies also used homologous sequences to identify conserved structures (32–36), but they predicted only structures for one genome and utilized sequence alignments to identify mutations. By contrast, LinearTurboFold is a global folding method without any limitations on sequence length or paring distance, and it jointly folds and aligns homologs to obtain conserved structures. Consequently, LinearTurboFold can capture long-range interactions even across the whole genome (the long arc in B and Fig. 3).Over a group of 25 SARS-CoV-2 and SARS-related homologous genomes, LinearTurboFold predictions are close to the canonical structures (40) and structures modeled with the aid of experimental data (32–34) for several well-studied regions. Due to global rather than local folding, LinearTurboFold discovers a long-range interaction involving 5′ and 3′ untranslated regions (UTRs) (∼29,800 nt apart), which is consistent with recent purely experimental work (39) and yet is out of reach for local folding methods used by existing studies (Fig. 1 B and C). In short, our in silico method of folding multiple homologs can achieve results similar to, and sometimes more accurate than, those of experimentally guided models for one genome. Moreover, LinearTurboFold identifies conserved structures supported by compensatory mutations, which are potential targets for small-molecule drugs (41) and antisense oligonucleotides (ASOs) (36). We further identify regions that are 1) sequence-level conserved; 2) at least 15 nt long; and 3) accessible (i.e., likely to be completely unpaired) as potential targets for ASOs (42), small interfering RNA (siRNA) (43), CRISPR-Cas13 guide RNA (gRNA) (44), and RT-PCR primers (45). LinearTurboFold is a general technique that can also be applied to other RNA viruses (e.g., influenza, Ebola, HIV, Zika, etc.) and full-length genome studies.     

Diffusion and distal linkages govern interchromosomal dynamics during meiotic prophase

The pairing of homologous chromosomes (homologs) in meiosis is essential for distributing the correct numbers of chromosomes into haploid gametes. In budding yeast, pairing depends on the formation of 150 to 200 Spo11-mediated double-strand breaks (DSBs) that are distributed among 16 homolog pairs, but it is not known if all, or only a subset, of these DSBs contribute to the close juxtaposition of homologs. Having established a system to measure the position of fluorescently tagged chromosomal loci in three-dimensional space over time, we analyzed locus trajectories to determine how frequently and how long loci spend colocalized or apart. Continuous imaging revealed highly heterogeneous cell-to-cell behavior of foci, with the majority of cells exhibiting a “mixed” phenotype where foci move into and out of proximity, even at late stages of prophase, suggesting that the axial structures of the synaptonemal complex may be more dynamic than anticipated. The observed plateaus of the mean-square change in distance (MSCD) between foci informed the development of a biophysical model of two diffusing polymers that captures the loss of centromere linkages as cells enter meiosis, nuclear confinement, and the formation of Spo11-dependent linkages. The predicted number of linkages per chromosome in our theoretical model closely approximates the small number (approximately two to four) of estimated synapsis-initiation sites, suggesting that excess DSBs have negligible effects on the overall juxtaposition of homologs. These insights into the dynamic interchromosomal behavior displayed during homolog pairing demonstrate the power of combining time-resolved in vivo analysis with modeling at the granular level.

During meiosis prophase I, homologous chromosomes undergo pairing, synapsis, and crossing over to ensure their proper segregation at meiosis I. An overarching question is how each chromosome identifies and pairs with its homolog partner within the complex nuclear environment that includes nonhomologous chromosomes (1–4). The general view is that pairing is achieved through homology-based mechanisms that can bring the axes of chromosome pairs into close juxtaposition such that discrete pairing interactions, in conjunction with the establishment of synapsis, are sufficient to align homologs end to end (5). While the intermediate steps leading to pairing are not well understood, the process itself is thought to be stochastic with heterogeneity from cell to cell.The budding yeast, Saccharomyces cerevisiae, is an important model for the study of homolog pairing as it has been used extensively for characterizing many of the other dynamic events that occur over the course of meiotic prophase I that are now known to be conserved across phyla. These include the transition from the Rabl (centromeres-clustered) to bouquet (telomeres-clustered) configurations; telomere-led chromosome movement driven by cytoskeletal motor proteins via the linker of nucleoskeleton and cytoskeleton complex; the formation and repair of Spo11-induced DNA double-strand breaks (DSBs); and the assembly and disassembly of the synaptonemal complex (SC), which is a ribbon-like structure that joins homologs together along their lengths (Fig. 1) (5–11). Several theoretical models of pairing in yeast have been developed that take into account chromosome size, linkage numbers, and the attachment and motion of telomeres at the nuclear envelope (12–16), yet no study to date has combined biophysical modeling together with empirical measurements of meiotic “pairing” dynamics in live cells.Open in a separate windowFig. 1.Overview of chromosome conformations in premeiotic cells (TM=T0) and in meiotic cells in midprophase (∼T3,T4) and late prophase (∼T5,T6). At T₀, cells are in the G0 stage prior to DNA replication, and chromosomes are arranged in the Rabl configuration with centromeres clustered at the nuclear periphery (59). Following transfer to sporulation media, the meiotic program begins with cells entering S phase, over which time the centromeres are dispersed and telomeres start to cluster in the bouquet (59–61, 74, 115). At early to midprophase, Spo11 initiates the formation of DSBs (116), shown as stars, of which the majority are repaired using the homologous chromosome as a substrate (9) (homologs are red and orange lines; note that each line in mid- and late prophase represents the pair of newly replicated sister chromatids). DSBs that go on to form class I or interfering COs, shown as the large stars, assemble the SIC (33, 41, 42), where the new SC is shown as blue lines. Concomitantly, telomeres are subject to motion driven by cytoskeletal motor proteins shown as gray arrows (7, 117). By late prophase, homologs are synapsed end to end and with CO intermediates maturing into CO products as shown.Homolog pairing in yeast has been studied using a number of different assays, including fluorescence in situ hybridization applied to spread chromosome preparations (17, 18), a “collision” assay based on Cre/loxP recombination measuring the relative position and accessibility of pairs of homologous loci (19), and a fluorescence reporter operator system (FROS) that enables specific chromosomal loci to be tagged and followed microscopically in live cells. When allelic loci on homologous chromosomes are tagged, this “one-spot, two-spot” assay has been used as a proxy for local homolog juxtaposition (Fig. 2) (2, 4, 20–26). However, with only a static snapshot, it is not possible to know if colocalization represents a true homolog pairing interaction: that is, if the foci remain colocalized until homologs are segregated at anaphase. While it has been proposed that homologs may undergo many transient interactions that become progressively stabilized throughout prophase (27), this has not yet been investigated.Open in a separate windowFig. 2.(A) A typical field of cells, highlighting example cells showing either two spots (Left) or one spot (Right). (B and C) Maximum intensity projections of the relative positions of fluorescent foci at 30-s intervals. In B, the vertical axis corresponds to a z stack (with step size 0.25 μm). For each x and y coordinate, the maximum value over all time points for that z stack is shown. In C, the vertical axis represents time (t; in seconds), and the projection is instead performed over z stacks. The positions of the loci and the distance between them are highlighted for select time points. (D–F) Kymographs showing the distance between the loci in a single cell over the 25-min imaging period. Each horizontal slice in the kymograph shows the fluorescence intensity along the line joining the centers of the two loci in a single frame. Example of cells where the loci are separated (D) or colocalize (F) for every frame. The mixed cell shown in E undergoes several transitions between the two states. (G) Fraction of cells in the mixed state vs. hours in SPM through meiosis for the URA3 and LYS2 loci in wild-type (WT) and spo11Δ cells. The plot was made from aggregating all available data for each meiotic stage. The error is the SEM with the sample count set to the number of trajectories. (H) Schematic representation of the genomic positions of the URA3 and LYS2 loci on chromosomes V and II, respectively.Although the mechanisms promoting homolog colocalization are not well understood, in yeast interhomolog linkages depend on the formation and repair of DSBs created by Spo11 and its partners during prophase I (9). For any given cell in meiosis, any sequence has the “potential” (albeit not all equally) to experience a DSB. While 150 to 200 DSBs are formed per cell, only ∼90 to 94 DSBs go on to form crossovers (COs). Another ∼66 are repaired using the homologous chromosome but do not lead to CO formation, called noncrossovers (NCOs), and the remaining ones are repaired with the sister chromatid (28–32). COs are divided into class I and class II. Class I COs account for ∼70% of total COs; their position and number are specified in midprophase by the ZMM proteins that make up the synapsis initiation complex (SIC), which functions to couple homologous recombination with the establishment of the SC (8, 33–43). Class II COs arise from an alternative repair process that does not involve the SIC and are “interference independent” (44–46). Thus, the following question arises. Are the excess DSBs necessary to mediate pairing, or is the smaller number that goes on to form COs (class I and/or class II) sufficient?Rather than the homolog pairing process being independent for each “paired” locus, several models relating meiotic homolog pairing to polymer theory predict that pairing at one locus will increase the probability that pairing at an adjacent site will occur (14–16, 47, 48). That is, a molecular linkage at one site on the chromosome is expected to restrict the diffusive properties of adjacent sites along that chromosome (49–51). However, this has not been explicitly evaluated experimentally in the case of meiotic homolog pairing. Furthermore, it is not known if the repair of Spo11 DSBs leads to any directed motion that could aid in bringing homolog axes into close juxtaposition, similar to the observed DSB-dependent directed motion that brings telomeres into proximity seen in ALT (alternative lengthening of telomeres) cells (52). For instance, it has been proposed that single-stranded DNA filaments, formed by resection of DSBs, might capture a locus of the homologous chromosome and processively “reel” the axis into alignment (53–55).To address these gaps in knowledge, we observed the behavior of FROS-tagged loci in three-dimensional (3D) space over time and show the highly dynamic behavior between loci on homologous chromosomes during meiosis prophase I. In contrast to static snapshots, continuous imaging revealed that the majority of cells show a “mixed” phenotype in which foci alternate between colocalized and separated states, indicating that once paired, homologous loci need not remain paired until anaphase. We then used our experimental measurements of the dynamic changes in distance between homologous loci to develop a theoretical model of interhomolog dynamics based on the presence of linkages and polymer diffusion in the viscoelastic medium of the nucleus. This modeling suggests that as chromosomes transition from an unlinked to a linked state, the chromosomes are subject to random fluctuations and not an active mechanism that progressively pulls or pushes them together. Moreover, the addition of a small number of linkages (between two and four) per chromosome pair, closely approximating the number of class I COs, accounts for the observed level of confinement, while the position of linkages and other factors account for the heterogeneous cell-to-cell behavior. These insights illustrate the utility of combining live imaging with biophysical modeling for the study of dynamic processes in living cells.     

Diffusion of a disordered protein on its folded ligand

Intrinsically disordered proteins often form dynamic complexes with their ligands. Yet, the speed and amplitude of these motions are hidden in classical binding kinetics. Here, we directly measure the dynamics in an exceptionally mobile, high-affinity complex. We show that the disordered tail of the cell adhesion protein E-cadherin dynamically samples a large surface area of the protooncogene β-catenin. Single-molecule experiments and molecular simulations resolve these motions with high resolution in space and time. Contacts break and form within hundreds of microseconds without a dissociation of the complex. The energy landscape of this complex is rugged with many small barriers (3 to 4 k_BT) and reconciles specificity, high affinity, and extreme disorder. A few persistent contacts provide specificity, whereas unspecific interactions boost affinity.

Specific molecular interactions orchestrate a multitude of simultaneous cellular processes. The discovery of intrinsically disordered proteins (IDPs) (1, 2) has substantially aided our understanding of such interactions. More than two decades of research revealed a plethora of functions and mechanisms (2–6) that complemented the prevalent structure-based view on protein interactions. Even the idea that IDPs always ought to fold upon binding has largely been dismantled by recent discoveries of high-affine–disordered complexes (7, 8). Classical shape complementary is indeed superfluous in the complex between prothymosin-α and histone H1, in which charge complementary is the main driving force for binding (7). However, complexes between IDPs and folded proteins can also be highly dynamic [e.g., Sic1 and Cdc4 (9), the Na⁺/H⁺ exchanger tail and ERK2 (10), nucleoporin tails, and nuclear transport receptors (11)]. Yet timescales of motions and their spatial amplitudes are often elusive, such that it is unclear how precisely the surfaces of folded proteins alter the dynamics of bound IDPs. Answering this question is a key step in understanding how specificity, affinity, and flexibility can be simultaneously realized in such complexes.To address this question, we focused on the dynamics of the cell adhesion complex between E-cadherin (E-cad) and β-catenin (β-cat), which is involved in growth pathologies and cancer (12). E-cad is a transmembrane protein that mediates cell–cell adhesions by linking actin filaments of adjacent epithelial cells (Fig. 1A). Previous NMR results showed that the cytoplasmic tail of E-cad is intrinsically disordered (13). E-cad binds β-cat, which establishes a connection to the actin-associated protein α-catenin (14–16). β-cat, on the other hand, is a multifunctional repeat protein (17–20) that mediates cadherin-based cell adhesions (21) and governs cell fate decisions during embryogenesis (22). It contains three domains: an N-terminal domain (130 amino acids [aa]), a central repeat domain (550 aa), and a C-terminal domain (100 aa). Whereas the N- and C-terminal domains of β-cat are in large parts unstructured (17), with little effect on the affinity of the E-cad/β-cat complex (23), the 12 repeats of the central domain arrange in a superhelix (24). The X-ray structure showed that the E-cad wraps around this central domain of β-cat (24) (Fig. 1B). However, not only is half of the electron density of E-cad missing, the X-ray unit cell also comprises two structures with different resolved parts of E-cad (Fig. 1B). In fact, only 45% of all resolved E-cad residues are found in both structures (Fig. 1C). Although this ambiguity together with the large portion of missing residues (25) suggests that E-cad is highly dynamic in the complex with β-cat, the timescales and amplitudes of these dynamics are unknown.Open in a separate windowFig. 1.Complex between the cytoplasmic tail of E-cad and β-cat. (A) Schematics of cell–cell junctions mediated by E-cad and β-cat. (B) The two X-ray structures of the complex between the tail of E-cad (red) and the central repeat domain of β-cat (white) resolve different parts of E-cad (Protein Data Bank: 1i7x), indicating the flexibility of E-cad in the complex. (Bottom) Cartoon representation of the resolved E-cad parts. (C) Scheme showing the resolved parts of E-cad (red).Here, we integrated single-molecule Förster resonance energy transfer (smFRET) experiments with molecular simulations to directly measure the dynamics of E-cad on β-cat with high spatial and temporal resolution. In our bottom-up strategy, we first probed intramolecular interactions within E-cad using smFRET to parameterize a coarse-grained (CG) model. In a second step, we monitored E-cad on β-cat, integrated this information into the CG model, and obtained a dynamic picture of the complex. We found that all segments of E-cad diffuse on the surface of β-cat at submillisecond timescales and obtained a residue-resolved understanding of these motions: A small number of persistent interactions provide specificity, whereas many weak multivalent contacts boost affinity, which confirms the idea that regulatory enzymes access their recognition motifs on E-cad and β-cat without requiring the complex to dissociate (24).     

Microscopic origins of the crystallographically preferred growth in evaporation-induced colloidal crystals

Unlike crystalline atomic and ionic solids, texture development due to crystallographically preferred growth in colloidal crystals is less studied. Here we investigate the underlying mechanisms of the texture evolution in an evaporation-induced colloidal assembly process through experiments, modeling, and theoretical analysis. In this widely used approach to obtain large-area colloidal crystals, the colloidal particles are driven to the meniscus via the evaporation of a solvent or matrix precursor solution where they close-pack to form a face-centered cubic colloidal assembly. Via two-dimensional large-area crystallographic mapping, we show that the initial crystal orientation is dominated by the interaction of particles with the meniscus, resulting in the expected coalignment of the close-packed direction with the local meniscus geometry. By combining with crystal structure analysis at a single-particle level, we further reveal that, at the later stage of self-assembly, however, the colloidal crystal undergoes a gradual rotation facilitated by geometrically necessary dislocations (GNDs) and achieves a large-area uniform crystallographic orientation with the close-packed direction perpendicular to the meniscus and parallel to the growth direction. Classical slip analysis, finite element-based mechanical simulation, computational colloidal assembly modeling, and continuum theory unequivocally show that these GNDs result from the tensile stress field along the meniscus direction due to the constrained shrinkage of the colloidal crystal during drying. The generation of GNDs with specific slip systems within individual grains leads to crystallographic rotation to accommodate the mechanical stress. The mechanistic understanding reported here can be utilized to control crystallographic features of colloidal assemblies, and may provide further insights into crystallographically preferred growth in synthetic, biological, and geological crystals.

As an analogy to atomic crystals, colloidal crystals are highly ordered structures formed by colloidal particles with sizes ranging from 100 nm to several micrometers (1–6). In addition to engineering applications such as photonics, sensing, and catalysis (4, 5, 7, 8), colloidal crystals have also been used as model systems to study some fundamental processes in statistical mechanics and mechanical behavior of crystalline solids (9–14). Depending on the nature of interparticle interactions, many equilibrium and nonequilibrium colloidal self-assembly processes have been explored and developed (1, 4). Among them, the evaporation-induced colloidal self-assembly presents a number of advantages, such as large-size fabrication, versatility, and cost and time efficiency (3–5, 15–18). In a typical synthesis where a substrate is immersed vertically or at an angle into a colloidal suspension, the colloidal particles are driven to the meniscus by the evaporation-induced fluid flow and subsequently self-assemble to form a colloidal crystal with the face-centered cubic (fcc) lattice structure and the close-packed {111} plane parallel to the substrate (2, 3, 19–23) (see Fig. 1A for a schematic diagram of the synthetic setup).Open in a separate windowFig. 1.Evaporation-induced coassembly of colloidal crystals. (A) Schematic diagram of the evaporation-induced colloidal coassembly process. “G”, “M”, and “N” refer to “growth,” “meniscus,” and “normal” directions, respectively. The reaction solution contains silica matrix precursor (tetraethyl orthosilicate, TEOS) in addition to colloids. (B) Schematic diagram of the crystallographic system and orientations used in this work. (C and D) Optical image (Top Left) and scanning electron micrograph (SEM) (Bottom Left) of a typical large-area colloidal crystal film before (C) and after (D) calcination. (Right) SEM images of select areas (yellow rectangles) at different magnifications. Corresponding fast-Fourier transform (see Inset in Middle in C) shows the single-crystalline nature of the assembled structure. (E) The 3D reconstruction of the colloidal crystal (left) based on FIB tomography data and (right) after particle detection. (F) Top-view SEM image of the colloidal crystal with crystallographic orientations indicated.While previous research has focused on utilizing the assembled colloidal structures for different applications (4, 5, 7, 8), considerably less effort is directed to understand the self-assembly mechanism itself in this process (17, 24). In particular, despite using the term “colloidal crystals” to highlight the microstructures’ long-range order, an analogy to atomic crystals, little is known regarding the crystallographic evolution of colloidal crystals in relation to the self-assembly process (3, 22, 25). The underlying mechanisms for the puzzling—yet commonly observed—phenomenon of the preferred growth along the close-packed <110> direction in evaporation-induced colloidal crystals are currently not understood (3, 25–29). The <110> growth direction has been observed in a number of processes with a variety of particle chemistries, evaporation rates, and matrix materials (3, 25–28, 30), hinting at a universal underlying mechanism. This behavior is particularly intriguing as the colloidal particles are expected to close-pack parallel to the meniscus, which should lead to the growth along the <112> direction and perpendicular to the <110> direction (16, 26, 31)^*.Preferred growth along specific crystallographic orientations, also known as texture development, is commonly observed in crystalline atomic solids in synthetic systems, biominerals, and geological crystals. While current knowledge recognizes mechanisms such as the oriented nucleation that defines the future crystallographic orientation of the growing crystals and competitive growth in atomic crystals (32–34), the underlying principles for texture development in colloidal crystals remain elusive. Previous hypotheses based on orientation-dependent growth speed and solvent flow resistance are inadequate to provide a universal explanation for different evaporation-induced colloidal self-assembly processes (3, 25–29). A better understanding of the crystallographically preferred growth in colloidal self-assembly processes may shed new light on the crystal growth in atomic, ionic, and molecular systems (35–37). Moreover, mechanistic understanding of the self-assembly processes will allow more precise control of the lattice types, crystallography, and defects to improve the performance and functionality of colloidal assembly structures (38–40).     

Kin selection for cooperation in natural bacterial populations

Bacteria produce a range of molecules that are secreted from the cell and can provide a benefit to the local population of cells. Laboratory experiments have suggested that these “public goods” molecules represent a form of cooperation, favored because they benefit closely related cells (kin selection). However, there is a relative lack of data demonstrating kin selection for cooperation in natural populations of bacteria. We used molecular population genetics to test for signatures of kin selection at the genomic level in natural populations of the opportunistic pathogen Pseudomonas aeruginosa. We found consistent evidence from multiple traits that genes controlling putatively cooperative traits have higher polymorphism and greater divergence and are more likely to harbor deleterious mutations relative to genes controlling putatively private traits, which are expressed at similar rates. These patterns suggest that cooperative traits are controlled by kin selection, and we estimate that the relatedness for social interactions in P. aeruginosa is r = 0.84. More generally, our results demonstrate how molecular population genetics can be used to study the evolution of cooperation in natural populations.

The growth and success of many bacteria appear to depend upon a stunning array of cooperative behaviors (1–3). Cells produce and secrete a range of factors that benefit the local group of cells and so, act as cooperative “public goods.” Examples include molecules to scavenge iron (siderophores) (4), enzymes that break down proteins (proteases) (5), and molecules to aid cell movement (rhamnolipids) (6).The potential problem with such cooperation is that it can be exploited by noncooperators (“cheats”) that do not produce public goods but can still benefit from those produced by others (7). A likely solution to this problem in bacteria is that clonal growth keeps close relatives together, and limited diffusion keeps public goods close to producers (8). Consequently, the benefits of cooperation tend to be shared with related cells that share the gene for cooperation, and so, cooperation is favored by kin selection (9).However, most evidence for cooperation and kin selection in bacteria has come from laboratory experiments (10–18). To what extent are test-tube cultures, often utilizing extreme gene knockouts, representative of natural populations (1, 12). A problem here is that while bacteria and other microorganisms offer many advantages for laboratory experiments, they can be very difficult to study in their natural environment.Population genetics offers a way to study natural populations because kin selection can leave signatures (“footprints”) of selection at the genomic level (10–12, 15, 19–28). In a clonal population, where the relatedness (r) between interacting cells is r = 1, the benefits of cooperating will always be passed onto other individuals that carry the gene for cooperation. In contrast, as relatedness decreases (r < 1), the benefits of cooperation will increasingly be passed onto individuals that do not carry the gene for cooperation (Fig. 1A). This reduces (dilutes) the kin-selected benefit of cooperation, making beneficial mutations less likely to fix and deleterious mutations more likely to fix (Fig. 1B) (9, 25).Open in a separate windowFig. 1.Population genetic theory for cooperative traits (15, 25). (A) Representation of how traits are categorized as private or cooperative. Cooperative traits are those involving the production and secretion of molecules where the fitness benefits can potentially be shared with other cells in the local group. Private traits are those where the fitness benefits are only felt by the individual expressing the gene (e.g., internal metabolism). (B) Probability of fixation for deleterious or beneficial mutations of varying effect (x axis) for mutations influencing private (black line) and cooperative (social; all lines) traits. In clonal populations, where the relatedness (r) between interacting individuals r = 1, the prediction is the same for mutations influencing private and cooperative traits (black line). As relatedness decreases, the prediction changes for mutations influencing cooperation, with beneficial mutations becoming less likely to fix and deleterious mutations becoming more likely to fix. Consequently, in nonclonal populations, there is relaxed selection on genes controlling cooperative traits relative to those controlling private traits. Adapted from ref. 15. (C) Prediction for relative polymorphism and divergence for cooperative (blue) relative to private (yellow) genes assuming a fixed r < 1. Due to the increased fixation likelihood of deleterious mutations and decreased fixation likelihood of beneficial mutation, genes for cooperative traits should have relatively greater levels of polymorphism and divergence. (D) Predicted polymorphism of private (yellow) and cooperative (blue) genes as relatedness varies for a trait where cooperation is favored when r>0.25. For private traits, polymorphism is independent of relatedness. For cooperative traits, expected polymorphism relative to a private trait is inversely proportional to r when cooperation is favored. When r = 1, there is no difference in polymorphism between cooperative and private traits. When r < 0.25, cooperation is not favored, so relatedness no longer predicts the level of polymorphism observed.Population genetic theory, therefore, predicts that, in nonclonal populations (r < 1), cooperative traits favored by kin selection will show increased polymorphism and divergence relative to traits that provide private benefits (Fig. 1 C and D) (15, 23, 25). Nonclonal populations appear to be very common in bacteria. At the scale of the social interaction, groups often contain multiple species, let alone multiple lineages of the same species (17, 29, 30). In addition, molecular and genomic studies have demonstrated selection for noncooperative cheats that exploit the cooperation of others as well as a diversity of mechanisms for attacking nonrelatives (14, 16, 31). Clonal interactions seem to be limited to extreme cases, such as cyanobacteria filaments (30).We tested for genomic signatures of kin selection for cooperation in the opportunistic pathogen Pseudomonas aeruginosa. Laboratory experiments have suggested that P. aeruginosa produces a range of cooperative public goods that facilitate both growth and virulence (4, 32, 33). A potential problem with genomic analyses is that they can be confounded by conditional gene expression. If a gene is only occasionally expressed, in certain conditions, this can also lead to relaxed selection, making beneficial mutations less likely to fix and deleterious mutations more likely to fix (10, 22). We controlled for this influence of conditional gene expression by making targeted comparisons between cooperative and private traits that are likely to be expressed at similar rates.     

Inside-out regulation of E-cadherin conformation and adhesion

Cadherin cell–cell adhesion proteins play key roles in tissue morphogenesis and wound healing. Cadherin ectodomains bind in two conformations, X-dimers and strand-swap dimers, with different adhesive properties. However, the mechanisms by which cells regulate ectodomain conformation are unknown. Cadherin intracellular regions associate with several actin-binding proteins including vinculin, which are believed to tune cell–cell adhesion by remodeling the actin cytoskeleton. Here, we show at the single-molecule level, that vinculin association with the cadherin cytoplasmic region allosterically converts weak X-dimers into strong strand-swap dimers and that this process is mediated by myosin II–dependent changes in cytoskeletal tension. We also show that in epithelial cells, ∼70% of apical cadherins exist as strand-swap dimers while the remaining form X-dimers, providing two cadherin pools with different adhesive properties. Our results demonstrate the inside-out regulation of cadherin conformation and establish a mechanistic role for vinculin in this process.

E-cadherins (Ecads) are essential, calcium-dependent cell–cell adhesion proteins that play key roles in the formation of epithelial tissue and in the maintenance of tissue integrity. Ecad adhesion is highly plastic and carefully regulated to orchestrate complex movement of epithelial cells, and dysregulation of adhesion is a hallmark of numerous cancers (1). However, little is known about how cells dynamically regulate the biophysical properties of individual Ecads.The extracellular region of Ecads from opposing cells bind in two distinct trans orientations: strand-swap dimers and X-dimers (Fig. 1 A and B). Strand-swap dimers are the stronger cadherin adhesive conformation and are formed by the exchange of conserved tryptophan (Trp) residues between the outermost domains of opposing Ecads (2–4). In contrast, X-dimers, which are formed by extensive surface interactions between opposing Ecads, are a weaker adhesive structure and serve as an intermediate during the formation and rupture of strand-swap dimers (5–7). Using cell-free, single-molecule experiments we previously showed that X-dimers and strand-swap dimers can be distinguished based on their distinctly different response to mechanical force. When a strand-swap dimer is pulled, its lifetime decreases with increasing force, resulting in the formation of a slip bond (8, 9) (Fig. 1B). In contrast, an X-dimer responds to pulling force by forming a catch bond, where bond lifetime initially increases up to a threshold force and then subsequently decreases (8, 10) (Fig. 1B). It has also been shown that wild-type Ecad ectodomains in solution can interconvert between X-dimer and strand-swap dimer conformations (9, 11). However, the biophysical mechanisms by which Ecad conformations (and adhesion) are regulated on the cell surface are unknown.Open in a separate windowFig. 1.Overview of experiment. (A) The extracellular region of Ecad from opposing cells mediates adhesion. The cytoplasmic region of Ecad associates either directly or indirectly with p120 catenin, β-catenin, α-catenin, vinculin, and F-actin. (B) Strand-swap dimers form slip bonds (blue) and X-dimers form catch bonds (red). Ecads interconvert between these two dimer conformations. Structures were generated from the crystal structure of mouse Ecad (PDB ID code 3Q2V); the X-dimer was formed by alignment to an X-dimer crystal structure (PDB ID code 3LNH). (C) Graphics showing the cell lines used in experiments and Western blot analysis of corresponding cell lysates.The cytoplasmic region of Ecad associates with the catenin family of proteins, namely, p120-catenin, β-catenin, and α-catenin. The Ecad–catenin complex, in turn, links to filamentous actin (F-actin) either by the direct binding of α-catenin and F-actin or by the indirect association of α-catenin and F-actin via vinculin (12) (Fig. 1A). Adhesive forces transmitted across intercellular junctions by Ecad induce conformational changes in α-catenin (13, 14), strengthen F-actin binding (15), and recruit vinculin to the sites of force application (16, 17). However, vinculin and α-catenin do not merely serve as passive cytoskeletal linkers; they also dynamically modulate cytoskeletal rearrangement and recruit myosin to cell–cell junctions (13, 18–20). Studies show that α-catenin and vinculin play important roles in strengthening and stabilizing Ecad adhesion: bead-twisting experiments show force-induced stiffening of Ecad-based junctions and cell doublet stretching experiments demonstrate reinforcement of cell–cell adhesion in vinculin- and α-catenin–dependent manners (18, 19, 21).Currently, actin anchorage and cytoskeletal remodeling are assumed to be the exclusive mechanisms by which α-catenin and vinculin strengthen Ecad adhesion (22–24). Here, we directly map the allosteric effects of cytoplasmic proteins on Ecad ectodomain conformation and demonstrate, at the single-molecule level, that vinculin association with the Ecad cytoplasmic region switches X-dimers to strand-swap dimers. We show that cytoskeletal tension, due to vinculin-mediated recruitment of myosin II, regulates Ecad ectodomain structure and adhesion. Finally, we demonstrate that only ∼50% of Ecads are linked to the underlying cytoskeleton and that while about 70% of Ecads form strand-swap dimers the remaining form X-dimers, which provides cells with two Ecad pools with different adhesive properties.     

An image-computable model of how endogenous and exogenous attention differentially alter visual perception

Attention alters perception across the visual field. Typically, endogenous (voluntary) and exogenous (involuntary) attention similarly improve performance in many visual tasks, but they have differential effects in some tasks. Extant models of visual attention assume that the effects of these two types of attention are identical and consequently do not explain differences between them. Here, we develop a model of spatial resolution and attention that distinguishes between endogenous and exogenous attention. We focus on texture-based segmentation as a model system because it has revealed a clear dissociation between both attention types. For a texture for which performance peaks at parafoveal locations, endogenous attention improves performance across eccentricity, whereas exogenous attention improves performance where the resolution is low (peripheral locations) but impairs it where the resolution is high (foveal locations) for the scale of the texture. Our model emulates sensory encoding to segment figures from their background and predict behavioral performance. To explain attentional effects, endogenous and exogenous attention require separate operating regimes across visual detail (spatial frequency). Our model reproduces behavioral performance across several experiments and simultaneously resolves three unexplained phenomena: 1) the parafoveal advantage in segmentation, 2) the uniform improvements across eccentricity by endogenous attention, and 3) the peripheral improvements and foveal impairments by exogenous attention. Overall, we unveil a computational dissociation between each attention type and provide a generalizable framework for predicting their effects on perception across the visual field.

Endogenous and exogenous spatial attention prioritize subsets of visual information and facilitate their processing without concurrent eye movements (1–3). Selection by endogenous attention is goal-driven and adapts to task demands, whereas exogenous attention transiently and automatically orients to salient stimuli (1–3). In most visual tasks, both types of attention typically improve visual perception similarly [e.g., acuity (4–6), visual search (7, 8), perceived contrast (9–11)]. Consequently, models of visual attention do not distinguish between endogenous and exogenous attention (e.g., refs. 12–19). However, stark differences also exist. Each attention type differentially modulates neural responses (20, 21) and fundamental properties of visual processing, including temporal resolution (22, 23), texture sensitivity (24), sensory tuning (25), contrast sensitivity (26), and spatial resolution (27–34).The effects of endogenous and exogenous attention are dissociable during texture segmentation, a visual task constrained by spatial resolution [reviews (1–3)]. Whereas endogenous attention optimizes spatial resolution to improve the detection of an attended texture (32–34), exogenous attention reflexively enhances resolution even when detrimental to perception (27–31, 34). Extant models of attention do not explain these well-established effects.Two main hypotheses have been proposed to explain how attention alters spatial resolution. Psychophysical studies ascribe attentional effects to modulations of spatial frequency (SF) sensitivity (30, 33). Neurophysiological (13, 35, 36) and neuroimaging (37, 38) studies bolster the idea that attention modifies spatial profiles of neural receptive fields (RFs) (2). Both hypotheses provide qualitative predictions of attentional effects but do not specify their underlying neural computations.Differences between endogenous and exogenous attention are well established in segmentation tasks and thus provide an ideal model system to uncover their separate roles in altering perception. Texture-based segmentation is a fundamental process of midlevel vision that isolates regions of local structure to extract figures from their background (39–41). Successful segmentation hinges on the overlap between the visual system’s spatial resolution and the levels of detail (i.e., SF) encompassed by the texture (39, 41, 42). Consequently, the ability to distinguish between adjacent textures varies as resolution declines toward the periphery (43–46). Each attention type differentially alters texture segmentation, demonstrating that their effects shape spatial resolution [reviews (1–3)].Current models of texture segmentation do not explain performance across eccentricity and the distinct modulations by attention. Conventional models treat segmentation as a feedforward process that encodes the elementary features of an image (e.g., SF and orientation), transforms them to reflect the local structure (e.g., regions of similarly oriented bars), and then pools across space to emphasize texture-defined contours (39, 41, 47). Few of these models account for variations in resolution across eccentricity (46, 48, 49) or endogenous (but not exogenous) attentional modulations (18, 50). All others postulate that segmentation is a “preattentive” (42) operation whose underlying neural processing is impervious to attention (39, 41, 46–49).Here, we develop a computational model in which feedforward processing and attentional gain contribute to segmentation performance. We augment a conventional model of texture processing (39, 41, 47). Our model varies with eccentricity and includes contextual modulation within local regions in the stimulus via normalization (51), a canonical neural computation (52). The defining characteristic of normalization is that an individual neuron is (divisively) suppressed by the summed activity of neighboring neurons responsive to different aspects of a stimulus. We model attention as multiplicative gains [attentional gain factors (15)] that vary with eccentricity and SF. Attention shifts sensitivity toward fine or coarse spatial scales depending on the range of SFs enhanced.Our model is image-computable, which allowed us to reproduce behavior directly from grayscale images used in psychophysical experiments (6, 26, 27, 29–33). The model explains three signatures of texture segmentation hitherto unexplained within a single computational framework (Fig. 1): 1) the central performance drop (CPD) (27–34, 43–46) (Fig. 1A), that is, the parafoveal advantage of segmentation over the fovea; 2) the improvements in the periphery and impairments at foveal locations induced by exogenous attention (27–32, 34) (Fig. 1B); and 3) the equivalent improvements across eccentricity by endogenous attention (32–34) (Fig. 1C).Open in a separate windowFig. 1.Signatures of texture segmentation. (A) CPD. Shaded region depicts the magnitude of the CPD. Identical axis labels are omitted in B and C. (B) Exogenous attention modulation. Exogenous attention improves segmentation performance in the periphery and impairs it near the fovea. (C) Endogenous attention modulation. Endogenous attention improves segmentation performance across eccentricity.Whereas our analyses focused on texture segmentation, our model is general and can be applied to other visual phenomena. We show that the model predicts the effects of attention on contrast sensitivity and acuity, i.e., in tasks in which both endogenous and exogenous attention have similar or differential effects on performance. To preview our results, model comparisons revealed that normalization is necessary to elicit the CPD and that separate profiles of gain enhancement across SF (26) generate the effects of exogenous and endogenous attention on texture segmentation. A preferential high-SF enhancement reproduces the impairments by exogenous attention due to a shift in visual sensitivity toward details too fine to distinguish the target at foveal locations. The transition from impairments to improvements in the periphery results from exogenous attentional gain gradually shifting to lower SFs that are more amenable for target detection. Improvements by endogenous attention result from a uniform enhancement of SFs that encompass the target, optimizing visual sensitivity for the attended stimulus across eccentricity.