Computational Insights into SARS-CoV-2 Main Protease Mutations and Nirmatrelvir Efficacy: The Effects of P132H and P132H-A173V.

Nirmatrelvir, a pivotal component of the oral antiviral Paxlovid for COVID-19, targets the SARS-CoV-2 main protease (Mpro) as a covalent inhibitor. Here, we employed combined computational methods to explore how the prevalent Omicron variant mutation P132H, alone and in combination with A173V (P132H-A173V), affects nirmatrelvir's efficacy. Our findings suggest that P132H enhances the noncovalent binding affinity of Mpro for nirmatrelvir, whereas P132H-A173V diminishes it. Although both mutants catalyze the rate-limiting step more efficiently than the wild-type (WT) Mpro, P132H slows the overall rate of covalent bond formation, whereas P132H-A173V accelerates it. Comprehensive analysis of noncovalent and covalent contributions to the overall binding free energy of the covalent complex suggests that P132H likely enhances Mpro sensitivity to nirmatrelvir, while P132H-A173V may confer resistance. Per-residue decompositions of the binding and activation free energies pinpoint key residues that significantly affect the binding affinity and reaction rates, revealing how the mutations modulate these effects. The mutation-induced conformational perturbations alter drug-protein local contact intensities and the electrostatic preorganization of the protein, affecting noncovalent binding affinity and the stability of key reaction states, respectively. Our findings inform the mechanisms of nirmatrelvir resistance and sensitivity, facilitating improved drug design and the detection of resistant strains.

Identification of and Mechanistic Insights into SARS-CoV-2 Main Protease Non-Covalent Inhibitors: An In-Silico Study.

The indispensable role of the SARS-CoV-2 main protease (Mpro) in the viral replication cycle and its dissimilarity to human proteases make Mpro a promising drug target. In order to identify the non-covalent Mpro inhibitors, we performed a comprehensive study using a combined computational strategy. We first screened the ZINC purchasable compound database using the pharmacophore model generated from the reference crystal structure of Mpro complexed with the inhibitor ML188. The hit compounds were then filtered by molecular docking and predicted parameters of drug-likeness and pharmacokinetics. The final molecular dynamics (MD) simulations identified three effective candidate inhibitors (ECIs) capable of maintaining binding within the substrate-binding cavity of Mpro. We further performed comparative analyses of the reference and effective complexes in terms of dynamics, thermodynamics, binding free energy (BFE), and interaction energies and modes. The results reveal that, when compared to the inter-molecular electrostatic forces/interactions, the inter-molecular van der Waals (vdW) forces/interactions are far more important in maintaining the association and determining the high affinity. Given the un-favorable effects of the inter-molecular electrostatic interactions-association destabilization by the competitive hydrogen bond (HB) interactions and the reduced binding affinity arising from the un-compensable increase in the electrostatic desolvation penalty-we suggest that enhancing the inter-molecular vdW interactions while avoiding introducing the deeply buried HBs may be a promising strategy in future inhibitor optimization.

Variances and covariances of linear summary statistics of segregating sites.

Each mutation in a population sample of DNA sequences can be classified by the number of sequences that inherit the mutant nucleotide, the resulting frequencies are known as mutations of different sizes or site frequency spectrum. Many summary statistics can be defined as a linear function of these frequencies. A flexible class of such linear summary statistics is explored analytically in this paper which include several well-known quantities, such as the number of segregating sizes and the mean number of nucleotide differences between two sequences. Some asymptotic variances and covariances are obtained while the analytical formulas for the variances and covariances of nine such linear summary statistics are derived, most of which are unknown to date. This study not only provides some theoretical foundations for exploring linear summary statistics, but also provides some newlinear summary statistics that may be utilized for analyzing sample polymorphism. Furthermore it is showed that a newly developed linear summary statistics has a smaller variance almost uniformly than Watterson's estimator, and that a class of linear summary statistics given too heavy weights on mutations of smaller sizes result in asymptotically non-zero variance.

