Understanding three-dimensional chromatin organization in diploid genomes.

The three-dimensional (3D) organization of chromatin in the nucleus of diploid eukaryotic organisms has fascinated biologists for many years. Despite major progress in chromatin conformation studies, current knowledge regarding the spatial organization of diploid (maternal and paternal) genomes is still limited. Recent advances in Hi-C technology and data processing approaches have enabled construction of diploid Hi-C contact maps. These maps greatly accelerated the pace of novel discoveries in haplotype-resolved 3D genome studies, revealing the role of allele biased chromatin conformation in transcriptional regulation. Here, we review emerging concepts and haplotype phasing strategies of Hi-C data in 3D diploid genome studies. We discuss new insights on homologous chromosomal organization and the interplay between allelic biased chromatin architecture and several nuclear functions, explaining how haplotype-resolved Hi-C technologies have been used to resolve important biological questions.

Introduction

Chromatin has a linear length spanning several meters but is compacted into the micrometer-sized nucleus of eukaryotic cells via an intricate and non-random folding process. High-throughput chromosome conformation capture (Hi-C) [1] is an essential technology providing genome-wide three dimensional (3D) chromatin organization by coupling proximity-based ligation with massively parallel sequencing [2]. The development of Hi-C and its derivatives [3], [4], [5] has revealed a complex genomic landscape of chromosomal structure across different cells and tissues of different species. Multilevel hierarchical chromatin structures are found in interphase cells of most higher eukaryotes, including chromosome territories (CTs), compartments of active (euchromatin, A-type) and inactive (heterochromatin, B-type) chromatin; chromosomal domains variably known as contact domains or topologically associating domains (TADs); chromatin loops; and promoter-enhancer interactions (PEIs) [1], [3], [6], [7], [8], [9], [10], [11], [12], [13]. The conformation of chromatin plays important roles in the regulation of gene expression, recombination, and DNA replication and repair [11], [14], [15], [16], [17], [18]. Hence, understanding 3D chromatin folding is essential to determining the mechanisms behind gene regulation in multicellular organisms across different cell types developing from the same genetic blueprint, as well as to reveal the genetic basis underlying distinct phenotypes and physical disorders.

Despite recent advances, at least two limitations affect the majority of 3D genome studies. The first is masking the differences in chromatin conformation that are inherent to individual cells. The other is assembling an averaged structure of the homologous chromosomes. The former has essentially been overcome by the development of single cell Hi-C technologies [19], [20], [21], whereas the latter remains a challenge. As a typical Hi-C experiment does not directly produce haplotype-resolved (phased) data, most Hi-C studies lack sufficient allele-specific information to distinguish the two parental haploids [22]. Apart from some genomic loci [23], [24], [25], little is known about the variability in higher order structures between homologous chromosomes in diploid genomes. In fact, substantial allelic differences in transcription, DNA methylation, and epigenetic chromatin states have already been found in diploid organisms [26], [27], [28], [29], [30], [31]. Hence, it has long been hypothesized that two homologous sets of chromosomes may differ in how they spatially fold and function.

Recently, significant research efforts have attempted to decipher 3D haplotype structures in diverse diploid organisms (Supplementary Table S1). In this review, we discuss recent methodological and biological advances in the search for 3D haploid-level genome architectures. We focus on the spatial organization of homologous chromosomes within the nucleus as well as their associations with different genomic functions. We also provide an outlook on the future trends and challenges for this field, as well as a glimpse into possible developments in 3D diploid genome research.

Reconstruction of the haplotype-resolved 3D genome

Conventional Hi-C data does not typically distinguish between different allelic chromosome copies. Hence, the main challenge in 3D haplotype-resolved genome studies is accurate allele assignment of each chromatin contact (pair of genomic loci that are joined by proximity ligation) [3], [32], [33], [34]. Observing a single Hi-C contact between locus i and j corresponds to one of four possible events, where either copy of locus i comes in contact with either copy of locus j [35]. When i and j come from the same chromosome, the contact represents a cis interaction. In contrast, when the two loci come from different chromosomes, it is considered a trans interaction. In particular, if i and j respectively belong to the homologous maternal and paternal chromosomes, the contact represents a trans-homolog (t-homo) interaction that reflects the spatial architecture of homologous chromosomes. Analogously, if i and j come from heterologs, the contact is classified as trans-heterolog (t-heter). Importantly, every 3D genome study at a true diploid scale must accurately account for this allelic uncertainty. In this section, we discuss the methodology currently used for reconstructing haplotype-resolved 3D genomes.

Phasing Hi-C data based on heterozygous variations

