Next Article in Journal
PCR-Based Equine Gene Doping Test for the Australian Horseracing Industry
Previous Article in Journal
Identification of IGF-1 Effects on White Adipose Tissue and Hippocampus in Alzheimer’s Disease Mice via Transcriptomic and Cellular Analysis
Previous Article in Special Issue
Transcriptomic Evidence of a Link between Cell Wall Biogenesis, Pathogenesis, and Vigor in Walnut Root and Trunk Diseases
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
IJMS
Volume 25
Issue 5
10.3390/ijms25052568
ijms-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Integrating Dynamic 3D Chromatin Architecture and Gene Expression Alterations Reveal Heterosis in Brassica rapa

by
Liu E
1,†,
Shanwu Lyu
2,†,
Yaolong Wang
1,
Dong Xiao
1,
Tongkun Liu
1,
Xilin Hou
1,
Ying Li
1 and
Changwei Zhang
1,*
1
National Key Laboratory of Crop Genetics & Germplasm Enhancement and Utilization, Nanjing Agricultural University, Nanjing 210095, China
2
Key Laboratory of South China Agricultural Plant Molecular Analysis and Genetic Improvement, Guangdong Provincial Key Laboratory of Applied Botany, South China Botanical Garden, Chinese Academy of Sciences, Guangzhou 510650, China
*
Author to whom correspondence should be addressed.
†
These authors contributed equally to this work.
Int. J. Mol. Sci. 2024, 25(5), 2568; https://doi.org/10.3390/ijms25052568
Submission received: 16 December 2023 / Revised: 25 January 2024 / Accepted: 9 February 2024 / Published: 22 February 2024
(This article belongs to the Special Issue Recent Analysis and Applications of Comparative and Functional Genomic on Horticultural Plants)
Browse Figures
Review Reports Versions Notes

Abstract

:
Heterosis plays a significant role in enhancing variety, boosting yield, and raising economic value in crops, but the molecular mechanism is still unclear. We analyzed the transcriptomes and 3D genomes of a hybrid (F1) and its parents (w30 and 082). The analysis of the expression revealed a total of 485 specially expressed genes (SEGs), 173 differentially expressed genes (DEGs) above the parental expression level, more actively expressed genes, and up-regulated DEGs in the F1. Further study revealed that the DEGs detected in the F1 and its parents were mainly involved in the response to auxin, plant hormone signal transduction, DNA metabolic process, purine metabolism, starch, and sucrose metabolism, which suggested that these biological processes may play a crucial role in the heterosis of Brassica rapa. The analysis of 3D genome data revealed that hybrid F1 plants tend to contain more transcriptionally active A chromatin compartments after hybridization. Supplementaryly, the F1 had a smaller TAD (topologically associated domain) genome length, but the number was the highest, and the expression change in activated TAD was higher than that of repressed TAD. More specific TAD boundaries were detected between the parents and F1. Subsequently, 140 DEGs with genomic structural variants were selected as potential candidate genes. We found two DEGs with consistent expression changes in A/B compartments and TADs. Our findings suggested that genomic structural variants, such as TADs and A/B chromatin compartments, may affect gene expression and contribute to heterosis in Brassica rapa. This study provides further insight into the molecular mechanism of heterosis in Brassica rapa.

1. Introduction

Heterosis is a widely observed biological phenomenon whereby hybrid offspring exhibit superior traits compared to their parents, including growth potential, viability, yield, and quality [1,2,3]. Given the increasing global population, the successful application of heterosis in agricultural production dramatically contributed to augmenting crop yields, mitigating food crises, and guaranteeing food security [4]. The concept of heterosis has been applied for more than 100 years [5]. At present, three hypotheses, namely dominance [6,7,8,9], over-dominance [9,10], and epistasis [11,12], were proposed to explain the genetic basis of heterosis in most crops. However, these hypotheses remain theoretical and cannot fully explain the phenomenon of heterosis. The QTL mapping methods developed in the late 1980s laid a foundation for exploring heterosis, and extensive genetic analyses of hybrids in maize, rice, rape, and other plants led to the identification of numerous quantitative trait loci [13]. With the development of high-throughput sequencing technology, a series of dominant genes related to heterosis traits were identified by integrating transcriptomics and other omics [14], thereby offering new insights into the molecular mechanisms underlying heterosis.
Genes encode genetic information in a one-dimensional linear sequence, which is expressed by forming a three-dimensional chromatin architecture. In addition, chromatin states and genome activity regulators also affect gene expression. Chromatin conformation plays a crucial role in executing biological activities via many regulatory elements [15]. Global and local chromatin rearrangements may occur upon sensing environmental and developmental cues, along with gene transcription changes [16,17]. Gene expression and regulatory network changes are entangled with heterosis [3]. In the interspecific hybridization between Arabidopsis thaliana and Arabidopsis lyata, the chromosomes from Arabidopsis thaliana in the hybrid are more compact, occupying a smaller nuclear volume and interacting more within the chromosomes [18]. These changes in the intensity of chromosomal interactions can affect gene expression. Studies have shown that the dynamic three-dimensional chromatin structure of Brassica napus also helps express plant hormone-related genes in hybrid F1 plants [19]. Systematically studying the frequency of interactions between genomic regions has led to the discovery of multiple principles of three-dimensional genome organization. Chromosome organization can be dissected into multiple functional and structural domains [20]. The units from high to low are chromosomal compartments (CTs), active (euchromatin, A-type) and inactive (heterochromatin, B-type) chromatin compartments, topologically associated domain (TADs), chromatin loops (CLs) and chromatin fibers. TAD regulates gene replication and epigenetic modification locally. The boundary region of TAD is genetically rich and conserved. The disruption of the TAD boundary can severely affect the regulation of gene expression and even lead to the occurrence of diseases [21]. The 3D genome and gene expression have complex relationships, with rapid dynamic changes in genome structure and gene expression accompanying hybridization. Nevertheless, current research on heterosis is mostly limited to one-dimensional linear sequences. Little is known about the changes in 3D structure for heterosis traits in Brassica rapa.
From appearance, Chinese cabbage (AA genome, 2n = 2x = 20) can be roughly divided into heading and non-heading Chinese cabbage. They are both highly nutritious, commercially viable, and easy to cultivate. They also have many differences, such as nutrient contentions and stress resistance. To exploit their advantages and explore the hybrid heterosis mechanism, we used non-heading Chinese cabbage 082 (P1) as the female parent and w30 (P2) as the male parent to generate the hybrid F1. We conducted a comparative analysis of Hi-C (high throughput chromatin conformation capture) and RNA-seq (RNA-sequencing) data, identifying several genes with specifically and significantly different expression levels in hybrids, and we discussed their potential functions and relationships. We also compared the 3D spatial structure of TADs and A/B compartments between the parents and hybrid F1 and observed the conserved and altered chromatin structures. We further investigated the associations between gene expression and the 3D chromatin structure. This study will contribute to the understanding of heterosis in Chinese cabbage, providing insight into achieving higher yields and better stress resistance through hybrid breeding.

2. Results

2.1. Heterosis Is Remarkable in the Hybrid F1

A remarkable biomass heterosis was observed in the hybrid F1 (see Figure 1A). To facilitate a quantitative comparison of their differences, we determined their fresh and dry weights at 21 days after sowing (Figure 1B). The values of MPH (mid-parent heterosis) and HPH (high-parent heterosis) both exceeded 50% (Figure 1C,D), and, in terms of fresh weight, MPH was as high as 159.28%. In comparison, HPH was as high as 150.11% (see Figure 1C). The F1 exhibited a remarkable degree of heterosis. The obvious phenotypic differences suited the subsequent integrative genomic analysis.