Amitosis as a strategy of cell division-Insight from the proliferation of Tetrahymena thermophila macronuclei.

Cell division is a necessity of life which can be either mitotic or amitotic. While both are fundamental, amitosis is sometimes considered a relic of little importance in biology. Nevertheless, eukaryotes often have polyploid cells, including cancer cells, which may divide amitotically. To understand how amitosis ensures the completion of cell division, we turn to the macronuclei of ciliates. The grand scheme governing the proliferation of the macronuclei of ciliate cells, which involves chromosomal replication and amitosis, is currently unknown, which is crucial for developing population genetics model of ciliate populations. Using a novel model that encompasses a wide range of mechanisms together with experimental data of the composition of mating types at different stages derived from a single karyonide of Tetrahymena thermophila, we show that the chromosomal replication of the macronucleus has a strong head-start effect, with only about five copies of chromosomes replicated at a time and persistent reuse of the chromosomes involved in the early replication. Furthermore the fission of a fully grown macronucleus is non-random with regard to chromosome composition, with a strong tendency to push chromosomes and their replications to the same daughter cell.

Exploring the Cold-Adaptation Mechanism of Serine Hydroxymethyltransferase by Comparative Molecular Dynamics Simulations.

Cold-adapted enzymes feature a lower thermostability and higher catalytic activity compared to their warm-active homologues, which are considered as a consequence of increased flexibility of their molecular structures. The complexity of the (thermo)stability-flexibility-activity relationship makes it difficult to define the strategies and formulate a general theory for enzyme cold adaptation. Here, the psychrophilic serine hydroxymethyltransferase (pSHMT) from Psychromonas ingrahamii and its mesophilic counterpart, mSHMT from Escherichia coli, were subjected to µs-scale multiple-replica molecular dynamics (MD) simulations to explore the cold-adaptation mechanism of the dimeric SHMT. The comparative analyses of MD trajectories reveal that pSHMT exhibits larger structural fluctuations and inter-monomer positional movements, a higher global flexibility, and considerably enhanced local flexibility involving the surface loops and active sites. The largest-amplitude motion mode of pSHMT describes the trends of inter-monomer dissociation and enlargement of the active-site cavity, whereas that of mSHMT characterizes the opposite trends. Based on the comparison of the calculated structural parameters and constructed free energy landscapes (FELs) between the two enzymes, we discuss in-depth the physicochemical principles underlying the stability-flexibility-activity relationships and conclude that (i) pSHMT adopts the global-flexibility mechanism to adapt to the cold environment and, (ii) optimizing the protein-solvent interactions and loosening the inter-monomer association are the main strategies for pSHMT to enhance its flexibility.

Molecular dynamics simulations reveal distinct differences in conformational dynamics and thermodynamics between the unliganded and CD4-bound states of HIV-1 gp120.

The entry of human immunodeficiency virus type I (HIV-1) into host cells is initiated by binding to the cell-surface receptor CD4, which induces a conformational transition of the envelope (Env) glycoprotein gp120 from the closed, unliganded state to the open, CD4-bound state. Despite many available structures in these two states, detailed aspects on the dynamics and thermodynamics of gp120 remain elusive. Here, we performed microsecond-scale (µs-scale) multiple-replica molecular dynamics (MD) simulations to explore the differences in the conformational dynamics, protein motions, and thermodynamics between the unliganded and CD4-bound/complexed forms of gp120. Comparative analyses of MD trajectories reveal that CD4 binding promotes the structural deviations/changes and conformational flexibility, loosens the structural packing, and complicates the molecular motions of gp120. Comparison of the constructed free energy landscapes (FELs) reveals that the CD4-complexed gp120 has more conformational substates, larger conformational entropy, and lower thermostability than the unliganded form. Therefore, the unliganded conformation represents a structurally and energetically stable "ground state" for the full-length gp120. The observed great increase in the mobility of V1/V2 and V3 along with their more versatile movement directions in the CD4-bound gp120 compared to the unliganded form suggests that their orientations with respect to each other and to the structural core determine the differences in the conformational dynamics and thermodynamics between the two gp120 forms. The results presented here provide a basis by which to better understand the functional and immunological properties of gp120 and, furthermore, to deploy appropriate strategies for the development of anti-HIV-1 drugs or vaccines.

