Global analysis of primary mesenchyme cell cis-regulatory modules by chromatin accessibility profiling

Abstract

Background: The developmental gene regulatory network (GRN) that underlies skeletogenesis in sea urchins and other echinoderms is a paradigm of GRN structure, function, and evolution. This transcriptional network is deployed selectively in skeleton-forming primary mesenchyme cells (PMCs) of the early embryo. To advance our understanding of this model developmental GRN, we used genome-wide chromatin accessibility profiling to identify and characterize PMC cis-regulatory modules (CRMs).

Results: ATAC-seq (Assay for Transposase-Accessible Chromatin using sequencing) analysis of purified PMCs provided a global picture of chromatin accessibility in these cells. We used both ATAC-seq and DNase-seq (DNase I hypersensitive site sequencing) to identify > 3000 sites that exhibited increased accessibility in PMCs relative to other embryonic cell lineages, and provide both computational and experimental evidence that a large fraction of these sites represent bona fide skeletogenic CRMs. Putative PMC CRMs were preferentially located near genes differentially expressed by PMCs and consensus binding sites for two key transcription factors in the PMC GRN, Alx1 and Ets1, were enriched in these CRMs. Moreover, a high proportion of candidate CRMs drove reporter gene expression specifically in PMCs in transgenic embryos. Surprisingly, we found that PMC CRMs were partially open in other embryonic lineages and exhibited hyperaccessibility as early as the 128-cell stage.

Conclusions: Our work provides a comprehensive picture of chromatin accessibility in an early embryonic cell lineage. By identifying thousands of candidate PMC CRMs, we significantly enhance the utility of the sea urchin skeletogenic network as a general model of GRN architecture and evolution. Our work also shows that differential chromatin accessibility, which has been used for the high-throughput identification of enhancers in differentiated cell types, is a powerful approach for the identification of CRMs in early embryonic cells. Lastly, we conclude that in the sea urchin embryo, CRMs that control the cell type-specific expression of effector genes are hyperaccessible several hours in advance of gene activation.

Keywords: ATAC-seq; DNase-seq; Gene regulatory network; Primary mesenchyme cells; Sea urchins; Skeletogenesis.

Fig. 1

ATAC-seq sample preparation and sequence analysis. a S. purpuratus embryos were cultured for 24 h at 15 °C in triplicate. PMCs and other cells were isolated and ATAC-seq libraries were generated and sequenced. Sequence reads were analyzed by the bioinformatics pipeline shown in Additional file 1: Fig. S1A. b Examples of ATAC-seq differential peaks. These differential peaks (yellow rectangles) are located near the Sp-kirrelL gene, which encodes a transmembrane protein required for PMC fusion [27]. The aligned reads for each replicate are visualized, and the difference in peak magnitudes can be seen when comparing differential peaks in the isolated PMC replicates (light green peak trace) to the other cell replicates (dark green trace). Nominal p-values for differential peaks are indicated. c The distribution of ATAC-seq peaks in the RPS with respect to the closest gene. See Methods for definitions of peak locations. d Distribution of ATAC-seq differential peaks with respect to the closest gene

Fig. 2

DNase-seq sample preparation and sequence analysis. a S. purpuratus embryos were treated with U0126 at the 2-cell stage to obtain PMC (−) embryos. Control and U0126-treated embryos were cultured for 28 h at 15 °C in triplicate. Nuclei were isolated and DNase-seq was carried out. Sequence reads were analyzed by the bioinformatics pipeline shown in Additional file 1: Figure S1A. b An example of DNase-seq differential peaks. The differential peaks (yellow rectangles) are located near the WHL22.245306 transcript. The aligned reads for each replicate are visualized as traces, and the differences in peak magnitude are clear when comparing control whole embryos (violet peak trace) to PMC(−) embryos (dark purple trace). Nominal p-values for differential peaks are indicated. c Distribution of DNase-seq peaks in the RPS with respect to the closest gene. See Methods for definitions of peak locations. d Distribution of DNase-seq differential peaks with respect to the closest gene

Fig. 3

Previously analyzed PMC-specific cis-regulatory modules. Previously analyzed cis-regulatory modules (CRMs) (orange and yellow rectangles) that control the spatio-temporal expression of four PMC DE genes are represented in both ATAC-seq (light green trace: isolated PMCs; dark green trace: other cells) and DNase-seq (violet trace: control whole embryos, purple trace: U0126-treated embryos) datasets. a Sp-alx1 cis-regulatory modules are not represented as differential peaks in the ATAC-seq or DNase-seq datasets. b The Sp-sm50 enhancer (orange rectangle) and the minimal element (yellow rectangle) required for correct spatio-temporal expression of Sp-sm50 are encompassed within a differential peak identified in the DNase-seq dataset (violet rectangle), but not identified as differential in the ATAC-seq dataset. c The Sp-sm30a enhancer overlaps a differential peak identified in the ATAC- seq dataset (light green), but is not identified as differential in the DNase-seq dataset. d Two of four previously studied Sp-tbr cis-regulatory modules overlap 2 differential peaks in the ATAC-seq (light green) dataset and 1 differential peak in the DNase-seq (violet) dataset

Fig. 4

Experimental validation of putative PMC CRMs. a The EpGFPII reporter construct: Of a total of 3073 PMC-enriched differential peaks identified using DNase-seq and ATAC-seq, 31 were cloned into the EpGFPII plasmid, upstream of the GFP coding sequence and the Sp-endo16 promoter, and injected into S. purpuratus eggs. b Representative images of S. purpuratus embryos injected with 7 reporter constructs, showing PMC-specific GFP expression (green fluorescence) at 48 hpf. DIC and Sp-kirrelL images show the same embryo; arrows indicate PMCs

Fig. 5

Computational analysis of high-confidence PMC CRMs. a ATAC-seq differential peaks (green rectangles) and DNase-seq differential peaks (violet rectangles) near Sp-p16 and Sp-mitf, both PMC DE genes. Aligned reads averaged across replicates, from isolated PMCs (light green trace) and non-PMC cells (dark green trace) using ATAC-seq, and control 28 hpf embryos (violet trace) and PMC (−) embryos (dark purple trace) using DNase-seq, are shown. b Temporal expression profiles (Tu et al., 2012) of 420 PMC DE genes identified previously (Rafiq et al., 2014). Each gene is represented by a single row. The color scale ranges from deep red (2.5-fold higher than mean expression) to deep blue (2.5-fold lower than mean expression). White indicates mean expression. Four clusters are delineated, corresponding to maximal gene expression at 0–10, 40–72, 24–40 and 18–24 hpf, respectively. c Temporal expression of the 62 PMC DE genes within 10 kb of overlapping differential peaks: these PMC DE genes were classified into four clusters, delineated in Fig. 3b. d PMC DE genes were classified into categories based on levels of gene expression in isolated PMCs (data obtained from (Rafiq et al., 2014). “High” expression genes: FPKM between 2512 and 100 (top 17% of all 420 DE genes); “very low” expression genes: FPKM between 14 and 0 (bottom 25% of all 420 DE genes)

Fig. 6

Examples of overlapping differential peaks accessible at the 128-cell stage. Overlapping differential peaks (yellow rectangles) around the Sp-msp130r gene are accessible at the 128-cell stage (red rectangles represent peaks called at the 128- cell stage). Hypersensitivity corresponding to the overlapping differential peaks is seen at the 128-cell stage (red trace), the 24 hpf stage isolated PMCs (light green trace) and other non-PMC cells (dark green trace), and control 28 hpf embryos (violet trace) and PMC (−) embryos (dark purple trace)