Analyzing the Geometry and Dynamics of Viral Structures: A Review of Computational Approaches Based on Alpha Shape Theory, Normal Mode Analysis, and Poisson-Boltzmann Theories.

The current SARS-CoV-2 pandemic highlights our fragility when we are exposed to emergent viruses either directly or through zoonotic diseases. Fortunately, our knowledge of the biology of those viruses is improving. In particular, we have more and more structural information on virions, i.e., the infective form of a virus that includes its genomic material and surrounding protective capsid, and on their gene products. It is important to have methods that enable the analyses of structural information on such large macromolecular systems. We review some of those methods in this paper. We focus on understanding the geometry of virions and viral structural proteins, their dynamics, and their energetics, with the ambition that this understanding can help design antiviral agents. We discuss those methods in light of the specificities of those structures, mainly that they are huge. We focus on three of our own methods based on the alpha shape theory for computing geometry, normal mode analyses to study dynamics, and modified Poisson-Boltzmann theories to study the organization of ions and co-solvent and solvent molecules around biomacromolecules. The corresponding software has computing times that are compatible with the use of regular desktop computers. We show examples of their applications on some outer shells and structural proteins of the West Nile Virus.

Conductance of concentrated electrolytes: Multivalency and the Wien effect.

The electric conductivity of ionic solutions is well understood at low ionic concentrations of up to a few millimolar but becomes difficult to unravel at higher concentrations that are still common in nature and technological applications. A model for the conductivity at high concentrations was recently put forth for monovalent electrolytes at low electric fields. The model relies on applying a stochastic density-functional theory and using a modified electrostatic pair-potential that suppresses unphysical, short-range electrostatic interactions. Here, we extend the theory to multivalent ions as well as to high electric fields where a deviation from Ohm's law known as the Wien effect occurs. Our results are in good agreement with experiments and recent simulations.

Sampling constrained stochastic trajectories using Brownian bridges.

We present a new method to sample conditioned trajectories of a system evolving under Langevin dynamics based on Brownian bridges. The trajectories are conditioned to end at a certain point (or in a certain region) in space. The bridge equations can be recast exactly in the form of a non-linear stochastic integro-differential equation. This equation can be very well approximated when the trajectories are closely bundled together in space, i.e., at low temperature, or for transition paths. The approximate equation can be solved iteratively using a fixed point method. We discuss how to choose the initial trajectories and show some examples of the performance of this method on some simple problems. This method allows us to generate conditioned trajectories with a high accuracy.

Conductivity of Concentrated Electrolytes.

The conductivity of ionic solutions is arguably their most important trait, being widely used in electrochemical, biochemical, and environmental applications. The Debye-Hückel-Onsager theory successfully predicts the conductivity at very low ionic concentrations of up to a few millimolars, but there is no well-established theory applicable at higher concentrations. We study the conductivity of ionic solutions using a stochastic density functional theory, paired with a modified Coulomb interaction that accounts for the hard-core repulsion between the ions. The modified potential suppresses unphysical, short-range electrostatic interactions, which are present in the Debye-Hückel-Onsager theory. Our results for the conductivity show very good agreement with experimental data up to 3 molars, without any fit parameters. We provide a compact expression for the conductivity, accompanied by a simple analytical approximation.

Parameterizing elastic network models to capture the dynamics of proteins.

Coarse-grained normal mode analyses of protein dynamics rely on the idea that the geometry of a protein structure contains enough information for computing its fluctuations around its equilibrium conformation. This geometry is captured in the form of an elastic network (EN), namely a network of edges between its residues. The normal modes of a protein are then identified with the normal modes of its EN. Different approaches have been proposed to construct ENs, focusing on the choice of the edges that they are comprised of, and on their parameterizations by the force constants associated with those edges. Here we propose new tools to guide choices on these two facets of EN. We study first different geometric models for ENs. We compare cutoff-based ENs, whose edges have lengths that are smaller than a cutoff distance, with Delaunay-based ENs and find that the latter provide better representations of the geometry of protein structures. We then derive an analytical method for the parameterization of the EN such that its dynamics leads to atomic fluctuations that agree with experimental B-factors. To limit overfitting, we attach a parameter referred to as flexibility constant to each atom instead of to each edge in the EN. The parameterization is expressed as a non-linear optimization problem whose parameters describe both rigid-body and internal motions. We show that this parameterization leads to improved ENs, whose dynamics mimic MD simulations better than ENs with uniform force constants, and reduces the number of normal modes needed to reproduce functional conformational changes.

