Inter- and intraspecific constraint of accessible chromatins maps regulatory loci involved in maize speciation and domestication

doi:10.21203/rs.3.rs-5009597/v1

Download PDF

Article

Inter- and intraspecific constraint of accessible chromatins maps regulatory loci involved in maize speciation and domestication

https://doi.org/10.21203/rs.3.rs-5009597/v1

This work is licensed under a CC BY 4.0 License

Version 1

posted

You are reading this latest preprint version

Comparative genomic studies have the capability to find genes that are either under evolutionary constraint for organismal living or specialized for trait innovations. Meanwhile, there is a growing body of evidence suggesting that evolutionary constraint can also happen on non-coding regulatory sequences, exerting profound effects on fitness traits. However, the majority of these evidences originate from primate studies and has yet to be thoroughly explored in plants. In this study, we adopt assay for transposase-accessible chromatin using sequencing (ATAC-seq) to profile the inter- and intraspecific constraint of over 80,000 maize accessible chromatin regions (ACRs). Our findings reveal that ACRs exhibit accelerated evolution over coding genes, resulting in approximately one-third of ACRs that specifically existed in maize lineage and exerting prominent transcriptional effects on genes associated with speciation. We highlight the role of transposable elements (TEs) in driving intraspecific innovation of ACRs, and uncover hundreds of candidate ACRs that likely involved in the transcriptional rewiring of genes during maize domestication. Furthermore, we demonstrate the importance of accessible chromatins in maintaining maize subgenome dominance and controlling complex trait variations. Our work establishes a new framework for analyzing the evolutionary trajectory of extensive plant regulatory sequences, and offers some candidate loci for downstream exploration and application during maize breeding.

Biological sciences/Genetics/Genomics/Comparative genomics

Biological sciences/Plant sciences/Plant genetics

Biological sciences/Genetics/Genomics/Epigenomics

Intricate regulation of gene expression is fundamentally important for the execution of complex biological processes in higher organisms. Maize, a classical genetic model in the field of gene regulation, has been extensively studied since Dr. Barbara McClintock’s groundbreaking discovery of transposable elements (TEs), commonly known as ‘jumping genes’ that are responsible for the mosaic pigmentation phenotype observed on maize kernels¹. In recent years, the assembly and analysis of maize genomes have revealed a highly sophisticated and hierarchical manner of gene regulation in maize, characterized by the prevalence of distal cis-regulatory elements (CREs)², particularly enhancers, that likely attributed to the proliferated maize genome (~ 2.3 Gb) with over 85% DNAs composed by TEs. Meanwhile, genetic mapping efforts in maize unveiled an enrichment of quantitative trait nucleotides (QTNs) within non-genic regulatory regions³, including well-known loci such as tb1^4,5, vgt1⁶, KRN4⁷, ZmCCT9⁸, ZmCCT10⁹, and UPA1¹⁰. In higher-dimensional scale, the widespread distal CREs in maize can orchestrate with transcription factors (TFs) and other trans-acting elements, and simultaneously interact with promoters or other gene-proximal regions to regulate gene expression through chromatin interactions². These findings underscore the crucial role of regulatory sequences in controlling development and trait variations in maize.

The progress of high-throughput chromatin profiling assays, such as ATAC-seq, CUT&Tag and DAP-seq technologies, has revolutionized genome-wide, base-resolution identification of CREs in plants¹¹. ATAC-seq utilizes engineered hyper-sensitive Tn5 transposons preloaded with sequencing adaptors to selectively target DNAs within accessible chromatin regions (ACRs)¹². ACRs offer biomacromolecules the physical access to DNAs, and are usually colocalized with CREs and transcription factor binding sites (TFBSs), therefore play a critical role in transcriptional regulation¹³. In maize, a series of studies have employed ATAC-seq to profile ACRs, which typically encompassed several to tens of megabases (Mbs) across major maize tissues including seedlings^14,15, roots¹⁶, leaf^2,17, tassel¹⁸ and ear^2,18. Except these progresses in bulk tissues, ATAC-seq has also been continually optimized to satisfy the increasingly smaller sample input, even down to a single cell, and a recent study utilized single-cell ATAC-seq to investigate the chromatin dynamics associated with cellular heterogeneity in a few maize tissues¹⁹. However, as a model species with substantially high tissue complexity, the existing ATAC-seq data in maize is undoubtedly unsaturated, and in-depth characterization of maize accessible chromatins is awaiting to be proceeded.

In comparative genomics, it is well established that highly constrained DNA sequences shared among related species are essential for organismal living and may have experienced more purifying selection, while those specifically existed in a single or a small group of species are usually under positive selection and responsible for innovative traits^19,20. In the past few years, the body of evidence is rapidly growing that a considerable proportion of genome sequences including not only coding regions but also massive CREs in non-coding regions can be evolutionary constrained^20–23. Since CREs are typically associated with some epigenetic features such as accessible chromatins, active histone modifications, low DNA methylation, as well as unique binding of TFs¹¹, an emerging research field termed ‘evolutionary epigenomics’ has combined both epigenomic and comparative genomic approaches to study the constraint of DNAs associated with CREs^20,22,24. In mammals, the Zoonomia project estimated that approximately half of the human CREs were evolutionary constrained in other 240 mammal genomes²¹, and another work reported that approximately 80% and 50% of the evolutionary constrained DNA bases in primates were outside of protein-coding exons and had no annotation in human ENCODE data²². Relevant scenarios are also found in plants that more than half of the highly constrained DNAs between diploid potato and other Solanaceae genomes happened in non-coding regions²³. With technical progress of CRE annotation methods and the availability of multiple high-quality genome sequences in diverse Poaceae species, we are now able to generate a comprehensive map of regulatory loci in maize and to trace their evolutionary trajectory during approximately 100 million years of Poaceae genome evolution²⁵, as well as to offer some new insights into maize trait innovations as compared with its wild ancestor and other grass species.

Here we portrayed a pan-tissue map of accessible chromatins in over ten major maize tissues and explored their interspecific and intraspecific constraint in multiple Poaceae and maize genomes. Unlike in primates, our analysis revealed strikingly low constraint (~ 1%) of ACRs in Poaceae, and highlighted the major role of TEs in generating maize specific ACRs that primarily happened in regions distal to genes. By comparison, approximately half of ACRs were constrained in diverse maize genomes, which exhibited a higher degree of enrichment in trait-associated SNPs over other highly variable ACRs. We also uncovered a set of ACRs that likely associated with maize domestication, demonstrating their entanglement with protein-coding genes that either artificially selected or differentially expressed between maize and wild teosintes. Our work established a new framework for the evolutionary analysis of massive plant regulatory elements, and unveiled some potentially important regulatory loci that involved in maize specification and domestication.

A pan-tissue map of accessible chromatins in maize

We profiled the accessible chromatins for 12 major tissues or cell types in maize (Supplementary Table 1). Extending from the earlier data of ear^2,18, tassel¹⁸, leaf², seedling^14,15, root¹⁶ and mesophyll¹⁷, we newly generated the ATAC-seq data of embryo, endosperm, pericarp, shoot apical meristem (SAM), silk and stem (Fig. 1a). Each data had two biological replicates that exhibited high reproducibility, canonical enrichment near gene TSS (transcription start sites), and clustered well in both principal component and heatmap analysis (Supplementary Fig. 1 and Supplementary Fig. 2). Therefore, the data of two biological replicates were combined for subsequent analysis. There were 26,743 to 39,561 ACRs (Fig. 1b) identified in each tissue, covering approximately 0.67–1.0% of the B73 reference genome (Fig. 1b). As expected, ACRs were mainly distributed in euchromatin regions (Fig. 1c), aggregating into some hotspots that either conserved or specifically existed in different tissues (Supplementary Fig. 3). We further classified ACRs into genic (36–52%), proximal (19–25%) and distal (29–46%) ACRs according to their intersection with genes (Fig. 1d). Due to the proliferation of maize genome, a considerable number of dACRs were located in “gene deserts” regions (100 kb away from the nearest gene, Supplementary Fig. 4a). Interestingly, these ultra-distal ACRs usually exhibited high conservation across the tissues we profiled, and could be independently verified by maize single-cell ATAC-seq data¹⁹ and ChIP-seq data² of various active histone modifications (Supplementary Fig. 4a and 4b). They were also found to be entwined with chromatin interactions, some pinpointed to genes in 3D genomic data² (Supplementary Fig. 4b), indicating the cis-regulatory role of maize ACRs in both local and distal manners.

