The Nonrandom and Geometric Role of Retroviral DNA in Genome

1. Introduction

The DNA structure in the human genome is the outcome of evolutionary processes governed by the synchronization of biological and environmental components driving the system near or far from equilibrium (Bundschuh, 2006; Wong, 2020; Basu, 2021). The redundancy of information stored in the DNA structure is reflected in its sequence segments following q-Gaussian distributions (Pavlos, 2015; Karakatsanis, 2018; 2021, Correia, 2022, Tsallis, 2022). Changes in the underlying information landscape are detectable via Machine Learning (ML) models with supervised or unsupervised learning (Manogaran et al., 2018; Washburn et al., 2019; Varma et al., 2019; Frey et al., 2019; Libbrecht and Noble, 2015; Karakatsanis et al., 2021).

The embodiment of retroviral DNA sequences into the human genome represents a crucial and significant research area. In this direction, the authors of Gonzalez-Cao et al. (2016) explored the correlation between Human Endogenous Retroviruses (ERVs) and human cancer types, including melanoma, breast cancer, germ cell tumors, renal cancer, and ovarian cancer, which express HERV proteins, primarily HERV-K (HML6) and HERV-K (HML2). In ref Alldredge et al. (Alldredge, 2023) investigated the expression of ERVs in cervical cancers, using publicly available RNA-seq data from 63 cervical cancer patients. ERV expression signatures in tumor biopsies may therefore be useful to help identify patients at greater risk of recurrence. Calero-Layana et al. (Calero-Layana, 2022) investigated the evolutionary history of endogenous retroviruses associated with five human genes (INPP5B, DET1, PSMA1, USH2A, and MACROD2), which are located within intron sections. They proposed that these elements play a relevant role in gene expression regulation for tumorigenesis control. Chuong (Chuong, 2018) found that ERVs implicate an extensive yet understudied role for retroviruses in shaping the evolution of placental gene regulatory networks. Kojima et al. (Kojima, 2021) applied machine learning to analyze endogenous RNA virus sequence signatures, enabling the identification of viruses in the human genome that are either undetected or extinct. Jansz and Faulkner (Jansz and Faulkner, 2021) explored the ERV co-option in development and innate immunity, the aberrant contribution of ERVs to tumorigenesis, and the wider biomedical potential of therapies directed at ERVs. Vargiu et al. (Vargiu, 2016) conducted a systematic study on the human ERV classification and characterization. Moreover, Russ and Iordanskiy studied the implications of HERVs in relation to innate immunity and their association with various pathological disease states. In a very recent study Ivancevic et al. (Ivancevic, 2024) showed that endogenous retroviruses mediate transcriptional rewiring in response to oncogenic signaling in colorectal cancer.

All previous studies in retroviral DNA means that the spatial organization which produced the spatial information of DNA can be characterized by complex character and is affected by the sub-set of retroviral DNA. The observed spatial organization indicates structured, nonrandom distributions of ERV elements within genomic DNA. The complexity metrics used reveal patterns exhibiting fractal and multifractal characteristics. These features align with the concept of strange attractors from nonlinear dynamics, suggesting complex, self-organizing genomic arrangements. This opens the scientific field to study the contribution of the sub-set retroviral DNA in the human genome.

In this study we measured the spatial organization and the non-extensive statistical characteristics in the retroviral DNA Human Endogenous Retro Viruses (HERV) database of the human genome studying the distribution of the position distance and the length of the ERV elements. Spatial organization refers explicitly to the positional arrangements and distribution patterns of ERV elements across chromosomes. The HERV database is a comprehensive, curated repository containing detailed information about retroviral elements integrated into the human genome, including their genomic location, structural composition, and sequence characteristics.

Following the ref Pavlos et all (Pavlos, 2015), we consider a DNA walker, showing that the step size distribution encodes all relevant information. We introduced a complexity factor (COFA) (Karakatsanis, 2021), incorporated the Hurst exponent (Weron, 2002) and the q-entropic index, which underlies nonextensive statistical mechanics (Tsallis, 1988; 2004; 2009), and explored, using machine learning techniques, whether it can summarize common patterns across chromosomes.