Fast computation of exact solutions of generic and degenerate assignment problems.

The linear assignment problem is a fundamental problem in combinatorial optimization with a wide range of applications, from operational research to data science. It consists of assigning "agents" to "tasks" on a one-to-one basis, while minimizing the total cost associated with the assignment. While many exact algorithms have been developed to identify such an optimal assignment, most of these methods are computationally prohibitive for large size problems. In this paper, we propose an alternative approach to solving the assignment problem using techniques adapted from statistical physics. Our first contribution is to fully describe this formalism, including all the proofs of its main claims. In particular we derive a strongly concave effective free-energy function that captures the constraints of the assignment problem at a finite temperature. We prove that this free energy decreases monotonically as a function of ß, the inverse of temperature, to the optimal assignment cost, providing a robust framework for temperature annealing. We prove also that for large enough ß values the exact solution to the generic assignment problem can be derived using simple roundoff to the nearest integer of the elements of the computed assignment matrix. Our second contribution is to derive a provably convergent method to handle degenerate assignment problems, with a characterization of those problems. We describe computer implementations of our framework that are optimized for parallel architectures, one based on CPU, the other based on GPU. We show that the latter enables solving large assignment problems (of the orders of a few 10 000s) in computing clock times of the orders of minutes.

Simultaneous Identification of Multiple Binding Sites in Proteins: A Statistical Mechanics Approach.

We present an extension of the Poisson-Boltzmann model in which the solute of interest is immersed in an assembly of self-orienting Langevin water dipoles, anions, cations, and hydrophobic molecules, all of variable densities. Interactions between charges are controlled by electrostatics, while hydrophobic interactions are modeled with a Yukawa potential. We impose steric constraints by assuming that the system is represented on a cubic lattice. We also assume incompressibility; i.e., all sites of the lattice are occupied. This model, which we refer to as the Hydrophobic Dipolar Poisson-Boltzmann Langevin (HDPBL) model, leads to a system of two equations whose solutions give the water dipole, salt, and hydrophobic molecule densities, all of them in the presence of the others in a self-consistent way. We use those to study the organization of the ions, cosolvent, and solvent molecules around proteins. In particular, peaks of densities are expected to reveal, simultaneously, the presence of compatible binding sites of different kinds on a protein. We have tested and validated the ability of HDPBL to detect pockets in proteins that bind to hydrophobic ligands, polar ligands, and charged small probes as well as to characterize the binding sites of lipids for membrane proteins.

Physics approach to the variable-mass optimal-transport problem.

Optimal transport (OT) has become a discipline by itself that offers solutions to a wide range of theoretical problems in probability and mathematics with applications in several applied fields such as imaging sciences, machine learning, and in data sciences in general. The traditional OT problem suffers from a severe limitation: its balance condition imposes that the two distributions to be compared be normalized and have the same total mass. However, it is important for many applications to be able to relax this constraint and allow for mass creation and/or destruction. This is true, for example, in all problems requiring partial matching. In this paper, we propose an approach to solving a generalized version of the OT problem, which we refer to as the discrete variable-mass optimal-transport (VMOT) problem, using techniques adapted from statistical physics. Our first contribution is to fully describe this formalism, including all the proofs of its main claims. In particular, we derive a strongly concave effective free-energy function that captures the constraints of the VMOT problem at a finite temperature. From its maximum we derive a weak distance (i.e., a divergence) between possibly unbalanced distribution functions. The temperature-dependent OT distance decreases monotonically to the standard variable-mass OT distance, providing a robust framework for temperature annealing. Our second contribution is to show that the implementation of this formalism has the same properties as the regularized OT algorithms in time complexity, making it a competitive approach to solving the VMOT problem. We illustrate applications of the framework to the problem of partial two- and three-dimensional shape-matching problems.

Numerical Simulation of Finite-Temperature Field Theory for Interacting Bosons.

We report a stable and efficient complex Langevin sampling scheme for performing approximation-free numerical simulations directly on the path-integral coherent-states field theory for an assembly of interacting bosons. We apply the method to generate the λ line of critical phase transitions associated with Bose-Einstein condensation in a model ϕ^{4} scalar field theory. The new approach enjoys near-linear scaling in the resolved (d+1) spatial and imaginary-time dimensions and should be particularly efficient for the study of dense systems at low temperature.

