Predictive screening of metal–organic framework (MOF) materials for their gas uptake properties has been previously limited by using data from a range of simulated sources, meaning the final predictions are dependent on the performance of these original models. In this work, experimental gas uptake data has been used to create a Gradient Boosted Tree model for the prediction of H2, CH4, and CO2 uptake over a range of temperatures and pressures in MOF materials. The descriptors used in this database were obtained from the literature, with no computational modeling needed. This model was repeated 10 times, showing an average R2 of 0.86 and a mean absolute error (MAE) of ±2.88 wt % across the runs. This model will provide gas uptake predictions for a range of gases, temperatures, and pressures as a one-stop solution, with the data provided being based on previous experimental observations in the literature, rather than simulations, which may differ from their real-world results. The objective of this work is to create a machine learning model for the inference of gas uptake in MOFs. The basis of model development is experimental as opposed to simulated data to realize its applications by practitioners. The real-world nature of this research materializes in a focus on the application of algorithms as opposed to the detailed assessment of the algorithms.

Chemical equilibria of multimerizations in systems with small numbers of particles exhibit a behavior seemingly at odds with that observed macroscopically. In this paper, we apply the recently proposed expression of equilibrium constant for binding, which includes cross-correlations in reactants’ concentrations, to write an equilibrium constant for the formation of clusters larger than two (e.g., trimer, tetramer, and pentamer) as series of two-body reactions. Results obtained by molecular dynamics simulations demonstrate that the value of this expression is constant for all concentrations and system sizes, as well as at an onset of a phase transition to an aggregated state, where densities in the system change discontinuously. In contrast, the value of the commonly utilized expression of equilibrium constant, which ignores correlations, is not constant and its variations can reach few orders of magnitude. Considering different paths for the same multimer formation, with elementary reactions of any order, yields different expressions for the equilibrium constant, yet, with exactly the same value. This is also true for routes with essentially zero probability to occur. Existence of different expressions for the same equilibrium constant imposes equalities between averages of correlated, along with uncorrelated, concentrations of participating species. Moreover, a relation between an average particle number and relative fluctuations derived for two-body reactions is found to be obeyed here as well despite couplings to additional equilibrium reactions in the system. Analyses of transfer reactions, where association and dissociation events take place on both sides of the chemical equation, further indicate the necessity to include cross-correlations in the expression of the equilibrium constant. However, in this case, the magnitudes of discrepancies of the uncorrelated expression are smaller, likely because of partial cancellation of correlations, which exist on both the reactant and product sides.

Transmissibility of SARS-CoV-2 initially relies on its trimeric Spike-RBDs to tether the ACE-2 on host cells, and enhanced self-association of ACE-2 engaged with Spike facilitates the viral infection. Two primary packing modes of Spike-ACE2 heteroproteins exist potentially due to discrepant amounts of RBDs loading on ACE-2, but the resultant self-association difference is inherently unclear. We used extensive coarse-grained dynamic simulations to characterize the self-association efficiency, the conformation relevance, and the molecular mechanism of ACE-2 with different RBD amounts. It was revealed that the ACE-2 hanging two/full RBDs (Mode-A) rapidly dimerized into the heteroprotein complex in a compact “linear” conformation, while the bare ACE-2 showed weakened self-association and a protein complex. The RBD-tethered ectodomains of ACE-2 presented a more upright conformation relative to the membrane, and the intermolecular ectodomains were predominantly packed by the neck domains, which was obligated to the rapid protein self-association in a compact pattern. Noted is the fact that the ACE-2 tethered by a single RBD (Mode-B) retained considerable self-association efficiency and clustering capability, which unravels the interrelation of ACE-2 colocalization and protein cross-linkage. The molecular perspectives in this study expound the self-association potency of ACE-2 with different RBD amounts and the viral activity implications, which can greatly enhance our comprehension of SARS-CoV-2 infection details.

