1.
Introduction
The outbreak of COVID 19 pandemic has caused an urgent need to understand a wide range of issues, both basic and applied, from epidemiology to diverse multi-scale molecular structures of its constitutive elements of the novel coronavirus (CoV-2) and its interaction with surroundings [1]–[4]. The RNA genome of CoV-2 is believed to encode only four major structural, sixteen non-structural and eleven accessory proteins. The four structural proteins including spike (S) protein, nucleocapsid (N) protein, membrane (M) protein, and the envelope (E) protein are key components for viral transmission and reproduction [4]–[9]. The CoVE protein is the smallest with 76 amino acid residues (1M2Y...76V) (for sequence see https://www.uniprot.org/uniprot/P59637) [10] but plays an essential role in replication cycle of the virus such as formation of envelope as an integral membrane protein in viral assembly, budding and pathogenesis [4]–[16].
Very recently, Kuzmin et al. [17] have examined the conformational variability of CoVE protein in presence of a membrane environment using an all-atom MD along with a coarse-grained simulation. Because of ‘prohibitively long time’ with the atomistic simulation in ‘exhaustive sampling of the conformational space’ utility of coarse-grained approach implemented in this study is very well emphasized. They find that CoVE induces curvature in membrane and ‘it has a preferable position’ in already curved membranes with ‘the C-termini localize to the convex region’. Apart from locating specific segments of the protein from the membrane sites, there is no estimate of the overall size of the protein or the underlying dynamics. For example, how the segments of protein perform its stochastic moves and how the overall dynamics (diffusive, non-diffusive, drift etc.) of the protein (i.e. that of its center of mass) sets in as its segments settles in the curved region. How its conformation and dynamics depends on temperature. It may be useful to explore the structure and dynamics of a protein monomer free from membrane as a function of temperature first. Then incorporate the environment such as membrane and other constitutive elements such as solvent to assess its effect on the structure and dynamics. We would like to focus on the folding dynamics of a free CoVE protein, a purely theoretical exercise as an idealized model system, the characteristic of which may be useful in hypothesizing its role in interpreting the laboratory data. Such a basic investigation has not been performed to our knowledge as presented here.
Sequence of CoVE protein and its specific characteristics such as formation of membrane ion channel, virulence intensity, etc. and other proteins are studied in depth [4]–[16]. Three distinct segments of the CoVE protein start with the hydrophilic amino (N) terminal with about 10 residues (1M–10G), followed by the hydrophobic transmembrane (TM) segment of about 25–27 residues (11T–37L), and ends with a relatively large hydrophilic carboxyl terminal consisting of the remaining residues (38R–76V). The structural details of how these segments interact, compete and cooperate is not fully known. The objective of this article is to investigate the effect of temperature on the structural dynamics by analyzing some of the local and global physical properties by a coarse-grained (CG) Monte Carlo simulation. The model and methods section is described in brief next followed by results and discussion and a conclusion.
2.
Model and methods
2.1. Model
The model adopted here is based on an effective and efficient method involving bond-fluctuating mechanism on a discrete lattice with ample degrees of freedom to investigate complex issues in probing long time properties in polymer systems [18] that cannot be accessed by traditional MC or MD in a continuum space. This coarse-grained approach has been further developed over the years [19]–[21] with options to fine-grain and incorporate the specificity of amino acids to model proteins. The relevant detail has been already described to examine structural dynamics of a number of proteins [19]–[21]. In coarse-grained description, we constructed a CG protein model of CoVE containing a chain of 76 amino acid residues [10] in a specific sequence on a discrete cubic lattice with ample degrees of freedom to describe essential motion and fluctuation for each residue in the protein. A residue is modeled using a single CG node with a cubic node of size (2a)3 where a refers to the lattice constant. Consecutive nodes are tethered together by fluctuating covalent bonds of which the bond length between consecutive nodes varies between 2 and √(10) in unit of lattice constant (a). A generalized Lennard-Jones potential (Uij) is used to describe the interaction between each residue and surrounding residues within a range (rc) as:
where rij is the distance between the amino acid residues at site i and j; rc = √8 and σ = 1 in units of lattice constant. The potential strength eij in phenomenological interaction (1) is based on a knowledge-based [22] residue-residue contact matrix which has been developed over many years [23]–[27] from a huge (and rapidly growing) ensemble of protein structures in Protein Data Bank (PDB).
2.2. Methods
The protein chain is placed in a random configuration initially since one can reach any microstate from any other in the vast landscape of statistical ensemble in thermodynamic equilibrium—a fundamental principle of statistical physics. Therefore, starting configuration is irrelevant. Each residue undergoes its stochastic movement with the Metropolis algorithm. The protein configurations are thus generated by movements of each residue with a probability proportional to the Boltzmann factor exp(–ΔE/T) where ΔE is the change in energy between new and old configurations of the underlying residues subject to strictly imposed excluded volume constraints. Monte Carlo time step is defined by attempts to move each residue once. Simulation parameters i.e. length in unit of lattice constant are measured in arbitrary unit and the temperature T in reduced units of the Boltzmann constant.
It is important to point out that connecting the arbitrary units with the real sample size and temperature i.e. Kelvin, and time scale (micro-second) quantitatively is difficult due to phenomenological nature of the residue-residue interaction and the coarse-grained representation. Calibration of the parameters used in simulation requires measurements of some common physical quantities in both computer simulations as well as laboratory experiments which is not feasible. However, based on the approximate size of residue and its representation and response of the physical quantities with respect to changes in temperature, one can estimate the order of magnitude.
All simulations reported here are performed on a 1503 cubic lattice for a range of temperature T = 0.010–0.030 with 100 independent samples each with 10 million time steps; different sample sizes are also used to make sure that the qualitative results are independent of finite size effects. Both local and global physical quantities such as radius of gyration (Rg), root mean square displacement (RMSD) of the center of mass, structure factor (S), contact map, etc. are investigated as a function of reduced temperature. In the following ‘temperature’ refers to ‘reduced temperature’ and should not be confused with the absolute value.
3.
Results and discussion
Snapshots of the protein (CoVE) at representative temperatures (T = 0.020, 0.024, 0.030) presented in Figure 1 illustrate the conformational changes as it denatures on raising the temperature. Obviously the protein conforms to a globular (compact) structure at low temperature (T = 0.020) and expands as it denatures (T = 0.030). The segmental globularity in the transmembrane segment of the protein is retained at least in part even at the high temperature (see below). Note that snapshots do not represent the average conformation; it is an illustrative mean that needs to be quantified. Therefore, some quantitative measures of segmental globularity are highly desirable. Average number (Nn) of residues within the range of interaction of each along the contour of the protein may provide some insight into the local restructuring as a function of temperature.
The contact profiles at representative temperatures are shown in Figure 2. At low temperatures (T = 0.012–0.020, in native phase), the coagulating centers are distributed sporadically particularly in N- and C-terminal regions along the backbone (e.g. 1M2Y, 12L13I14V (in N-terminal), 51L52V53K, 56V57Y58V59Y, 74L75L76V (in C-terminus)) with relatively large Nn. The transmembrane segment (18L–44C), however, has the largest globular cluster with high number Nn of surrounding residues. Clearly, there is no significant change in contact profiles with the temperature in native phase. However, raising the temperature (T = 0.024) further leads to a substantial decrease in segmental globularity in both N- and C-terminals with low contact number Nn while retaining some globularity in TM segment with a relatively large Nn (see Figure 2). Thus, three distinct (N- and C-terminals separated by TM) segments can be identified easily from the evolution of contact profiles with the temperature.
The global structure and dynamics of the protein emerge from the stochastic movement of each residues and their configurational stabilities. Variation of the root mean square displacement of the center of mass (Rc) with the time step (t) at a range of temperature presented in Figure 3 show, how the asymptotic dynamics depends on the temperature. One can identify the type of dynamics by estimating the power-law exponent (ν) by fitting the data with a scaling, Rc ∝ tν where ν = ½ represents diffusion. Accordingly, the protein is found to move extremely slowly at low temperatures with ultra-sub-diffusive dynamics (i.e. ν = 0.05 at T = 0.020) and resumes diffusion at higher temperatures (i.e. ν = 0.47 at T = 0.024, see Figure 3). Dynamics do affect the conformation of the protein which can be quantified by analyzing radius of gyration and the structure factor (see below).
How do we know that the protein chain has reached its steady-state equilibrium in 10 million time steps? It is important to check the convergence of physical quantities to make sure that system has reached steady-state within the range of simulation time steps in all such simulations. We have done here as in all our previous studies. For an illustration, variation of the radius of gyration with the time steps for a range of temperature is presented in Figure 4. We see that the radius of gyration has reached its steady-state equilibrium at almost all temperatures. It is not feasible to reach as good steady-state at low temperatures as is the case in most such interacting systems and even in laboratories. However, the probability of reaching the same asymptotic magnitude of Rg at all temperatures (T below 0.020) in sub-native phase (i.e. T = 0.018, 0.016 in Figure 4) is extremely low. Since it is not feasible to extend this simulation for such a long asymptotic time regime, it is worth issuing a caution to very cautious readers. The data on Rg in this sub-native phase (T = 0.010–0.020) in the following Figures should be therefore considered with caution.
The radius of gyration (Rg) of the protein is a measure of its average size. Variation of the average (over configurations from 5 million to 10 million time steps each with 100 independent samples) radius of gyration (Rg) with the temperature (T) is presented in Figure 5. In contrast to denature phase (T = 0.020–0.030) where the radius of gyration (Rg) increases on raising the temperature, it decays on increasing the temperature in the native phase (T = 0.010–0.020). Why the protein contract on raising the temperature in sub-native phase? At very low temperatures, segmental interactions dominate over thermal noise. On raising the temperature, the frozen or frustrated trapped segments begin to explore and equilibrate better and therefore result into a more compact conformation. This process continues until a certain temperature beyond which thermal noise take over the segmental interactions with the temperature and protein begins to expand. Such an opposite thermal response in native and denature phases is also observed in other membrane proteins in recent years [20]. Based on a crude calibration of the measurements of the physical quantities with the same coarse-grained approach and an all-atom MD simulation [20], we believe that the range of temperature T = 0.020–0.030 may lie in the range of T = 300–400 K roughly. In order to claim such thermal response, a universal or sub-universal characteristics for the membrane proteins, more studies on many different membrane proteins are needed. The thermal response of the nucleocapsid protein (CoVN) of COVID 19 is also non-monotonic [21] but very different from that of CoVE. We hope this results will stimulate further investigations in this direction.
It is worth pointing out in passing that the validation of this model is based on the response of protein at extreme temperature where residue-residue interaction becomes irrelevant and conformation of a protein chain behaves similar to that of an athermal polymer chain. That is the case here as in all our previous work based on this coarse-grained model. Conformational dynamics of a range of proteins with various sizes are examined using this coarse-grained approach and observed different thermal response of proteins with comparable size. For example, proteins H3.1 and H2AX [19] are about the same size (136 and 140 residues) and exhibit very different thermal response to temperature: on raising the temperature, Rg of H3.1 increases continuously before saturation while Rg of H2AX increases to a maximum value followed by decay (a non-monotonic response different from the one observed here for CoVE). Thus, the specificity of residues and its distribution (i.e. sequence) in a protein pays a critical role in their collective response. Unfortunately, we do not have a data on a protein of the COVID related protein with the size comparable to CoVE to compare here. However, based on the analysis on many proteins (some cited in this study) with this approach, the conclusion on the thermal response of CoVE is reliable. A reader can always check it out by carry out simulation independently on an idealized model like the one presented here.
Analysis of the structure factor can also provide insight into the global conformation of the protein. The structure factor (S(q)) is defined as:
where rj is the position of each residue. The magnitude of the wave vector (q) of wavelength λ is given by |q| = 2π/λ. Using a power-law scaling for the structure factor with the wave vector, S(q) ∝ q–1/γ, one can evaluate the exponent γ. Variations of the structure factor with the wavelength (comparable to radius of gyration) is presented in Figure 5 for representative temperatures.
In general, the radius of gyration (Rg) of a polymer chain shows a power-law scaling with number (N) of monomer (residues), Rg ∝ Nx with a well-defined exponent x. Since a protein is a hetero-polymer chain of residues, one may use the same scaling law here for the protein with the radius of gyration comparable to wavelength (Rg ~ λ). Then the power-law exponent x = γ and one may be able to estimate the globularity index g which is the effective dimension (D) for the spread of the protein (g ≅ D). Since, N ∝ RgD, D = 1/γ. Higher magnitude of D represents more compact structure with higher degree of globularity. Note that this scaling analysis is an estimate and should not be taken as a precise measurement. Obviously D cannot be higher than 3, however high value indicates more compact conformations. Clearly, the protein remains globular (three dimensional) at low temperatures in native phase and it denatures into a ramified random coil (D ~ 2) structures at high temperatures (see Figure 6). It is important to point out that the such a vast changes in conformation of the protein in sub-native i.e. low temperature (T = 0.010–0.020) to native-to-denature phase i.e. at higher temperature regime (T = 0.020–0.030) may pose the possibility of COVE to be an intrinsically disordered protein at least in some range of temperature and time scales. For example, the segmental interactions in the native phase may be somewhat similar to interaction of a protein with its underlying environment e.g. a specific surface [28] even in the absence of explicit external surface.
4.
Conclusions
Monte Carlo simulations are performed with a coarse-grained model of envelope membrane protein CoVE of the human coronavirus to assess its structural variability as a function of temperature. The evolution of the protein's conformation as a function of temperature is clearly seen in visualization via snapshots and animations. Segmental organizing of residues is quantified by examining local physical quantities such as average contact and mobility profiles as a function of temperature. From the variations in contact profiles as a function of temperature, it is feasible to identify three distinct segments (N-terminal and C-terminal separated by a transmembrane segment) which is consistent with the previous studies on CoVE [4],[5],[11]–[16]. The radius of gyration exhibits a non-monotonic thermal response: the protein CoVE contracts on heating in native phase and expands continuously on raising the temperature before reaching a saturation in denatured phase. Such a unique non-monotonic thermal response has been recently observed in other membrane proteins (unrelated to coronavirus). This may lead to speculate that thermal response may be a unique characteristic to identify different classes of proteins—obviously more studies are needed to confirm this observation. Thermal response of CoVE is different from another structural protein CoVN (a nucleocapsid protein) of the coronavirus [21]. Variation of the structure factor with the wave vector is studied in detail for a wide range of temperature. Scaling analysis of the structure factor provides a quantitative estimate of the overall globularity with an effective dimension (D). CoVE is found to remain globular (D ~ 3) in its native phase and unfolds continuously to a random coil (D ~ 2) conformation in denature phase on raising the temperature.
DownLoad:
