Synchronized attachment and the Darwinian evolution of coronaviruses CoV-1 and CoV-2.

CoV2019 has evolved to be much more dangerous than CoV2003. Experiments suggest that structural rearrangements dramatically enhance CoV2019 activity. We identify a new first stage of infection that precedes structural rearrangements by using biomolecular evolutionary theory to identify sequence differences enhancing viral attachment rates. We find a small cluster of mutations which show that CoV-2 has a new feature that promotes much stronger viral attachment and enhances contagiousness. The extremely dangerous dynamics of human coronavirus infection is a dramatic example of evolutionary approach of self-organized networks to criticality. It may favor a very successful vaccine. The identified mutations can be used to test the present theory experimentally.

The current coronavirus belongs to a virus family that also causes common colds. CoV-1 and CoV-2 are much more dangerous. Worse still, many numbers now show that CoV-2 (2019) spikes have evolved to be much more dangerous than CoV-1(2003) spikes. Much of this “nearly perfect” 2019 viral action [1] can be explained in terms of two well-studied cleavage–reassembly structural differences of coronal spikes [2]. The spikes contain > 1200 amino acids, and BLAST shows that 300 of these are mutated. Can only a few of these 300 mutations, far from the two cleavage sites, also be contributing to the extremely strong viral interaction of CoV-2? Maybe, but how do we identify these possibly critical distal and so far unobserved sites?

Here we will show that the exceptional viral strength of CoV-2 [1] can be used to find such sites. Our search confirms the importance of the two known cleavage sites [2]. It selects from the nearly 300 remaining mutated sites one new small mutated cluster, which can also contribute to the very strong viral activity of CoV-2 and precede cleavage. Our search is based on the concept of self-organized criticality, which goes back to van der Waals (1873) and has been discussed extensively by physicists for decades, but with only a few applications to proteins [3].

The abstract of [3] is an excellent introduction to this paper’s content. For the readers’​ convenience, here is how it begins: “A celebrated and controversial hypothesis suggests that some biological systems-parts, aspects, or groups of them-may extract important functional benefits from operating at the edge of instability, halfway between order and disorder, i.e., in the vicinity of the critical point of a phase transition. Criticality has been argued to provide biological systems with an optimal balance between robustness against perturbations and flexibility to adapt to changing conditions as well as to confer on them optimal computational capabilities, large dynamical repertoires, unparalleled sensitivity to stimuli, etc. Criticality, with its concomitant scale invariance, can be conjectured to emerge in living systems as the result of adaptive and evolutionary processes that, for reasons to be fully elucidated, select for it as a template upon which further layers of complexity can rest. This hypothesis is suggestive as it proposes that criticality could constitute a general and common organizing strategy in biology stemming from the physics of phase transitions”.

The present application uses a technique called hydropathic scaling, which has successfully explained small sequence differences in the evolutions of many protein families The thermodynamic principles of self-organized criticality [3] are universal, and have successfully produced strong signals from the evolution of many proteins [4], [5], [6], [7], [8], [9], [10], [11]. The shape of a globular protein depends strongly on its ins and outs relative to its center. The limits of these ins and outs are associated with hydrophobic and hydrophilic extrema. The shape of the protein changes during phase transitions. A given protein can function mainly through either first- or second-order phase transitions, and second order is more common [4], [5], [6], [7], [8], [9], [10], [11]. Protein unfolding from water to air is a first-order transition, well-described by the most popular hydropathic scale (KD) [12]. Second-order transitions involve much smaller structural (long-range or allosteric) transitions, such as may occur in viral cellular attachment and penetration.

Second-order phase transitions are characterized in principle by fractal (non-integral) exponents, which are difficult to measure. The discovery of 20 universal amino-acid specific fractals in the solvent accessible surface areas (SASA) of > 5000 protein segments [13] is by far the most important application of phase transition theory in the 21st century and the latter part of the 20th century. The existence of these hydrofractals proves that most proteins function near critical points of second-order phase transitions, with human proteins often closest to their critical points [5] with level sets [10]. As the reader has guessed, the fractal hydropathic SASA scale (aa) (MZ) is the tool that enables us here to find the new distal site cluster which probably contributes to the very strong viral activity of CoV-2.