2.2. More Genes Were Actively Expressed in the Hybrid F1

High-quality transcriptome reference genomes of 124.57 Mb, 133.47 Mb, and 134.98 Mb were obtained for 082, w30, and F1 (Supplementary file 1: Tables S1 and S2). Genes with FPKM ≥ 1 were defined as actively expressed. The number of actively expressed genes in the hybrid F1 (24,601) was higher than those in the maternal line 082 (24,320) and paternal line w30 (24,524) (Supplementary file 2: Tables S3–S5 and Figure S2A). We speculate that the expression of more active genes drives the formation of heterosis. However, no distinct differences were observed based on the genome-wide FPKM expression profiles for the three varieties (Supplementary file 3: Table S6, Figures S1 and S2B). A pairwise analysis of these actively expressed genes further confirmed that many genes (more than 92%) were commonly expressed between the parents and that the male parent exhibited greater specificity than the female parent (Figure 2C). F1 had more genes that were specifically expressed when compared to either parent alone (Figure 2D,E). More genes between F1 and P2 (23,226) shared common expression profiles (Figure 2E). Concurrently, we compared the correlation of all expressed genes with nine samples and found that the correlation coefficient between F1 and P2 (R2 = 0.919) was higher than that between F1 and P1 (R2 = 0.888) (Figure 2G). Interestingly, when compared with both parents, F1 displayed the least specific expression of genes (485) compared to P1 (921) and P2 (891) (Figure 2F). We hypothesized that 485 genes with specific expressions in F1 may be involved in heterosis.

2.3. Transcriptome Analysis of DEGs

We analyzed the expression patterns and expression levels of differentially expressed genes (DEGs) to investigate the differential expression profiles of the genes. Firstly, we conducted a cluster analysis of the expression patterns. Some of the genes showed the same expression pattern between the parents and F1. In contrast, others showed significant differences (Figure 3A), suggesting that these genes may be responsible for the heterosis of the F1. Next, we analyzed DEG expression levels. A total of 2156 genes displayed significantly differential expression levels between the parents (p < 0.05, |log2FC| > 1) (Supplementary file 4: Table S8, Figure 3B, and Table 1), and the ratio between the maternal and paternal lines was similar. F1 and its parents differed significantly in 2568 genes, and most of those genes were expressed at higher levels in the hybrid F1 (Supplementary file 4: Tables S1, S8 and S9). Among these genes, we focused on 296 genes in F1 that had significantly different expression levels from both parents, and the expression level of most genes (173) in F1 was higher than that of its parents (Supplementary file 5: Tables S1 and S10). We used the volcano map to summarize the significant DEGs (Figure 3B). These results further confirmed that hybridization activated the expression of more genes in F1 and made the expression level of these genes in F1 significantly higher than in its parents (Table 1). In addition, a Venn diagram analysis was conducted on the DEGs of the three groups (082 vs. w30, F1 vs. 082, F1 vs. w30), and 102 shared DEGs were identified (Supplementary file 6: Table S11 and Figure S3C). Among them, ten homologous hyperparental genes were identified in Arabidopsis, which were enriched in developmental growth, stress and resistance, response to light, the nitrogen compound biosynthetic process, the carboxylic acid metabolic process, and the lipid biosynthetic process (Supplementary file 7: Table S12). These genes may play an essential role in the growth, development, resistance, and other advantages of hybrid F1 plants.

2.4. GO and KEGG Analysis of DEGs

We conducted a gene ontology (GO) and KEGG pathway enrichment analysis to annotate the DEGs [22,23,24]. Using an over-represented p-value threshold of less than 0.018, we identified significant changes in the GO terms between the F1 and its parents. The top GO terms for biological processes were as follows: “response to auxin”, “nucleoside transport”, and “nicotinamide riboside transport”, while the majority of DEGs were found in “DNA metabolic process” (78 DEGs). In most studies, the metabolic process was unanimously recognized as related to growth or biomass heterosis [25,26,27]. In terms of cell components, the most common GO term was “Golgi-associated vesicle”, while “Cytoplasmic membrane-bound vesicle” had the highest number of DEGs (11 DEGs). The most common GO term for molecular function was “ADP binding” with the most DEGs (Supplementary file 8: Tables S13 and S14 and Figure 4A). We also identified the top 20 pathways assigned to our DEGs (Figure 4B). Most of the DEGs were involved in “Plant hormone signal transduction” (17 DEGs), “Purine metabolism” (17 DEGs), and “Starch and sucrose metabolism” (12 DEGs) (Supplementary file 9: Tables S15 and S16). The “Plant hormone signal transduction” pathway helps the expression of plant hormone-related genes in F1, which is related to plant growth [19]. The main transport form of carbohydrates within plants is sucrose. The enrichment level of metabolic pathways in different tissues and development stages of different plants is different, and many studies suggested that carbohydrate metabolism is related to the formation of heterosis [24,28,29]. These findings suggest that hybridization results in a substantial proportion of transcriptome alterations, and more genes were activated than repressed (Figure S2).

2.5. Genome-Wide Interaction Matrices of the Parents and F1

To explore the dynamic 3D spatial structural changes of the P1, P2, and F1 genomes during B. rapa hybridization, we conducted Hi-C sequencing and obtained 117.8 G, 141.1 G, and 133.0 G of raw data, respectively (Supplementary file 10: Table S17). We then mapped the clean Hi-C data to the high-quality B. rapa (v2.5) reference genome [30] and filtered out invalid pairs for subsequent comparative 3D structural analysis (Supplementary file 10: Table S18). The genome-wide simulation images showed that the P1, P2, and F1 nuclei were nearly spherical, and the chromosomes were localized in a limited volume (Figure 5A, Supplementary file 15). The results confirmed previous findings that each chromosome occupies an exclusive region in the nucleus, a concept termed “chromosomal territory” [31]. At a resolution of 1 Mb, chromatin interaction analysis showed that the intra-chromosome interaction (cis interaction) of the three varieties was higher than the inter-chromosome interaction (trans interaction) (Figure S3). Chromosomal territory also revealed why intra-chromatin interactions were higher than inter-chromatin interactions. But F1 has a higher proportion of interchromosome/intrachromosome (trans/cis) contact than its parents (Supplementary file 10: Table S19 and Figure S4). Therefore, each chromosome in the F1 nucleus may have a different space. To further explore the interaction patterns within chromosomes, we enhanced the resolution of the genome-wide interaction matrix to 100 kb (Figure S8). For example, among the interactions within chromosome 8 (Figure 5B), the dark red diagonal indicated the strongest interaction. Horizontal distance decreased the occurrence of intra-chromosome interactions. Supplementaryly, more refined components, like TADs (Figure 5C), were available. Each diagonally distributed triangle corresponded to a topological association domain (TAD) (Figure 5E). Meanwhile, there was no significant difference in the distribution of compartments (Figure 5D) between the F1 and parents.

2.6. Identification of A/B Compartment Shifts