Large numbers of vertebrates began rapid population decline in the late 19th century.

Accelerated losses of biodiversity are a hallmark of the current era. Large declines of population size have been widely observed and currently 22,176 species are threatened by extinction. The time at which a threatened species began rapid population decline (RPD) and the rate of RPD provide important clues about the driving forces of population decline and anticipated extinction time. However, these parameters remain unknown for the vast majority of threatened species. Here we analyzed the genetic diversity data of nuclear and mitochondrial loci of 2,764 vertebrate species and found that the mean genetic diversity is lower in threatened species than in related nonthreatened species. Our coalescence-based modeling suggests that in many threatened species the RPD began ∼123 y ago (a 95% confidence interval of 20-260 y). This estimated date coincides with widespread industrialization and a profound change in global living ecosystems over the past two centuries. On average the population size declined by ∼25% every 10 y in a threatened species, and the population size was reduced to ∼5% of its ancestral size. Moreover, the ancestral size of threatened species was, on average, ∼22% smaller than that of nonthreatened species. Because the time period of RPD is short, the cumulative effect of RPD on genetic diversity is still not strong, so that the smaller ancestral size of threatened species may be the major cause of their reduced genetic diversity; RPD explains 24.1-37.5% of the difference in genetic diversity between threatened and nonthreatened species.

Age-dependent transition from cell-level to population-level control in murine intestinal homeostasis revealed by coalescence analysis.

In multi-cellular organisms, tissue homeostasis is maintained by an exquisite balance between stem cell proliferation and differentiation. This equilibrium can be achieved either at the single cell level (a.k.a. cell asymmetry), where stem cells follow strict asymmetric divisions, or the population level (a.k.a. population asymmetry), where gains and losses in individual stem cell lineages are randomly distributed, but the net effect is homeostasis. In the mature mouse intestinal crypt, previous evidence has revealed a pattern of population asymmetry through predominantly symmetric divisions of stem cells. In this work, using population genetic theory together with previously published crypt single-cell data obtained at different mouse life stages, we reveal a strikingly dynamic pattern of stem cell homeostatic control. We find that single-cell asymmetric divisions are gradually replaced by stochastic population-level asymmetry as the mouse matures to adulthood. This lifelong process has important developmental and evolutionary implications in understanding how adult tissues maintain their homeostasis integrating the trade-off between intrinsic and extrinsic regulations.

Effect of the Solvent Temperatures on Dynamics of Serine Protease Proteinase K.

To obtain detailed information about the effect of the solvent temperatures on protein dynamics, multiple long molecular dynamics (MD) simulations of serine protease proteinase K with the solute and solvent coupled to different temperatures (either 300 or 180 K) have been performed. Comparative analyses demonstrate that the internal flexibility and mobility of proteinase K are strongly dependent on the solvent temperatures but weakly on the protein temperatures. The constructed free energy landscapes (FELs) at the high solvent temperatures exhibit a more rugged surface, broader spanning range, and higher minimum free energy level than do those at the low solvent temperatures. Comparison between the dynamic hydrogen bond (HB) numbers reveals that the high solvent temperatures intensify the competitive HB interactions between water molecules and protein surface atoms, and this in turn exacerbates the competitive HB interactions between protein internal atoms, thus enhancing the conformational flexibility and facilitating the collective motions of the protein. A refined FEL model was proposed to explain the role of the solvent mobility in facilitating the cascade amplification of microscopic motions of atoms and atomic groups into the global collective motions of the protein.

Highly variable recessive lethal or nearly lethal mutation rates during germ-line development of male Drosophila melanogaster.