Multicanonical Monte Carlo ensemble growth algorithm.

We present an ensemble Monte Carlo growth method to sample the equilibrium thermodynamic properties of random chains. The method is based on the multicanonical technique of computing the density of states in the energy space. Such a quantity is temperature independent, and therefore microcanonical and canonical thermodynamic quantities, including the free energy, entropy, and thermal averages, can be obtained by reweighting with a Boltzmann factor. The algorithm we present combines two approaches: The first is the Monte Carlo ensemble growth method, where a "population" of samples in the state space is considered, as opposed to traditional sampling by long random walks, or iterative single-chain growth. The second is the flat-histogram Monte Carlo, similar to the popular Wang-Landau sampling, or to multicanonical chain-growth sampling. We discuss the performance and relative simplicity of the proposed algorithm, and we apply it to known test cases.

Statistical Physics Approach to the Optimal Transport Problem.

Originally defined for the optimal allocation of resources, optimal transport (OT) has found many theoretical and practical applications in multiple domains of science and physics. In this Letter we develop a new method for solving the discrete version of this problem using techniques derived from statistical physics. We derive a strongly concave free energy function that captures the constraints of the OT problem at a finite temperature. Its maximum defines an optimal transport plan, or registration between the two discrete probability measures that are compared, as well as a pseudodistance between those measures that satisfies the triangular inequalities. The computation of this pseudodistance is fast and numerically stable. The temperature dependent OT pseudodistance is shown to decrease monotonically with respect to the inverse of the temperature and to converge to the standard OT distance at zero temperature, providing a robust framework for temperature annealing. We illustrate applications of this framework to the problem of image comparison.

Optimal transport at finite temperature.

Optimal transport (OT) has become a discipline by itself that offers solutions to a wide range of theoretical problems in probability and mathematics. Despite its appealing theoretical properties, solving the OT problem involves the resolution of a linear program whose computational cost can quickly become prohibitive whenever the size of the problem exceeds a few hundred points. The recent introduction of entropy regularization, however, has led to the development of fast algorithms for solving an approximate OT problem. The successes of those algorithms have resulted in a popularization of the applications of OT in several applied fields such as imaging sciences and machine learning, and in data sciences in general. Problems remain, however, as to the numerical convergence of those regularized approximations towards the actual OT solution. In addition, the physical meaning of this regularization is unclear. In this paper, we propose an approach to solving the discrete OT problem using techniques adapted from statistical physics. Our first contribution is to fully describe this formalism, including all the proofs of its main claims. In particular we derive a strongly concave effective free energy function that captures the constraints of the optimal transport problem at a finite temperature. Its maximum defines a pseudo distance between the two set of weighted points that are compared, which satisfies the triangular inequalities. The temperature dependent OT pseudo distance decreases monotonically to the standard OT distance, providing a robust framework for temperature annealing. Our second contribution is to show that the implementation of this formalism has the same properties as the regularized OT algorithms in time complexity, making it a competitive approach to solving the OT problem. We illustrate applications of the framework to the problem of protein fold recognition based on sequence information only.

A Polymer Model for the Quantitative Reconstruction of Chromosome Architecture from HiC and GAM Data.

It is widely believed that the folding of the chromosome in the nucleus has a major effect on genetic expression. For example, coregulated genes in several species have been shown to colocalize in space despite being far away on the DNA sequence. In this manuscript, we present a new, to our knowledge, method to model the three-dimensional structure of the chromosome in live cells based on DNA-DNA interactions measured in high-throughput chromosome conformation capture experiments and genome architecture mapping. Our approach incorporates a polymer model and directly uses the contact probabilities measured in high-throughput chromosome conformation capture experiments and genome architecture mapping experiments rather than estimates of average distances between genomic loci. Specifically, we model the chromosome as a Gaussian polymer with harmonic interactions and extract the coupling coefficients best reproducing the experimental contact probabilities. In contrast to existing methods, we give an exact expression of the contact probabilities at thermodynamic equilibrium. The Gaussian effective model reconstructed with our method reproduces experimental contacts with high accuracy. We also show how Brownian dynamics simulations of our reconstructed Gaussian effective model can be used to study chromatin organization and possibly give some clue about its dynamics.

Dielectric constant of ionic solutions: Combined effects of correlations and excluded volume.

