Global Genetic Response in a Cancer Cell: Self-Organized Coherent Expression Dynamics

Understanding the basic mechanism of the spatio-temporal self-control of genome-wide gene expression engaged with the complex epigenetic molecular assembly is one of major challenges in current biological science. In this study, the genome-wide dynamical profile of gene expression was analyzed for MCF-7 breast cancer cells induced by two distinct ErbB receptor ligands: epidermal growth factor (EGF) and heregulin (HRG), which drive cell proliferation and differentiation, respectively. We focused our attention to elucidate how global genetic responses emerge and to decipher what is an underlying principle for dynamic self-control of genome-wide gene expression. The whole mRNA expression was classified into about a hundred groups according to the root mean square fluctuation (rmsf). These expression groups showed characteristic time-dependent correlations, indicating the existence of collective behaviors on the ensemble of genes with respect to mRNA expression and also to temporal changes in expression. All-or-none responses were observed for HRG and EGF (biphasic statistics) at around 10–20 min. The emergence of time-dependent collective behaviors of expression occurred through bifurcation of a coherent expression state (CES). In the ensemble of mRNA expression, the self-organized CESs reveals distinct characteristic expression domains for biphasic statistics, which exhibits notably the presence of criticality in the expression profile as a route for genomic transition. In time-dependent changes in the expression domains, the dynamics of CES reveals that the temporal development of the characteristic domains is characterized as autonomous bistable switch, which exhibits dynamic criticality (the temporal development of criticality) in the genome-wide coherent expression dynamics. It is expected that elucidation of the biophysical origin for such critical behavior sheds light on the underlying mechanism of the control of whole genome.

contains biological and experimental noise [1][2][3]. In our prior works, for LPS innate responses of macrophages and naive CD4+ T cell differentiation with a few experimental time points, groups are made according to changes in expression with respect to initial (at time t 0 ) wild-type genome, whereas for HL-60 cell (human promyelocytic leukemia cell) differentiation with enough experimental time points like MCF-7 cell in this report, grouping is based on an expression variation, where average value of each gene expression in time is subtracted from gene expression, and genes are sorted by the standard deviation of expression in time. In any case, groups are based on pure expression level basis with no reference to the physiological role of the genes.
We found that the average value ( ) of the i th group between different time points or different genotypes ( ) exhibited nonlinear correlation (curve) between and as the group size n (the number of probes in each group) increased: ! = ! [4][5][6].
Interestingly, the emergence of global trend (manifold) from scattered points to (asymptotic) time-dependent correlation was noticed in the distribution of the whole groups when the grouping size is about n = 30-50. Especially, as the group size (n) increased, the ordinary square root distance (Euclidian distance) from the correlation decreases based on the inverse square root law (α/ n: α is a constant coefficient); the average value of a group approaches to an asymptotic value, which indicates a distinct xi y i xi y i asymptotic frequency distribution for temporal expression of genes within the group [5].
We progressively approach to the 'true' (i.e., asymptotic) distribution relative to the reference population, in which each gene expression in a group is stochastically fluctuating around the asymptotic average value. Such global trend on the whole gene expression imply the presence of a hidden governing principle, which rationalizes groups of genes to guide all the genes allowing for such a convergent behavior when the entire Notably, as the biological and biophysical significance, i) for the innate immune grouped genes that can be linearly superposed to decipher the global gene regulatory differential control of the transcriptional and mRNA decay machineries between wildtype and mutant genomes [4]; ii) for T-cell differentiation, we observed global switching behavior of gene expression between help T1 and T2 differentiations [5]; iii) for HL-60 cell differentiation, gene expression dynamics exhibits coherent oscillating behaviors during differentiation process [6]. Furthermore, we demonstrated, in a statistically rigorous manner, the existence of specific gene ensembles (collectively named "genome vehicle") responsible for the formation of a neutrophil attractor; the collective motion of lowly and moderately variable genes within a genome vehicle guides the whole gene expression toward the neutrophil attractor.
These results revealed the emergent hidden global regulations in genome-wide expression, which drives orchestrated diverse regulations of biological processes. Thus, there must be organized principles to explain such highly coordinated averaging responses of biological processes, which might imply the existence of simple rules that control and regulate complex biological processes. Next, we address a basic underlying mechanism to derive the global genetic response.

Dynamic Criticality
We demonstrated that a change from a unimodal to a bimodal frequency to DEAB (Figures 4 and 6) anticipates the occurrence of criticality, which encompasses genome-wide coherent expression. In genetic criticality, the expression possesses a characteristic of an order parameter such as density, and the parameter rmsf, which describes the standard deviation of temporal fluctuation in mRNA expression, depicts the degree of symmetry breaking such as temperature to categorize critical phenomena as above-criticality (supercritical), near-criticality, and below-criticality (subcritical) ( Figure   8); a local minimum of energy potential corresponds to a hilltop of a CES and a peak of the frequency distribution of the expression (Figures 6A-B).
We considered gene expression near criticality, where the potential energy (or energy like function) could be expanded by the lowest order of the order parameter At the symmetric-symmetry breaking, the whole system can move between two new states with an equal transition probability to produce two equivalent Boltzmann distributions. To make asymmetric-symmetry breaking ( ≠ 0), the linear term aη such as a role of electric field contributes to left (a > 0) or right side (a < 0) asymmetric local minima ( Figure S1 C).
To understand a non-equilibrium process in genetic activity, we point to an example that the bimodal frequency distribution of the expression in the static domain ( Figure 6) does not show two independent Boltzmann distributions as expected in a equilibrium system, where two peaks should correspond to peaks of a bimodal profile, but rather the mixing of the expression, the broken Boltzmann distribution, appears (Figures 6   and 8). We expect that the mixing expression dynamics stems from nonlinear interactions among gene expression in a non-equilibrium system. In addition, the two high-expression states in the unimodal frequency distribution in the dynamic domain further support a nonequilibrium genetic activity ( Figure 6A).

Figure S1: Criticality and Symmetry Breaking
Dynamic energy-potential, the time-dependent change in the energy profile is expected by the energy flow (due to a non-equilibrium process) in coherent expression dynamics ( Figure 9). Self-organized CESs exhibit criticality in DEAB of the expression at Combining the pictures of Figures 8 and 9, the time-dependent of criticality for the three domains in the HRG response with the symmetry argument is summarized: iii) In below-criticality: the static domains, a local energy minimum bifurcates from a high energy minimum through symmetry breaking (a negative quadratic and a positive quartic term) of the energy potential, which results in the creation of a potential barrier that separates coherent expression into high and low. Furthermore, at the right side a deeper well (valley) of the genetic energy profile indicated a key effect of a linear external term (a > 0). During the period 10-20min, the energy profile for DEAB of the expression (rmsf < 0.21 with a fixed change in the expression) should be temporally invariant (due to overlapping of thee DEABs at 10, 15, and 20 min: Figure 1A). On the other hand, DEAB of the expression change (i.e., a fixed expression) for HRG (rmsf < 0.21) changed from up-to down-regulation ( Figure 1B) as a part of the pendulum oscillation of LES1 between 10-15min and 15-20min (Figure 9), which suggests the occurrence of a single-well shift of the energy potential during the period 10-15min ( Figure 9).
Upon elucidation of coherent interactions of genomic DNA folding/unfolding dynamics, a basic mechanism of the formation of ABS by a pair of CESs could be deciphered, and the biophysical significance of the opposite pendulum oscillations of CES between the HRG dynamic and static domains (Figure 9) might be revealed.