The principal finding of this study is the identification of structured complexity and significant nonrandom patterns in the positions and lengths of ERV elements in the genome. This structured complexity is measurable using complexity metrics (Hurst exponent, Tsallis q-index, and COFA index). Such complexity potentially impacts genomic stability, functional adaptability, and evolutionary processes, highlighting the biological significance of ERV distribution patterns.

2. Theoretical Framework

2.1. Nonrandom DNA Walker

Let $a_i=\mathrm{0,1},\dots ,N\in \mathrm{A}\complement \mathrm{R}$^d be a collection of ordered coordinate centers representing molecular buildings blocks of a DNA molecule. d is the dimension of the embedding space and may be adjusted according to the observer’s scopes. Ordering means that the sequence i=0,1,...,N is the outcome of some optimization process running across an evolutionary time scale T.

The interevent vector, i.e., the vector connecting any pair of ordered A coordinates reads (Pavlos, 2015):

$$\overrightarrow{r_i}=a_{i+1}-a_i$$

(1)

$‖{\overrightarrow{r}}_i=r_i‖$, where $‖.‖:R^n\to R$ is some arbitrary norm operation, is intuitively understood as the size of the ’step’ that a DNA walker undertakes in the abstract configuration DNA space.

Let now $p(r,T)$ denote the probability distribution function of the step size, r, at a given T. Note that for clarity and ease of reading, the subscript i is omitted in this instance and thereafter.

To define a dynamic equation for p, we base our approach on the following two assumptions.

First, we assume the influence of a linear constraint,

$$F\left(r\right)=-\gamma r,\gamma >0,$$

(2)

biasing the DNA walker towards shorter ’jumps’. Shorter steps restrict the DNA walker’s ability to traverse the genomic landscape swiftly, thereby limiting its exploration of distant configurations. While this reduced mobility may slow adaptation to dynamic environments, it simultaneously promotes the stability of existing, well-adapted functions by curbing excessive exploration and minimizing the risk of harmful changes. Consequently, γ reflects a delicate balance between evolutionary adaptability and the maintenance of functional robustness in the DNA system.

Second, we assume an isotropic diffusive force that ’disorients’ the DNA walker. Diffusion acts here as a moderator of the DNA walker’s motion by counteracting the bias introduced by (2). This ensures a nuanced exploration of DNA configuration space, where both localized and distant regions are probabilistically accessible, preventing the walker from becoming trapped in specific regions of the configuration space.

Given the pair of assumptions outlined above, the dynamics of p across T can be described by the following Fokker-Planck equation:

$${\displaystyle \frac{\partial p(r,T)}{\partial T}}=-{\displaystyle \frac{\partial \left[F\left(r\right)p\left(r,T\right)\right]}{\partial r}}+D{\displaystyle \frac{{\partial }^2{\left[p\left(r,T\right)\right]}^{2-q}}{{\partial r}^2}},p\left(r_0,0\right)=\delta (r_0)$$

(3)

where D>0 is the diffusion coefficient, i.e., the ’disorientation strength’. Note that (3) describes a Uhlenbeck-Ornstein process initiated for δ(r₀) (Tsallis & Bukman, 1995, Tsallis, 2009). q is the so-called q-entropic index introduced by C. Tsallis (Tsallis, 1998).

Its biophysical role is discussed below.

Solving 3 for p, yields (Tsallis & Bukman, 1995, Tsallis, 2009):

$$p\left(r,T\right)={\displaystyle \frac{{{\{1-\left(1-q\right)\beta [r-r_M]}^2\}}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$1-q$}\right.}}{Z}},$$

(4)

with,

• $\beta ≔\beta \left(T\right)\ne 0$,
• $Z≔Z\left(T\right)>0$
• $r_M≔r_M\left(T\right)=r_M,0^{\exp \left(-\gamma T\right)},r_{M,0}=r_M\left(0\right)\in R$