In contrast to the bloated maize genome that primarily composed by TEs (~ 85%), only ~ 25% of ACRs overlapped with TEs (Fig. 1e), suggesting the enrichment of ACRs for low-copy sequences in maize. Specifically, dACRs (~ 40%) had higher TE composition over pACRs (~ 23%) and gACRs (~ 15%). Regarding to TE superfamilies, Gypsy accounted for over half of TEs in dACRs, while Helitrons enhanced substantially in pACRs and gACRs which was consistent with their strong preference for gene-rich regions in maize²⁶ (Supplementary Fig. 5a). It was worth noting that the insertion time of LTRs associated with ACRs seemed to be more ancient than the remaining non-association LTRs, suggesting the purification of recent TE insertions that may disrupt ACRs’ function (Supplementary Fig. 5b). Unlike TEs, simple repeats were reported to be enriched in CREs in humans²⁷ and also found to be overrepresented in maize ACRs (Supplementary Fig. 5c). Specifically, they exhibited a notable valley in the core region of ACRs while elevated in flanking regions, in line with their presumed role to assist the binding of TFs²⁷.

By integrating the data of PlantTFDB²⁸, JASPAR²⁹ and DAP-seq^2,30, the DNA binding motifs of 306 maize TFs from 36 families were curated (Supplementary Table 2). We discovered the motifs of 243 TFs within 50 bp sequences flanking the summit of ACRs, with each tissue having 157 to 191 TFs (Supplementary Table 3). Similar to a report in cotton³¹, the motifs of BRR-BPC, C2H2, Dof, ERF and TCP families, which were fundamentally important for plant growth and development, were highly enriched in all profiled tissues (Fig. 1f and Supplementary Fig. 6). Specifically, the enrichment of SBP families was high in the ACRs of ear, tassel and SAM, in consistent with one of their major roles in establishing meristem boundaries within inflorescence tissues³². Additionally, we also found the motifs of some WRKY family genes that overrepresented in roots and mesophyll, a phenomenon that was accorded with their significantly enhanced transcription in similar tissues^33,34.

We further merged the ACRs of all tissues into 80,365 pan-tissue ACRs (~ 61.4 Mb, ~ 2.9% of the maize genome), accounting for 2.8 to 4.4-fold of the ACR length in a single tissue. Our estimation approached the maize single-cell ATAC-seq assays (~ 4%) that had integrated the data of 56,573 nuclei from 6 major tissues¹⁹. However, our ACRs were apparently not saturated, representing approximately 85% of saturation according to Michaelis-Menten equation analysis (Fig. 1g), indicating the further necessity to incorporate more specialized tissues or cell-types. Pan-tissue ACRs were further classified into core, variable and tissue-specific ones (Fig. 1g). The contrasting numbers of variable (n = 46,483) and core ACRs (n = 6,475) again suggested the high tissue specificity of accessible chromatins in maize. To further reveal the selective pressure on ACRs, we estimated the SNP density within and flanking ACRs using data of maize HapMap3³⁵. As expected, the core ACRs had the lowest SNP density in their central region (Supplementary Fig. 7), followed by variable and private ACRs that were accorded with their conservation degree. Interestingly, the flanking regions of core and variable ACRs shown significantly higher SNP density than private ACRs, possibly due to their correlation with previously reported higher genetic diversity of genic and proximal regions in maize³⁶. In summary, we generated a pan-tissue map of over 80k maize ACRs that could allow their evolutionary characterization in subsequent sections.

Evolutionary constraint of maize ACRs in other grass genomes

Comparative genomic studies in plants have commonly characterized the conservation or innovation of protein-coding genes, yet the evolutionary constraint of more widespread regulatory sequences remains less investigated. Based on the pan-tissue ACRs, we employed a Cactus-based comparative approach (Supplementary Fig. 8) to trace their evolutionary trajectory in 34 Poaceae species and 4 outgroups (Supplementary Fig. 9 and Supplementary Table 4). To demonstrate the robutness of our method, we identified the regions in sorghum genome that highly constrained (alignment coverage > = 90%) with maize ACRs, which could successfully reconstruct the two maize subgenomes (Supplementary Fig. 8b). When expanding to other Poaceae species, the number of highly constrained ACRs was negatively correlated with their pairwise divergence (Fig. 2a). For example, there were 7,117, 3,592 and 1,421 ACRs that highly constrained in sorghum (diverged ~ 11.1 Mya), foxtail millet (~ 22.8 Mya) and rice (~ 49.1 Mya) respectively, which exhibited a good consistency with their overall phylogeny (Fig. 2a and Supplementary Table 4). In general, the constraint of ACRs situated between that of protein-coding sequences and randomly selected DNA sequences (Fig. 2a), suggesting an accelerated evolutionary rate of ACRs over protein-coding sequences in crops.

We then projected ACRs into a planar map (N1-N2, Fig. 2b) according to the number of highly constrained species (N1) and the number of species without constraint (N2, alignment coverage < = 10%), which exhibited a similar triangle as the map generated by comparing human CREs against hundreds of mammal genomes²¹. It was worth noting that only 1% (n = 808) of the maize ACRs were highly constrained in most Poaceae species (N1 > = 20) (Fig. 2b), contrasting with approximately half of human CREs that were constrained in primates. In addition, there were a substantial number of ACRs (~ 34.6%) that specifically existed in maize (N2 > = 34) (Fig. 2b). These data implied unprecedented high genetic diversity of regulatory sequences in grasses, and further prompted us to explore the driving forces of maize specific ACRs. We found a profound proportion (~ 70%) of maize specific ACRs located in distal regions (Fig. 2c), while both highly and moderately constrained ACRs were exclusively located (> 90%) in genic regions. Meanwhile, nearly half of the sequences in maize specific ACRs were annotated as TEs (Fig. 2d), particularly Gypsy elements, contrasting with very rare TE composition in highly and moderately constrained ACRs. These data highlighted the major role of TEs, mostly resided in regions distal to genes, that involved in the innovation of maize specific ACRs.

To further characterize the interspecific constraint of ACRs, we applied the Uniform Manifold Approximation and Projection (UMAP) algorithm on all ACRs based on their alignment coverage in each of the non-maize genomes (Fig. 2e). This analysis resulted in a UMAP resembling ‘a bunch of flowers’ with constrained (N1 > = 10) ACRs enriched in stalks, while the maize specific ACRs formed two distinct groups, one (class I, n = 7,897) dispersed within the ‘petal’ region, and the other (class II, n = 19,780) formed an ‘isolated island’. We found almost no sequence homology of class II ACRs in other Poaceae genomes as compared with sparse homology of class I ACRs. Although class II ACRs exhibited higher affinity with TEs over class I ACRs, there were still more than half of both class II and class I ACRs that not originated from TEs (Fig. 2d). These non-TE maize specific ACRs might result from sequence mutations that accumulated due to genome redundancy after whole genome duplication of maize.