The standard tool for comparing two or more proteins is BLAST, which compares different protein sequences site-by-site. While there are many refinements of BLAST, all utilize site-by-site comparisons. It is not enough to profile (aa) site by site for two proteins, as this produces hundreds of oscillations, with no clear differences between the two proteins. Instead we plot (aa,W), where W is a rectangular box of length W = 35, centered on aa, over which (aa) is averaged. This W value is chosen to maximize the hydropathic shape differences between CoV-1 and CoV-2, as measured by their variance ratio, which is dominated by in/out extrema [10], [11]. Smoothing the profile over W enhances dynamical resolution by reducing the number of in/out extrema to the level which corresponds to critical domain motion differences. Choosing the best value of W to enhance resolution of dynamical evolutionary differences is similar to adjusting the focal plane of a microscope.

Our main result is shown in Fig. 1, which displays the 400–800 central region of the spikes. This region is part of the N-terminal domain (S1) responsible for cellular attachment [14]. Fig. 1 shows the (aa,W) hydropathic profiles of for CoV-1 and CoV-2, using the (aa) values of MZ [4]. While the two (aa,35) profiles are similar, there are important differences in their extrema. The maxima are most hydrophobic and are located in the globular interior, while the minima are most hydrophilic and are located on the globular surface. The two cleavage sites S1/S2 and  [2] of CoV-1 have moved lower (hydrophilically, further outside) in CoV-2 (Fig. 1), consistent with the very accurate MZ scale. When a cleavage segment is further outside, there is more space for cleavage and reassembly, which will occur more rapidly. The insertion PRRA in CoV-2 was identified [2] with BLAST (W = 1) as unique to CoV-2, but with BLAST alone one cannot show that this change has made CoV-2 more dangerous.

Fig. 1

The hydropathic profiles of CoV-1 and CoV-2 reveal a hidden symmetry when plotted using the MZ scale (second-order phase transitions). The three hydrophilic minima of CoV-2, labeled 1–3, are nearly equal at 140, whereas the similar minima of CoV-1 range from 131 to 147. Note especially the very deep minimum of CoV-1 at 559. The new sequence PRRA in CoV-2, inserted at 681 in CoV-1 (the S1/S2 cleavage interface) [2], has an average MZ hydropathicity 108.5. This lowers (aa,35) to 140.9 at the 3 minimum, aligning it with 140 minima 1 and 2. The three minima span 250 amino acids sites, which makes their water-driven synchronization for CoV-2, but not CoV-1, outside the range of most simulation or modeling methods.

A characteristic feature of second order phase transitions is that they break what physicists often call a “hidden” symmetry, here the key symmetry of extrema which lead hydropathic fronts [15]. Proteins are typically composed of domains of 100 amino acids, and these domains can rotate from the resting state to the functional state, and then return reversibly to the resting state. This phase transition requires a balance between stability and flexibility [16], [17]. To maximize the reversibility and extend protein life one can synchronize this domain motion by leveling the pivotal hydrophobic (inside) extrema forming level sets [10], [11], [18], [19], [20], [21]. Viruses must act rapidly before being destroyed by antibodies, and they could do this through synchronized motion differently, by leveling their hydrophilic (outside) extrema. As shown in Fig. 1, such a leveling of minima 1–3 occurs in CoV-2, while it is absent from CoV-1. The change in minimum 2 is especially striking: it is caused by a cluster of four key mutations from CoV-1 to CoV-2. These are (CoV-1 site numbering from Uniprot P59594): 546Gln to Leu; 556 and 561Ser to Ala; and 568Ser to Leu. The differences associated with each of these key mutations are hydropathically large (50-100 in the MZ scale [4]; all 20 amino acids span a range from most hydrophilic to most hydrophobic of 170).

How important is the hydrofractal MZ scale? In Fig. 2 the MZ panoramic profiles of CoV-1 and CoV-2 are compared across the entire 1255 amino acid spike. We see immediately a strong hydrophobic peak above 1200, which is responsible for anchoring the spike to the central cushion; it is almost unchanged from CoV-1 to CoV-2. The central hydrophilic level set, absent from CoV-1 and present in CoV-2, is our main result (Fig. 1). Similar plots with the first order KD scale, shown in Fig. 3, Fig. 4, do not exhibit synchronization by a hydrophilic level set (tilted arrow). One can see a shorter and more nearly level hydrophobic set, but only in CoV-1 and not in CoV-2.