To explore the changes in A/B compartments between the F1 and its parents, each chromosome’s A/B compartment distribution of the F1 and its parents was displayed (Figure 6A). The F1 contained more A compartments, while the parents tended to contain more B compartments (Supplementary file 11: Table S20 and Figure 6B). A total of 473 genes from P1 to F1 were identified as A-to-B shifts, and 950 genes were identified as B-to-A shifts. A combined analysis with transcriptomics revealed 15 DEGs in A-to-B shifts (2 DEGs down-regulated in F1) and 20 DEGs in B-to-A shifts (16 DEGs up-regulated in F1). Meanwhile, a total of 281 genes from P2 to F1 were identified as A-to-B shifts, and 115 genes were B-to-A shifts. Moreover, 7 DEGs were in A-to-B shifts (1 DEGs down-regulated in F1) and 68 DEGs in B-to-A shifts (48 DEGs up-regulated in F1) (Supplementary file 11: Table S20). The inconsistent relationship resembled the results in Drosophila [32]. We found that hybridization can make it more inclined to transform into A compartments. In general, A compartments are in the transcriptional active region of euchromatin, and B compartments are in the transcriptional repressed region of heterochromatin. Therefore, the B-to-A compartment shift is considered beneficial to the formation of heterosis.
On this basis, we correlated the gene densities of 082-/w30-/F1 with A/B compartments. Results from F1 and its parents showed that the A compartment had a higher gene density than the B compartment (Figure S5). Further analysis of the gene expression of the A/B compartments found that there was no significant difference between the F1 and its parents. Nevertheless, the gene expression in the A compartment was significantly higher than that in the B compartment (Figure 6C). This was consistent with previous studies in Arabidopsis [33]. In summary, the A compartment, with its higher gene expression and higher gene density, may be a key factor in heterosis.

2.7. Changes in Compartments during Hybridization

During hybridization, 90.31% (A-to-A-to-A shifts are 72.96%, B-to-B-to-B shifts are 20.35%) of compartments did not change, suggesting that compartments in the nucleus are relatively stable. Furthermore, the F1 showed a higher proportion of compartment A, which increased by 4.04% (including 3.14% for A-to-B-to-A shifts and 0.90% for B-to-A-to-A shifts), while compartment B decreased by 1.44%, compared to the parents. Moreover, very few compartments underwent completely different transformations from their parent; of those that did, A-to-A-to-B shifts accounted for 0.06%, and B-to-B-to-A shifts accounted for 1.15% (Figure 2D). Heterochromatin remodeling is critical for a variety of cellular processes [34]. Therefore, the B-to-A compartment shift is considered beneficial to the formation of heterosis. Interestingly, 485 SEGs in the F1 showed similar results, but A-to-A-to-B shifts did not occur. When the parents’ compartments were not aligned, the compartments in the F1 were more inclined towards compartment A (Figure S6). The A/B compartment transformation analysis of 296 DEGs with significantly different expression levels in the parents and F1 was similar to that of the SEGs (Figure S6). These results suggested that different compartment activities between parents may influence compartment transitions in hybrids, ultimately impacting gene expression.

2.8. Identification of Different Kinds of Topologically Associating Domains

TADs, as relatively static physical and regulatory domains, facilitate specific and intentional gene expression programs in various ways [35,36]. To explore the 3D genomic differences between the parents and F1 more closely, we compared the change in TADs. The F1 shared similarities and differences with its parents regarding TADs (Figure S7). Further statistical analyses of the number and length of TADs showed that the number of TADs in the F1 was greater than that of its parents, and the average length of each TAD was smaller, as was the length of the whole genome (Supplementary file 12: Table S21 and Figure S7A). TADs are usually considered a relatively conservative structure, but the number and size of TADs change during cell differentiation. Therefore, we believe that changes in the number and size of TADs in hybrid F1 plants benefit heterosis.
Compared with P1, F1 exhibited 46 specific TADs and 105 conserved TADs (12 activated TADs, 10 repressed TADs, and 83 other TADs) (Figure 7B). Similarly, compared with P2, there were 91 conserved TADs (10 activated TADs, 9 repressed TADs, and 72 other TADs) and 53 specific TADs in F1 (Figure 7B). A total of 141 DEGs were recognized in activated TADs (53 up-regulated DEGs in P1/P2), and 43 DEGs were identified in repressed TADs (42 down-regulated DEGs in P1/P2) (Supplementary file 12: Table S22). Together with transcriptome analysis, some DEGs were not expressed consistently with TAD types, and some DEGs were not detected in active or repressed TADs. These observations suggested that heterosis formation may be largely due to trans-regulatory mechanisms. Figure 7C showed that activated TADs had a percentage of fold change (P1 or P2/F1) >1 than repressed TADs. Moreover, the expression change in activated TADs was higher than that of repressed TADs (Figure 7D). The correlation analysis of different boundary regions was carried out based on the distribution of DI values in the upstream and downstream windows. We found more specific TAD boundaries between the parents and F1 and less specific TAD boundaries between the parents (Figure S9). This may be due to the large-scale rearrangement of chromatin during hybridization. Indeed, when the TAD boundary is disrupted or reconstructed, it can promote the formation of new promoter–enhancer interactions, thereby altering gene expression. These results indicated that the changes in 3D chromatin structure between the parents and F1 were greater than those between the parents.

2.9. Identification of Candidate Genes for Heterosis and Verification of qRT-PCR

Based on the above transcriptome and 3D genome joint analysis results, we identified 140 DEGs associated with genomic structural alterations, including 56 DEGs in A/B compartments and 87 DEGs in TADs (Supplementary file 14: Table S24). Among them, three photosynthetic genes (BraA02003217, BraA07001020, and BraA07001021), eight plants cell size/division/cycle-related genes (BraA06003576, BraA05002844, BraA07001159, BraA03006134, BraA05001791, BraA02002099, BraA02003420, BraA02002337), four carbohydrate metabolism genes (BraA07001318, BraA06002248, BraA07001149, BraA01003773), nine resistance- and stress-related genes (BraA04000638, BraA07000080, BraA09003537, BraA09003725, BraA02003318, BraA07001184, BraA06002291, BraA07001349, BraA06002292), three genes related to the response to auxin (BraA02003420, BraA06000631, BraA07001182), one development and cell death gene (BraA05002866), and one senescence-associated gene (BraA08002232) were identified. For example, the P1 vs. F1 group used F1 as a reference. If the TAD type of P1 relative to F1 was a repressed TAD, then P1 relative to F1 should be down-regulated. Notably, two candidate DEGs (BraA02003296 and BraA02003297) were shared between these categories. To verify the accuracy of RNA-seq data, we selected 10 DEGs between the parents and F1 in 173 DEGs mentioned above by real-time quantitative PCR (qRT-RCR) (Supplementary file 13: Table S23). The expression patterns of these 10 DEGs were highly consistent with the data obtained by RNA-seq (Figure S10).

3. Discussion