Genes regulated by highly constrained or species-specific ACRs are likely to encode for traits of conservation or innovation, respectively²¹. Using our previously generated Chromatin Interaction Analysis by Paired-End Tag Sequencing (ChIA-PETs) data³⁷, we identified these two different sets of genes (n = 2,041 and 10,779), many of which that unable to be identified by the simple principle of proximity (Supplementary Fig. 10). Genes targeted by highly constrained ACRs were found to be enriched in photosynthesis pathways including NADPH: quinone reductase activity and response to light stimulus (Supplementary Fig. 11a). They were also participated in the transmembrane transport of some polyols (Supplementary Fig. 11a) among which mannitol served as a major product of photosynthesis³⁸, as well as Xylose which could play a crucial role during the formation of plant cell walls^39,40 to mediate plant response to abiotic stresses (Supplementary Fig. 11a). On the other hand, genes (n = 10,779) regulated by maize specific ACRs were primarily enriched in plant immune or defense response pathways as were commonly observed for PAV genes^41,42 (Supplementary Fig. 11b). They include jasmonic acid and salicylic acid metabolic process (Supplementary Fig. 11b), which could regulate plant resistance against pathogen invasion⁴³ during the continuous evolution of plant immune system⁴⁴. Meanwhile, genes regulated by maize specific ACRs were also enriched in flower development and morphogenesis pathways (Supplementary Fig. 11b). One special example was ramosa1 (ra1), a pivotal gene that could regulate maize inflorescence branching architecture⁴⁵ and specifically expressed in immature tassel and ear in our profiled tissues (Supplementary Fig. 12). It was notable that only one ACR was identified in both the gene body and 20-kb flanking regions of ra1 (Fig. 2f), which was very likely to be a decisive regulatory locus of ra1. This ACR was maize-specific since it had no constraint in even sorghum and other Poaceae species, therefore it might function through ra1 to contribute to the substantial changes of inflorescence architecture in maize compared to other grasses.

We further investigate the expression profiles of aforementioned two set of genes to find more clues regarding their roles in maize speciation. Although genes regulated by highly constrained ACRs exhibited overall higher expression level (Supplementary Fig. 13a), those regulated by maize specific ACRs, by contrast, exhibited significantly higher tissue specificity (Fig. 2g). This result was consistent with previous findings that tissue-specific expressed genes usually evolved faster⁴⁶ and were more likely to be affected by speciation and gene duplication⁴⁷. We further utilized an independent dataset that applied the phylogenetic approach to estimate the age of genes in Poaceae⁴⁸, and found a higher proportion of young genes among genes that regulated by maize specific ACRs (Supplementary Fig. 13b). Collectively, we demonstrated the massive number of maize specific ACRs that primarily derived from the activity of TEs, and genes they regulated might involve in the innovation of some maize specific traits such as inflorescence development and immune response.

Intraspecific constraint of ACRs and their association with maize domestication

Recent release of diverse maize reference genomes^14,49,50 offers us an opportunity to further explore the intraspecific constraint of maize ACRs. By adopting the same comparative approach, we identified highly constrained ACRs in a total of 50 non-B73 maize genomes (Fig. 3a and Supplementary Table 5), finding approximately half of ACRs were fixed in nearly all maize genomes (Supplementary Fig. 14a). This ratio was comparable to the estimation of intraspecific constraint of CREs (~ 50%) in primates²², and accorded with a previous estimation of genetic diversity in maize that was comparable or even higher than the divergence between human and chimpanzee⁵¹. On the other hand, only ~ 7% (n = 5,239) of ACRs were constrained in less than 10 maize genomes (Fig. 3a), among which 342 were specifically found in B73. It seemed that the conservation of ACRs we estimated here was stronger than genes, since a pan-genome analysis in maize NAM population had reported a higher degree of gene variability¹⁴ (~ 70% were dispensable or private genes). To make clarification, we simultaneously estimated the constraint of ACRs and genes in our panel. We found apparently stronger constraint of genes over ACRs, and their constraint were highly correlated (r = 0.98) with only a few exceptions (e.g., LH244 and PH207) that the constraint of ACRs was slightly stronger than genes (Supplementary Fig. 14b). Interestingly, ACRs that constrained intraspecifically characterized quite different genomic properties with interspecifically constrained ACRs, for example, TEs were common in intraspecifically constrained ACRs while rare in interspecifically constrained ones (Fig. 3b). In addition, a considerable proportion of intraspecifically constrained ACRs were found in distal regions (Fig. 3c), while almost no distal ACRs was found in interspecifically constrained ones. These findings underscored the importance of historical TEs that shared in the common maize ancestors, which might reside in regulatory regions and subsequently be fixed after maize diversification.

Genome comparison between wild teosinte and maize had uncovered hundreds to thousands of candidate genes associated with maize domestication^5,7, yet there was still relatively few information regarding the regulatory sequences involved in this process. Based on the classical case of tb1, which was a ‘gain-of-function’ Hopscotch insertion that fixed in maize while absent in teosintes⁵, we hypothesized that ACRs which were highly constrained in most maize genomes while absent in two wild teosintes (Zea mays ssp. parviglumis and Zea mays ssp. mexicana)⁵² should play a pivotal role in rewiring transcriptional networks during maize domestication. Under this scenario, we could successfully validate the ACR that colocalized with tb1⁵³ (Supplementary Fig. 15a). On genome-wide scale, we totally identified 497 ACRs that potentially involved in maize domestication. Interestingly, as compared with all highly constrained ACRs in maize, this small group of highly-constrained and domestication-associated ACRs shown uncanonically higher TE composition (Fig. 3c), as well as a profound ratio within distal non-genic regions. These results again highlighted the prevalence and importance of TEs, particularly those located distal to genes during domestication of maize. We were interested in the TFs that involved in regulating these domestication associated ACRs and further predicted their TFBSs. We found significant motif enrichment of TCP family (Fig. 3d), which was accorded with the master regulator of maize domestication gene tb1. In addition, we also discovered the enriched motifs of C2H2 family including ramosa1 gene, as well as ERF and bZIP families that were fundamentally important during plant development.

To offer more functional clues into these domestication-associated ACRs, we again used the ChIA-PETs data to find a total of 262 genes that entangled with these ACRs through chromatin interaction. Notably, we noticed significant enrichment of these genes in another independent domestication gene list (Supplementary Fig. 15b) that identified by scanning SNPs between teosinte and maize populations⁵⁴. Moreover, among the differentially expressed genes (DEGs) identified in three major tissues between teosintes and maize⁵⁵, we found a significantly higher ratio (17.6%, Supplementary Table 6) of DEGs within these interacting genes over randomly selected genes without interactions (12.7%, Fig. 3e). As an example, we identified an ACR in ~ 51 kb upstream of bHLH68, also known as phytochrome-interacting factor1 or ZmPIF1 that could play versatile roles in regulating plant growth, development and drought tolerance in maize⁵⁶. Chromatin interactions between this ACR and bHLH68 can be found not only in our ChIA-PETs data but also in other 3D genomic datasets (Fig. 3f ), and this ACR was fixed in multiple maize lines while absent in teosintes, which correlated well with its up-regulation in maize over teosintes (Supplementary Fig. 13). A similar case was also found for a domestication-associated ACR and its interacting gene cmu2 (Supplementary Fig. 15c), which encoded a chorismate mutase with elevated expression in maize ear and leaf and might be involved in plant pathogen infection resistance⁵⁷ (Supplementary Fig. 13). In summary, our results demonstrated the great promise of analyzing evolutionary constraint of ACRs to identify candidate regulatory loci that associated with crop domestication and subsequent adaptation processes.