The most commonly used approach to assign alleles of each Hi-C contact is employing phased single nucleotide variations (SNVs) to define haplotypes from the Hi-C reads (Fig. 1). Parent-specific SNVs make it possible to resolve the haplotypes of Hi-C data by mapping reads containing allele-specific genetic variations (i.e., informative Hi-C reads) to their parental chromosome. Because many Hi-C reads overlap SNVs, this strategy efficiently assigns Hi-C contacts to specific chromosomal homologs in order to finally construct a ‘diploid’ Hi-C map comprising both maternal and paternal contact matrices. However, only in rare instances is it possible to obtain phased SNV data, as is the case of the International Hapmap Project [36] and the 1000 Genomes Project [37], [38], [39]. Such information is not available for the majority of studies. Genotyping parent–child trios can directly determine whole-genome haplotype information for the F1 hybrid based on heterozygous sites. However, in order to obtain a high density of heterozygous sites, and consequently high resolution of haplotype phasing, it is usually necessary to establish hybrids with divergent breeds or lineages. The hybrids should be sufficiently divergent to allow for straightforward sequence-level discrimination of homologs but maintain a nearly complete synteny in order to guarantee the mapping rates. Hi-C reads can be aligned to reference genome N – masking variable sites and then performing haplotype imputation based on phased variants. Alternatively, Hi-C reads can be mapped to the diploid F1 hybrid genome, which is then rebuilt with the reference genome and the phased variants (replacing the corresponding bases) to directly obtain allele assignments of the Hi-C reads.

Fig. 1

Schematics of haplotype phasing strategies for Hi-C read pairs. (A) The strategy of haplotype phasing for paired-end Hi-C reads based on heterozygous variants. (B) Two local imputation phasing methods for Hi-C reads. Left: The haplotype imputation algorithm in Dip-C developed by Tan et al. [21]. An unknown haplotype of chromatin contact (circled blue dot) can be imputed based on the statistical properties of interchromosomal and long-range intrachromosomal contacts (blue dots) in the super-elliptical neighborhood (light blue shaded area) of the unphased contact (circled blue dot). The plot is modified from Tan et al. [21]. Right: A novel phasing technique for Hi-C data, named ‘HaploHiC’. A paired-end Hi-C read with no allele-specific SNV or short InDels (insertions and deletions) has its origin imputed using information from nearby reads. The ratios of paternally and maternally mapped reads in the flanking regions of the unknown-haplotype reads are used to determine the likelihood of a given parental origin of unknown-haplotype reads.