The underlying molecular mechanism of heterosis remains elusive, and the use of three-dimensional structural dynamics on heterosis is still in its infancy. Therefore, integrating 3D genomics and transcriptomics is imperative for gaining a more comprehensive understanding of the role of chromatin in 3D space in heterosis.
Here, we present a study on the effects of hybridization on the 3D chromatin architecture and gene expression in Brassica rapa. The transcriptome data confirmed that hybridization leads to the activation of gene expression in F1 plants. A total of 485 SEGs, 173 DEGs, above parental expression level, more actively expressed genes, and up-regulated DEGs in F1 confirmed this conclusion. However, only a few genes were specifically expressed or displayed significantly different expression levels. These indicated that the genome-wide expression profiles of the hybrid F1 and parents were similar. Previous studies on heterosis produced by intraspecific hybridization also confirmed that most genes in hybrids and parents showed identical expression profiles, and the expression level of those genes was close to MPV [25,28,29]. In Brassica napus, Hu [19] also found that these hybrids had more up-regulated genes compared to the mid-parental value of gene expression, suggesting that altered gene expression patterns may influence their heterotic phenotypes. The DEGs detected in the F1 and its parents were mainly involved in the response to auxin, plant hormone signal transduction, the DNA metabolic process, purine metabolism, and starch and sucrose metabolism, suggesting that these biological processes may play a crucial role in the heterosis of B. rapa. However, the specific genes and pathways involved in heterosis may differ from species to species. In Brassica napus, the study found that hybrids with superior heterosis (better performance than both parents) had larger leaf sizes, which was attributed to increases in both cell size and cell number [19]. The differences in heterosis between species may be attributed to the specific genes and pathways involved, as well as the extent to which chromatin architecture influences gene expression. Further research is needed to fully elucidate the molecular mechanisms underlying heterosis and how they may vary between different Brassica species.
Although Hi-C technology has been widely used in many crops [37,38,39], little is known about genome architecture and its effect on B. rapa in heterosis. A preformed and stable topology (TAD) organizes the physical proximity between enhancers and their target genes. We compared the specific TAD and TAD boundary and the length and quantity of the TAD genomes. The results showed that the number of TADs in the F1 was greater than that of its parents, and the average length of each TAD was smaller, as was the length of the whole genome. This change results in the chromosomes of the hybrid being more compact than those of the parents, with higher interactions within the chromatin. Furthermore, we found that the F1 was more active than that between the chromosomes of both parents by calculating the interaction ratio of inter/intra. This may also be one of the reasons why the F1 had more active genes and up-regulated genes. Parents and hybrids exhibit significant differences in terms of TAD. More specific TAD boundaries between parents and the F1 were found. Changes in TAD boundaries can lead to gene rearrangement, thereby altering gene expression. It is precisely these changes in gene expression that give hybrid varieties heterosis traits. The chromatin structure of TADs is widely present in many crops, including rice and cotton [38,39]. However, there is no classical TAD interaction pattern in rice and Arabidopsis. Many species in the Brassica genus have been identified to have chromatin structures containing TADs, including B. napus, B. rapa, and B. oleracea. This study utilized three-dimensional spatial structural data to explore the impact of chromatin structural variation on heterosis traits, providing us with a new research perspective.
Hybridization allows F1 plants more A compartments, resulting in more high-density and highly expressed genes. In Brassica napus, F1 hybrids with superior heterosis had more transcriptionally active A compartments in their chromatin compared to those with inferior heterosis [19]. Although there are no obvious differences in A/B compartments, the change in TADs may be one of the reasons why hybridization activates more genes in F1 plants. Subsequently, 140 DEGs with genomic structural variants were selected as potential candidate genes. Their functions were related to various processes, such as photosynthesis, plant cell size/division/cycle, resistance and stress, response to auxin, and carbohydrate metabolism. However, further functional verification is required for these candidate DEGs. In the future, we will obtain a mutant of the candidate DEGs through gene editing or other methods, observe its corresponding phenotype, and then obtain overexpressed plants to further determine the phenotype. In addition, the expression pattern analysis of the candidate DEGs is also needed (for example, analyzing expression levels in different tissues at different stages or expression sites in cells). If possible, we will further investigate interacting proteins and the regulatory relationships between upstream and downstream.
Gene expression is precisely regulated by the multi-layered three-dimensional structure of chromatin [36]. Different layers of the 3D genome have various levels of regulatory control [40]. This study specifically focused on the beneficial aspects of the 3D genome. However, further exploration is needed to understand the specific regulations. With the advancement of technology, it is believed that the impact of the spatial structure changes of chromatin on heterosis will be able to be further explored and interpreted in the future.

4. Materials and Methods

4.1. Plant Materials, Growth Conditions, and Sample Collection

The male parent w30 (P2), female parent 082 (P1), and their hybrid F1 used in this study were all provided by Nanjing Agricultural University (Nanjing, China). The F1 was generated by hand pollination. To reduce bias, all plant materials used for Hi-C and RNA-seq were grown in the same environment (24 °C, 16 h of light/8 h of darkness, 75% RH). Samples for Hi-C and RNA-seq were collected using fully expanded parental and hybrid F1 leaves (3rd upper leaf of a plant) with three replicates per sample. After being quick-frozen in liquid nitrogen, they were placed in a −80 °C refrigerator for later use.

4.2. Evaluation of Heterosis

To measure the phenotypic data, including the fresh weight (FW) and dry weight (DW), 30 similar plants (10 plants per replicate) were collected at 21 days. The collected material was washed with distilled water immediately and then baked at 105 °C for 30 min before further drying at 80 °C for 60 h [19]. MPH (mid-parent heterosis) was calculated as follows: MPH = (the mean value of hybrid F1—the average value of both parents)/the average value of both parents × 100%. HPH (high-parent heterosis) was calculated as HPH = (the mean value of hybrid F1—the optimal parental value)/the optimal parental value × 100%.

4.3. RNA-Seq Experiment and Sequencing

Per the kit’s instructions, the total RNA of nine samples was extracted using the RNAprep Pure Plant Kit (Tiangen, Beijing, China). In addition, the construction and sequencing of cDNA libraries were also completed by the Novogene Company (Beijing, China). The construction of cDNA libraries included the synthesis, purification, and end repair of double-stranded cDNA, adding A-tail, ligation of sequencing adapter, PCR enrichment, and so on. After the completed libraries were qualified, different libraries were pooled according to the requirements of an effective concentration and target data volume, and finally, RNA sequencing was performed.

4.4. Hi-C Libraries Construction and Sequencing

The plant leaves were ground in liquid nitrogen, as described in Section 4.1 of this study. Then, 2 g of uniform powder was taken to establish a Hi-C library. The establishment of Hi-C libraries and completion of Illumina sequencing were conducted by the Novogene Bioinformatics technology company (Beijing, China). The libraries were built according to the previous standard protocol, with some modifications [41]. After fixing it with polyformaldehyde and DNA, restriction endonuclease was used to make a gap at the cross-linking point, and a biotin marker was used simultaneously. The adjacent DNA fragments were treated with nucleic acid ligase, then the protein at the junction was digested with protease. Finally, the fragments with a length of 350 bp were broken by the Covaris crusher for recovery. After the libraries were constructed, Qubit 2.0 and Q-PCR were used for preliminary and precise quantification, respectively, while Agilent 2100 (Agilent Technology Co., Ltd.; USA) was used to detect the insertion size of the library. Finally, the different libraries were pooled and sequenced once they had been qualified.

4.5. RNA-Seq Data Analysis

Spliced reads were effectively compared to RNA-Seq data using HISAT 2 (version 2.1.0) [42]. In the RNA-Seq analysis, to make different genes and sequencing data comparable, we introduced the concept of FPKM. FPKM represents the number of fragments per million from a gene per kilobase length [43] and utilizes DESeq2 R (version 1.18.0) software for differences in gene analysis (DEG), selection criteria for |log2 fold change| > 1, or error rate (FDR) < 0.05. The corrected p-value was obtained by multiple hypothesis testing based on Benjamini and Hochberg’s method, and the square of the Pearson correlation coefficient between biological repetitions (R2 > 0.8) was used to test the experimental reliability. The Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases (www.kegg.jp/kegg/kegg1.html. accessed on 4 August 2023.) [22,23,24] were used to identify possible biological functions and pathways of DEGs. Among them, the GOseq (Release 2.12) software package was used to conduct GO enrichment analysis of the DEGs in R [44]. Cluster Profiler R package [45] was used to study the enrichment of DEGs in the KEGG pathway.