Constraint and dominance of ACRs between two maize subgenomes

It was well documented that genes in one maize subgenome always dominated expression over genes in the other one⁵⁸, a phenomenon known as subgenome dominance which was still less investigated for chromatin accessibility. Using 3,719 high-quality genes pairs with exact 2:1 synteny to sorghum genes (Supplementary Fig. 16), we reconstructed two maize subgenomes (m1, 1,181 Mb; m2, 726 Mb) and classified ACRs into m1 (n = 47,666) and m2 (n = 28,442), respectively (Supplementary Table 7). It was notable that the numbers of highly constrained ACRs in m1 (n = 4,647, 9.7%) and m2 (n = 2,326, 8.1%) when compared to sorghum were higher than the highly constrained ACRs when mutually compared between m1 and m2 (~ 1,200, Fig. 4a and 4b), a phenomenon accorded with the well-known reciprocal subgenome fractionation observed for genes⁵⁹. We then characterized key metrics of ACRs between two subgenomes, finding both the number and density of ACRs in m1 were higher than m2, although their overall differences were subtle (Fig. 4c and Supplementary Fig. 17a). Accordingly, the transcriptional differentiation was also weak when profiled across two subgenomes (Supplementary Fig. 18a). However, when we focused on these 3,719 pairs of syntenic genes, the dominance of both chromatin accessibility (Fig. 4d and Supplementary Fig. 17b-c) and gene transcription was significantly enhanced (Supplementary Fig. 18b). Meanwhile, the chromatin accessibility of these conserved regions was also much stronger than the remaining less-conserved regions (Fig. 4c and 4d). These results suggested the importance of chromatin accessibility that could function to differentiate duplicated gene loci to introduce subsequent gene expression differentiation⁶⁰.

To further explore accessible chromatin in specific duplicated loci, we focused on teosinte branched1 (tb1, m1) and its m2 paralog tcp6, both of which encoded the TCP family transcription factors. Although gene density and chromatin accessibility were fairly high in up- and downstream regions of tcp6, no ACR identified in its gene body and 10-kb flanking regions (Fig. 4f), which was consistent with its silenced expression in major maize tissues (Supplementary Fig. 19a). In contrast, tb1 was exclusively activated in ear, which was likely resulted from the functional Hopscotch insertion that shown colocalization with two independent ACRs in our data (Fig. 4f). We also investigated two other syntenic genes pairs, ub2 (m1) and ub3 (m2), gt1 (m1) and Zm00001d048172 (m2), both were crucial in regulating plant architecture^61,62. Cases of ub3 and gt1 were similar to tb1 that their causative insertions, named KRN4⁷ and prol1.1⁶², exhibited accessible chromatin in most of our profiled tissues (Supplementary Fig. 19b and 19c). However, no parallel functional insertion was found near ub2 and Zm00001d048172, and consequently their transcription was weakened or completely silenced (Supplementary Fig. 19a). These data suggested that subgenome-specific functional sequence insertions, mostly TEs, might be one of the driving forces underlying the subgenome differentiation of chromatin accessibility and gene expression.

We next sought to decipher the divergence of TFBSs between two maize subgenomes. Unlike genes and ACRs that exhibited strong fractionation⁶³ (Fig. 4a), approximately 80% of the TFBSs (69.1–86.3% across our profiled tissues) can be identified in both two subgenomes (Fig. 4e and Supplementary Fig. 20a). However, the divergence of TFBSs enhanced in promoter and distal ACRs (Fig. 4e and Supplementary Fig. 20b and 20c), e.g., only one-third of the TFBSs in distal ACRs of seedling were shared between two subgenomes. Therefore, the strong conservation of TFBSs between two subgenomes might come from the complementation of ACRs from different functional regions. Additionally, we found a higher number of subgenome-specific TFBS in m1 than m2, which was accorded with the general principal of m1 dominance over m2 (Fig. 4e and Supplementary Fig. 20a-c). For instance, the TFBS of GLK family, which could regulate some earliest photosynthesis steps by activating genes encoding chloroplast-localized and photosynthesis-related proteins⁶⁴, were exclusively found in m1 (Supplementary Fig. 21). Interestingly, genes that encoding these GLK TFs were also existed in m1, indicating their regulatory functions that preferentially exerted within the same subgenome. These results indicated that the regulomes of two maize subgenome were largely shared except for a few TF families that specifically functioned within a single subgenome to precisely regulate the development of particular tissues.

Trait-associated SNPs are enriched in ACRs that highly constrained in maize

In species with giant and complex genomes, for example in maize⁶⁵ and humans⁶⁶, trait-associated SNPs (TASs) identified by Genome Wide Association Studies (GWAS) often fall into non-coding regulatory regions⁶⁷. We anticipated a better annotation of maize TASs using our pan-tissue map of ACRs. Therefore, we firstly performed GWAS for 15 complex traits in a maize association mapping panel (n = 260) that genotyped across HapMap3 SNPs, resulting in 110 to 363 TASs (minor allele frequency > 0.05, P-value < 10^− 5) for each trait (Supplementary Fig. 22 and Supplementary Table 8). We then integrated all TASs (n = 2,984) to improve their statistical power when assessed with ACRs. As expected, TASs were overrepresented in ACRs with their density being 1.64 (cob-weight) to 6.49-fold (leaf-length) higher over the remaining non-ACR regions (Fig. 5a). In addition, the density of TASs in m1 ACRs was overall higher than m2 ACRs (Fig. 5b), indicating a predominant role of m1 subgenome in the genetic regulation of maize complex traits. Interestingly, although all the traits we investigated here were believed to be classical quantitative traits, six out of them (for example, upper-leaf-angle and days-to-tassel) exhibited non-canonical results that the TASs density was higher in ACRs of m2 subgenome (Fig. 5b). The existence of hotspots with an excessive of TASs (Supplementary Fig. 23a and b), which usually entangled with strong LD and can be mapped within m1 or m2 subgenomes for different traits, might be one of the appropriate explanations to these observations.

Researches in humans have revealed that TASs of diseases or other complex traits are more likely to be enriched in interspecific constrained CREs^21,22. However, due to low and limited number of interspecifically constrained ACRs in Poaceae (Fig. 2b), the power to investigate maize TASs in ACRs with a varying degree of interspecific constraint was statistically weak. However, we observed that TASs were overrepresented in ACRs with higher intraspecific constraint (Fig. 5c). Specifically, the density of TASs in ACRs that highly constrained in most maize genomes exhibited 1.4-fold increase compared with ACRs with weak constraint in maize. We further verified this result (Supplementary Fig. 24) by using another independent TASs dataset (n = 2,360) that also mapped on maize B73 V4 genome and tested in an ultra-large maize population with more than 1,000 lines⁶⁸. These findings suggested that complex traits in maize are primarily driven by genetic variants in ACRs, especially for those that fixed in diverse maize population.