Paired-end Hi-C reads can be classified into three different groups based on SNVs: reads containing SNVs at a) one or b) both ends, and c) reads containing no SNVs, respectively termed partially informative, informative, and non-informative reads (Fig. 1A). Haplotype phasing for reads with allele-specific SNVs at both ends is simple and accurate: it is possible to unambiguously map these reads to their respective parental genome based on the SNVs present (Fig. 1A). All proposed methods can process this group of Hi-C reads, including the algorithm Haplotype-Assisted Read Parsing (HARP, https://github.com/dvera/harp) [40]. Notably, HARP performs allele-specific alignment of Hi-C reads to the two parental reference genomes, thus yielding a higher number of usable reads for phasing due to the separate alignment of each haploid. Rivera-Mulia et al. [40] took advantage of HARP methodology and the high SNV density between the genomes of mouse subspecies to explore the allele-specific control of genome organization and replication timing during development. This kind of approach has been widely used to obtain 3D diploid genomes [3], [41], [42].

With partially informative Hi-C reads, it is possible to determine parental origin by arbitrarily using the SNVs at one end of the read, or by combining information from local imputation phasing based on empirical probabilities or statistical modeling in order to increase accuracy (Fig. 1A). Some studies have adopted the former strategy, including Giorgetti et al. [34], Du et al. [41], Bonora et al. [43], and Greenwald et al. [44], while others have opted for the latter, such as the studies proposing Dip-C [21] and HaploHiC [45]. Dip-C, an algorithm originally designed for single-cell Hi-C data by Tan et al. [21], can impute unknown haplotypes from ‘neighboring’ contacts (measured as linear distances) by assuming that the two homologs would typically contact different partners (Fig. 1B). Briefly, considering the statistical properties of interchromosomal and long-range intrachromosomal contacts, Tan et al. [21] defined a contact neighborhood as a super-ellipse with an exponent of 0.5 and a radius of 10-Mb, where the haplotypes of nearby contacts were weighted to impute the origin of the unphased contacts. Then, the assignment of allele-ambiguous contacts was refined by the 3D structures. In this case, the imputation accuracy for each haplotype was estimated as ~ 96% by cross-validation.

Notably, a substantial proportion of total Hi-C reads are non-informative, ranging from 21.1% (patski) to 94.1% (GM12878) as reported in [46], and are typically discarded in most phasing methods. Even for hybrids, maternal and paternal haploids are often identical in most parts of the genome, and most phasing approaches can only infer the parental origin at a small subset of genomic loci [45]. This means that the segregated maternal and paternal contact matrices are extremely sparse. To tackle this problem, a few methods [21], [45], [46], [47] have attempted to reintroduce unphased reads to both parental contact maps. The Poisson-Gamma model developed by Deng et al. [47] can impute both partially informative and non-informative contacts using an iterative expectation–maximization (EM) algorithm. However, this model might not work robustly in fine-resolution analyses due to lack of diploid 3D structures to assist the assignment of allele-ambiguous contacts.

HaploHiC [45], a novel phasing algorithm, also allows for more Hi-C data to be imputed using local imputation combined with SNVs/InDels-based phasing, including Hi-C reads with only one or no informative read end, thus reducing the sparsity of the allele-specific Hi-C contact matrices (Fig. 1B). With HaploHiC, assigned reads containing informative SNVs/InDels can be employed to predict the origin of nearby unknown reads using a data-derived ratio. The haplotype-predicting accuracy of HaploHiC (97.66%) was comparable to Dip-C (~96%) for reads of unknown origin using simulation data. Notably, rigorous evaluations are required to further assess its efficiency and accuracy on more Hi-C data, given that HaploHiC is currently a preprint version.

More recently, Ye and Ma [46] proposed ASHIC methods, a hierarchical Bayesian framework to provide a probabilistic model of all possible allele origins, to impute allele-specific chromatin organizations and simultaneously infer allelic 3D structures in diploid genomes using both informative and non-informative Hi-C contacts. The ASHIC methods include the Poisson-multinomial (ASHIC-PM) model and the zero-inflated Poisson-multinomial (ASHIC-ZIPM) model. Both models outperformed the above-mentioned allele-certain methods [3], [42] as well as mate-rescue approaches in which allele-ambiguous read-ends are directly assigned to the same allele as with allele-certain mate-ends [34], [41], [43], especially under low coverage and low SNV density conditions using both simulated and real Hi-C data.

There at least two kinds of errors that may be introduced to allele assignment of Hi-C data based on heterozygous variants: homolog misassignment and systematic bias caused by variable genetic variant density in the genome. Erceg et al. [48] developed a computational approach, called Ohm (oversight of homolog misassignment), which can minimize homolog misassignment and thus robustly distinguish t-homo from cis contacts for haplotype-resolved Hi-C studies. This approach can be extended to any number of methods in order to improve the assignment of parental origin to Hi-C reads. As for the bias in diploid Hi-C contact maps caused by asymmetric genetic variant density along the genome, a novel strategy has been proposed by Luo et al. [49] to facilitate whole-genome identification of the 3D organization of diploid chromatin. This approach is included as an integrated tool called HiCHap.

Inferring the 3D structure of diploid genomes

Hi-C data allows for reconstruction of 3D genome structures with the aid of sophisticated computational algorithms. Even though Hi-C contact matrices yield valuable insights, modeling and visualizing 3D genome structures can unveil higher-order structural patterns and relationships that are not apparent from the raw data [11], [50], [51]. In particular, this might provide a 3D architecture that can be interpreted by humans, orienting genomic regions relative to various nuclear landmarks, and serving as a framework for integrating other data types [52]. Most inference methods have modeled a single structure without differentiating the two alleles, but many fewer algorithms have so far been developed to account for diploidy in Hi-C 3D modeling [6], [21], [46], [53], [54], [55], [56]. Hence, 3D structural inference of diploid genomes remains challenging.

Two separate research groups have developed modeling methods that can be exclusively applied to single cell Hi-C data to generate 3D physical models of diploid genomes [21], [54]. One approach [21] adopted an iterative modeling strategy, starting with coarse resolution modeling and gradually refining the structures to higher resolutions. For more widely available bulk Hi-C data, some genome modeling approaches [6], [50], [55] have been developed to generate a diploid 3D genome structure using a simulated annealing/molecular dynamics sampling method. Ay et al. [53] and Cauer et al. [35] proposed an extension of the PASTIS software to infer haploid 3D genome structures [57]. One published algorithm [35] can infer the 3D chromatin structure at a true diploid level from population Hi-C data. Specifically, it incorporates two constraints into structural inference – it enforces even spacing of the beads along the chromosome (the homolog separation constraint), and spatially separates the homologs (the bead connectivity constraint). Notably, the weight assigned to the homolog-separating constraint should be tuned for each organism based on prior knowledge of genome architecture, especially the relationship between the number of contacts and pairwise distances. Despite its good performance with real Hi-C data, this diploid PASTIS method can be further improved via joint estimation of the biophysical model alongside the 3D structure, or by incorporating a multiscale optimization strategy, in which a high-resolution structure is deduced in a stepwise fashion with gradually increasing resolution.

Coincidentally, Han et al. [22] have built an iteratively weighted adjusting model to fit each haploid chromosome into the 3D nuclear space from the allele-specific interaction matrices. The iterative process continuously minimizes the sum of errors between coordinate-based distance and the 'real' distance that was converted from the allele-specific interaction matrix. Briefly, Han et al. first initialized random 3D coordinates for each chromosome, and then iteratively adjusted the coordinates based on the distance errors between a specific chromosome and the remaining chromosomes. The model is thereby more likely to reach local optima and a steady state with finite iterations.

In addition, Ye and Ma [46] designed an iterative expectation–maximization (EM) algorithm to infer 3D structures of homologs. They developed an inter-homologous optimization strategy, that is, they estimated homologous chromosomes separately, then predicted the relative position between the two homologs to improve the final 3D structures. However, this method cannot be used to build an accurate genome-wide 3D model, due to lack of an additional estimation step to model the relative position of multiple haploid chromosomes.

Features of the spatial organization of homologous chromosomes

Dissection of 3D genome structures in diploid genomes has revealed homolog-specific features. In this section, we discuss the biophysical characteristics and links between these structural features.

Similar Hi-C contact patterns between homologous chromosomes

In interphase cells, homologous chromosomes have the most similar Hi-C contact patterns (Fig. 2A), a prominent characteristic of diploid genome architecture. Both inter- and intra-chromosomal contact profiles of a chromosome resemble those of its respective homolog during interphase in human and mouse cells, and this similarity is strongly correlated with allelic co-expression levels [3], [22]. Principal component analysis (PCA) on diploid Hi-C interaction matrices showed that homologous chromosomes are in close proximity on the PCA projection and that PCA distances between homologous pairs are significantly shorter compared to non-homologs [22].

Fig. 2

Schematics of highly structured homolog pairing in diploid genomes. (A) Haplotype-resolved Hi-C map of a diploid F1 hybrid, which is assumed to have three pairs of chromosomes. The three boxes on the map show the chromatin interactions within the haploid chromosome (cis contacts, box 1), between homologous chromosomes (trans homolog contacts, box 2), and between heterologous chromosomes (trans heterolog contacts, box 3). The t-homo diagonals along the main cis diagonal of the map indicate the homolog pairing. (B) Spatial organization of haploid-level chromosomes in the 3D nucleus of diploid F1 hybrids obtained by 3D modeling of the haploid-resolved Hi-C map. Cis (box 1), t-homo (box 2), and t-heter (box 3) contacts are indicated on the map. (C) Zoomed in Hi-C maps of boxes 1, 2, and 3 shown in (A). (D) Local contact map and 3D organization of variable homolog pairing. Homolog pairing encompasses tightly and loosely paired regions, and displays highly structured trans-domains and trans-interaction peaks.

Global folding patterns of homologous chromosomes are highly and consistently similar during interphase (Fig. 2B) [32]. Studies on 3D nucleus reconstruction have shown that homologs are nearly equidistant from the center of the 3D nucleus, and are usually located next to each other, confirming the existence of homolog pairing based on Hi-C data. The distance of a chromosome to the nuclear center is highly correlated with that of its homologous counterpart, and the distance between any chromosome and a homologous pair is highly correlated, supporting similarity in the spatial organization of homologs [22]. This is in agreement with traditional chromatin conformation studies showing that loci with similar genomic composition and chromatin status tend to be located close to each other [1], [58], [59], which is expected given that homologous chromosomes are nearly identical in size, sequence composition, and associated proteins. However, there are exceptions. For example, maternal and paternal X chromosomes in mammalian cells are located far away from each other due to X chromosome inactivation (XCI) [34]. In a study of 3D diploid nuclei [22], larger chromosomes tended to locate in one pole, with smaller chromosomes in the opposite pole, which was consistent with observations in non-haplotype-resolved Hi-C studies.

Homolog pairing in diploid genomes

One of the most intriguing biological phenomena in nuclear architecture of diploid genomes involves the colocalization of homologous chromosomes, known as ‘homolog pairing’ (Fig. 2C, 2D). Homolog pairing was first found in Drosophila somatic cells nearly a century ago, occurring from embryogenesis to adulthood, and various degrees of somatic homolog pairing have since been observed during interphase in diploid yeasts and mammals. Such pairing has been implicated in DNA repair, V(D)J recombination, imprinting, XCI, and cell fate commitment across species [60], [61]. As attention towards this phenomenon grows, it is becoming increasingly clear that homolog pairing can facilitate regulatory interactions between different alleles of the same gene, a process known as transvection or intragenic complementation [61], [62]. In particular, allele-specific Hi-C technologies shed light on the 3D organization and formation of homolog pairing, and partially explain whether pairing is uniform across the genome or how homolog proximity correlates with alignment precision and genome function.

3D organizational characteristics and their association with function

One of the main characteristics of homolog pairing is the remarkable variation in pairing intensity, which has been observed in diploid yeasts [63] and fruit flies [48], [64]. After achieving a high-resolution map of homolog pairing in the diploid Saccharomyces genome, Kim and colleagues found that the strength of pairing fluctuates across the genome and as a function of growth conditions [63]. Their results on condition-specific pairing also suggested a role for the nuclear pore complex in homolog pairing of Saccharomyces. Hence, the degree of pairing in Saccharomyces depends on both the genomic location and the external environment, and sometimes on an interaction between these factors.

Homolog pairing in Drosophila embryos is highly structured in different forms, including the well-aligned form ‘railroad track’, the less-aligned form ‘laissez-faire’, and highly disordered pairing [48]. By introducing pairing scores (PS) to quantify the variation of homolog colocalization, AlHaj Abed et al. [64] identified at least two distinct forms of homolog pairing in the somatic cells of Drosophila: tight and loose pairing (analogous to the ‘railroad track’ and ‘laissez-faire’ forms, respectively) (Fig. 2D), which have quite different 3D structural features. Tightly paired regions, usually spanning contiguous small TADs, are extensively interspersed with loosely paired segments typically consisting of single large domains. The two types of pairing also differ in the decay of cis and t-homo contact frequencies. Moreover, tighter pairing is more likely to be located at domain boundaries, suggesting that some insulator proteins and elements can contribute to homolog pairing and promote transvection, since the domain boundaries are enriched for insulators and architectural proteins [8], [65], [66].

Most recently, Hoencamp and colleagues [67] have determined that homologs tend to be separated or paired across species depending on whether they carry functional condensin II subunits. Species lacking all or partial condensin II subunits (e.g., fruit flies and mosquitoes), exhibit Rabl-like architectural features (centromere clustering, telomere clustering, and a telomere-to-centromere axis), and tend to be homolog pairing. In contrast, species with full sets of condensin II subunits (e.g., humans and mice), display typical chromosome territories and are more likely to show homolog separation or highly disordered pairing.

Besides the variation in pairing, significant genome-wide correlations were found between pairing, compartments, and gene expression in Drosophila [48], [64], [68]. Focusing on trans interactions, Erceg et al. [48] discovered that pairing correlates with the opening of the chromatin mediated by the transcription factor Zelda during zygotic genome activation. Taking this analysis further, AlHaj Abed and colleagues [64] uncovered complex associations between pairing and chromatin status, compartments, gene expression, and regulatory epigenetic states. First, the active compartment A or highly expressed regions were usually tightly paired, while compartment B showed a bimodal pairing distribution – inactive and repressed chromatin contained regions of both tight and loose pairing. Second, the pairing types were differentially associated with gene expression: highly expressed genes tended to be tightly paired, whereas genes expressed at low levels were associated with both high and low pairing. These observations demonstrate the complex relationship between homolog pairing and chromatin activity or, alternatively, the imperfect classification of pairing/chromatin status. Nevertheless, this suggests that pairing has functional implications.

Importantly, the close proximity of homologs may influence the way in which loops and domains are formed in Drosophila. As observed by AlHaj Abed et al. [64], the boundaries of loose pairings are tightly paired, supporting a model that integrates pairing, loop formation, and chromosome compaction via a mechanism in which the chromosomes are buckled out to form loops by extrusion, anti-pairing, or a combination of both [64], [69]. In the model, tightly paired regions form barriers to extrusion. A series of RNAi experiments showed that the disruption of pairing-promoting factors (Slmb and TopII) resulted in global changes in pairing and disruption of some interaction peaks, suggesting that the close proximity of homologs may influence the formation of loops and domains in Drosophila [64].

Potential somatic pairing in human and mouse, albeit rare, was captured on 3D models by Tan et al. [21]. Nonetheless, homolog pairing has not yet been systematically investigated in mammals. Joyce et al. [61] suggest in a review that pairing in mammals has so far remained unnoticed because of unknown mechanisms to prevent it. Thus, conducting studies on somatic paring of mammals can not only provide an overview of homolog positioning in these animals, but also clarify the inhibitory mechanisms underlying anti-pairing.

Mechanical models of pairing formation

Two main models explain the mechanisms driving colocalization of homologous chromosomes: the ‘zipper’ and ‘button’ models. In the ‘zipper’ model, the different regions of the genome have an equal pairing affinity based on sequence homology, and homolog pairing begins at the centromere and proceeds distally to the telomeres [62]. In contrast, the ‘button’ model proposes that regions of high pairing ability are interspersed across the chromosomes and come together through a random walk that launches pairing [70], [71]. Viets et al. [68] identified TADs (the ‘buttons’) interspersed across the fly genome that brought homologous chromosomes together in order to facilitate cell-type-specific inter-chromosomal gene regulation. These specialized TADs may take on unique chromatin conformations or uniquely bind different protein combinations to create nuclear microcompartments that enable homologous TAD association and pairing. In addition, using a button spanning the spineless gene as a paradigm, they revealed that pairing was necessary but not sufficient to enable transvection, indicating that these two processes were mechanistically independent. The ‘button’ model was also supported by the distribution and common features of tightly paired regions found by AlHaj Abed et al. [64].

However, there are many uncharacterized features of homolog pairing. For example, it is still unclear whether homologs become structurally similar because of pairing, or instead have comparable conformations prior to pairing [48]. Similarly, the existence of different pairing mechanisms remains unknown, in light of the contrary correlations between tight pairing and gene activity found by Li et al. [72] and AlHaj Abed et al. [64]. Present studies have not provided definitive evidence, meaning that these issues need to be addressed in the future.

Fundamental 3D architecture of homologous chromosomes

Multi-layered organizational structures

Allele-specific Hi-C research has not only validated the central cis diagonal (representing short-distance intra-chromosomal interactions) and hierarchical chromatin structures within chromosomes (Fig. 2C), but also revealed more detailed features of the multi-layered organization and the rich structure of t-homo interactions (Fig. 2D).

Haplotype-resolved chromosomes occupy distinct regions within the simulated 3D nucleus of most animals. These regions are known as chromosome territories (CTs) and their nuclear positioning mainly correlates with genomic properties and functions [21], [35]. Haploid-level Hi-C maps demonstrate the characteristic chromatin compartmentalization in diploid genomes [21]. As expected, broad similarity has been observed between the patterning of compartments and TADs in maternal and paternal genomes [22], [32], [45]. Moreover, a diploid Hi-C study revealed that A compartments have higher heterogeneity compared to B compartments in interphase cells of mice [22]. For example, the local boundary score (LBS) – a quantitative value with peaks representing TAD boundaries – within A compartments exhibit stronger fluctuations and higher variability than that in B compartments. The multiple TADs found in A compartments contrast with the often single TAD found in B compartments, demonstrating the complexity of A compartments in diploid genomes.

Furthermore, allele-specific Hi-C research on Drosophila has revealed the existence of t-homo diagonals in the diploid Hi-C maps, and showed that paired homologs can form t-homo compartments, TADs, and loops (Fig. 2D) [48], [64]. The hierarchical structures of t-homo interactions largely coincided with analogous features in cis contacts. In particular, the observed concordance between the t-homo, cis-maternal, and cis-paternal Hi-C maps in terms of the sizes and positions of domains and loops supports the idea that cis and trans interactions are structurally coordinated and may be formed by similar mechanisms [48], and demonstrates a high level of agreement between paired homologs [64].

Similar to what was found in haploid mouse embryonic stem cells (mESCs) [73], chromosomes in diploid mESCs preferred a more parallel conformation (i.e. Rabl configuration), whereas terminally differentiated cells displayed a more centromere-facing-out organization [21]. The differentiated cell type-specific chromatin organization has been confirmed by Han et al. [22], who showed that inter-chromosomal interactions increased along the chromosomes from the centromeres to the telomeres. These observations suggested that chromatin conformation may function in pluripotency and the characterization of cell identity.

Allelic-biased chromatin conformation

Rare genomic regions show changes in compartment status between alleles. Despite the high similarities in the patterns of compartments between homologs, ~0.6–2.3% of the genome exhibited different compartments between homologs in human embryonic stem cells (hESCs) and four hESC-derived lineages [32]. Furthermore, allele switching between active and inactive compartments was observed genome-wide in the cell cycle of human B-lymphocytes, and the switched alleles were not the same in the different cell cycle phases [74]. Moreover, the vast majority of AB switching occurred at the borders between compartments, instead of within entire compartments [45], [74]. This indicated that the boundaries may be more unstable compared to the inner parts of the compartments. The compartments identified in two alleles of a single cell can vary in an almost independent manner, with an average Spearman correlation close to zero [21]. Han et al. [22] found that the A/B compartment status of haploids in hybrid mice was highly correlated with its parent-of-origin due to cis-effects, and parental divergent compartments basically transited into the same status in hybrid mice due to a shared cellular environment, i.e., the trans-effect. Approximately 58% of the divergent A/B compartment bins between two parents retained the same status as the parent-of-origin in hybrid mice (cis-effect), while 95% of the remaining parentally divergent compartments converged into the same compartment status (trans-effect). These observations confirmed the essential function of micro-environments in epigenetic regulation of gene expression.

Even though TAD is considered to have a relatively conservative organizational structure among different species and across tissues, allele-specific TADs have been detected in diploid genomes. TADs with the largest allelic differences are smaller in size and contain higher gene density compared to randomly selected TADs [74]. Lindsly et al. [45] carefully identified sub-TADs and local structural changes, and determined that sub-domains were not always consistent between alleles or cell cycle phases, and were even more variable between cell types.

Several studies focusing on chromatin loop variation across haplotypes have provided conflicting conclusions. Rao et al. [3] studied allelic differences in looping with diploid human in situ Hi-C contact matrices and observed few disparities between the maternal and paternal 3D genome maps except for imprinted regions. In contrast, studies using dilution Hi-C and 4C-seq [32], ChIA-PET [15], and H3K27ac Hi-ChIP [75] technologies reported that allelic imbalance in chromatin looping occurred throughout the genome. Analyzing phased human multi-omics data including high resolution haplotype-resolved Hi-C maps, Greenwald et al. [44] quantitatively investigated chromatin looping across different haplotypes and cell types. Aided by the increased power of combined cell type data, the authors found subtle allelic differences in the majority of loops, and identified infrequent, large allelic imbalances in looping in the genome (only 114 haplotype-associated chromatin loops across all individuals [FDR < 0.05]), which were mainly driven by imprinting and/or copy number variants (CNVs). These contradictory results are likely a result of distinct experimental designs, detection power, and the types of effects examined in these studies.

From homologous chromatin conformation to genome function

Allelic differences have been detected in gene expression, protein binding, and 3D chromatin organization [76], [77], [78], [79], [80]. Given the basic biological principle that structure (form) determines function, it follows that 3D chromatin architecture contributes to allele-biased (or allele-specific) gene expression and regulation. At present, both direct and indirect evidence support an intimate relationship between allele-specific genome structures and nuclear functions.

Allele biased expression (ABE)

Diploid 3D genome studies showed that ABE can, to some extent, be attributed to allele-biased chromatin structures. By integrating haplotype-resolved 3D genome maps with epigenomic and transcriptomic data, Dixon et al. [32] observed that genomic regions containing allele biased or imprinted genes show a subtle but statistically significant increase in the variability of A/B compartment scores between alleles. However, A/B compartment switches between alleles constituted rare regions of the genome and were not enriched for either allele-biased or known imprinted genes. These results indicate that despite most allele-biased genes having no differential compartment status between alleles, subtle differences may exist in higher-order chromatin structure between homologs at ABE regions. In addition, they found that allele-biased enhancer activity possibly underlies ABE through either stable DNA looping between alleles or potential allele-specific PEIs.

Interaction domains also function in ABE. Focusing on the top 10% of bins with the largest degree of Hi-C change to fully characterize Hi-C allelic differences, Chen et al. [74] discovered a subset of ABE genes that were preferentially localized near TAD boundaries in lymphoblastoid cells. These genes were enriched with chromatin organization transcription factor binding sites and contained a higher number of SNVs/InDels in the CCCTC-binding factor (CTCF) binding sites, providing a view of how the two haploid genomes differ in form and function. Parental-specific domains participate in parentally biased gene expression in early murine embryogenesis, including parentally pre-formed domains associated with a lower average expression on the structured allele and a higher frequency of strongly biased genes (including some transiently imprinted loci) [81]. These findings offer important clues to understand how dynamic allele-specific chromatin organization relates to gene expression during embryogenesis.

More recently, a genome-wide proportional relationship between differential contact propensity and differential expression across haplotypes and cell types has been quantitatively investigated at the chromatin loop level [44]. Greenwald et al. [44] determined that haplotype-associated differences in chromatin looping were considerably smaller than cell-type-associated differences, and both were smaller than the range at which gene expression and H3K27ac fold changes occurred, indicating that subtle changes in chromatin contact propensity were associated with large functional effects. Additionally, haplotype-associated loops that were highly enriched for imprinted regions, were related to CNVs but not eQTLs, suggesting that regulatory genetic variation does not exert large effects on chromatin interaction but might affect gene expression through small modifications of contact frequency at chromatin loops. Hence, contact propensity is a mechanism likely involved in gene regulation in a similar way to enhancer activity or transcription factor binding strength. Notably, contact propensity may either affect, or be affected by, gene expression and/or regulation.

Interestingly, allelic imbalance in Hi-C contacts is more relevant to nascent gene transcripts than to mature mRNAs [74]. Hence, the form-function relationship can be more accurately captured when comparing Hi-C data with nascent RNA Bru-seq (bromouridine labeling and sequencing) data, instead of steady-state RNA-seq (RNA sequencing) data. However, the majority of the aforementioned studies employed mRNA-seq data to analyze the interplay between chromatin contact imbalance and gene expression. This may lead to inaccurate conclusions and requires further in-depth investigation using Bru-seq data.

Imprinting

Haplotype-resolved 3D genome assays enable an understanding of the genetic basis of imprinting, which is essential for mammal development. At imprinted loci, the two alleles can differ drastically in transcriptional activity [82]. Therefore, imprinted regions dictate allele-specific gene activity [83], [84], which can only be studied at the haploid chromosome level.

Rao et al. [3] observed many homolog-specific features within diploid Hi-C maps, including imprinting-specific loops. The 3D structural differences between alleles at imprinted loci were also captured by Tan et al. [21]. Through visualizing the maternally transcribed H19 gene and the paternally transcribed IGF2 gene in simulated 3D nuclei of single cells, Tan et al. found that, despite cell-to-cell heterogeneity, the maternal allele more frequently separated IGF2 from both H19 and the nearby HIDAD site, and disrupted the IGF2-HIDAD CTCF loop; in contrast, the paternal allele more frequently stayed fully intermingled, which was consistent with the observed allelic differences found in contact profiles of bulk Hi-C [3]. More recently, Collombet et al. [81] found that maternally pre-formed domains in early mice embryos encompassed most transiently maternally imprinted genes (19 out of 27 genes). They also observed a similar pattern of shifting, from maternal imprinted domains at the two-cell stage to TADs at the blastocyst stage, for some imprinted genes, including Enc1, Mbnl2, Tle3, and Xist.

X chromosome inactivation (XCI)

Similar to imprinting, XCI means that one of the X chromosomes in female mammals is expressed exclusively from a single allele in somatic cells, which represents a striking case of homolog differences [82]. In the 3D space, the active X chromosome tends to exhibit an extended chromatin morphology, whereas the inactive X chromosome is generally compacted [21]. This pair of homologs is characterized by distinct patterning of compartments defined based on chromatin contacts: the active X chromosome features clear compartmentalization of euchromatin and heterochromatin; by contrast, the silent X chromosome contains a more uniform compartment patterning with multiple ‘super-domains’ or ‘super-loops’ [3], [21], [34], [42]. In mice, the inactive X chromosome shows a global loss of the TAD structure, except in genes that are being expressed [34], [81], [85], which corroborates the association between chromatin architecture and nuclear activity. Furthermore, a study on DNA organization dynamics related to XCI during early development suggested that loss of the TAD structure on the silenced X chromosome follows gene silencing, rather than preceding it [81].

Additionally, allele-biased 3D genome organization, which is dynamically coordinated with transcriptional activity, can facilitate the programming of replication-timing (RT) during development. Rivera-Mulia et al. [40] showed that RT asynchrony in mESCs was strongly correlated with changes in Hi-C compartments between alleles but not with sequence variation, gene expression, or chromatin accessibility. Long-range PEIs in RT asynchronous domains were restricted to the early replicating allele. Furthermore, the RT asynchrony and allelic discordance observed in the genome organization of mESCs were gradually lost during cell fate commitment. Dynamic allelic imbalances in chromatin organization also function directly in early mammalian development, which has also been investigated in several 3D diploid genome studies [41], [81], [86], [87], [88]. These studies depicted notably asymmetric chromatin architectures between alleles as well as significant reorganization during early development, partially explaining how haploid chromosomal conformation flexibly facilitates genomic functions during different biological processes.

Summary and outlook

Many fundamental and long-standing biological questions are linked to 3D genome architecture [59]. Haplotype-resolved Hi-C studies provide new ways to answer these questions. In this review, we described some of the emerging strategies and algorithms that are used to phase genome-scale Hi-C data, and summarized the highly structured chromatin organization of the diploid genome, especially homolog pairing. We also discussed the functional implications of the recently reported allele-biased chromosomal conformation, which attempts to explain how differences in chromosome architecture between alleles affects differences in genome function.

One of the main challenges of haplotype-resolved Hi-C studies is to improve the efficiency of data generation and processing. In the coming years, the integration of comprehensive studies on diploid genome conformation, multi-omics data, and live-cell imaging will enable the generation of a panorama of the genome. This trend necessitates derivation of comprehensive models that allow for integrating haplotype imputation and 4D Nucleome (combining 3D genome with time series) with gene expression and epigenetic signals to investigate chromatin form-function dynamics under different contexts (e.g., cold stress or elevated glucose) or during important biological processes (e.g., cell reprogramming and early development). The dissection of the functional significance of haploid chromosomal conformation, especially higher-order structural elements, with genetic perturbations will likely be predominant in the near future.

CRediT authorship contribution statement

Jing Li: Visualization, Funding acquisition. Yu Lin: Methodology. Qianzi Tang: Methodology, Writing - original draft. Mingzhou Li: Conceptualization, Writing - review & editing, Funding acquisition.

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.