4.6. Hi-C Read Mapping

The quality of the obtained Hi-C sequencing data was controlled before output. Paired reads, including adapter-contaminated sequences, unknown nucleotide “N” ratios >10%, and more than 50% base Q < 5 were filtered out. Next, the qualified reads obtained were processed using the HiCUP pipeline (version 0.57) [46]. For data comparison, Bowtie2 software (version 2.2.3) [47] was used to compare the obtained reads with the B. rapa genome (version 2.5) [30]. The observed interaction matrices were constructed for the final influential contacts according to the statistical interaction matrices with a certain resolution interaction matrix. The maximum likelihood method was used to standardize the observation interaction matrix. The MDS algorithm of the PASTIS software (Release 0.3.3) [48] was employed to imitate the 3D position of chromatin. On this basis, the constructed heat map of chromatin interactions was divided into A/B compartments with a resolution of 100 kb by principal component analysis (PCA). TadLib (hitad 0.1.1-r1) software was applied to estimate the TAD topology, with default parameters at a 40 Kb resolution. RNA-seq was performed on the genes in TADs, and we classified the TADs into three categories based on the proportion of genes with positive and negative FC in each TAD of P1/P2 vs. F1: activated TADs, repressed TADs, and other TADs. We sorted the foldchange from largest to smallest, with top 10% being considered activated TADs, bottom 10% being expressed TADs, and the rest being considered others TADs.

4.7. Quantitative Reverse-Transcription PCR

QRT-PCR was used to perform RNA-seq verification on 10 randomly selected DEGs. In the qRT-PCR, cDNA was derived from the reverse transcription of three biologically repetitive RNAs of the parents and hybrids. Total RNA extraction was carried out according to the instructions of the RNAprep Plant Kit (DP419; Tiangen, Beijing, China). For reverse transcription, an Evo M-MLV RT Mix Kit with gDNA Clean [AG11728, Accurate Biotechnology Co., Ltd. Changsha, China] was used, the QuantStudio Q3 instrument was used for quantitative analysis, and Excel was used for mapping and data analysis. The primers of corresponding candidate DEGs are listed in Supplementary File S13. BrActin (Bra028615) was used as a quantitative reference [49]. The 2−ΔΔCt approach was applied to quantify the relative gene expression levels [50].

5. Conclusions

Our study revealed the gene expression profiles and changes in the 3D chromatin structure of B. rapa hybrids and their parents. Several gene regulatory networks may play a crucial role in the heterosis of B. rapa, including those involved in “plant hormone signal transduction” and “starch and sucrose metabolism”. The identification of 173 genes with higher expression in the F1 and 485 genes with specific expression in the F1 may contribute to the generation of heterosis traits. The differences in heterosis between species may be attributed to the specific genes and pathways involved, as well as the extent to which chromatin architecture influences gene expression. Further research is needed to fully elucidate the molecular mechanisms underlying heterosis and how they may vary between different Brassica species.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ijms25052568/s1.

Author Contributions

L.E., S.L., C.Z., Y.W., X.H., Y.L., T.L. and D.X. designed this study. L.E. and S.L. performed the experiments. L.E., S.L. and C.Z. prepared the article. All authors have read and agreed to the published version of the manuscript.

Funding

This study was financially supported by the National Key R&D plan (2022YFD1200502) and the Natural Science Foundation of Jiangsu Province (BK20200560).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The datasets supporting the results of the present study are included within this article (and its Supplementary files). All sequencing data generated for this study have been submitted to the NCBI Sequence Read Archive under accession number (PRJNA970959).

Conflicts of Interest

The authors declare no conflicts of interest.

Abbreviations

DEGsDifferentially expressed genes;
FPKMFragments per kilobase of transcript sequence per million base pairs;
GOGene ontology;
MPVMid-parent value;
qRT-PCRQuantitative real-time PCR;
QTLQuantitative trait locus;
Hi-CHigh throughput chromatin conformation capture;
RNA-seqRNA-sequencing data.