Fig. 2

With the MZ scale (second-order phase transitions) panoramic spike profiles of CoV-1 and CoV-2 reveal a set of hydrophilic level extrema in CoV-2, but not in CoV-1. The new level set was identified with a cluster of four single mutations in 546–568. Addition of a new mutation, D796H, brings the hydrophilic minimum S2’at 796, into excellent alignment with minima 1–3 in Fig. 1.

Fig. 3

Even on the panoramic picture, it is clear that with the KD scale (first-order phase transitions) the three minima that were aligned in CoV-2 are not level. This is shown in more detail in Fig. 4.

Fig. 4

The enlarged central region 400–800 of the hydropathic profile using the KD scale (first-order phase transitions) is similar to what was shown in Fig. 1 using the MZ profile (second-order phase transitions), but there is an important difference. With the MZ scale all three minima 1–3 are aligned in CoV-2, but here only the 2 and 3 are aligned near 150. The 1 minimum has improved from 135 in Cov-1 to 142 in CoV-2, but this is not alignment, compared to minima 1–3 with the MZ scale (Fig. 1).

There is an excellent review of the principles of self-organized criticality in living matter [3]. It mentions examples of synchronization (also used here for the hydrodynamics of domain motion of molecular level sets [15], [18]). Another example is melting phase transitions in DNA strands. These can be compared to relaxation of homogeneous glass alloys, whose composition has been adjusted to bring the glass network to a critical point. There the glass (commercial Gorilla glass) is strong and flexible – not brittle [21]. Another broad and deep review includes fractals and their many implications for quantifying protein dynamic criticality [22]. In the genomic age abstract tools have many applications to analyzing the evolution of protein functions.

Some readers may be interested in the connections between hydropathic scaling theory of proteins and the more general synchronization of complex networks. The analogy is closest for networks that are scale-free because they are near a second-order phase transition. Social networks contain many such examples, usually described by one or two fractals [23]. Generally it is found that synchronizability is robust against random removal of nodes, but it is fragile to specific removal of the most highly connected nodes. In proteins the most highly connected nodes are hydrophobic extrema (deeply interior), where the number of van der Waals contacts is largest. In coronavirus spikes, it is the hydrophilic extrema that are most exposed to water and most likely to attach to protein targets. Thus the symmetry of extrema is broken between hydrophobic (proteins) and hydrophilic (CoV-2). One can also regard Darwinian selection as a kind of percolation [24], [25]. Further extensions to multiphase hydrodynamic level sets require a special theory [15].

Our primary focus here has been on enhanced attachment caused by the synchronization of CoV-2 hydrophilic minima, but the mutations of CoV-1 are of some interest (P59594). Only a fragment (318–510) of CoV-1 accounts for its binding, and smaller fragments (327–510) and 318–490) do not bind [26]. The successful fragment is centered on the hydrophilic minimum 1 in Fig. 1, and it spans the region on both sides up to the two nearest hydrophobic maxima in Fig. 2. It was found that hydroneutralizing (alanining, Ala) six extremely hydrophobic (Cys) or hydrophilic (Asp, Glu) sites caused five of them to lose their binding: 452, 454, 463, 467, 474. These grouped sites are all close to the hydrophilic minimum 1 in Fig. 1 (MZ scale) and Fig. 3 (KD scale). Fragmental binding is a feature common to both scales, and thus not sensitive to whether the phase transition is first- or second-order. Instead, by synchronizing minima 2 and 3 with minimum 1, CoV-2 is able to greatly increase its attachment to ACE2 (angiotensin-converting enzyme 2) and probably other receptors as well; synchronization is specifically sensitive to allosteric interactions and is definitely a second-order effect.