1/β is the effective temperature of the DNA walk. Note that negative effective temperatures are also possible. This is rationalized on the basis that the DNA configuration space may be constrained within specific boundaries due to evolutionary pressures and/or thermodynamic limitations associated with DNA crumbling.

• Z is the strictly positive normalization function (partition sum) guaranteeing that,
• $\int drp\left(r,T\right)=1,\forall T\ge 0.$
• β and A are co-determined by the following equation:

$${\displaystyle \frac{{\beta }_0}{\beta \left({\rm T}\right)}}={\left({\displaystyle \frac{{\rm Z}\left({\rm T}\right)}{{{\rm Z}}_0}}\right)}^2={\left[\left(1-{\displaystyle \frac{1}{k}}\right)exp\left({\displaystyle \raisebox{1ex}{$-{\rm T}$}\!\left/ \!\raisebox{-1ex}{$\tau $}\right.}\right)+{\displaystyle \frac{1}{k}}\right]}^{{\displaystyle \frac{2}{3-q}}},{\beta }_0=\beta \left(0\right),{{\rm Z}}_0={\rm Z}\left(0\right)$$

(5)

with

• $k≔{\displaystyle \frac{1}{C_{0}2\left(2-q\right)D}}{\displaystyle \frac{\gamma }{D}}>0,C_0={\beta }_0{Z_0}^{q-1}\ne 0,q\ne 2$
• $\tau ={\displaystyle \frac{1}{\gamma (3-q)}}>0,q<3.$

To guarantee that $k>0$, the sign of C₀ is adjusted by requiring that β₀ < 0 and β₀ > 0 for 2 < q < 3 and 1<q <2, respectively.

For a sufficiently long evolutionary trajectory (i.e., for $T\to \infty$), (4) can be written as:

$$p_{\infty }=p\left(r,T\to \infty \right)={\displaystyle \frac{{\left\{1-(1-q){\beta }_{\infty }{r}^2\right\}}^{\frac{1}{1-q}}}{Z_{\infty }}},$$

(6)

with

• ${\beta }_{\mathrm{\infty }}=\beta \left(T\to \mathrm{\infty }\right)={\displaystyle \frac{{\beta }_0}{\lambda }}\ne 0,\lambda ≔{\left({\displaystyle \raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$k$}\right.}\right)}^{{\displaystyle \raisebox{1ex}{$2$}\!\left/ \!\raisebox{-1ex}{$3-q$}\right.}}>0$
• $Z_{\mathrm{\infty }}=Z\left(T\to \mathrm{\infty }\right)=Z_0\sqrt{{\displaystyle \frac{{\beta }_0}{{\beta }_{\mathrm{\infty }}}}}$.

We notice that the amplitude of the absolute effective temperature is determined is proportional to ${\displaystyle \frac{\gamma }{D}}$. This implies that ’hot’ regions of the DNA configuration space are those having a high likelihood to be revisited, since a large value of γ (and/or a small value of D) discourages configuration exploration.

Also, (6) can be obtained as the solution of an optimization process where the following entropy functional:

$$S_q\left[p_{\infty }\right]={\displaystyle \frac{1-\int dr{\left[p_{\infty }\right]}^q}{q-1}}$$

(7)

is maximized (Tsallis, 2009). In the light of (7), it is clarified that q is indeed an entropic index.

• Generally, q can be understood as the degree of nonrandomness of the ordering process under scrutiny since for q → 1 the Boltzmann-Gibbs entropy functional is obtained:

  $$S_{q\to 1}\left[p_{\infty }\right]=-\int dr{p}_{\infty ,q\to 1}ln\left(p_{\infty ,q\to 1}\right)$$

  (8)

with

• $p_{\mathrm{\infty },q\to 1}=exp\left(-{\beta }_{\mathrm{\infty },q\to 1}{r}^2\right)/Z_{\mathrm{\infty },q\to 1}$,
• ${\beta }_{\mathrm{\infty },q\to 1}={\displaystyle \frac{\gamma }{2D}}$,
• $Z_{\mathrm{\infty },q\to 1}=Z_0\sqrt{{\displaystyle \frac{{\beta }_0}{{\beta }_{\mathrm{\infty },q\to 1}}}}$

indicating that as tend to be randomly ordered or equivalently that the DNA walker explores the configuration space in a nearly random fashion.

2.2. Data Acquisition

2.2.1. HERV Database

The Human Endogenous Retro Viruses database (HERVd) (see https://herv.img.cas.cz/) provides complex information on and analysis of retroviral elements found in the human genome. It can be used for searches of individual ERV families, identification of ERV parts, graphical output of ERV structures, comparison of ERVs and identification of retrovirus integration sites (Paces 2002; 2004).

2.2.2. Construction of Data

For each chromosome, we created some sequences from ERV elements, such as the ‘position distance’, and the ‘length’. The ‘position distance’ raw data (Figure 1a) corresponds to the distance between the start point of two consecutive repeats and the ‘length’ raw data (Figure 2a) corresponds to the length of each repeat. The idea here is to assess whether the location of each element and its length follow any pattern identifiable by the complexity metrics. Evaluating the patterns in positions and lengths of ERV elements through complexity metrics provides potential roles these elements may have in gene regulation and chromosomal stability.

4. Results

4.1. Hurst Exponent

The Hurst exponent was estimated, for both ‘position distance’ and ‘length’ sequences and for each chromosome. Figure 2 presents the estimated values of the Hurst exponent. The blue line corresponds to the ‘position distance’ while the orange one to the ‘length’ sequence.

A small value of Hurst exponent shows a higher fractal dimension and a rougher surface. Biologically, fractal dimension characterizes the complexity or irregularity of genomic element arrangements. A lower fractal dimension indicates a smoother, structured genomic landscape, likely associated with regions of high regulatory or functional importance. In contrast, higher fractal dimensions suggest complex, irregular arrangements, reflecting dynamic genomic interactions.

A larger Hurst exponent shows a smaller fractional dimension and a smoother surface. The values of the Hurst exponent range between 0 and 1. A value of 0.5 indicates a true random process (a Brownian raw data). A Hurst exponent value, between 0.5 and 1 indicates “persistent behavior”. A Hurst exponent value between 0 and 0.5 indicates “anti-persistent behavior”. As we observe from Figure 2, all chromosomes characterized by “persistent behavior”. Persistent behavior describes the sustained correlation of sequential genomic events or positions across genomic distances, indicating nonrandom, structured patterns. This persistence biologically implies evolutionary constraints and genomic functional stability, essential for maintaining genomic integrity and regulatory mechanisms.

For the estimation of the Hurst exponent in this study we use Rescaled Range Analysis (R/S) (Weron, 2002). There is a clear discrimination profile of hurst in positions distance data and length data. The independent two-sample t-test was employed to statistically compare complexity metrics between 'position distance' and 'length' sequences. In t-test calculator, the t-value is 8.06309. The p-value is < .00001. The result is significant at p < .05. The p-values obtained (< 0.05) indicate statistically significant differences, validating distinct complexity profiles between these two types of genomic sequences.

Figure 3. Example calculation of the Hurst exponent for chromosome 6, demonstrating the method used (Rescaled Range Analysis, R/S) and its interpretation for genomic complexity.

Figure 3. Example calculation of the Hurst exponent for chromosome 6, demonstrating the method used (Rescaled Range Analysis, R/S) and its interpretation for genomic complexity.

4.2. q Stationary (Tsallis)

In Figure 4, we present the estimation of qstat index for both sequences and for all chromosomes. Concerning the qstat index, as one can see, the value in all chromosomes in both sequences are higher than 1 and suggests the presence of long-range correlations, a distinctive property of open nonequilibrium systems, with underlying spatiotemporal organization characterized by non-Gaussian (q Gaussian) distributions. The emphasis on non-equilibrium open systems illustrates that genomic organization is dynamically maintained, shaped by continuous external and internal genomic interactions, reflecting the ongoing adaptive and evolutionary processes.

As we observe, the qstat index in ‘position distance’ sequence is higher than the ‘length’ sequence for all chromosomes. This means that the non-extensive character of the spatial organization is higher in ‘position distance’ than the ‘length’ sequence and presents stronger long-range correlations in ‘position distance’ sequence. Long-range correlations describe dependencies and interactions among genomic positions widely separated within chromosomes. Biologically, these correlations imply that distant genomic regions influence each other's functional and structural characteristics, potentially affecting gene regulation, chromosomal stability, and adaptive genomic responses.

Moreover, in some chromosomes, we observe a significant differentiation of the qstat index between the two sequences (i.e. Chr 15, Chr21, Chr 22), while there is a high variation of the qstat index in Chromosomes 15-22 between the two sequences. There is a clear discrimination profile of qstat in positions distance data and length data. In t-test calculator, the t-value is 10.38241. The p-value is < .00001. The result is significant at p < .05. In Figure 5 we present a sample of the qstat index calculation for chromosome 6.

4.3. Complexity Factor (COFA)

The Complexity Factor (COFA) synthesizes diverse complexity metrics, offering a comprehensive measure of genomic organization connecting the real space with the phase space of the DNA structure. In Figure 6, the estimation of the technical term COFA per chromosome for both sequences is presented. There is a clear discrimination profile of COFA in positions distance data and length data. In t-test calculator, the t-value is 2.4071. The p-value is < .00001. The result is significant at p < .05.

4.4. K-Means Clustering

In this section we applied the unsupervised k-means clustering as an input in ML algorithms with the thought to see if the variation of the metrics that correspond to each HERVd entity for all chromosomes can be identified as a common dynamical feature which is characterizing these geometrical indices positions distance and lengths. We prepared the model using a distinct set of complexity metrics every time and we run the clustering process. To evaluate each clustering process, we used the Davies-Bouldin (DB) index (Davies, 1979). The DB index provide an internal evaluation schema (the score is based on the cluster itself and not on external knowledge such as labels) and is bounded from 0 to 1, where a lower score is better.

In Figure 7, the DB performance of clustering process for the Complexity Factor (COFA) position distance and lengths for different values of k parameter is presented. In COFA positions we calculate the best (lower) DB index performance for k = 8 parameter. Similarly, in COFA lengths we calculate the best (lower) DB index for k=6.

In Table 1 the clusters for the best DB index performance are presented. Each cluster included a set of different chromosomes with a common geometrical center of the variations of the COFA index. With this method of clustering based on the COFA index we discriminated sets of chromosomes, which appear to have similar complex behavior on position distance or lengths of retroviral elements as we observe in Figure 8.

The K-means clustering based on COFA values reveals groups of chromosomes with shared complexity profiles, reflecting potential structural, regulatory, or evolutionary constraints. These biologically meaningful clusters offer valuable insights into genomic architecture, suggest common functional or evolutionary patterns, and support hypotheses about coordinated gene regulation and evolutionary dynamics.

5. Discussion

In this study we analyze the statistical characteristics of the Human Endogenous Retro Viruses database (HERVd) focusing on the subcase of positions distances and lengths of ERVs elements, based on complexity metrics (Hurst exponent, q stationery index of Tsallis statistics and COFA technical index) for the purpose of understanding the degree of complexity behavior and internal organization of chromosomes in relation to the embodiment of retroviral DNA into the human genome. The internal organization refers to structured interactions and dependencies within the genome and the complexity behavior indicates measurable patterns of nonrandomness and structured organization that reflect underlying regulatory or adaptive genomic processes.

The analysis was based upon complexity metrics to phase or physical space with the estimation of Hurst exponent, q-stationary of Tsallis and COFA index and presented variations in the degree of complexity behavior per chromosomes in relation with the positions and lengths of ERVs elements. This analysis shows memory effects and long-range correlations in the distributions of ERVs elements in all chromosomes regarding their positions and lengths within the chromosomes. "Memory" in genomic sequences refers to the sustained positional correlations observed over long genomic distances, indicating functional links between distant regions potentially mediated by regulatory elements or structural constraints. This phenomenon, along with the "multiplicity" of roles that these positions can assume across various regulatory and structural contexts, plays a crucial role in maintaining genomic stability, adaptability, and overall functionality.

Moreover, the findings of this study reveal that the geometry of ERVs elements (positions and lengths) create a complex environment that communicates informationally with the spatial information of DNA as a separate subsystem that affects the biological functionality of the genome. Structured genomic complexity, characterized by ERV element positions and lengths, influences genomic functionality, potentially modulating gene expression and other cellular regulatory processes.

In this direction, complexity theory and computational tools can lead to further decoding of hidden information within the DNA. In addition, the Tsallis theory were used in this study showed the existence of the non-Gaussian character regarding the positions and lengths of ERVs elements and the embodiment of retroviral 261 DNA into the human genome.

The results of the Hurst exponent reveal that the size distributions of positions distance and lengths of retroviral elements in the genome are characterized by memory character or persistent behavior in all chromosomes. Specifically, this memory character has a differential profile so much between positions and lengths of retroviral elements among all chromosomes as well. It is observed that position distance retroviral elements maintain a higher Hurst exponent in all chromosomes suggesting that the distribution of position distance retroviral elements possess an enriched multiplicity character with a high degree of organization, as opposed to lengths of retroviral elements that maintain a lower degree of multiplicity and therefore a lower degree of organization. The above conclusion does not apply for the chromosomes 1,12,17,19,22 where the hurst exponent has approximately the same value. The above, in biological terms, may suggest that the positions of retroviral elements are engaged in multiple structural or functional roles, while lengths of retroviral elements are more limited. The multiplicity and the positional complexity refer to the capability of genomic elements, such as ERVs, to participate in diverse functional or structural roles within genomic contexts, reflecting versatility in regulatory and evolutionary dynamics.

The results of the q stationary reveal that the size distribution of the position distance and lengths of retroviral elements in all chromosomes is characterized by long range correlations. Non-extensive statistical mechanics refers to systems where traditional (extensive) thermodynamic relationships do not scale linearly with system size. In genomic terms, this implies that the observed complexity and correlations extend beyond local genomic segments, indicating widespread genomic interactions.

This non-extensive behavior is stronger in position distance of retroviral elements as compared to lengths of retroviral elements with some degree of variations per chromosome. Similarly, the variations are also significant in position distance data, reflecting long range correlations within chromosomes and specifically in chromosomes 15-22. Both position distance and lengths of retroviral elements distributions are independent of the chromosomal size except chromosomes 15-22. It seems that the long-range correlations are getting stronger in position distance data from the lengths data and specifically in chromosomes 15-22. These results would suggest that the positions distance of retroviral elements is coordinated with the position distance located in other distant regions of the same single chromosome and that all chromosomes have similar interactive structural relationships dictated by the same principles as suggested by the Tsallis q stationary index. The q stationary index clearly demonstrates that there is a coordination of the distributions of the position distance and lengths of retroviral elements within chromosomes characterized by specific profiles.

With COFA index and ML models we identify sets of position distance and lengths of retroviral elements among all chromosomes which present similar complex behavior. These new sets may contain interactions of information among chromosomes based on internal laws and geometrical symmetries that affect the biological functions of the cell. The reference to cellular processes emphasizes that structural genomic features, influenced by ERV distribution, likely affect gene expression and cellular functionality.

In future studies, we plan to analyze genomic information using sliding windows within the framework of Tsallis non-extensive entropy. By applying this lens with varying zoom levels (both zooming in and out), we aim to observe the distribution of complexity metrics and examine how these metrics vary across different regions of DNA chromosomes. Additionally, we will explore how these variations relate to the thermodynamics of the cell, providing insights into the relationship between genetic information and cellular energy.