References

  1. Birchler, J.A.; Yao, H.; Chudalayandi, S.; Vaiman, D.; Veitia, R.A. Heterosis. Plant Cell 2010, 22, 2105–2112. [Google Scholar] [CrossRef]
  2. Hochholdinger, F.; Hoecker, N. Towards the molecular basis of heterosis. Trends Plant Sci. 2007, 12, 427–432. [Google Scholar] [CrossRef]
  3. Birchler, J.A.; Auger, D.L.; Riddle, N.C. In search of the molecular basis of heterosis. Plant Cell 2003, 15, 2236–2239. [Google Scholar] [CrossRef]
  4. Hochholdinger, F.; Baldauf, J.A. Heterosis in plants. Curr. Biol. 2018, 28, R1089–R1092. [Google Scholar] [CrossRef] [PubMed]
  5. Liu, J.; Li, M.; Zhang, Q.; Wei, X.; Huang, X. Exploring the molecular basis of heterosis for plant breeding. J. Integr. Plant Biol. 2020, 62, 287–298. [Google Scholar] [CrossRef]
  6. Jones, D.F. Dominance of Linked Factors as a Means of Accounting for Heterosis. Proc. Natl. Acad. Sci. USA 1917, 3, 310–312. [Google Scholar] [CrossRef] [PubMed]
  7. Bruce, A.B. The Mendelian theory of heredity and the augmentation of vigor. Science 1910, 32, 627–628. [Google Scholar] [CrossRef]
  8. Davenport, C.B. Degeneration, albinism, and inbreeding. Science 1908, 28, 454–455. [Google Scholar] [CrossRef]
  9. Crow, J.F. 90 years ago: The beginning of hybrid maize. Genetics 1998, 148, 923–928. [Google Scholar] [CrossRef] [PubMed]
  10. East, E.M. Heterosis. Genetics 1936, 21, 375–397. [Google Scholar] [CrossRef] [PubMed]
  11. Williams, W. Heterosis and the genetics of complex characters. Nature 1959, 184, 527–530. [Google Scholar] [CrossRef]
  12. Richey, F.D. Mock-dominance and hybrid vigor. Science 1942, 96, 280–281. [Google Scholar] [CrossRef]
  13. Fujimoto, R.; Uezone, K.; Ishikura, S.; Osabe, K.; Peacoc, W.J.; Dennis, E.S. Recent research on the mechanism of heterosis is important for crop and vegetable breeding systems. Breed. Sci. 2018, 68, 145–158. [Google Scholar] [CrossRef] [PubMed]
  14. Chen, Z.J. Genomic and epigenetic insights into the molecular bases of heterosis. Nat. Rev. Genet. 2013, 14, 471–482. [Google Scholar] [CrossRef] [PubMed]
  15. Dixon, J.R.; Gorkin, D.U.; Ren, B. Chromatin Domains: The Unit of Chromosome Organization. Mol. Cell 2016, 62, 668–680. [Google Scholar] [CrossRef]
  16. Probst, A.V.; Scheid, O.M. Stress-induced structural changes in plant chromatin. Curr. Opin. Plant Biol. 2015, 27, 8–16. [Google Scholar] [CrossRef] [PubMed]
  17. Rosa, S.; De Lucia, F.; Mylne, J.S.; Zhu, D.; Ohmido, N.; Pendle, A.; Kato, N.; Shaw, P.; Dean, C. Physical clustering of FLC alleles during Polycomb-mediated epigenetic silencing in vernalization. Genes Dev. 2013, 27, 1845–1850. [Google Scholar] [CrossRef]
  18. Ouyang, W.; Xiong, D.; Li, G.; Li, X. Unraveling the 3D Genome Architecture in Plants: Present and Future. Mol. Plant 2020, 13, 1676–1693. [Google Scholar] [CrossRef]
  19. Hu, Y.; Xiong, J.; Shalby, N.; Zhuo, C.; Jia, Y.; Yang, Q.-Y.; Tu, J. Comparison of dynamic 3D chromatin architecture uncovers heterosis for leaf size in Brassica napus. J. Adv. Res. 2022, 42, 289–301. [Google Scholar] [CrossRef]
  20. Sexton, T.; Cavalli, G. The Role of Chromosome Domains in Shaping the Functional Genome. Cell 2015, 160, 1049–1059. [Google Scholar] [CrossRef]
  21. Lupianez, D.G.; Kraft, K.; Heinrich, V.; Krawitz, P.; Brancati, F.; Klopocki, E.; Hom, D.; Kayserili, H.; Opitz, J.M.; Laxova, R.; et al. Disruptions of Topological Chromatin Domains Cause Pathogenic Rewiring of Gene-Enhancer Interactions. Cell 2015, 161, 1012–1025. [Google Scholar] [CrossRef] [PubMed]
  22. Kanehisa, M.; Furumichi, M.; Sato, Y.; Kawashima, M.; Ishiguro-Watanabe, M. KEGG for taxonomy-based analysis of pathways and genomes. Nucleic Acids Res. 2022, 51, D587–D592. [Google Scholar] [CrossRef] [PubMed]
  23. Kanehisa, M.; Goto, S. KEGG: Kyoto Encyclopedia of Genes and Genomes. Nucleic Acids Res. 2000, 28, 27–30. [Google Scholar] [CrossRef] [PubMed]
  24. Kanehisa, M. Toward understanding the origin and evolution of cellular organisms. Protein Sci. 2019, 28, 1947–1951. [Google Scholar] [CrossRef] [PubMed]
  25. Song, G.-S.; Zhai, H.-L.; Peng, Y.-G.; Zhang, L.; Wei, G.; Chen, X.-Y.; Xiao, Y.-G.; Wang, L.; Chen, Y.-J.; Wu, B.; et al. Comparative Transcriptional Profiling and Preliminary Study on Heterosis Mechanism of Super-Hybrid Rice. Mol. Plant 2010, 3, 1012–1025. [Google Scholar] [CrossRef] [PubMed]
  26. Gaertner, T.; Steinfath, M.; Andorf, S.; Lisec, J.; Meyer, R.C.; Altmann, T.; Willmitzer, L.; Selbig, J. Improved Heterosis Prediction by Combining Information on DNA- and Metabolic Markers. PLoS ONE 2009, 4, e5220. [Google Scholar] [CrossRef]
  27. Lisec, J.; Steinfath, M.; Meyer, R.C.; Selbig, J.; Melchinger, A.E.; Willmitzer, L.; Altmann, T. Identification of heterotic metabolite QTL in Arabidopsis thaliana RIL and IL populations. Plant J. 2009, 59, 777–788. [Google Scholar] [CrossRef]
  28. Wei, G.; Tao, Y.; Liu, G.; Chen, C.; Luo, R.; Xia, H.; Gan, Q.; Zeng, H.; Lu, Z.; Han, Y.; et al. A transcriptomic analysis of superhybrid rice LYP9 and its parents. Proc. Natl. Acad. Sci. USA 2009, 106, 7695–7701. [Google Scholar] [CrossRef]
  29. Hu, X.; Wang, H.; Diao, X.; Liu, Z.; Li, K.; Wu, Y.; Liang, Q.; Wang, H.; Huang, C. Transcriptome profiling and comparison of maize ear heterosis during the spikelet and floret differentiation stages. BMC Genom. 2016, 17, 959. [Google Scholar] [CrossRef]
  30. Cai, C.; Wang, X.; Liu, B.; Wu, J.; Liang, J.; Cui, Y.; Cheng, F.; Wang, X. Brassica rapa Genome 2.0: A Reference Upgrade through Sequence Re-assembly and Gene Re-annotation. Mol. Plant 2017, 10, 649–651. [Google Scholar] [CrossRef]
  31. Tanabe, H.; Muller, S.; Neusser, M.; von Hase, J.; Calcagno, E.; Cremer, M.; Solovei, I.; Cremer, C.; Cremer, T. Evolutionary conservation of chromosome territory arrangements in cell nuclei from higher primates. Proc. Natl. Acad. Sci. USA 2002, 99, 4424–4429. [Google Scholar] [CrossRef] [PubMed]
  32. Ghavi-Helm, Y.; Jankowski, A.; Meiers, S.; Viales, R.R.; Korbel, J.O.; Furlong, E.E.M. Highly rearranged chromosomes reveal uncoupling between genome topology and gene expression. Nat. Genet. 2019, 51, 1272–1282. [Google Scholar] [CrossRef] [PubMed]
  33. Grob, S.; Schmid, M.W.; Grossniklaus, U. Hi-C Analysis in Arabidopsis Identifies the KNOT, a Structure with Similarities to the flamenco Locus of Drosophila. Mol. Cell 2014, 55, 678–693. [Google Scholar] [CrossRef]
  34. Zhang, X.; Liu, X.; Du, Z.; Wei, L.; Fang, H.; Dong, Q.; Niu, J.; Li, Y.; Gao, J.; Zhang, M.Q.; et al. The loss of heterochromatin is associated with multiscale three-dimensional genome reorganization and aberrant transcription during cellular senescence. Genome Res. 2021, 31, 1121–1135. [Google Scholar] [CrossRef]
  35. Gibcus, J.H.; Dekker, J. The Hierarchy of the 3D Genome. Mol. Cell 2013, 49, 773–782. [Google Scholar] [CrossRef] [PubMed]
  36. Bonev, B.; Cavalli, G. Organization and function of the 3D genome. Nat. Rev. Genet. 2016, 17, 661. [Google Scholar] [CrossRef]
  37. Xie, T.; Zhang, F.-G.; Zhang, H.-Y.; Wang, X.-T.; Hu, J.-H.; Wu, X.-M. Biased gene retention during diploidization in Brassica linked to three-dimensional genome organization. Nat. Plants 2019, 5, 822–832. [Google Scholar] [CrossRef]
  38. Wang, M.; Wang, P.; Lin, M.; Ye, Z.; Li, G.; Tu, L.; Shen, C.; Li, J.; Yang, Q.; Zhang, X. Evolutionary dynamics of 3D genome architecture following polyploidization in cotton. Nat. Plants 2018, 4, 90–97. [Google Scholar] [CrossRef]
  39. Dong, Q.; Li, N.; Li, X.; Yuan, Z.; Xie, D.; Wang, X.; Li, J.; Yu, Y.; Wang, J.; Ding, B.; et al. Genome-wide Hi-C analysis reveals extensive hierarchical chromatin interactions in rice. Plant J. 2018, 94, 1141–1156. [Google Scholar] [CrossRef]
  40. Kong, S.; Zhang, Y. Deciphering Hi-C: From 3D genome to function. Cell Biol. Toxicol. 2019, 35, 15–32. [Google Scholar] [CrossRef]
  41. Van Berkum, N.L.; Lieberman-Aiden, E.; Williams, L.; Imakaev, M.; Gnirke, A.; Mirny, L.A.; Dekker, J.; Lander, E.S. Hi-C: A method to study the three-dimensional architecture of genomes. J. Vis. Exp. 2010, 39, e1869. [Google Scholar]
  42. Kim, D.; Langmead, B.; Salzberg, S.L. HISAT: A fast spliced aligner with low memory requirements. Nat. Methods 2015, 12, 357–360. [Google Scholar] [CrossRef] [PubMed]
  43. Trapnell, C.; Roberts, A.; Goff, L.; Pertea, G.; Kim, D.; Kelley, D.R.; Pimentel, H.; Salzberg, S.L.; Rinn, J.L.; Pachter, L. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat. Protoc. 2012, 7, 562–578. [Google Scholar] [CrossRef] [PubMed]
  44. Young, M.D.; Wakefield, M.J.; Smyth, G.K.; Oshlack, A. Gene ontology analysis for RNA-seq: Accounting for selection bias. Genome Biol. 2010, 11, R14. [Google Scholar] [CrossRef] [PubMed]
  45. Yu, G.; Wang, L.-G.; Han, Y.; He, Q.-Y. clusterProfiler: An R Package for Comparing Biological Themes Among Gene Clusters. Omics A J. Integr. Biol. 2012, 16, 284–287. [Google Scholar] [CrossRef]
  46. Wingett, S.; Ewels, P.; Furlan-Magaril, M.; Nagano, T.; Schoenfelder, S.; Fraser, P.; Andrews, S. HiCUP: Pipeline for mapping and processing Hi-C data. F1000Research 2015, 4, 1310. [Google Scholar] [CrossRef]
  47. Langmead, B.; Salzberg, S.L. Fast gapped-read alignment with Bowtie 2. Nat. Methods 2012, 9, 357–359. [Google Scholar] [CrossRef]
  48. Szalaj, P.; Tang, Z.; Michalski, P.; Pietal, M.J.; Luo, O.J.; Sadowski, M.; Li, X.; Radew, K.; Ruan, Y.; Plewczynski, D. An integrated 3-Dimensional Genome Modeling Engine for data-driven simulation of spatial genome organization. Genome Res. 2016, 26, 1697–1709. [Google Scholar] [CrossRef]
  49. Dheda, K.; Huggett, J.F.; Bustin, S.A.; Johnson, M.A.; Rook, G.; Zumla, A. Validation of housekeeping genes for normalizing RNA expression in real-time PCR. Biotechniques 2004, 37, 112–119. [Google Scholar] [CrossRef]
  50. Schmittgen, T.D.; Livak, K.J. Analyzing real-time PCR data by the comparative C-T method. Nat. Protoc. 2008, 3, 1101–1108. [Google Scholar] [CrossRef]