Each cell of higher organism adults is derived from a fertilized egg through a series of divisions, during which mutations can occur. Both the rate and timing of mutations can have profound impacts on both the individual and the population, because mutations that occur at early cell divisions will affect more tissues and are more likely to be transferred to the next generation. Using large-scale multigeneration screening experiments for recessive lethal or nearly lethal mutations of Drosophila melanogaster and recently developed statistical analysis, we show for male D. melanogaster that (i) mutation rates (for recessive lethal or nearly lethal) are highly variable during germ cell development; (ii) first cell cleavage has the highest mutation rate, which drops substantially in the second cleavage or the next few cleavages; (iii) the intermediate stages, after a few cleavages to right before spermatogenesis, have at least an order of magnitude smaller mutation rate; and (iv) spermatogenesis also harbors a fairly high mutation rate. Because germ-line lineage shares some (early) cell divisions with somatic cell lineage, the first conclusion is readily extended to a somatic cell lineage. It is conceivable that the first conclusion is true for most (if not all) higher organisms, whereas the other three conclusions are widely applicable, although the extent may differ from species to species. Therefore, conclusions or analyses that are based on equal mutation rates during development should be taken with caution. Furthermore, the statistical approach developed can be adopted for studying other organisms, including the human germ-line or somatic mutational patterns.

Estimating DNA polymorphism from next generation sequencing data with high error rate by dual sequencing applications.

BACKGROUND: As the error rate is high and the distribution of errors across sites is non-uniform in next generation sequencing (NGS) data, it has been a challenge to estimate DNA polymorphism (θ) accurately from NGS data. RESULTS: By computer simulations, we compare the two methods of data acquisition - sequencing each diploid individual separately and sequencing the pooled sample. Under the current NGS error rate, sequencing each individual separately offers little advantage unless the coverage per individual is high (>20X). We hence propose a new method for estimating θ from pooled samples that have been subjected to two separate rounds of DNA sequencing. Since errors from the two sequencing applications are usually non-overlapping, it is possible to separate low frequency polymorphisms from sequencing errors. Simulation results show that the dual applications method is reliable even when the error rate is high and θ is low. CONCLUSIONS: In studies of natural populations where the sequencing coverage is usually modest (~2X per individual), the dual applications method on pooled samples should be a reasonable choice.

Estimating population genetic parameters and comparing model goodness-of-fit using DNA sequences with error.

It is known that sequencing error can bias estimation of evolutionary or population genetic parameters. This problem is more prominent in deep resequencing studies because of their large sample size n, and a higher probability of error at each nucleotide site. We propose a new method based on the composite likelihood of the observed SNP configurations to infer population mutation rate theta = 4N(e)micro, population exponential growth rate R, and error rate epsilon, simultaneously. Using simulation, we show the combined effects of the parameters, theta, n, epsilon, and R on the accuracy of parameter estimation. We compared our maximum composite likelihood estimator (MCLE) of theta with other theta estimators that take into account the error. The results show the MCLE performs well when the sample size is large or the error rate is high. Using parametric bootstrap, composite likelihood can also be used as a statistic for testing the model goodness-of-fit of the observed DNA sequences. The MCLE method is applied to sequence data on the ANGPTL4 gene in 1832 African American and 1045 European American individuals.

Significantly fewer protein functional changing variants for lipid metabolism in Africans than in Europeans.

BACKGROUND: The disorders in metabolism of energy substances are usually related to some diseases, such as obesity, diabetes and cancer, etc. However, the genetic background for these disorders has not been well understood. In this study, we explored the genetic risk differences among human populations in metabolism (catabolism and biosynthesis) of energy substances, including lipids, carbohydrates and amino acids. RESULTS: Two genotype datasets (Hapmap and 1000 Genome) were used for this study. The genetic risks of protein functional changing variants (PFCVs) on genes involved in lipid, carbohydrate and amino acid metabolism were calculated using two genetic risk indices: the total number of PFCVs (Num) and the total possibly harmful score of PFCVs (R). Observations in these two genotype datasets consistently showed that Africans had lower genetic risk in lipid metabolism (both catabolic and biosynthetic processes) compared to Europeans. However this relationship was not observed in carbohydrate and amino acid metabolism. CONCLUSIONS: Our results suggested that Africans had higher efficiency of utilizing lipids as energy substances than Europeans. In other words, lipids might be more preferred as energy substances in Africans than in Europeans.