As a versatile polymer in many applications, synthesized polyethylenimine (PEI) is polydisperse with diverse branched structures that attain pH-dependent protonation states. Understanding the structure–function relationship of PEI is necessary for enhancing its efficacy in various applications. Coarse-grained (CG) simulations can be performed at length and time scales directly comparable with experimental data while maintaining the molecular perspective. However, manually developing CG forcefields for complex PEI structures is time-consuming and prone to human errors. This article presents a fully automated algorithm that can coarse-grain any branched architecture of PEI from its all-atom (AA) simulation trajectories and topology. The algorithm is demonstrated by coarse-graining a branched 2 kDa PEI, which can replicate the AA diffusion coefficient, radius of gyration, and end-to-end distance of the longest linear chain. Commercially available 25 and 2 kDa Millipore-Sigma PEIs are used for experimental validation. Specifically, branched PEI architectures are proposed, coarse-grained using the automated algorithm, and then simulated at different mass concentrations. The CG PEIs can reproduce existing experimental data on PEI’s diffusion coefficient and Stokes–Einstein radius at infinite dilution as well as its intrinsic viscosity. This suggests a strategy where probable chemical structures of synthetic PEIs can be inferred computationally using the developed algorithm. The coarse-graining methodology presented here can also be extended to other polymers.

G-protein coupled receptors (GPCRs) are the most prominent family of membrane proteins that serve as major targets for one-third of the drugs produced. A detailed understanding of the molecular mechanism of drug-induced activation and inhibition of GPCRs is crucial for the rational design of novel therapeutics. The binding of the neurotransmitter adrenaline to the β2-adrenergic receptor (β2AR) is known to induce a flight or fight cellular response, but much remains to be understood about binding-induced dynamical changes in β2AR and adrenaline. In this article, we examine the potential of mean force (PMF) for the unbinding of adrenaline from the orthosteric binding site of β2AR and the associated dynamics using umbrella sampling and molecular dynamics (MD) simulations. The calculated PMF reveals a global energy minimum, which corresponds to the crystal structure of β2AR–adrenaline complex, and a meta-stable state in which the adrenaline is moved slightly deeper into the binding pocket with a different orientation compared to that in the crystal structure. The orientational and conformational changes in adrenaline during the transition between these two states and the underlying driving forces of this transition are also explored. Based on the clustering of MD configurations and machine learning-based statistical analyses of time series of relevant collective variables, the structures and stabilizing interactions of these two states of the β2AR–adrenaline complex are also investigated.

Knowledge of critical properties, such as critical temperature, pressure, density, as well as acentric factor, is essential to calculate thermo-physical properties of chemical compounds. Experiments to determine critical properties and acentric factors are expensive and time intensive; therefore, we developed a machine learning (ML) model that can predict these molecular properties given the SMILES representation of a chemical species. We explored directed message passing neural network (D-MPNN) and graph attention network as ML architecture choices. Additionally, we investigated featurization with additional atomic and molecular features, multitask training, and pretraining using estimated data to optimize model performance. Our final model utilizes a D-MPNN layer to learn the molecular representation and is supplemented by Abraham parameters. A multitask training scheme was used to train a single model to predict all the critical properties and acentric factors along with boiling point, melting point, enthalpy of vaporization, and enthalpy of fusion. The model was evaluated on both random and scaffold splits where it shows state-of-the-art accuracies. The extensive data set of critical properties and acentric factors contains 1144 chemical compounds and is made available in the public domain together with the source code that can be used for further exploration.

Tyrosine kinases are a subfamily of kinases with critical roles in cellular machinery. Dysregulation of their active or inactive forms is associated with diseases like cancer. This study aimed to holistically understand their flexibility–activity relationships, focusing on pockets and fluctuations. We studied 43 different tyrosine kinases by collecting 120 μs of molecular dynamics simulations, pocket and residue fluctuation analysis, and a complementary machine learning approach. We found that the inactive forms often have increased flexibility, particularly at the DFG motif level. Noteworthy, thanks to these long simulations combined with a decision tree, we identified a semiquantitative fluctuation threshold of the DGF+3 residue over which the kinase has a higher probability to be in the inactive form.