Ijms 25 02568 g001
Figure 1. Comparisons of heterosis in w30(P2)-/082(P1)-/F1 traid. (A). Phenotypes of the hybrid F1 and its parents in Brassica rapa. (B). W30(P2)-/082(P1)-/F1 traid at 21 DAS. (C). The total fresh weight of the w30(P2)-/082(P1)-/F1 traid. (D). The total dry weight of the w30(P2)-/082(P1)-/F1 traid. MPH indicated the mid-parent heterosis in the w30(P2)-/082(P1)-/F1 traid, while HPH indicated over-parent heterosis. Values with different letters differ significantly (p < 0.05).
Figure 1. Comparisons of heterosis in w30(P2)-/082(P1)-/F1 traid. (A). Phenotypes of the hybrid F1 and its parents in Brassica rapa. (B). W30(P2)-/082(P1)-/F1 traid at 21 DAS. (C). The total fresh weight of the w30(P2)-/082(P1)-/F1 traid. (D). The total dry weight of the w30(P2)-/082(P1)-/F1 traid. MPH indicated the mid-parent heterosis in the w30(P2)-/082(P1)-/F1 traid, while HPH indicated over-parent heterosis. Values with different letters differ significantly (p < 0.05).
Ijms 25 02568 g001
Ijms 25 02568 g002
Figure 2. Comparative analysis of the actively expressed genes detected in w30(P2)-/082(P1)-/F1. (A). Number of actively expressed genes in F1 and its parents. (B). Genome-wide expression profile of w30(P2)-/082(P1)-/F1 traid. (C–F). Venn diagrams of the actively expressed genes drawn by comparative analysis. (G). The expression profiles of all detected genes were used to calculate the correlation of nine samples.
Figure 2. Comparative analysis of the actively expressed genes detected in w30(P2)-/082(P1)-/F1. (A). Number of actively expressed genes in F1 and its parents. (B). Genome-wide expression profile of w30(P2)-/082(P1)-/F1 traid. (C–F). Venn diagrams of the actively expressed genes drawn by comparative analysis. (G). The expression profiles of all detected genes were used to calculate the correlation of nine samples.
Ijms 25 02568 g002
Ijms 25 02568 g003
Figure 3. Overview of the DEGs of the hybrid F1 and its parents. (A). Cluster analysis of differentially expressed genes for w30(P2)−/082(P1)−/F1. The heat maps were constructed based on the log2 (fold−change in the normalized expression levels) of two arbitrary samples in w30(P2)−/082(P1)−/F1. The color scale represents the log2 (fold−change in the normalized expression levels) of two arbitrary samples, with blue denoting low expression and red denoting high expression. (B). Volcano map of DEGs. Red points denote up−regulated genes; green points denote down−regulated genes; blue points denote non−differentiated genes. Each group is referenced by the one after ‘vs.’, such as P1 vs. P2 is referenced by P2. (C). Venn diagrams of the differentially expressed genes were drawn by comparative analysis of the A, B, and C groups. Group A represents P1 vs. P2, group B represents P1 vs. F1, and group C represents P1 vs. P2.
Figure 3. Overview of the DEGs of the hybrid F1 and its parents. (A). Cluster analysis of differentially expressed genes for w30(P2)−/082(P1)−/F1. The heat maps were constructed based on the log2 (fold−change in the normalized expression levels) of two arbitrary samples in w30(P2)−/082(P1)−/F1. The color scale represents the log2 (fold−change in the normalized expression levels) of two arbitrary samples, with blue denoting low expression and red denoting high expression. (B). Volcano map of DEGs. Red points denote up−regulated genes; green points denote down−regulated genes; blue points denote non−differentiated genes. Each group is referenced by the one after ‘vs.’, such as P1 vs. P2 is referenced by P2. (C). Venn diagrams of the differentially expressed genes were drawn by comparative analysis of the A, B, and C groups. Group A represents P1 vs. P2, group B represents P1 vs. F1, and group C represents P1 vs. P2.
Ijms 25 02568 g003
Ijms 25 02568 g004
Figure 4. Analysis of diferentially expressed genes (DEGs) between the F1 and its parents. (A). The most enriched GO terms with all DEGs. The ordinate is the enriched GO term, and the abscissa is the number of DEGs. (B). The top 20 enriched KEGG pathways of DEGs. The pathway label is on the vertical axis, the Rich factor is on the horizontal axis, the size of the dots represents the number of DEGs in the pathway, and the color of the dots represents the distinct Q-value levels. Each group is referenced by the one after ‘vs.’; for example, F1 references P1 vs. F1. Kyoto Encyclopedia of Genes and Genomes (KEGG) databases (www.kegg.jp/kegg/kegg1.html. accessed on 4 August 2023.) were used.
Figure 4. Analysis of diferentially expressed genes (DEGs) between the F1 and its parents. (A). The most enriched GO terms with all DEGs. The ordinate is the enriched GO term, and the abscissa is the number of DEGs. (B). The top 20 enriched KEGG pathways of DEGs. The pathway label is on the vertical axis, the Rich factor is on the horizontal axis, the size of the dots represents the number of DEGs in the pathway, and the color of the dots represents the distinct Q-value levels. Each group is referenced by the one after ‘vs.’; for example, F1 references P1 vs. F1. Kyoto Encyclopedia of Genes and Genomes (KEGG) databases (www.kegg.jp/kegg/kegg1.html. accessed on 4 August 2023.) were used.
Ijms 25 02568 g004
Ijms 25 02568 g005
Figure 5. Hi-C analysis of chromatin contacts in F1 and parents on chromatin 8. The x-axis represents the chromosome. (A). Three-dimensional model of whole chromosomes. A different hue represents each chromosome. (B). Intrachromosomal interactions of the chromosome at 100 Kb resolution. The Y-axis represents log2 (interaction matrix). (C). Each triangle, distributed diagonally, is represented as a topologically associated domain (TAD). The y-axis represents Triangular. (D). The first principal component values show A/B compartment status. (E). Distribution of TADs. The y-axis represents domains.
Figure 5. Hi-C analysis of chromatin contacts in F1 and parents on chromatin 8. The x-axis represents the chromosome. (A). Three-dimensional model of whole chromosomes. A different hue represents each chromosome. (B). Intrachromosomal interactions of the chromosome at 100 Kb resolution. The Y-axis represents log2 (interaction matrix). (C). Each triangle, distributed diagonally, is represented as a topologically associated domain (TAD). The y-axis represents Triangular. (D). The first principal component values show A/B compartment status. (E). Distribution of TADs. The y-axis represents domains.
Ijms 25 02568 g005
Ijms 25 02568 g006
Figure 6. The analysis of the A/B compartments in the F1 and its parents. (A). The distributions of the A/B compartment on each chromosome. With the y−axis 0 scale line as the reference, above is Compartment A, and below is Compartment B. Different colors correspond to different chromosomes. (B). The number of genes contained in the A and B compartments. (C). Box plots representing the gene expression of the A and B compartments. The Y−axis represents log2 (FPKM + 1 of genes). (D). The percentage of the A/B compartment shifts between the F1 and its parents.
Figure 6. The analysis of the A/B compartments in the F1 and its parents. (A). The distributions of the A/B compartment on each chromosome. With the y−axis 0 scale line as the reference, above is Compartment A, and below is Compartment B. Different colors correspond to different chromosomes. (B). The number of genes contained in the A and B compartments. (C). Box plots representing the gene expression of the A and B compartments. The Y−axis represents log2 (FPKM + 1 of genes). (D). The percentage of the A/B compartment shifts between the F1 and its parents.
Ijms 25 02568 g006
Ijms 25 02568 g007
Figure 7. The analysis of the TADs in the F1 and parents. (A). The TAD length genomic distribution of the F1 and parents. From left to right: 082, F1, and w30. The Y−axis represents size in Mb. The X−axis represents the chromosome. (B). P1 vs. F1 Venn diagrams of the TAD. The overlapping part represents conservative TADs, and vice versa represents specific TADs. (C). Proportional stacking maps of genes with different expression levels in different TADs. For each gene contained in the TADs, the change in expression level was calculated for both sets of samples, and the percentage of all genes corresponding to that TAD with a fold change greater than 1. The percentages were divided into four categories, including 0–25%, 25–50%, 50–75%, and 75–100%. The results for the different categories of TADs were presented as bar charts. (D). Gene expression of the three types of TADs (repressed TADs/other TADs/activated TADs). The Y−axis represents expression change (log2).
Figure 7. The analysis of the TADs in the F1 and parents. (A). The TAD length genomic distribution of the F1 and parents. From left to right: 082, F1, and w30. The Y−axis represents size in Mb. The X−axis represents the chromosome. (B). P1 vs. F1 Venn diagrams of the TAD. The overlapping part represents conservative TADs, and vice versa represents specific TADs. (C). Proportional stacking maps of genes with different expression levels in different TADs. For each gene contained in the TADs, the change in expression level was calculated for both sets of samples, and the percentage of all genes corresponding to that TAD with a fold change greater than 1. The percentages were divided into four categories, including 0–25%, 25–50%, 50–75%, and 75–100%. The results for the different categories of TADs were presented as bar charts. (D). Gene expression of the three types of TADs (repressed TADs/other TADs/activated TADs). The Y−axis represents expression change (log2).
Ijms 25 02568 g007
Table 1. Differentially expressed genes detected in the hybrid F1 and their parents.
Table 1. Differentially expressed genes detected in the hybrid F1 and their parents.
Hybrid Set
Samples		UpDownB2PTotal
082 vs. w30888 (41.19%)1268 (58.81%) 2156
F1 vs. 082885 (80.38%)216 (19.62%) 1101
F1 vs. w301115 (63.24%)648 (36.76%) 1763
F1 vs. (082 and w30)173 (58.45%)32 (10.81%)91 (30.74%)296
Note: DGPP indicated the differentially expressed genes (DEGs) between two parents. DGHP indicated the DEGs between the hybrid and parents. B2P indicated the expression levels of the DEGs between the two parental lines. Each group is referenced by the one after ‘vs.’; for example, 082 vs. w30 is referenced by w30.
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