Genomic inference of a severe human bottleneck during the Early to Middle Pleistocene transition.

Population size history is essential for studying human evolution. However, ancient population size history during the Pleistocene is notoriously difficult to unravel. In this study, we developed a fast infinitesimal time coalescent process (FitCoal) to circumvent this difficulty and calculated the composite likelihood for present-day human genomic sequences of 3154 individuals. Results showed that human ancestors went through a severe population bottleneck with about 1280 breeding individuals between around 930,000 and 813,000 years ago. The bottleneck lasted for about 117,000 years and brought human ancestors close to extinction. This bottleneck is congruent with a substantial chronological gap in the available African and Eurasian fossil record. Our results provide new insights into our ancestry and suggest a coincident speciation event.

Mechanistic Origin of Different Binding Affinities of SARS-CoV and SARS-CoV-2 Spike RBDs to Human ACE2.

The receptor-binding domain (RBD) of the SARS-CoV-2 spike protein (RBDCoV2) has a higher binding affinity to the human receptor angiotensin-converting enzyme 2 (ACE2) than the SARS-CoV RBD (RBDCoV). Here, we performed molecular dynamics (MD) simulations, binding free energy (BFE) calculations, and interface residue contact network (IRCN) analysis to explore the mechanistic origin of different ACE2-binding affinities of the two RBDs. The results demonstrate that, when compared to the RBDCoV2-ACE2 complex, RBDCoV-ACE2 features enhanced dynamicsand inter-protein positional movements and increased conformational entropy and conformational diversity. Although the inter-protein electrostatic attractive interactions are the primary determinant for the high ACE2-binding affinities of both RBDs, the significantly enhanced electrostatic attractive interactions between ACE2 and RBDCoV2 determine the higher ACE2-binding affinity of RBDCoV2 than of RBDCoV. Comprehensive comparative analyses of the residue BFE components and IRCNs between the two complexes reveal that it is the residue changes at the RBD interface that lead to the overall stronger inter-protein electrostatic attractive force in RBDCoV2-ACE2, which not only tightens the interface packing and suppresses the dynamics of RBDCoV2-ACE2, but also enhances the ACE2-binding affinity of RBDCoV2. Since the RBD residue changes involving gain/loss of the positively/negatively charged residues can greatly enhance the binding affinity, special attention should be paid to the SARS-CoV-2 variants carrying such mutations, particularly those near or at the binding interfaces with the potential to form hydrogen bonds and/or salt bridges with ACE2.

A Deep Learning Approach for Predicting Antigenic Variation of Influenza A H3N2.

Modeling antigenic variation in influenza (flu) virus A H3N2 using amino acid sequences is a promising approach for improving the prediction accuracy of immune efficacy of vaccines and increasing the efficiency of vaccine screening. Antigenic drift and antigenic jump/shift, which arise from the accumulation of mutations with small or moderate effects and from a major, abrupt change with large effects on the surface antigen hemagglutinin (HA), respectively, are two types of antigenic variation that facilitate immune evasion of flu virus A and make it challenging to predict the antigenic properties of new viral strains. Despite considerable progress in modeling antigenic variation based on the amino acid sequences, few studies focus on the deep learning framework which could be most suitable to be applied to this task. Here, we propose a novel deep learning approach that incorporates a convolutional neural network (CNN) and bidirectional long-short-term memory (BLSTM) neural network to predict antigenic variation. In this approach, CNN extracts the complex local contexts of amino acids while the BLSTM neural network captures the long-distance sequence information. When compared to the existing methods, our deep learning approach achieves the overall highest prediction performance on the validation dataset, and more encouragingly, it achieves prediction agreements of 99.20% and 96.46% for the strains in the forthcoming year and in the next two years included in an existing set of chronological amino acid sequences, respectively. These results indicate that our deep learning approach is promising to be applied to antigenic variation prediction of flu virus A H3N2.

Inferring population mutation rate and sequencing error rate using the SNP frequency spectrum in a sample of DNA sequences.

One challenge of analyzing samples of DNA sequences is to account for the nonnegligible polymorphisms produced by error when the sequencing error rate is high or the sample size is large. Specifically, those artificial sequence variations will bias the observed single nucleotide polymorphism (SNP) frequency spectrum, which in turn may further bias the estimators of the population mutation rate theta =4N mu for diploids. In this paper, we propose a new approach based on the generalized least squares (GLS) method to estimate theta, given a SNP frequency spectrum in a random sample of DNA sequences from a population. With this approach, error rate epsilon can be either known or unknown. In the latter case, epsilon can be estimated given an estimation of theta. Using coalescent simulation, we compared our estimators with other estimators of theta. The results showed that the GLS estimators are more efficient than other theta estimators with error, and the estimation of epsilon is usable in practice when the theta per bp is small. We demonstrate the application of the estimators with 10-kb noncoding region sequence sampled from a human population and provide suggestions for choosing theta estimators with error.

Recombination facilitates neofunctionalization of duplicate genes via originalization.

BACKGROUND: Recently originalization was proposed to be an effective way of duplicate-gene preservation, in which recombination provokes the high frequency of original (or wild-type) allele on both duplicated loci. Because the high frequency of wild-type allele might drive the arising and accumulating of advantageous mutation, it is hypothesized that recombination might enlarge the probability of neofunctionalization (Pneo) of duplicate genes. In this article this hypothesis has been tested theoretically. RESULTS: Results show that through originalization recombination might not only shorten mean time to neofunctionalizaiton, but also enlarge Pneo. CONCLUSIONS: Therefore, recombination might facilitate neofunctionalization via originalization. Several extensive applications of these results on genomic evolution have been discussed: 1. Time to nonfunctionalization can be much longer than a few million generations expected before; 2. Homogenization on duplicated loci results from not only gene conversion, but also originalization; 3. Although the rate of advantageous mutation is much small compared with that of degenerative mutation, Pneo cannot be expected to be small.

Stairway Plot 2: demographic history inference with folded SNP frequency spectra.

Inferring the demographic histories of populations has wide applications in population, ecological, and conservation genomics. We present Stairway Plot 2, a cross-platform program package for this task using SNP frequency spectra. It is based on a nonparametric method with the capability of handling folded SNP frequency spectra (that is, when the ancestral alleles of the SNPs are unknown) of thousands of samples produced with genotyping-by-sequencing technologies; therefore, it is particularly suitable for nonmodel organisms.

The Energetic Origin of Different Catalytic Activities in Temperature-Adapted Trypsins.

Psychrophilic enzymes were always observed to have higher catalytic activity (k cat) than their mesophilic homologs at room temperature, while the origin of this phenomenon remains obscure. Here, we used two different temperature-adapted trypsins, the psychrophilic Atlantic cod trypsin (ACT) and the mesophilic bovine trypsin (BT), as a model system to explore the energetic origin of their different catalytic activities using computational methods. The results reproduce the characteristic changing trends in the activation free energy, activation enthalpy, and activation entropy between the psychrophilic and mesophilic enzymes, where, in particular, the slightly decreased activation free energy of ACT is determined by its considerably reduced activation enthalpy rather than by its more negative activation entropy compared to BT. The calculated electrostatic contributions to the solvation free energies in the reactant state/ground sate (RS/GS) and transition state (TS) show that, going from BT to ACT, the TS stabilization has a predominant effect over the RS stabilization on lowering the activation enthalpy of ACT. Comparison between the solvation energy components reveals a more optimized electrostatic preorganization to the TS in ACT, which provides a larger stabilization to the TS through reducing the reorganization energy, thus resulting in the lower activation enthalpy and hence lower activation free energy of ACT. Thus, it can be concluded that it is the difference in the protein electrostatic environment, and hence its different stabilizing effects on the TS, that brings about the different catalytic activities of different temperature-adapted trypsins.