Single-cell RNA sequencing (scRNA-seq) is widely used to reveal heterogeneity in cells, which has given us insights into cell–cell communication, cell differentiation, and differential gene expression. However, analyzing scRNA-seq data is a challenge due to sparsity and the large number of genes involved. Therefore, dimensionality reduction and feature selection are important for removing spurious signals and enhancing the downstream analysis. We present Correlated Clustering and Projection (CCP), a new data-domain dimensionality reduction method, for the first time. CCP projects each cluster of similar genes into a supergene defined as the accumulated pairwise nonlinear gene–gene correlations among all cells. Using 14 benchmark data sets, we demonstrate that CCP has significant advantages over classical principal component analysis (PCA) for clustering and/or classification problems with intrinsically high dimensionality. In addition, we introduce the Residue-Similarity index (RSI) as a novel metric for clustering and classification and the R-S plot as a new visualization tool. We show that the RSI correlates with accuracy without requiring the knowledge of the true labels. The R-S plot provides a unique alternative to the uniform manifold approximation and projection (UMAP) and t-distributed stochastic neighbor embedding (t-SNE) for data with a large number of cell types.

The relative energy gradient (REG) method is paired with the topological energy partitioning method interacting quantum atoms (IQA), as REG-IQA, to provide detailed and unbiased knowledge on the intra- and interatomic interactions. REG operates on a sequence of geometries representing a dynamical change of a system. Its recent application to peptide hydrolysis of the human immunodeficiency virus-1 (HIV-1) protease (PDB code: 4HVP) has demonstrated its full potential in recovering reaction mechanisms and through-space electrostatic and exchange–correlation effects, making it a compelling tool for analyzing enzymatic reactions. In this study, the computational efficiency of the REG-IQA method for the 133-atom HIV-1 protease quantum mechanical system is analyzed in every detail and substantially improved by means of three different approaches. The first approach of smaller integration grids for IQA integrations reduces the computational overhead by about a factor of 3. The second approach uses the line-simplification Ramer–Douglas–Peucker (RDP) algorithm, which outputs the minimal number of geometries necessary for the REG-IQA analysis for a predetermined root mean squared error (RMSE) tolerance. This cuts the computational time of the whole REG analysis by a factor of 2 if an RMSE of 0.5 kJ/mol is considered. The third approach consists of a “biased” or “unbiased” selection of a specific subset of atoms of the whole initial quantum mechanical model wave-function, which results in more than a 10-fold speed-up per geometry for the IQA calculation, without deterioration of the outcome of the REG-IQA analysis. Finally, to show the capability of these approaches, the findings gathered from the HIV-1 protease system are also applied to a different system named haloalcohol dehalogenase (HheC). In summary, this study takes the REG-IQA method to a computationally feasible and highly accurate level, making it viable for the analysis of a multitude of enzymatic systems.

Machine-learning and deep-learning models have been extensively used in cheminformatics to predict molecular properties, to reduce the need for direct measurements, and to accelerate compound prioritization. However, different setups and frameworks and the large number of molecular representations make it difficult to properly evaluate, reproduce, and compare them. Here we present a new PREdictive modeling FramEwoRk for molecular discovery (PREFER), written in Python (version 3.7.7) and based on AutoSklearn (version 0.14.7), that allows comparison between different molecular representations and common machine-learning models. We provide an overview of the design of our framework and show exemplary use cases and results of several representation–model combinations on diverse data sets, both public and in-house. Finally, we discuss the use of PREFER on small data sets. The code of the framework is freely available on GitHub.

The contradictory behaviors in light harvesting and non-photochemical quenching make xanthophyll lutein the most attractive functional molecule in photosynthesis. Despite several theoretical simulations on the spectral properties and excited-state dynamics, the atomic-level photophysical mechanisms need to be further studied and established, especially for an accurate description of geometric and electronic structures of conical intersections for the lowest several electronic states of lutein. In the present work, semiempirical OM2/MRCI and multi-configurational restricted active space self-consistent field methods were performed to optimize the minima and conical intersections in and between the 1Ag–, 2Ag–, 1Bu+, and 1Bu– states. Meanwhile, the relative energies were refined by MS-CASPT2(10,8)/6-31G*, which can reproduce correct electronic state properties as those in the spectroscopic experiments. Based on the above calculation results, we proposed a possible excited-state relaxation mechanism for lutein from its initially populated 1Bu+ state. Once excited to the optically bright 1Bu+ state, the system will propagate along the key reaction coordinate, i.e., the stretching vibration of the conjugated carbon chain. During this period of time, the 1Bu– state will participate in and forms a resonance state between the 1Bu– and 1Bu+ states. Later, the system will rapidly hop to the 2Ag– state via the 1Bu+/2Ag– conical intersection. Finally, the lutein molecule will survive in the 2Ag– state for a relatively long time before it internally converts to the ground state directly or via a twisted S1/S0 conical intersection. Notably, though the photophysical picture may be very different in solvents and proteins, the current theoretical study proposed a promising calculation protocol and also provided many valuable mechanistic insights for lutein and similar carotenoids.

The folding/misfolding of membrane-permiable Amyloid beta (Aβ) peptides is likely associated with the advancing stage of Alzheimer’s disease (AD) by disrupting Ca2+ homeostasis. In this context, the aggregation of four transmembrane Aβ17–42 peptides was investigated using temperature replica-exchange molecular dynamics (REMD) simulations. The obtained results indicated that the secondary structure of transmembrane Aβ peptides tends to have different propensities compared to those in solution. Interestingly, the residues favorably forming β-structure were interleaved by residues rigidly adopting turn-structure. A combination of β and turn regions likely forms a pore structure. Six morphologies of 4Aβ were found over the free energy landscape and clustering analyses. Among these, the morphologies include (1) Aβ binding onto the membrane surface and three transmembrane Aβ; (2) three helical and coil transmembrane Aβ; (3) four helical transmembrane Aβ; (4) three helical and one β-hairpin transmembrane Aβ; (5) two helical and two β-strand transmembrane Aβ; and (6) three β-strand and one helical transmembrane Aβ. Although the formation of the β-barrel structure was not observed during the 0.28 ms─long MD simulation, the structure is likely to form when the simulation time is further extended.

Traditional Chinese medicine (TCM) not only maintains the health of Asian people but also provides a great resource of active natural products for modern drug development. Herein, we developed a Database of Constituents Absorbed into the Blood and Metabolites of TCM (DCABM-TCM), the first database systematically collecting blood constituents of TCM prescriptions and herbs, including prototypes and metabolites experimentally detected in the blood, together with the corresponding detailed detection conditions through manual literature mining. The DCABM-TCM has collected 1816 blood constituents with chemical structures of 192 prescriptions and 194 herbs and integrated their related annotations, including physicochemical, absorption, distribution, metabolism, excretion, and toxicity properties, and associated targets, pathways, and diseases. Furthermore, the DCABM-TCM supported two blood constituent-based analysis functions, the network pharmacology analysis for TCM molecular mechanism elucidation, and the target/pathway/disease-based screening of candidate blood constituents, herbs, or prescriptions for TCM-based drug discovery. The DCABM-TCM is freely accessible at http://bionet.ncpsb.org.cn/dcabm-tcm/. The DCABM-TCM will contribute to the elucidation of effective constituents and molecular mechanism of TCMs and the discovery of TCM-derived drug-like compounds that are both bioactive and bioavailable.

Due to the complex nature of noncovalent interactions and their long-range effects, analyzing protein conformations using network theory can be enlightening. Protein Structure Networks (PSNs) provide a convenient formalism to study protein structures in relation to essential properties such as key residues for structural stability, allosteric communication, and the effects of modifications of the protein. PSNs can be defined according to very different principles, and the available tools have limitations in input formats, supported models, and version control. Other outstanding problems are related to the definition of network cutoffs and the assessment of the stability of the network properties. The protein science community could benefit from a common framework to carry out these analyses and make them easier to reproduce, reuse, and evaluate. We here provide two open-source software packages, PyInteraph2 and PyInKnife2, to implement and analyze PSNs in a reproducible and documented manner. PyInteraph2 interfaces with multiple formats for protein ensembles and incorporates different network models with the possibility of integrating them into a macronetwork and performing various downstream analyses, including hubs, connected components, and several other centrality measures, and visualizes the networks or further analyzes them thanks to compatibility with Cytoscape.PyInKnife2 that supports the network models implemented in PyInteraph2. It employs a jackknife resampling approach to estimate the convergence of network properties and streamline the selection of distance cutoffs. We foresee that the modular structure of the code and the supported version control system will promote the transition to a community-driven effort, boost reproducibility, and establish common protocols in the PSN field. As developers, we will guarantee the introduction of new functionalities and maintenance, assistance, and training of new contributors.

Two-photon photodynamic therapy (TP-PDT), as a treatment technology with deep penetration and less damage, provides a broad prospect for cancer treatment. Nowadays, the development of TP-PDT suffers from the low two-photon absorption (TPA) intensity and short triplet state lifetime of photosensitizers (PSs) used in TP-PDT. Herein, we propose some novel modification strategies based on the thionated NpImidazole (the combination of naphthalimide and imidazole) derivatives to make efforts on those issues and obtain corresponding fluorescent probes for detecting ClO– and excellent PSs for TP-PDT. Density functional theory (DFT) and time-dependent DFT (TD-DFT) are used to help us characterize the photophysical properties and TP-PDT process of the newly designed compounds. Our results show that the introduction of different electron-donating groups at the position 4 of NpImidazole can effectively improve their TPA and emission properties. Specifically, 3s with a N,N-dimethylamino group has a large triplet state lifetime (τ = 699 μs) and TPA cross section value (δTPA = 314 GM), which can effectively achieve TP-PDT; additionally, 4s (with electron-donating group 2-oxa-6-azaspiro[3.3]heptane in NpImidazole) effectively realizes the dual-function of a PS for TP-PDT (τ = 25,122 μs, δTPA = 351 GM) and a fluorescent probe for detecting ClO– (Φf = 29% of the product 4o). Moreover, an important problem is clarified from a microscopic perspective, that is, why the transition property of 3s and 4s (1π–π*) from S1 to S0 is different from that of 1s and 2s (1n−π*). It is hoped that our work can provides valuable theoretical clues for the design and synthesis of heavy-atom-free NpImidazole-based PSs and fluorescent probes for the detection of hypochlorite.

Water network rearrangement from the ligand-unbound state to the ligand-bound state is known to have significant effects on the protein–ligand binding interactions, but most of the current machine learning-based scoring functions overlook these effects. In this study, we endeavor to construct a comprehensive and realistic deep learning model by incorporating water network information into both ligand-unbound and -bound states. In particular, extended connectivity interaction features were integrated into graph representation, and graph transformer operator was employed to extract features of the ligand-unbound and -bound states. Through these efforts, we developed a water network-augmented two-state model called ECIFGraph::HM-Holo-Apo. Our new model exhibits satisfactory performance in terms of scoring, ranking, docking, screening, and reverse screening power tests on the CASF-2016 benchmark. In addition, it can achieve superior performance in large-scale docking-based virtual screening tests on the DEKOIS2.0 data set. Our study highlights that the use of a water network-augmented two-state model can be an effective strategy to bolster the robustness and applicability of machine learning-based scoring functions, particularly for targets with hydrophilic or solvent-exposed binding pockets.

Protein tyrosine phosphatase 1B (PTP1B) is a negative regulator of the insulin and leptin signaling pathways, making it a highly attractive target for the treatment of type II diabetes. For PTP1B to perform its enzymatic function, a loop referred to as the “WPD loop” must transition between open (catalytically incompetent) and closed (catalytically competent) conformations, which have both been resolved by X-ray crystallography. Although prior studies have established this transition as the rate-limiting step for catalysis, the transition mechanism for PTP1B and other PTPs has been unclear. Here we present an atomically detailed model of WPD loop transitions in PTP1B based on unbiased, long-timescale molecular dynamics simulations and weighted ensemble simulations. We found that a specific WPD loop region─the PDFG motif─acted as the key conformational switch, with structural changes to the motif being necessary and sufficient for transitions between long-lived open and closed states of the loop. Simulations starting from the closed state repeatedly visited open states of the loop that quickly closed again unless the infrequent conformational switching of the motif stabilized the open state. The functional importance of the PDFG motif is supported by the fact that it is well conserved across PTPs. Bioinformatic analysis shows that the PDFG motif is also conserved, and adopts two distinct conformations, in deiminases, and the related DFG motif is known to function as a conformational switch in many kinases, suggesting that PDFG-like motifs may control transitions between structurally distinct, long-lived conformational states in multiple protein families.

Our skin constitutes an effective permeability barrier that protects the body from exogenous substances but concomitantly severely limits the number of pharmaceutical drugs that can be delivered transdermally. In topical formulation design, chemical permeation enhancers (PEs) are used to increase drug skin permeability. In vitro skin permeability experiments can measure net effects of PEs on transdermal drug transport, but they cannot explain the molecular mechanisms of interactions between drugs, permeation enhancers, and skin structure, which limits the possibility to rationally design better new drug formulations. Here we investigate the effect of the PEs water, lauric acid, geraniol, stearic acid, thymol, ethanol, oleic acid, and eucalyptol on the transdermal transport of metronidazole, caffeine, and naproxen. We use atomistic molecular dynamics (MD) simulations in combination with developed molecular models to calculate the free energy difference between 11 PE-containing formulations and the skin’s barrier structure. We then utilize the results to calculate the final concentration of PEs in skin. We obtain an RMSE of 0.58 log units for calculated partition coefficients from water into the barrier structure. We then use the modified PE-containing barrier structure to calculate the PEs’ permeability enhancement ratios (ERs) on transdermal metronidazole, caffeine, and naproxen transport and compare with the results obtained from in vitro experiments. We show that MD simulations are able to reproduce rankings based on ERs. However, strict quantitative correlation with experimental data needs further refinement, which is complicated by significant deviations between different measurements. Finally, we propose a model for how to use calculations of the potential of mean force of drugs across the skin’s barrier structure in a topical formulation design.

Molecular mechanics rely on existing experimental and theoretical inputs to confidently calculate the trajectories of molecular systems. These calculations, however, are often hindered by missing force field parameters. A notable subject of this problem is metal centers of proteins. This study parameterized, through an adaptive force matching (AFM) workflow, the copper cofactor of plastocyanin in its two oxidation states. New 12-6 Lennard-Jones (LJ) parameters and atomic partial charges were generated to complete the non-bonded description of the copper site. Our models show uniform distorted tetrahedral structures for reduced plastocyanin, Cu(I), and oxidized plastocyanin, Cu(II). These structures align with the QM/MM MD results and existing crystallography studies. TD-DFT calculations, meanwhile, showed that conformations with elongated axial Cu–SMet and shortened equatorial Cu–SCys bonds retain the experimental UV-Vis profile of blue copper (BC) proteins, thus signifying the importance of Cu–S interactions on BC proteins’ unique spectroscopic properties.

Pf2001 esterase (Pf2001) from Pyrococcus furiosus has hyperthermophilic properties and exerts a biocatalytic function in a dimeric state. Crystal structures revealed that the structural rearrangement of the cap domain is responsible for the Pf2001 dimer formation. However, the details of the cap domain remodeling and the effects of temperature on the dimerization process remain elusive at the molecular level, taking into account that experimental methods are difficult to capture the dynamic process of dimerization to some extent. Herein, four dimer models based on the monomeric crystal structure (PDB ID: 5G59) were constructed to investigate the conformational transition details and temperature effects in the dimerization by conventional molecular dynamics and accelerated molecular dynamics simulations. Our simulation results indicate that the monomer undergoes a conformational change into a “preparatory state” at high temperatures, which is more favorable for its transformation into a stable dimer. The subsequent free energy landscape analysis further identifies four intermediate states (from separated state to dimeric state) and discloses that a more accessible α-helix driven by stronger hydrophobic interactions induces a rearrangement of the cap domain, displaying a “tic-tac-toe” activation feature that is important for stabilizing the dimer interface and facilitating the formation of hydrophobic pockets. In addition, the electrostatic potential surface analysis illustrates that the weaker electrostatic repulsion (Lys and Arg) in the dimer interface at high temperatures is also a key factor for dimer stabilization. Altogether, our results can provide molecular-level insight into the dimer formation process of hyperthermophilic esterase and would be useful to understand the enzymatic specificity of α/β-hydrolase.