E, L.; Lyu, S.; Wang, Y.; Xiao, D.; Liu, T.; Hou, X.; Li, Y.; Zhang, C. Integrating Dynamic 3D Chromatin Architecture and Gene Expression Alterations Reveal Heterosis in Brassica rapa. Int. J. Mol. Sci. 2024, 25, 2568. https://doi.org/10.3390/ijms25052568

AMA Style

E L, Lyu S, Wang Y, Xiao D, Liu T, Hou X, Li Y, Zhang C. Integrating Dynamic 3D Chromatin Architecture and Gene Expression Alterations Reveal Heterosis in Brassica rapa. International Journal of Molecular Sciences. 2024; 25(5):2568. https://doi.org/10.3390/ijms25052568

Chicago/Turabian Style

E, Liu, Shanwu Lyu, Yaolong Wang, Dong Xiao, Tongkun Liu, Xilin Hou, Ying Li, and Changwei Zhang. 2024. "Integrating Dynamic 3D Chromatin Architecture and Gene Expression Alterations Reveal Heterosis in Brassica rapa" International Journal of Molecular Sciences 25, no. 5: 2568. https://doi.org/10.3390/ijms25052568

APA Style

E, L., Lyu, S., Wang, Y., Xiao, D., Liu, T., Hou, X., Li, Y., & Zhang, C. (2024). Integrating Dynamic 3D Chromatin Architecture and Gene Expression Alterations Reveal Heterosis in Brassica rapa. International Journal of Molecular Sciences, 25(5), 2568. https://doi.org/10.3390/ijms25052568

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop