Patterns of Transposable Element Distribution around Chromatin Ligation Points Revealed by Micro-C Data Analysis

2.1. Public Chromatin Conformation Capture Datasets

We analyzed chromatin interaction data from four independent datasets to identify ligation points (LPs). Dataset 1 (DS1) and Dataset 2 (DS2) were produced by Micro-C: SRR12625672 and SRR12625674 (biological replicates from HUDEP cell line, ~43M paired-end reads each, 150nt; https://www.ncbi.nlm.nih.gov /sra/SRR1 2625672). The Micro-C protocol employs micrococcal nuclease digestion, achieving nucleosome-level (~200bp) resolution.

Datasets 3 and 4 (DS3 and DS4) were produced by Hi-C: SRR27906244 and SRR27906243 from the GM13977 cell line (~75M paired-end reads each). Hi-C libraries were constructed using a standard protocol with restriction enzyme digestion (MboI, MseI, or NlaIII), which typically produces 2-6 kilobase fragments, followed by biotin incorporation and proximity ligation. The Hi-C data were obtained from NCBI Bioproject PRJNA1074296 (https://www.ncbi.nlm.nih.gov /bioproject /1074296) and processed using JuicerTools v1.14.08. Despite the methodological differences between Micro-C and Hi-C, all datasets were analyzed using identical parameters for TE density patterns and correlations.

2.2. Identification of Ligation Points

Although contact points in the chromatin conformation capture are called points, they are imprecisely defined and are actually represented by short fragments. Yet, it is possible to define precisely to a single nucleotide which sequences were ligated in the assay. The first step in the analysis was to identify the exact positions of the ligation points on the genome. The method identifies ligation points with single-nucleotide precision by using the alignment of paired-end reads to the genome. In Micro-C, most paired reads overlap, forming a continuous sequence. This sequence corresponds to the real physical fragment (ligated_fragment) produced by ligation in chromatin conformation capture. For each paired read, we reconstructed the ligated_fragment sequence by aligning paired reads. Next, we aligned the ligated_fragment against the genome. Only those paired reads were retained where two parts of the ligated_fragments (called here arms) aligned to the same chromosome and were positioned at least 1Kb apart. Once the arms were aligned, the exact position of coordinates of the chromosomes that were ligated to each other in the assay were identified (LP1 and LP2). The pipeline for the identification of the ligation points was called LigP_finder_v2 (the Python code for LigP_finder_v2 can be downloaded from: https://github.com/maxrempel /DRRF/tree/ main/LigP_ finder-main).

All datasets were processed using whole-genome alignment in LigP_finder_v2.

When applied to Hi-C datasets, the program produced fewer ligation points since the reads in Hi-C did not overlap, but the number of found ligation points was sufficient to produce TE density pattern plots and correlations. In general, the amount of data from the experiments was in excess since the plots did not visually deteriorate even when a single small chromosome was used for pattern plotting. The pattern was consistent across chromosomes, although the density of ligation points varied widely, about 10-fold between chromosomes, reflecting variable compaction of chromosomes in given cell lines.

To test whether the program was correct, we ran tests using synthetic data and performed a spot-check of the results. Also, the final program was described in English algorithm, and reprogrammed from scratch by another remote programmer. The recreated program produced identical results, demonstrating that both programs were correct.

The ligation points, LPs, were obtained with hg38, unmasked version from: https://hgdownload.soe.ucsc.edu /goldenPath /hg38/bigZips /hg38.fa.gz.

We used the hg38 transposable element (TE) annotation from UCSC (https://hgdownload.soe.ucsc.edu /goldenPath /hg38/bigZips /hg38.fa.out.gz).

2.3. Plotting Density Patterns

For each identified LP, we examined the genomic regions extending 50kb in both directions, which we termed “exbors” (from EXtended harBORS). Within these regions, we cataloged the density of transposable elements using bedtools and function intersect. In that, the rigth and left-oriented elements were analysed separately. We focused on the most abundant TE families and subfamilies based on their genomic copy numbers.

To quantify TE density patterns, we calculated TE density in 1kb bins across the 100kb window centered on each LP. Importantly, we conducted separate analyses for elements on the plus and minus strands to investigate potential strand-specific (same as orientation-specific) patterns. To visualize these density patterns, we smoothed the density curves with the Kernel Density Estimation (KDE) smoothing method.

For negative controls, we performed random control analyses. This involved generating random LP positions across the genome in excess and selecting those random coordinates that were located in non-repeat-masked (unique) regions across the whole genome. Then, real and random coordinates of LPs were compared using density and correlation plots. To control for programming errors, three programmers (DO, AV, and MM) reproduced the prlots independently, and the plots came out identical. Correlations between density pattern plots were calculated using the Pearson correlation and visualized as diagonal heatmaps. The Python code for LigP_finder_v2 also contains TE density pattern plotting with KDE smoothing (https:// github.com /maxrempel /DRRF/ tree/main /LigP_ finder- main).

Some plots are shown in the results section, and others are shown in the Supplementary Materials pdf file.

2.4. Analysis of Harbor Homology

To measure sequence similarity between harbor pairs, we extracted ±5 Kb sequences around each ligation point (10kb total harbor size). Sequence alignment was performed using minimap2 with the 'sr' preset and minimum 70% identity threshold. To compare within-pair versus between-pair homology, we analyzed sequence similarity between harbors from the same ligation point pair (REAL) versus harbors from different ligation point pairs (CONTROL). Both unique/low-copy sequences (regions remaining after repeat masking) and total sequences (unmasked) were analyzed. Pairs with overlapping harbors (distance < 10kb) were excluded from the analysis. The Python code for homology analysis is available at https:// github.com /maxrempel /DRRF/ tree/main /LigP_ finder- main.

3.1. Additional Transposable Elements

Our analysis of transposable element (TE) distribution around chromatin ligation points (LPs) revealed that each TE family and subfamily has a distinct density pattern around ligation points. We examined the density of six major TE families: Alu, L1, L2, MER, MIR, and MLT elements, as well as several subfamilies within these groups, within a 100 kb window centered on LPs identified through Micro-C experiments (Fig. 1 for L1, Fig. 2 for Alu; Supplementary Figures for all tested TE families and subfamilies).

The distribution patterns observed in every examined TE family had very similar trends to the trends observed for L1 elements (Fig. 1). However, the density patterns around the ligation points were unique for each tested TE family and subfamily.

Moreover, the patterns for all the TE families and subfamilies were asymmetric relative to the ligation points, oriented with respect to chromosome direction (with the short arm (p) positioned towards the left, following standard genomic convention). This asymmetry indicates that transposon directionality relative to chromosome orientation participates in chromatin folding. The consistency of these patterns suggests that the evolutionary selection of transposon placement produced long-range chromosome-wide orientations of transposons. Such large-scale genomic patterns may represent fundamental principles of chromosome organization that persist through evolution.

Since every tested TE family and subfamily produced consistent density patterns around the ligation points, this suggests the fundamental role of transposable elements in sequence-specific chromatin folding that transcends the specific characteristics of individual TE types. The observed strand biases and dataset-specific variations appear to be general features of how TEs are distributed relative to chromatin contact sites.

To exclude potential methodological artifacts, we utilized Hi-C data (Datasets 3 and 4) to validate our findings. While both Hi-C and Micro-C capture chromatin conformation, they differ fundamentally in their digestion methods; Hi-C uses restriction enzymes producing fragments of several kilobases, while Micro-C employs micrococcal nuclease digestion, achieving nucleosome-level (~200bp) resolution. Despite this resolution difference, Hi-C data reproduced the key patterns of transposable element distribution around contact points, with correlation coefficients between biological replicates reaching 0.7-0.8 for same-strand comparisons and remaining below 0.2 for opposite strands (Supplementary Figures). This cross-method validation confirms the strand-specific organization of transposable elements around chromatin contact points. While absolute correlation values between datasets from different laboratories and cell lines were not expected due to the tissue-specific nature of chromatin organization, the fundamental patterns were consistently reproduced: strong correlations between biological replicates and between members of ligation point pairs and strong strand asymmetry. These patterns were reproduced despite using cell lines from different blood lineages (HUDEP erythroid progenitors (datasets 1 and 2) vs. GM13977 lymphoblastoid cells (datasets 3 and 4)), suggesting that the results reflect a common principle in chromatin organization.

3.2. Asymmetric Distribution Suggests Homological Adhesion

We noticed that the asymmetry of transposable element density patterns relative to chromosome orientation suggests tandem-like organization. However, we found no periodic repetition of these patterns, suggesting an aperiodic repetition of patterns that maintains directional consistency. Coincidentally, this aligns with Erwin Schrödinger's description of the hereditary material as aperiodic crystal [13].

The figure shows how patterns of identical transposon sequences (arbitrarily placed colored arrows) can serve as contact points through sequence-specific homological adhesion. The 100 Kb window demonstrates the alignment of identical transposons creating interactions between distant chromatin regions.

We propose that this aperiodic but directional repetition of TE patterns functions in chromatin architecture through homological adhesion, where identical sequences, particularly transposable elements, can form contact points through sequence-specific adhesion of identical sequences. Here, we propose for the first time that two parallel double helices of DNA can adhere to each other when their sequences are identical, forming the molecular basis for homological adhesion. The density plots reveal that the patterns are asymmetric and face in one direction on chromosomes. Since they face in one direction on chromosomes, they could serve as specific anchors for homological adhesion and the formation of large-scale helices (Fig. 4).

We next tested if the sequences that come together in chromatin display homology. As a negative control, we measured homology between unpaired unligated harbors from different harbor pairs. Comparing harbors (ligated fragments) from the same harbor pair versus different pairs revealed strong enrichment of homology in ligated regions. In unique and low-copy sequences, within-pair homology reached 70%, while between-pair homology remained near zero (mean 0.47% vs 0.02%). The signal persisted in an unmasked sequence, with within-pair homology reaching 100% compared to a maximum of 40% between pairs. This provides quantitative evidence for sequence-specific homological dsDNA-dsDNA adhesion in chromatin contacts.

4.1. Homological Adhesion

The discovery of the density patterns of transposable elements around the chromatin contact sites and of the homology between contacting DNA sequences offers support for homological adhesion. This is referred to as “homadhesion” for brevity. Homadhesion is the adhesion of two DNA duplexes (double helixes) to each other if they have similar sequences.

Several studies provide support for sequence-dependent attraction between DNA duplexes. Yoo et al. [15] demonstrated sequence-dependent attraction between double-stranded DNA molecules using molecular dynamics simulations and single-molecule FRET experiments. They showed that DNA duplexes can attract each other over distances up to 2-3 nm in the presence of polyamines like spermine. Importantly, AT-rich sequences showed stronger attraction than GC-rich sequences, and DNA methylation enhanced these interactions.

Lee et al. [16] provided theoretical support for sequence-dependent attraction between intact DNA duplexes. Their model incorporated electrostatic forces and sequence-dependent DNA shape variations, predicting lower interaction energies for pairs of DNA fragments with parallel homologous sequences compared to those with uncorrelated sequences. This interaction was modeled without strand separation, considering the aqueous environment through parameters that account for electrostatic screening in solution.

Barzel and Kupiec reviewed the evidence for the pairing of similar DNA sequences across different organisms [17]. They noted that in yeast, matching DNA sequences can locate each other and recombine efficiently even when in different genomic locations despite the large amount of genomic DNA present. Their review suggested that similar sequences are paired as part of the genome's basic organization.

For the adhesion mechanism itself, several forces could be considered: electrostatic (including ionic), hydrophilic-hydrophobic, hydrogen bonding, van der Waals, and other weak forces. The dynamic nature of chromatin condensation and decondensation suggests that DNA duplexes remain nearly intact during these interactions. While hydrated DNA structure may be primary in mediating adhesion, other nucleoplasmic components likely participate, including histones, other proteins, and low molecular weight chemicals. The highly negatively charged DNA duplexes must overcome electrostatic repulsion to adhere. This requires neutralization by positive ions present in the nucleoplasm, including protons (H⁺), hydronium ions (H3O⁺), positively charged histones and other proteins, and ions such as Na⁺, K⁺, Mg²⁺, and polyamines like spermine and spermidine.