Recapitulation-like developmental transitions of chromatin accessibility in vertebrates

Abstract

The relationship between development and evolution has been a central theme in evolutionary developmental biology. Across the vertebrates, the most highly conserved gene expression profiles are found at mid-embryonic, organogenesis stages, whereas those at earlier and later stages are more diverged. This hourglass-like pattern of divergence does not necessarily rule out the possibility that gene expression profiles that are more evolutionarily derived appear at later stages of development; however, no molecular-level evidence of such a phenomenon has been reported. To address this issue, we compared putative gene regulatory elements among different species within a phylum. We made a genome-wide assessment of accessible chromatin regions throughout embryogenesis in three vertebrate species (mouse, chicken, and medaka) and estimated the evolutionary ages of these regions to define their evolutionary origins on the phylogenetic tree. In all the three species, we found that genomic regions tend to become accessible in an order that parallels their phylogenetic history, with evolutionarily newer gene regulations activated at later developmental stages. This tendency was restricted only after the mid-embryonic, phylotypic periods. Our results imply a phylogenetic hierarchy of putative regulatory regions, in which their activation parallels the phylogenetic order of their appearance. One evolutionary mechanism that may explain this phenomenon is that newly introduced regulatory elements are more likely to survive if activated at later stages of embryogenesis. Possible relationships between this phenomenon and the so-called recapitulation are discussed.

Keywords: ATAC-seq; Development; Developmental hourglass model; Evolution; Gene regulatory evolution; Parallelism; Recapitulation.

Fig. 1

Strategy for assessing accessible chromatin landscapes in vertebrate embryos. a Genome browser views showing enrichment of whole-embryo ATAC-seq reads in representative regions of mouse, chicken, and medaka genomes, respectively. ATAC-seq read enrichment is presented as the mean of three biological replicates. Colors below the read enrichment represent the estimated evolutionary ages of genomic regions that correspond to the tracks of the evolutionary trajectories, which are shown as a phylogenetic tree on the right. b The representative ACRs (blue boxes) overlapping with annotated enhancers (red regions) acquired from the VISTA Enhancer Database [23]. For each enhancer, ATAC-seq read enrichment in E10.5 and E12.5 mice and in vivo enhancer activity in E11.5 mice are shown with the VISTA Enhancer ID, the flanking gene, and an embryo image from the VISTA Enhancer Database [23]. c Schematic diagram showing the three steps for estimating the relative ATAC-seq signal for each evolutionary age: (1) ACRs were identified by ATAC-seq signal intensity; (2) ACRs were categorized according to their estimated evolutionary ages; and (3) for each developmental stage, the percentage corresponding to the summed signal intensities of each evolutionary category divided by the total signal intensities for all evolutionary categories was calculated (relative ATAC-seq signals)

Fig. 2

Numbers of ACRs and expressed protein-coding genes categorized according to evolutionary ages. Stacked bar graphs show the numbers of evolutionarily categorized ACRs (a, d, g) and expressed (FPKM > 1) protein-coding genes (b, e, h) at each developmental stage in mouse (a, b), chicken (d, e), and medaka (g, h). Evolutionary ages of ACRs were estimated based on Method I (for details, see Methods and Additional file 2: Figure S4). The evolutionary ages of protein-coding genes were estimated according to the most recent common ancestors of all the species sharing the homologs; the expressed genes that were estimated to be lost secondarily in any of the compared species were excluded (see Methods for details). Colors in each stacked bar graph indicate the categories of the evolutionary ages of each element. Each evolutionary category includes ACRs or expressed protein-coding genes that originated during the correspondingly colored period in the phylogenetic trees shown in c, f, and i

Fig. 3

Transition of developmental stages with the maximum evolutionarily categorized chromatin accessibility during vertebrate embryogenesis. For each developmental stage in three vertebrate species (mouse, chicken, and medaka), the percentages on the vertical axis represent the summed signal intensity for each evolutionary category of ACRs divided by the total signal intensity for all categories (relative ATAC-seq signal). The color of each category indicates the estimated evolutionary age of the region (shown at right). In each graph, the developmental stages with the highest signal from the potential phylotypic period are highlighted in the corresponding colors, as is the range that showed a recapitulative pattern for unknown reasons. Error bars indicate the standard deviation of three biological replicates for each developmental stage. Changes in the relative ATAC-seq signals were statistically significant (Kruskal–Wallis rank sum test) in all cases, except for the vase tunicate category in medaka. Detailed statistical information is provided in Additional file 1: Table S7

Fig. 4

No recapitulative pattern was observed in developmental gene expression levels For each developmental stage, summed expression levels of evolutionarily categorized protein-coding genes were shown as the percentage relative to the total expression levels of all evolutionary categories. The color of each category corresponds to the estimated evolutionary age of the protein-coding genes (shown at right). The evolutionary age of each protein-coding gene was estimated according to the most recent common ancestors of all the species sharing the homologs. Genes that were estimated to be secondarily lost in any of the compared species were excluded (see Methods for details). In each graph, developmental stages with the highest value after the potential phylotypic period are highlighted in the corresponding color. Error bars indicate the standard deviation of biological replicates for each developmental stage. Statistical information for the Kruskal–Wallis rank sum test is given in Additional file 1: Table S8.