The dielectric constant of ionic solutions is known to reduce with increasing ionic concentrations. However, the origin of this effect has not been thoroughly explored. In this paper, we study two such possible sources: long-range Coulombic correlations and solvent excluded-volume. Correlations originate from fluctuations of the electrostatic potential beyond the mean-field Poisson-Boltzmann theory, evaluated by employing a field-theoretical loop expansion of the free energy. The solvent excluded-volume, on the other hand, stems from the finite ion size, accounted for via a lattice-gas model. We show that both correlations and excluded volume are required in order to capture the important features of the dielectric behavior. For highly polar solvents, such as water, the dielectric constant is given by the product of the solvent volume fraction and a concentration-dependent susceptibility per volume fraction. The available solvent volume decreases as a function of ionic strength due the increasing volume fraction of ions. A similar decrease occurs for the susceptibility due to the correlations between the ions and solvent, reducing the dielectric response even further. Our predictions for the dielectric constant fit well with experiments for a wide range of concentrations for different salts in different temperatures, using a single fit parameter related to the ion size.

Self consistent field theory of virus assembly.

The ground state dominance approximation (GSDA) has been extensively used to study the assembly of viral shells. In this work we employ the self-consistent field theory (SCFT) to investigate the adsorption of RNA onto positively charged spherical viral shells and examine the conditions when GSDA does not apply and SCFT has to be used to obtain a reliable solution. We find that there are two regimes in which GSDA does work. First, when the genomic RNA length is long enough compared to the capsid radius, and second, when the interaction between the genome and capsid is so strong that the genome is basically localized next to the wall. We find that for the case in which RNA is more or less distributed uniformly in the shell, regardless of the length of RNA, GSDA is not a good approximation. We observe that as the polymer-shell interaction becomes stronger, the energy gap between the ground state and first excited state increases and thus GSDA becomes a better approximation. We also present our results corresponding to the genome persistence length obtained through the tangent-tangent correlation length and show that it is zero in case of GSDA but is equal to the inverse of the energy gap when using SCFT.

Numerical Encodings of Amino Acids in Multivariate Gaussian Modeling of Protein Multiple Sequence Alignments.

Residues in proteins that are in close spatial proximity are more prone to covariate as their interactions are likely to be preserved due to structural and evolutionary constraints. If we can detect and quantify such covariation, physical contacts may then be predicted in the structure of a protein solely from the sequences that decorate it. To carry out such predictions, and following the work of others, we have implemented a multivariate Gaussian model to analyze correlation in multiple sequence alignments. We have explored and tested several numerical encodings of amino acids within this model. We have shown that 1D encodings based on amino acid biochemical and biophysical properties, as well as higher dimensional encodings computed from the principal components of experimentally derived mutation/substitution matrices, do not perform as well as a simple twenty dimensional encoding with each amino acid represented with a vector of one along its own dimension and zero elsewhere. The optimum obtained from representations based on substitution matrices is reached by using 10 to 12 principal components; the corresponding performance is less than the performance obtained with the 20-dimensional binary encoding. We highlight also the importance of the prior when constructing the multivariate Gaussian model of a multiple sequence alignment.

RNA Base Pairing Determines the Conformations of RNA Inside Spherical Viruses.

Many simple RNA viruses enclose their genetic material by a protein shell called the capsid. While the capsid structures are well characterized for most viruses, the structure of RNA inside the shells and the factors contributing to it remain poorly understood. We study the impact of base pairing on the conformations of RNA and find that it undergoes a swollen coil to globule continuous transition as a function of the strength of the pairing interaction. We also observe a first order transition and kink profile as a function of RNA length. All these transitions could explain the different RNA profiles observed inside viral shells.

Ab initio sampling of transition paths by conditioned Langevin dynamics.

We propose a novel stochastic method to generate Brownian paths conditioned to start at an initial point and end at a given final point during a fixed time tf under a given potential U(x). These paths are sampled with a probability given by the overdamped Langevin dynamics. We show that these paths can be exactly generated by a local stochastic partial differential equation. This equation cannot be solved in general but we present several approximations that are valid either in the low temperature regime or in the presence of barrier crossing. We show that this method warrants the generation of statistically independent transition paths. It is computationally very efficient. We illustrate the method first on two simple potentials, the two-dimensional Mueller potential and the Mexican hat potential, and then on the multi-dimensional problem of conformational transitions in proteins using the "Mixed Elastic Network Model" as a benchmark.