We finally tried to annotate some new GWAS loci that resided in ACRs. As an illustration, we identified a lead SNP (chr3: 215,791,183) that significantly associated with both days-to-silk and days-to-tassel, which fell into a non-coding region that was ~ 5.3 kb away from nearest gene. Notably, this lead SNP was colocalized with an intraspecifically highly constrained ACR that accessible in three reproductive tissues (ear, tassel, and SAM) (Fig. 5d), suggesting their potential role of regulating genes involved in flowering. Further incorporation of 3D chromatin interaction data, for example the HiChIP data (Fig. 5d), could further pinpoint to a list of candidate genes that were likely interacting with this ACR, although the exact function of these GWAS loci and interacting genes needed further experimental validations. To end up our work, we have integrated the pan-tissue ATAC-seq data with our recently established HEMU online database (https://shijunpenglab.com/HEMUdb/databrowse)^69–71 to facilitate easy visualization and analysis of any loci that recorded by B73 V4 genome coordinates.

With the rapid advancement of plant genomics, we are on a highway to enter the post-genomic era with the primary research goals progressively shifted from gene-centric studies to encompass a wider spectrum of genome regulatory sequences. In this work, we have generated a comprehensive map of accessible chromatins by integrating the ATAC-seq data of 12 major maize tissues. Although ATAC-seq has been commonly used by many maize studies to explore the tissue-specificity of accessible chromatins^2,14,16–19, most of them are carried out on a single or relatively few tissues. They include the recently published single-cell atlas of ACRs across six maize tissues, which has highlighted the strong association between GWAS SNPs and tissue/cell-type specificity of ACRs. To our knowledge, our map might be one of the most comprehensive ACR maps in plants to date⁷², although it is still far from reaching saturation that awaiting the incorporation of more ATAC-seq data from distinct tissues or cell types. We anticipate the growing of ATAC-seq data in maize and other plant species in the forthcoming years, especially by applying the single-cell ATAC-seq technology that has shown its unique advantages in decoding cell heterogenity and the flexibility to obtain chromatin accessibility of bulk tissues when data of tens of thousands of cells were integrated.

We performed large-scale comparative genomic analysis to estimate the inter- and intraspecific evolutionary constraint of our pan-tissue ACRs. Although it is well-documented that highly constrained DNA sequences in a group of related species are usually functional, there are also many circumstances that ‘loss-of-function’ mutations can happen in these highly constrained regions, for example, large-effect SNPs that can result in pseudogenes. Therefore, data resources of ATAC-seq are required for other crop species to unveil the true extent of chromatin accessibility for any DNA sequences that constrained with an ACR in maize. In addition, by analyzing the intraspecific constraint, we also identified several hundreds of candidate ACRs that fixed in most maize lines while absent in teosintes, which were likely associated with maize domestication. However, as compared with these ‘gain-of-function’ regulatory elements, the effect of genetic bottleneck happened during crop domestication are more usually to result in the ‘loss-of-function’ of beneficial alleles in their ancestors. Therefore, the generation of similar regulatory maps in wild teosintes will help to shed new light on these untapped regulatory loci that should also play an important role of maize domestication process. We also demonstrated the power of utilizing ATAC-seq data to annotate the genetic loci as discovered by GWAS in maize, although experimental tools, e.g., transient transcription activation assay or genome editing is required for downstream verification. In summary, our study extends our knowledge of the plant evolutionary epigenomic, and offers a wealthy data resources for the excavation and utilization to enhance the genetic improvement of maize.

Plant materials

The seeds of maize B73 inbred line were planted at the experimental station of Shandong Agricultural University, Tai’an, China in 2021. After 12 days after self-pollination, the embryo, endosperm and pericarp from roughly 50 seeds were separated, then freeze quickly in liquid nitrogen. SAM tissues were dissected from sixty two-week-old field-grown plants and immediately frozen in liquid nitrogen.

Library construction and ATAC sequencing

Plant tissues were ground into powder in liquid nitrogen, then incubated with lysis buffer at 4°C for 10 min. After filtration with 40 um cell strainer, cell solution was transferred to the top 60% percoll solution and centrifuged at 1800 g for 20 min. The intermediate phase of the centrifuged liquid was extracted and then diluted with 10 mL RSB buffer. Cell nuclei pellets were harvested by centrifuge at 1000 g for 10 min and washed with RSB buffer once time. Approximate 50,000 nucleus were added to transposition reaction solution to perform tagmentation. Tn5 transposed DNA were purified by AMPure DNA magnetic beads, then amplified by PCR. The final qualified ATAC-seq libraries were sequenced on the Illumina Novo-seq platform (San Diego, CA, United States) with PE150 mode.

ATAC-seq analysis

ATAC-seq raw data was processed using ATAC-pipe pipeline⁷³. The quality of mapping reads was evaluated using the “--MappingQC” module in ATAC-pipe. After trimming adapter sequences, reads were mapped to Z. mays AGPv4 reference genome using Bowtie2 (v.2.3.0)⁷⁴ and PCR duplicates were removed as described⁷³. Mapped reads were then shifted + 4/-5bp depending on the strand of the reads, so that the first base of each mapped read can represent the Tn5 cleavage position. All mapped reads were further extended to 50 bp centered by the cleavage position. Reads mapping to highly repeated regions, chromosome M and P, and unanchored contigs were removed. Two biological replicates from a same sample were merged using “--Merge” module in ATAC-pipe.

Samples clustering and principal component analysis

The clustering among all samples was performed by the “dba.plotHeatmap” function of DiffBind (v.3.4)⁷⁵ with default parameters. The PCA analysis of all samples was performed using the “dba.plotPCA” function of DiffBind (v.3.4) with default parameters.

Transcriptome data analysis

Firstly, FASTP (v.0.20.1)⁷⁶ was used to trim adaptors and filter out low quality reads. Trimmed reads were further mapped to the Z. mays AGPv4 reference genome by HISAT2 (v.2.2.1)⁷⁷ with default parameters. Alignments were subsequently sorted and converted into bam files using Samtools (v.1.12)⁷⁸. Gene expression was estimated by StringTie (v.2.1.5)⁷⁹ in the forms of both fragments per kilobase of transcript per million mapped reads (FPKM) and transcripts per kilobase million reads (TPM). Genes with FPKM or TPM > 1 were considered to be expressed in the corresponding samples.

ACR identification

MACS2⁸⁰ was used to call peaks by “--PeakCalling” module in ATAC-pipe with the parameters of “macs2 callpeak -f BED -g 2.1e + 9 -q 0.01 --nomodel --shift 0”. Tn5 integration site density (insertions per kb) was estimated for each peak, as well as the matched and randomized control regions (established using bedtools (v.2.26.0)⁸¹, shuffled specifically excluding ACRs from the randomized selection space). ACRs with Tn5 integration site density less than the 99.9% upper quantile of randomized control regions were removed according to the published method¹⁴.

Annotation of TEs and their association with ACRs

Repetitive elements were identified using Extensive de novo TE annotation (EDTA, v.1.9.6)⁸² pipeline with the parameters “--species Maize --anno 1”. Next, the unclassified transposons from EDTA was re-annotated using DeepTEs⁸³. Then, TE composition in ACRs was calculated by intersectBed (v.2.29.2)⁸¹ with the parameters “intersectBed -wo”.

The insertion time of intact LTRs associated with ACRs (50% TE sequence overlap with ACR) was calculated. Firstly, intact LTRs were extracted by two long terminal repeat sequences and mutually aligned by Muscle⁸⁴. The distance K with Kimura Two-Parameter approach between LTRs was calculated by distmat program in EMBOSS (v6.6.0)⁸⁵. The activity of each LTR was calculated by the formular: T = K/(2 × r), where “r” refers to a general substitution rate of 1.3 × 10^− 8 per site per year in the grass family⁸⁶. RepeatMasker was used to annotate simple repeats in the ACR with the parameters of “-no_is -norna”.

TF annotation and motif analysis

TF genes (n = 2,361) in the B73 V4 genome were annotated using Transcription Factor Prediction function of PlantTFDB²⁸ (http://planttfdb.gao-lab.org/). A total of 259 and 45 TF motif models were downloaded from PlantTFDB and JASPAR (https://jaspar.genereg.net/), respectively. Additional 12 TF motif models from De novo motif identification were detected by the discriminative motif discovery workflow of MEME-ChIP (v5.0.2) with default settings based on DAP-seq data^2,30. The motifs of 306 TF genes from 36 families were curated. Finally, the enriched motifs in the interest region of ACRs were identified by using the AME v.5.0.2 tool of MEME Suite with the parameters; “--scoring avg --method fisher --hit-lo-fraction 0.25 --evalue-report-threshold 10.0 --control”.

Evolutionary constraint analysis of ACRs

To measure the sequence constraint of ACRs, Progressive Cactus (v2.6.0)⁸⁷ was used to perform pairwise alignment of maize with other genomes, with the parameters “cactus ./js_grass ./B73.othergenome.txt ./ B73.othergenome.hal --consMemory 320G --binariesMode local --stats ”. After obtaining the hal file, the haltochains function of Cactus was utilized to convert into chain file with the parameters " cactus-hal2chains ./js_grass ./ B73.othergenome..hal chain/ --refGenome othergenome ". Then, we followed the published code to employed the LiftOver⁸⁸ software to determine the conservation of ACR sequences from B73 genome in other genomes²¹. Alignment > = 90% was considered as constrained ACRs (see Supplementary Fig. 7a). Phylogenetic relationships among the 38 Poaceae species and 4 outgroups were obtained from TimeTree⁸⁹.

UMAP analysis on ACRs

For each ACRs, LiftOver software was used to get the percentages of aligned positions in the other 34 genomes. The resulting matrix of n. ACRs by 34 genomes was then used as input to run UMAP⁹⁰ with default parameters.

GO analysis

GO enrichment analysis was performed using GOATOOLS⁹¹ with the Z. mays AGPv4 GO annotation and go-basic.obo (http://current.geneontology.org/ontology/go-basic.obo). The module “find_enrichment.py” with Fisher's exact test was used.

Analysis of differentially expressed genes

The published RNA-seq data between 9 teosintes and 6 maize inbred lines (generated by John Doebley’s lab⁵⁵, NCBI SRA accession number of SRP047528) from ear, leaf and stem tissues were downloaded. For the identification of differentially expressed genes between teosintes and maize inbred lines, we treated the bam files from teosintes (n = 18, 2 replicates from each of 9 lines) and maize inbred lines (n = 12, 2 replicates from each of 6 lines) as replicates, then used Cuffdiff⁹² to calculate the FPKM in teosinte and maize. Genes with the fold-change of FPKM > 2 and q-value < 0.05 were identified as DEGs between teosintes and maize.

Maize two subgenome identification

The two subgenome of Z. mays AGPv4 were identified according to our published pipeline⁹³. Briefly, we used the Synmap pipeline in CoGe (https://genomevolution.org/CoGe/SynMap.pl). We used last to blast⁹⁴ the CDS sequences between maize and sorghum, then detected syntenic blocks with DAGchainer⁹⁵ with options -D 20 -A 5. QuotaAlign⁹⁶ was further used to merge adjacent syntenic blocks. The syntenic depth was set to 2:1 for maize and sorghum comparison, and the overlapped distance was set to 40 to permit overlapped syntenic regions. Fractionation bias was applied to determine subgenome organization in maize compared with sorghum. For each pair of chromosomes, the copy with a greater number of unique genes was assigned to the maize1 subgenome, whereas the pair with fewer uniquely retained genes was assigned to the maize2 subgenome. The CodeML utility in the PAML⁹⁷ software package was used to calculate the non-synonymous (Ka) and synonymous (Ks) rates between orthologous genes.

High-quality syntenic gene pairs between two subgenomes

According to the above COGE output file of QuotaAlign results, the gene pairs in maize having 2:1 relationship with sorghum genes were extracted, and among those with Ks > 2 were further removed (potential pre-grass homeologs). Finally, remaining gene pairs were assigned to m1 and m2 subgenomes according to the result of subgenome classification.

Genome wide association studies (GWAS)

Plink¹⁰⁰ was used to convert vcf files to bed files, then filter out sites with missrate > 0.5 or MAF < 0.05. GCTA¹⁰¹ was used to perform principal component analysis. The kinship matrix was obtained using emmax-kin for calculating genetic relatedness. EMMAX¹⁰² was used to conduct GWAS analysis and correcte for population structure by including the top 20 principal components (to account for population structure) and the kinship matrix (to correct for genetic relatedness and other fixed effects).

Data Availability

All the published ATAC-seq data that used in this study were summarized in the Supplementary Table 1. All the newly generated ATAC-seq data in this study were deposited into NCBI SRA under BioProject ID PRJNA1048297.

Acknowledgement

We want to thank the colleagues from Frasergen company in Wuhan, China for their assistance in preparing the ATAC-seq sequencing libraries. We also appreciate the assistance and suggestions from Dr. Zhong, Silin from The Chinese University of HongKong. This work was supported by grants from the National Natural Science Foundation of China (32172014 and 32372115), the Shenzhen Outstanding Youth Science Fund project (RCYX20231211090225038), the Young Elite Scientists Sponsorship Program by CAST (2021QNRC001), and the Fundamental Research Funds for the Central Universities, Sun Yat-sen University (23lgbj016) to J.S.

Author Contributions

J.S. and X.G. conceived this research project. Y.L. completed most of the data analysis. X.G., Y.Z., L.H., Q.L., H.L., X.Y., F.X., F.Y., and J.L. participated in data analysis. X.G., H.L., X.L., F.X., and X.Y. involved in sample collection and sequencing. J.S., X.G. and Y.L. wrote and finalized this manuscript.

Competing interests

The authors declare no competing interests.

Wessler, S.R. Phenotypic diversity mediated by the maize transposable elements Ac and Spm. Science 242, 399-405 (1988).
Ricci, W.A. et al. Widespread long-range cis-regulatory elements in the maize genome (vol 5, pg 1237, 2019). Nature Plants 6, 328-328 (2020).
Liang, Y., Liu, H.J., Yan, J. & Tian, F. Natural Variation in Crops: Realized Understanding, Continuing Promise. Annu Rev Plant Biol 72, 357-385 (2021).
Clark, R.M., Wagler, T.N., Quijada, P. & Doebley, J. A distant upstream enhancer at the maize domestication gene tb1 has pleiotropic effects on plant and inflorescent architecture. Nature Genetics 38, 594-597 (2006).
Studer, A., Zhao, Q., Ross-Ibarra, J. & Doebley, J. Identification of a functional transposon insertion in the maize domestication gene tb1. Nature Genetics 43, 1160-1163 (2011).
Castelletti, S., Tuberosa, R., Pindo, M. & Salvi, S. A MITE transposon insertion is associated with differential methylation at the maize flowering time QTL Vgt1. G3 (Bethesda) 4, 805-12 (2014).
Liu, L. et al. KRN4 Controls Quantitative Variation in Maize Kernel Row Number. PLoS Genet 11, e1005670 (2015).
Huang, C. et al. ZmCCT9 enhances maize adaptation to higher latitudes. Proc Natl Acad Sci U S A 115, E334-e341 (2018).
Yang, Q. et al. CACTA-like transposable element in ZmCCT attenuated photoperiod sensitivity and accelerated the postdomestication spread of maize. Proc Natl Acad Sci U S A 110, 16969-74 (2013).
Tian, J. et al. Teosinte ligule allele narrows plant architecture and enhances high-density maize yields. Science 365, 658-664 (2019).
Marand, A.P., Eveland, A.L., Kaufmann, K. & Springer, N.M. cis-Regulatory Elements in Plant Development, Adaptation, and Evolution. Annu Rev Plant Biol (2022).
Buenrostro, J.D., Giresi, P.G., Zaba, L.C., Chang, H.Y. & Greenleaf, W.J. Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position. Nature Methods 10, 1213 (2013).
Klemm, S.L., Shipony, Z. & Greenleaf, W.J. Chromatin accessibility and the regulatory epigenome. Nature Reviews Genetics 20, 207-220 (2019).
Hufford, M.B. et al. De novo assembly, annotation, and comparative analysis of 26 diverse maize genomes. Science 373, 655-+ (2021).
Tu, X. et al. Reconstructing the maize leaf regulatory network using ChIP-seq data of 104 transcription factors. Nature Communications 11, 5089 (2020).
Crisp, P.A. et al. Stable unmethylated DNA demarcates expressed genes and their cis-regulatory space in plant genomes. Proc Natl Acad Sci U S A 117, 23991-24000 (2020).
Dong, P. et al. 3D Chromatin Architecture of Large Plant Genomes Determined by Local A/B Compartments. Mol Plant 10, 1497-1509 (2017).
Sun, Y. et al. 3D genome architecture coordinates trans and cis regulation of differentially expressed ear and tassel genes in maize. Genome Biol 21, 143 (2020).
Marand, A.P., Chen, Z., Gallavotti, A. & Schmitz, R.J. A cis-regulatory atlas in maize at single-cell resolution. Cell 184(2021).
Wirthlin, M.E. et al. Vocal learning-associated convergent evolution in mammalian proteins and regulatory elements. Science, eabn3263 (2024).
Andrews, G. et al. Mammalian evolution of human cis-regulatory elements and transcription factor binding sites. Science 380, eabn7930 (2023).
Kuderna, L.F.K. et al. Identification of constrained sequence elements across 239 primate genomes. Nature 625, 735-742 (2024).
Wu, Y. et al. Phylogenomic discovery of deleterious mutations facilitates hybrid potato breeding. Cell 186, 2313-2328.e15 (2023).
Hu, Y. et al. Comparative analysis reveals epigenomic evolution related to species traits and genomic imprinting in mammals. Innovation (Camb) 4, 100434 (2023).
Huang, W. et al. A well-supported nuclear phylogeny of Poaceae and implications for the evolution of C(4) photosynthesis. Mol Plant 15, 755-777 (2022).
Morgante, M. et al. Gene duplication and exon shuffling by helitron-like transposons generate intraspecies diversity in maize. Nature Genetics 37, 997-1002 (2005).
Horton, C.A. et al. Short tandem repeats bind transcription factors to tune eukaryotic gene expression. Science 381, eadd1250 (2023).
Jin, J. et al. PlantTFDB 4.0: toward a central hub for transcription factors and regulatory interactions in plants. Nucleic Acids Research 45, D1040-D1045 (2017).
Castro-Mondragon, J.A. et al. JASPAR 2022: the 9th release of the open-access database of transcription factor binding profiles. Nucleic Acids Res 50, D165-d173 (2022).
Galli, M. et al. The DNA binding landscape of the maize AUXIN RESPONSE FACTOR family. Nat Commun 9, 4526 (2018).
Han, J. et al. Genome-wide chromatin accessibility analysis unveils open chromatin convergent evolution during polyploidization in cotton. Proceedings of the National Academy of Sciences 119, e2209743119 (2022).
Chuck, G., Whipple, C., Jackson, D. & Hake, S. The maize SBP-box transcription factor encoded by tasselsheath4 regulates bract development and the establishment of meristem boundaries. Development 137, 1243-50 (2010).
Zhang, C., Li, X., Wang, Z., Zhang, Z. & Wu, Z. Identifying key regulatory genes of maize root growth and development by RNA sequencing. Genomics 112, 5157-5169 (2020).
Tao, S. et al. Single-Cell Transcriptome and Network Analyses Unveil Key Transcription Factors Regulating Mesophyll Cell Development in Maize. Genes 13, 374 (2022).
Bukowski, R. et al. Construction of the third-generation Zea mays haplotype map. Gigascience 7, 1-12 (2018).
Hufford, M.B. et al. Comparative population genomics of maize domestication and improvement. Nature Genetics 44, 808-811 (2012).
Renny-Byfield, S., Rodgers-Melnick, E. & Ross-Ibarra, J. Gene Fractionation and Function in the Ancient Subgenomes of Maize. Mol Biol Evol 34, 1825-1832 (2017).
Rumpho, M.E., Edwards, G.E. & Loescher, W.H. A pathway for photosynthetic carbon flow to mannitol in celery leaves : activity and localization of key enzymes. Plant Physiol 73, 869-73 (1983).
Hatfield, R.D., Rancour, D.M. & Marita, J.M. Grass Cell Walls: A Story of Cross-Linking. Front Plant Sci 7, 2056 (2016).
Faik, A. Xylan biosynthesis: news from the grass. Plant Physiol 153, 396-402 (2010).
Lai, J. et al. Genome-wide patterns of genetic variation among elite maize inbred lines. Nature Genetics 42, 1027-1030 (2010).
Gabur, I. et al. Gene presence-absence variation associates with quantitative Verticillium longisporum disease resistance in Brassica napus. Scientific Reports 10, 4131 (2020).
Ye, J., Zhang, L., Zhang, X., Wu, X. & Fang, R. Plant Defense Networks against Insect-Borne Pathogens. Trends Plant Sci 26, 272-287 (2021).
Han, G.Z. Origin and evolution of the plant immune system. New Phytol 222, 70-83 (2019).
Sigmon, B. & Vollbrecht, E. Evidence of selection at the ramosa1 locus during maize domestication. Mol Ecol 19, 1296-311 (2010).
Liao, B.Y. & Zhang, J. Low rates of expression profile divergence in highly expressed genes and tissue-specific genes during mammalian evolution. Mol Biol Evol 23, 1119-28 (2006).
Mantica, F. et al. Evolution of tissue-specific expression of ancestral genes across vertebrates and insects. Nature Ecology & Evolution 8, 1140-1153 (2024).
Tong, Y.-B. et al. GenOrigin: A comprehensive protein-coding gene origination database on the evolutionary timescale of life. Journal of Genetics and Genomics 48, 1122-1129 (2021).
Wang, B. et al. De novo genome assembly and analyses of 12 founder inbred lines provide insights into maize heterosis. Nature Genetics 55, 312-323 (2023).
Chen, J. et al. A complete telomere-to-telomere assembly of the maize genome. Nature Genetics 55, 1221-1231 (2023).
Buckler, E.S. et al. The Genetic Architecture of Maize Flowering Time. Science 325, 714-718 (2009).
Yang, N. et al. Two teosintes made modern maize. Science 382, eadg8940 (2023).
Liu, J., Fernie, A.R. & Yan, J. The Past, Present, and Future of Maize Improvement: Domestication, Genomics, and Functional Genomic Routes toward Crop Enhancement. Plant Communications 1, 100010 (2020).
Chen, W. et al. Convergent selection of a WD40 protein that enhances grain yield in maize and rice. Science 375, eabg7985 (2022).
Lemmon, Z.H., Bukowski, R., Sun, Q. & Doebley, J.F. The role of cis regulatory evolution in maize domestication. PLoS Genet 10, e1004745 (2014).
Gao, Y. et al. A maize phytochrome-interacting factors protein ZmPIF1 enhances drought tolerance by inducing stomatal closure and improves grain yield in Oryza sativa. Plant Biotechnol J 16, 1375-1387 (2018).
Jan, R. et al. Over-Expression of Chorismate Mutase Enhances the Accumulation of Salicylic Acid, Lignin, and Antioxidants in Response to the White-Backed Planthopper in Rice Plants. Antioxidants (Basel) 10(2021).
(!!! INVALID CITATION !!! 24,25).
(!!! INVALID CITATION !!! 24).
Fang, C. et al. Dynamics of cis-regulatory sequences and transcriptional divergence of duplicated genes in soybean. Proc Natl Acad Sci U S A 120, e2303836120 (2023).
Kong, D. et al. UB2/UB3/TSH4-anchored transcriptional networks regulate early maize inflorescence development in response to simulated shade. Plant Cell (2022).
Wills, D.M. et al. From many, one: genetic control of prolificacy during maize domestication. PLoS Genet 9, e1003604 (2013).
Schnable, J.C., Springer, N.M. & Freeling, M. Differentiation of the maize subgenomes by genome dominance and both ancient and ongoing gene loss. Proceedings of the National Academy of Sciences of the United States of America 108(2011).
Li, X. et al. Maize GOLDEN2-LIKE genes enhance biomass and grain yields in rice by improving photosynthesis and reducing photoinhibition. Communications Biology 3, 151 (2020).
Rodgers-Melnick, E., Vera, D.L., Bass, H.W. & Buckler, E.S. Open chromatin reveals the functional maize genome. Proceedings of the National Academy of Sciences 113, E3177-E3184 (2016).
Maurano, M.T. et al. Systematic localization of common disease-associated variation in regulatory DNA. Science 337, 1190-5 (2012).
(!!! INVALID CITATION !!! 9,36).
Li, C. et al. Genomic insights into historical improvement of heterotic groups during modern hybrid maize breeding. Nature Plants 8, 750-763 (2022).
Zhu, Y. et al. HEMU: an integrated comparative genomics database and analysis platform for Andropogoneae grasses. Plant Commun, 100786 (2023).
Zhu, Y., Wang, Z., Zhou, Z., Liu, Y. & Shi, J. HEMU: an integrated Andropogoneae comparative genomics database and analysis platform. (bioRxiv, 2023).
Yuzhi, Z., Zijie, W., Zanchen, Z., Yuting, L. & Junpeng, S. HEMU: an integrated Andropogoneae comparative genomics database and analysis platform, (2023).
Zhu, T. et al. Comprehensive mapping and modelling of the rice regulome landscape unveils the regulatory architecture underlying complex traits. Nature Communications 15, 6562 (2024).
Zuo, Z. et al. ATAC-pipe: general analysis of genome-wide chromatin accessibility. Brief Bioinform 20, 1934-1943 (2019).
Langmead, B., Trapnell, C., Pop, M. & Salzberg, S.L. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol 10, R25 (2009).
Stark, R. & Brown, G. DiffBind differential binding analysis of ChIP-Seq peak data. In R package version 100(2011).
Chen, S., Zhou, Y., Chen, Y. & Gu, J. fastp: an ultra-fast all-in-one FASTQ preprocessor. Bioinformatics 34, i884-i890 (2018).
Kim, D., Paggi, J.M., Park, C., Bennett, C. & Salzberg, S.L. Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype. Nature Biotechnology 37, 907-915 (2019).
Li, H. et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078-9 (2009).
Kovaka, S. et al. Transcriptome assembly from long-read RNA-seq alignments with StringTie2. Genome Biol 20, 278 (2019).
Zhang, Y. et al. Model-based analysis of ChIP-Seq (MACS). Genome Biol 9, R137 (2008).
Quinlan, A.R. & Hall, I.M. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics 26, 841-2 (2010).
Ou, S. et al. Benchmarking transposable element annotation methods for creation of a streamlined, comprehensive pipeline. Genome Biol 20, 275 (2019).
Yan, H., Bombarely, A. & Li, S. DeepTE: a computational method for de novo classification of transposons with convolutional neural network. Bioinformatics 36, 4269-4275 (2020).
Edgar, R.C. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Res 32, 1792-7 (2004).
Rice, P., Longden, I. & Bleasby, A. EMBOSS: the European Molecular Biology Open Software Suite. Trends Genet 16, 276-7 (2000).
Ma, J. & Bennetzen, J.L. Rapid recent growth and divergence of rice nuclear genomes. Proc Natl Acad Sci U S A 101, 12404-10 (2004).
Armstrong, J. et al. Progressive Cactus is a multiple-genome aligner for the thousand-genome era. Nature 587, 246-251 (2020).
Hinrichs, A.S. et al. The UCSC Genome Browser Database: update 2006. Nucleic Acids Research 34, D590-D598 (2006).
Hedges, S.B., Marin, J., Suleski, M., Paymer, M. & Kumar, S. Tree of Life Reveals Clock-Like Speciation and Diversification. Molecular Biology and Evolution 32, 835-845 (2015).
Mcinnes, L. & Healy, J. UMAP: Uniform Manifold Approximation and Projection for Dimension Reduction. The Journal of Open Source Software 3, 861 (2018).
Klopfenstein, D.V. et al. GOATOOLS: A Python library for Gene Ontology analyses. Scientific Reports 8, 10872 (2018).
Trapnell, C. et al. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc 7, 562-78 (2012).
Sun, S. et al. Extensive intraspecific gene order and gene structural variations between Mo17 and other maize genomes. Nat Genet 50, 1289-1295 (2018).
Kiełbasa, S.M., Wan, R., Sato, K., Horton, P. & Frith, M.C. Adaptive seeds tame genomic sequence comparison. Genome Res 21, 487-93 (2011).
Haas, B.J., Delcher, A.L., Wortman, J.R. & Salzberg, S.L. DAGchainer: a tool for mining segmental genome duplications and synteny. Bioinformatics 20, 3643-6 (2004).
Tang, H. et al. Screening synteny blocks in pairwise genome comparisons through integer programming. BMC Bioinformatics 12, 102 (2011).
Yang, Z. PAML 4: phylogenetic analysis by maximum likelihood. Mol Biol Evol 24, 1586-91 (2007).
Wang, Y. et al. MCScanX: a toolkit for detection and evolutionary analysis of gene synteny and collinearity. Nucleic Acids Res 40, e49 (2012).
Krzywinski, M.I. et al. Circos: An information aesthetic for comparative genomics. Genome Research (2009).
Purcell, S. et al. PLINK: a tool set for whole-genome association and population-based linkage analyses. Am J Hum Genet 81, 559-75 (2007).
Yang, J., Lee, S.H., Goddard, M.E. & Visscher, P.M. Genome-wide complex trait analysis (GCTA): methods, data analyses, and interpretations. Methods Mol Biol 1019, 215-36 (2013).
Kang, H.M. et al. Variance component model to account for sample structure in genome-wide association studies. Nat Genet 42, 348-54 (2010).

The supplementary figures and tables are not available with this version.

There is NO Competing Interest.

Download PDF

Version 1

posted

You are reading this latest preprint version

Inter- and intraspecific constraint of accessible chromatins maps regulatory loci involved in maize speciation and domestication

Status:

Version 1

Abstract

Figures

Introduction

Results

A pan-tissue map of accessible chromatins in maize

Evolutionary constraint of maize ACRs in other grass genomes

Intraspecific constraint of ACRs and their association with maize domestication

Constraint and dominance of ACRs between two maize subgenomes

Trait-associated SNPs are enriched in ACRs that highly constrained in maize

Discussion

Methods

Plant materials

Library construction and ATAC sequencing

ATAC-seq analysis

Samples clustering and principal component analysis

Transcriptome data analysis

ACR identification

Annotation of TEs and their association with ACRs

TF annotation and motif analysis

Evolutionary constraint analysis of ACRs

UMAP analysis on ACRs

GO analysis

Analysis of differentially expressed genes

Maize two subgenome identification

High-quality syntenic gene pairs between two subgenomes

Genome wide association studies (GWAS)

Declarations

References

Supplementary Figures and Tables

Additional Declarations

Status:

Version 1