The answer to the question posed in [1] is that, compared to CoV-1, CoV-2 has moved closer to its functional critical point [3]. CoV-2 has an added synchronized feature in its central attachment region that contributes to its extreme contagiousness in an asymptomatic phase [27]. The attachment of the S1 region precedes cleavage [2], and it is likely that both mechanisms enhance the spread of CoV-2. An early estimate of the basic contagious reproductive number (R0) of Coronavirus 2019 was 2.2–2.7, but this was later revised to 5.7 [28]. By comparison, the numbers for H1N1 flu (1918, 1.8; 2009, 1.5) were much smaller [29] The uncontrolled R0 for CoV-1 was 3.6, and this dropped to 0.7 after the rapid implementation of control measures [30].

The new feature explains the unexpected yet often asymptomatic behavior of early stages of CoV-2 infection, when the new feature is acting alone. It should be possible to desynchronize this added attachment feature and reduce R0 by attaching a small molecule in the region 530–575. In any case, the “hidden” symmetry of the CoV-2 sequence should be of interest to both biologists and physicists. The effects of the four critical mutations on CoV contagiousness could be tested on mice [31]. If the predictions of their importance should prove to be correct, then the general principle of Darwinian evolution of proteins would be confirmed.

The level hydrophilic extrema of the spikes of CoV-2 may be unique among viruses. These level extrema bring CoV-2 contagion very close to a critical point. There very small differences in many genes can cause large fluctuations in infection levels, not only between individuals, but also even in neighboring countries [32]. Because CoV-2 is very close to a critical point, its symptoms can vary widely between individuals. One of the symptoms is loss of smell and taste. A large statistical survey reported that this was often an “all-or-nothing” effect [33], characteristic of a critical point. Because of interactions with other proteins, the individual reaction can bifurcate on either side of the critical point.

The quasi-linear spike morphology, with its large surface/volume ratio, is ideal for strong protein–water interactions, explaining [1]. The Oxford vaccine is based on the CoV-2 spike protein hydrophobically anchored to an adenovirus vaccine vector [34], [35], [36]. It seems likely that this vaccine will be more effective for CoV-2 than for CoV-1 because the synchronized attachment mechanism is stronger for CoV-2 and thus more easily disrupted by induced antibodies. Thus vaccines based on the entire spike are expected to be more effective than vaccines based on only part of the spike [37]. Another way of thinking about this is by analogy with the van der Waals equation. There when one is at or near the P–V critical point on the isotherm, a small change in pressure (i.e., the vaccine) produces a large change in volume and drives the system away from the critical point.

The most dangerous type of flu virus is N2, and it does not exhibit level sets [38]. Flu viruses mutate significantly annually, while CoV-2 has shown only a few mutations [39]. The spike stability against mutations is consistent with maximized attachment. There are many examples of scale-free networks. Readers unfamiliar with the general properties of scale-free networks, such as (often only 1 or 2) fractals, preferential attachment, and synchronization, may find early reviews interesting [40], [41], [42]. The sequences used here for CoV-1 and CoV-2 are the ones used in [2], that is, AAP13441.1 and YP_009724390.1. Many other CoV sequences are aligned (W = 1) in their Supp. Table [2]. Similar W = 1 sequence alignments have been reported by others [43], [44], but the genetic origin of CoV-2 remains mysterious [1], [45].

Historically it has long been the case that the amount of information that could be obtained from hydropathic scales was limited, because the hydrofractals describing second-order phase transitions were not known. Now that we are in the genomic age, with a very large sequence data base available for many proteins and many species, the discovery of these 20 hydrofractals [5], [13] opens a new biophysics field of accurate thermodynamic analysis of small but medically important evolutionary differences.

The methods used here on CoV-1 and CoV-2 do not require elaborate calculation to obtain new results. The hydropathic scaling methods are easily implemented on a spread sheet, but are still not routine, and vary between proteins. The modular nature of protein domains responsible for protein assembly is well known [45], [46]. Here we maximized small and otherwise mysterious evolutionary differences by using hydropathic scaling to optimize W at 35. This enables us to focus on domain-scale motion responsible for stronger attachment in CoV-2. Evidence for large-scale motion has been found in the evolution of many proteins [4], [5], [6], [7], [9], including molecular motors [11], while it is not accessible to any of the traditional W = 1 BLAST-based alignment methods [2], [47], [48], [49]. From detailed genomic studies it has been suggested that “a more comprehensive theory of molecular evolution must be sought” [50]. Readers wishing to explore these novel scaling methods can obtain copies of sample EXCEL spread sheets from the author. All the details obtained by analyzing the (aa,W) matrices are straightforward, and should be easily reproduced with or without our sample spreadsheets by any physical scientist.

Many physicists have long suspected that proteins must function near a critical point very near thermodynamic equilibrium [3]. However, to implement this idea in practice with protein sequences, one needs to compress the existing vast array of three-dimensional structural data into a one-dimensional form. If you are one in a million, namely among the small community of co-authors of the 500+ references of [3], what comes next may not surprise you. As discussed above, this near-miracle compression has been achieved in the 21st century by two Brazilian physicists [13], located far from the main stream. They described the discovery of 20 (not just one) universal amino-acid specific hydrofractals in the solvent accessible surface areas (SASA) of proteins. This amazing 21st century discovery was achievable only for proteins. It has made possible the quantitative analysis of Darwinian evolution of many different proteins [4], [5], [6], [7], [8], [9], [10], [11], which brings us to the present application.

The last steps of evolution are usually the most informative. Recent molecular dynamics (MD) simulations [51] suggest “that the formation of the intermediate state is a prerequisite for the host cell recognition”, which is similar to our synchronized state that precedes structural rearrangements. How close is CoV19 to its critical point? Very close, because there are very few mutations. These are discussed at length in analyzing in vitro data on the ACE2 affinity of a CoV19 receptor fragment supported by yeast [21]. Comparing our Fig. 1 to the CoV19 mutations in Fig. 7B of [52], our thermodynamic analysis suggests that their most dangerous mutation is S477N. We have seen that CoV19 is much more contagious than CoV03, because the three minima in Fig. 1 are nearly level at 137.7, 138.4, and 140.1. The S477N mutation changes the first minimum to 138.0 (slightly more level). This change is 25 times smaller than the evolutionary leveling shift of min. 2 of CoV1 from 131.3 to CoV2 138.4. In general, single amino acid mutations are unlikely to make CoV19 both more flexible and more stable, and are unlikely to proliferate [53]. A current topic is the significance of the mutation D614G [54]. A physical explanation for the small success of this single mutation [4] is that the replacement of hydrophilic D with hydroneutral G should provide additional stability for CoV-1/2 in aerosol drops (see Table 1).

Table 1

The hydropathic first-(KD) and second-(MZ) order amino acid parameters. Both sets have been rescaled linearly in order to fit into the convenient range (50, 250) with almost identical average values and ranges.

	A	C	D	E	F	G	H	I
MZ	157	246	87	94	218	156	152	222
KD	200.5	214.3	96.1	96.1	220.2	157.1	102.0	253.7
	K	L	M	N	P	Q	R	S
MZ	69	197	221	113	121	105	78	100
KD	88.2	239.9	202.4	96.1	133.5	96.1	76.4	149.2
	T	V	W	Y
MZ	135	238	174	222
KD	151.2	247.7	147.3	139.4

At least two very sophisticated explanations of how Darwinian evolution involving multiple mutations could have led from CoV03 to CoV19 are now available [55], [56]. Note that each of our four new key mutations improves synchronization. We expect that vaccines based on slightly modified spikes should be very successful because any disturbance in the spike population will probably drive the virus far away from the critical point, rendering it harmless.

All the sequence data used in these calculations are available online, for instance, at Uniprot., as spike P59594 and P0DTC2. The very large success of the Oxford vaccine [34], [35] was predicted in the abstract of an early version of this paper in arXiv2008.12168 (Aug. 28, 2020). It was first reported in media news Nov. 20, 2020, and in [57]. The greater success with half an initial dose can be understood as the result of training immunogenicity. On the initial half, the immunogenic system chooses to target only the more dangerous spike virus, and ignores the adenovirus component. With a full first dose, it can do both, making the second dose less effective.

Postscript. In a later note this thermodynamic model successfully explains the enhanced spreading of several later strains of CoV19 (the UK, SA and Calif. strains). This note is available at bioRxiv 2021.02.16.431437 under the title “Phase transitions may explain why SARS-CoV-2 spreads so fast and why new variants are spreading faster”. Readers who are unfamiliar with fractals may find Gene Stanley’s books helpful [58].

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.