Combinatorial chromatin dynamics foster accurate cardiopharyngeal fate choices (original) (raw)

Abstract

During embryogenesis, chromatin accessibility profiles control lineage-specific gene expression by modulating transcription, thus impacting multipotent progenitor states and subsequent fate choices. Subsets of cardiac and pharyngeal/head muscles share a common origin in the cardiopharyngeal mesoderm, but the chromatin landscapes that govern multipotent progenitors competence and early fate choices remain largely elusive. Here, we leveraged the simplicity of the chordate model Ciona to profile chromatin accessibility through stereotyped transitions from naive _Mesp_+ mesoderm to distinct fate-restricted heart and pharyngeal muscle precursors. An FGF-Foxf pathway acts in multipotent progenitors to establish cardiopharyngeal-specific patterns of accessibility, which govern later heart vs. pharyngeal muscle-specific expression profiles, demonstrating extensive spatiotemporal decoupling between early cardiopharyngeal enhancer accessibility and late cell-type-specific activity. We found that multiple _cis_-regulatory elements, with distinct chromatin accessibility profiles and motif compositions, are required to activate Ebf and Tbx1/10, two key determinants of cardiopharyngeal fate choices. We propose that these ‘combined enhancers’ foster spatially and temporally accurate fate choices, by increasing the repertoire of regulatory inputs that control gene expression, through either accessibility and/or activity.

Research organism: C. intestinalis

Introduction

How a species’ genome encodes its diverse and specific biological features has fascinated generations of biologists, and answers regarding the genetic control of body plan, organ, tissue and cell type formation have emerged from steady progress in developmental biology. Cell types arise as cells divide and the progeny of pluripotent embryonic stem cells progress through multipotent and fate-restricted states. The ontogeny of diverse terminal cell identities involves differential expression of hundreds to thousands of genes. Their dynamic activities are orchestrated by complex gene regulatory networks, whereby DNA-binding proteins and co-factors act upon specific _cis_-regulatory elements to control gene expression (Davidson, 2010). Technical and conceptual revolutions in genome biology have extensively characterized the chromatin dynamics that govern the function of _cis_-regulatory elements (Klemm et al., 2019). Specifically, as the nuclear genome is packaged in nucleosomes, DNA-binding transcription factors compete with histones to interact with _cis_-regulatory elements and control gene expression. Thus, identifying changing landscapes of accessible chromatin governing the transition from multipotent to fate restricted progenitors offers privileged insights into the genomic code for progressive cell type specification.

Dynamic chromatin states underlying cardiomyocyte differentiation have been extensively profiled (Paige et al., 2012; Wamstad et al., 2012), and chromatin state regulation is essential for heart development (He et al., 2014; Rosa-Garrido et al., 2013; Zaidi et al., 2013). However, different parts of the heart originate from separate first and second fields of progenitor cells, including those referred to as cardiopharyngeal, which can also produce branchiomeric head muscles (Diogo et al., 2015; Lescroart et al., 2010). Bulk and single cell transcription profiling have begun to illuminate gene expression changes underlying cardiopharyngeal fate choices (Lescroart et al., 2018; Lescroart et al., 2014), but the corresponding chromatin dynamics remains largely elusive.

The tunicate Ciona emerged as a powerful chordate model to study early cardiopharyngeal development with high spatio-temporal resolution (Diogo et al., 2015; Kaplan et al., 2015). In Ciona, the cardiopharyngeal lineages arise from naive Mesp+ mesodermal progenitors that emerge at the onset of gastrulation, and divide into two multipotent cardiopharyngeal progenitors (aka trunk ventral cells, TVCs) and two anterior tail muscles (ATMs), on either side of the embryo (Figure 1A). Following induction by FGF-MAPK signaling, cardiopharyngeal progenitors migrate collectively, before dividing asymmetrically and medio-laterally to produce small median first heart precursors (FHPs), and large lateral second trunk ventral cells (STVCs) (Davidson et al., 2005; Stolfi et al., 2010; Wang et al., 2013). The latter are also multipotent cardiopharyngeal progenitors, which upregulate Tbx1/10 and then divide again to produce small median second heart precursors (SHPs), and large lateral atrial siphon muscle founder cells (ASMFs). ASMFs activate Ebf, which is necessary and sufficient to induce pharyngeal muscle specification (Razy-Krajka et al., 2014; Stolfi et al., 2014; Stolfi et al., 2010; Tolkin and Christiaen, 2016). Importantly, spatially and temporally accurate activation of Tbx1/10 and Ebf in the STVC and ASMF, respectively, is essential to permit the emergence of all cardiopharyngeal cell lineages, as their ectopic expression would inhibit proper heart fate specification (Figure 1A).

Figure 1. Profiling chromatin accessibility dynamics during early cardiopharyngeal cell development.

(A) Embryos, larvae and lineage diagram showing B7.5 blastomeres, their cardiopharyngeal progeny, and the main stages sampled for ATAC-seq. Anterior tail muscle (ATM, gray), trunk ventral cell (TVC, green), secondary TVC (STVC, yellow), first heart precursor (FHP, red), second heart precursor (SHP, orange), atrial siphon precursor cells (ASMF, blue). Stages (St.) according to Hotta et al. (2007) with hours post fertilization (hpf). (B) Spearman correlation of RPKM (reads per kb per million mapped reads) values in 14,178 regions changing accessibility over time or between B7.5 and B-line mesenchyme lineages. (C) Temporal changes in chromatin accessibility for 5,450 regions. ‘B7.5 6 > 10’: 3,691 regions more accessible at Mesp>LacZ 6 hpf than Mesp>LacZ 10 hpf. ‘B7.5 6 < 10’: 1,759 regions more accessible at Mesp>LacZ 10 than Mesp>LacZ 6 hpf. The accessibility of these regions is shown for Mesp>LacZ 6 hpf, Mesp>LacZ 10 hpf, and Hand-r>LacZ 18 hpf vs. the average (avg) accessibility in the control cells. Cell-type-specific chromatin accessibility is shown in the comparison of Mesp>LacZ and MyoD905>GFP at 10 and Hand-r>LacZ and MyoD905>GFP 18 hpf. (D) Gene Set Enrichment Analysis (GSEA) normalized enrichment score of defined gene sets in regions ranked by difference in accessibility between time points as indicated (see Materials and methods).

Figure 1.

Figure 1—figure supplement 1. General characterization of the accessome.

Figure 1—figure supplement 1.

(A) Fragment size for B7.5 control samples at developmental time indicated (first three panels from the left) showing consistent approximately 147 bp periodicity, likely corresponding to nucleosome-protected fragments. This periodicity is not present in the purified genomic DNA (gDNA) input control (far right panel). Dotted lines show the cumulative fraction of reads larger than a given size. (B) Proportion of genomic features covered by the accessome and proportion of the accessome covered by genomic features (in bp). (C) Two-tailed binomial test for enrichment of accessible elements overlapping a genomic feature. Accessible elements associated with any DE gene from any scRNA-seq or bulk RNA-seq experiment were considered DE peaks. Bars show the predicted (based on bp coverage of the genome by a feature) and observed probabilities that an accessible element will overlap a genomic feature. Error bars show the 99% confidence interval. (D) Comparison of GC content of genomic features and accessible regions overlapping these features. (E) A 6 kb region displaying MACS2-called peaks in seven ATAC-seq libraries (upper panel), gene model (black) and accessome (gray). (F) Accessible element size distribution.

Figure 1—figure supplement 2. Characterization of promoter regions.

Figure 1—figure supplement 2.

(A) Read density within 2 kb upstream or downstream of TSS shows global closing of TSS over time. (B) Accessibility of the TSS is maintained in response to FgfrDN perturbation at 10 hpf. (C) Intersections of overlapping promoter elements. ‘ATAC-seq promoter peaks’ are accessible elements overlapping the putative promoter (−1,107 to +107 bp from the TSS). ‘DE promoters’ are differentially expressed putative promoters in any bulk RNA-seq condition (see Materials and methods). (D) Two-tailed binomial test for enrichment of accessible elements overlapping TSS-seq sites by genomic feature. Bars show the observed probability that an accessible element overlapping a genomic feature will also overlap a TSS-seq element. Error bars show the 99% confidence interval. The expected probability was calculated from the percent of the accessome overlapping TSS-seq elements (3.8%). (E) One-tailed hypergeometric test for enrichment of promoter motifs in genomic features. Eukaryotic core promoter motifs were taken from Haberle and Stark (2018). (F) Spearman correlation of inferred motif enrichment (see Materials and methods). Motifs in each class of element were ranked based on a one-tailed hypergeometric test.

Figure 1—figure supplement 3. Annotation of the accessome.

Figure 1—figure supplement 3.

(A) Example diagram of accessible elements surrounding a gene locus. (B) Distribution of number of genes associated with each peak. (C) Distribution of number of peaks associated with each gene. (D) Gene size vs. number of associated peaks. Red and blue dots show all transcription factor (TF) and signaling molecule (SM) genes, respectively. (E) Density of peaks per kb of each gene. (F) Two-tailed binomial test for enrichment of accessible elements associated with TF or signaling molecule genes. Predicted probability of an element associating to a TF or signaling molecule is given by the proportion of bps in all gene bodies ± 10 kb covering TFs or signaling molecule gene bodies ± 10 kb. Error bars show the 99% confidence interval of the binomial test. Only genes that were differentially expressed in any scRNA-seq or bulkRNA-seq comparison (Wang et al., 2019) were considered.

Building on previous extensive transcription profiles (Christiaen et al., 2008; Razy-Krajka et al., 2014), including single cell RNA-seq (scRNA-seq) from multipotent cardiopharyngeal progenitors to first and second heart lineages and pharyngeal muscle precursors (Wang et al., 2019), here we characterized the genome-wide chromatin accessibility dynamics underlying cardiopharyngeal fate specification. We identified regulatory inputs that govern _cis_-regulatory element accessibility and activity, as well as cell-type-specific enhancers for key cardiopharyngeal determinants. We found that, in multipotent progenitors, an FGF-Foxf pathway controls cardiopharyngeal-specific patterns of accessibility, which govern later heart vs. pharyngeal muscle-specific expression profiles. We further characterized temporal patterns of chromatin accessibility during cardiopharyngeal development. In particular, activation of fate determinants Tbx1/10 and Ebf specifically in the STVCs and ASMF, respectively, require multiple _cis_-regulatory elements with distinct spatio-temporal patterns of accessibility, which precede gene expression. We propose that these elements function as ‘combined enhancers’, which mediate distinct inputs, including from determinants of chromatin accessibility, to regulate gene activation. The observation that _cis_-regulatory inputs from multiple elements control expression of a single gene is consistent with the ‘shadow-’ and ‘super-enhancer’ paradigms (Barolo, 2012; Hnisz et al., 2013; Hong et al., 2008; Kvon et al., 2014; Perry et al., 2011; Perry et al., 2010; Pott and Lieb, 2015; Whyte et al., 2013). However, while shadow enhancer promotes robust transcription through the actions of multiple elements mediating similar regulatory inputs (Frankel, 2012; Frankel et al., 2010; Lam et al., 2015; Zeitlinger et al., 2007), we propose that combined enhancers promote spatially and temporally accurate fate choices, by augmenting the repertoire of _trans_-acting inputs controlling gene activation through enhancer activity and/or chromatin accessibility.

Results

A reference accessome for cardiopharyngeal development

To characterize the chromatin landscape underlying early cardiopharyngeal development, we used the assay for transposon-accessible chromatin (ATAC-seq; Buenrostro et al., 2013) on lineage-specific samples isolated at successive time points, and following defined perturbations (Figure 1A; Supplementary file 1; Razy-Krajka et al., 2018; Wang et al., 2019). Using the B7.5 lineage-specific Mesp>tagRFP reporter (Wang et al., 2018), we used FACS to collect ~4,000 cells per biological replicate from embryos dissociated at five time points encompassing key transitions in cardiopharyngeal development (Figure 1A): naive _Mesp_+ mesoderm (aka founder cells; Cooley et al., 2011), ATMs, TVCs, STVCs as well as fate-restricted first and second heart precursors (FHPs and SHPs), and pharyngeal muscle precursors (aka atrial siphon muscle founder cells -ASMF- and their progeny, the ASM precursors -ASMP; Razy-Krajka et al., 2014). The latter fate-restricted progenitors were obtained from larvae dissociated at three time points (15, 18 and 20 hours post-fertilization, hpf; Figure 1A). From the same embryonic cell populations, we used a co-transfected MyoD905>GFP reporter to isolate B-line mesenchymal cells (Christiaen et al., 2008). In the present study, we predominantly focused our analysis on the cardiopharyngeal progenitors at 6, 10 and 18 hpf (Figure 1A).

We obtained ~500 million unique ATAC-seq reads, with fragment-size distributions showing the characteristic ~150 bp periodicity and patterns of mono-, di- and tri-nucleosomal fragments (Buenrostro et al., 2013), which were absent in the genomic DNA control (Figure 1—figure supplement 1A). We identified ATAC-seq peaks using MACS2 (Zhang et al., 2008), and generated a combined atlas of 56,090 unique and non-overlapping accessible regions covering 9.25% of the C. robusta genome, which we used as our reference ‘accessome’ (Figure 1—figure supplement 1B,E). General metrics including peak numbers, size, GC content and genomic distribution were comparable to consensus peaksets reported in other studies of chromatin accessibility in developmental contexts (Materials and methods; Figure 1—figure supplement 1; Daugherty et al., 2017; Hockman et al., 2019; Jänes et al., 2018; Li et al., 2007; Madgwick et al., 2019).

Next, we annotated the reference accessome by associating accessible regions with other genomic features, especially gene models. In Ciona, the transcripts of approximately half of the protein-coding genes undergo spliced-leader (SL) _trans_-splicing, which replaces the original 5’ end sequence of pre-mRNAs by a short non-coding RNA, causing the 5’ end of mRNAs to differ from the transcription start site (TSS) (Ganot et al., 2004; Hastings, 2005; Satou et al., 2006; Vandenberghe et al., 2001). Using annotated TSSs (Satou et al., 2006; Vandenberghe et al., 2001; Yokomori et al., 2016), RNA-seq datasets (Wang et al., 2019), and our ATAC-seq data (Figure 1—figure supplement 2C), we determined that promoter regions and 5’ untranslated regions (5’UTR) were over-represented in the accessome (p < 0.001, two-tailed binomial test; Figure 1—figure supplement 1C), indicating that promoter proximal regions tend to be accessible as observed in other systems (Mayran et al., 2018). We also detected nucleosome footprints immediately upstream of TSSs, consistent with a tendency for constitutive accessibility (Figure 1—figure supplement 2A,B; Mavrich et al., 2008a; Mavrich et al., 2008b). By contrast, intronic and intergenic regions were significantly under-represented in our reference accessome, compared to the whole genome, although they were the most abundant elements (32.8% and 20.8%, respectively; Figure 1—figure supplement 1B). This suggests that most of these elements are accessible in specific contexts, as expected for tissue-specific _cis_-regulatory elements (Long et al., 2016).

We associated annotated genes with ATAC-seq peaks located within 10 kb of the TSS or transcription termination site (TTS) (Figure 1—figure supplement 3A) (Brozovic et al., 2018), thus assigning median values of 11 peaks per gene, and three genes per peak, owing to the compact Ciona genome (Figure 1—figure supplement 3B,C). Notably, active regulatory genes encoding transcription factors and signaling molecules were associated with significantly more peaks than other expressed genes (p < 0.001, two-tailed binomial test; Figure 1—figure supplement 3F). This high peak density surrounding regulatory genes is reminiscent of previously described super-enhancers (Whyte et al., 2013) and Clusters of Open _Cis_-Regulatory Elements (COREs) surrounding developmental regulators (Gaulton et al., 2010; Khan et al., 2018; Pott and Lieb, 2015).

Cardiopharyngeal accessibility profiles are established in multipotent progenitors

Using this reference accessome, we investigated lineage-specific and dynamic patterns of chromatin accessibility during fate decisions. We observed the greatest contrast in accessibility between the B7.5 and B-line mesenchyme lineages, with biological replicates correlating most highly (Spearman’s ρ > 0.93), indicating reproducible detection of extensive lineage-specific accessibility (Figure 1B). Within the B7.5 lineage, correlation analysis suggested that most changes occur between 6 and 10 hpf, during the transition from naive _Mesp_+ mesoderm to multipotent cardiopharyngeal progenitors (TVCs). Higher correlation between multipotent progenitors and mixed heart and pharyngeal muscle precursors, obtained from 18 hpf larvae, suggested more stable accessibility profiles during and immediately following early cardiopharyngeal fate choices (Figure 1B). Consistent with correlation analyses, most significant temporal changes in accessibility occurred during the transition from naive Mesp+ mesoderm to multipotent progenitors (5,450 regions, FDR < 0.05; Figure 1C). Specifically, about two thirds (64.7%, 3,525/5,450) of these regions showed reduced accessibility at 10 hpf, in multipotent progenitors, compared to 6 hpf naive _Mesp_+ mesoderm (Figure 1C). Conversely, 1,252 regions become accessible between 6 and 10 hpf or later, and 38.8% (486/1,252) of these regions were more accessible in the B7.5-lineage compared to the mesenchyme (Figure 1C). Moreover, the subset of regions opening between 6 and 10 hpf or later was enriched in genomic elements associated with cardiopharyngeal markers, including primed pan-cardiac and pharyngeal muscle markers, while elements flanking tail muscle markers (ATMs) or multipotent progenitor-specific genes were predominantly closing between 6 and 10 to 18 hpf (Figure 1D; see Materials and methods for definition of gene sets). Taken together, these observations suggest that cardiopharyngeal accessibility profiles are established specifically in the B7.5 lineage, upon induction of multipotent progenitors, and persist in fate-restricted cells.

To further analyze changes in accessibility associated with multipotent progenitor induction, we performed ATAC-seq on B7.5 lineage cells isolated at 10 hpf following defined perturbations of FGF-MAPK signaling: a constitutively active form of Mek (MekS216D,S220E), which converts all B7.5 lineage cells into multipotent cardiopharyngeal progenitors or a dominant negative form of the Fgf receptor (FgfrDN), which blocks induction and transforms all B7.5-derived cells into ATMs (Figure 1A; Davidson et al., 2006; Razy-Krajka et al., 2018). We used DESeq2 (Love et al., 2014) to compute differential accessibility of the elements in the reference accessome, and identified 2,728 and 2,491 differentially accessible regions following either inhibition or activation of FGF-MAPK signaling, respectively (Figure 2A,B; Figure 2—figure supplement 1A). Using peak-to-gene annotations (Figure 1—figure supplement 3A), we cross-referenced ATAC-seq with expression microarray data obtained from B7.5 lineage cells expressing the same FgfrDN(Christiaen et al., 2008), and observed a positive correlation between changes in differential accessibility and differential gene expression at 10 hpf (Spearman’s ρ = 0.47; Figure 2A; Figure 2—figure supplement 2A). Specifically, 48% of FGF-MAPK-regulated genes were associated with at least one element showing consistent differential accessibility, including 260 candidate FGF-MAPK-activated TVC markers associated with 557 regions predicted to open specifically in multipotent cardiopharyngeal progenitors at 10 hpf (Supplementary file 2, Figure 2—figure supplement 2B). Conversely, the majority (603 regions) of differentially accessible ATAC-seq peaks associated with 263 FGF-MAPK-inhibited tail muscles markers, and were also more accessible upon inhibition of FGF signaling (Supplementary file 2; Figure 2—figure supplement 1A). Taken together, these observations indicate that cardiopharyngeal accessibility profiles are established in the multipotent progenitors by opening regions associated with genes upregulated upon induction by FGF-MAPK signaling.

Figure 2. Cardiopharyngeal accessibility profiles are established in multipotent progenitors.

(A–B) Correlations between differential gene expression (DE) and differential chromatin accessibility (DA) in response to FGF-MAPK perturbation in the multipotent progenitors (A) and between chromatin accessibility in response to FGF-MAPK perturbation and in multipotent progenitors (10 hpf) versus founder cells (6 hpf). (B). Colored dots are DA peaks associated with cell type-specific DE genes. ρ is the Spearman correlation of expression and accessibility for DA regions associated with DE genes (A) or of region response to MAPK perturbation with accessibility in founder cells versus multipotent progenitors (B). (C) Relationship between expression and accessibility of DE genes associated with DA regions for genes in the bottom 0.75% quantile of fold change between expression in FgfrDN and control (log2(FC) < −1.32). Microarray log2(fold change (FC)) values are shown on the left. The fold change for all time points is versus the average. (**D**) A 24 kb region on chromosome eight displaying expression (RNA-seq) and chromatin accessibility (ATAC-seq; normalized by total sequencing depth). Gray shaded boxes show validated ATM-specific promoters and a newly identified TVC-specific enhancer in _Nk4/Nkx2-5_ intron. (**E**) Enhancer-driven in vivo reporter expression (green) of tested ATAC-seq regions (_KhC8.2200_ and ._2201_). TVCs marked with _Mesp_>H2B::mCherry (red). Numbers indicate observed/total of half-embryos scored. (F) Endogenous expression of Nk4/Nkx2-5 visualized by in situ (green) in TyrosinaseCRISPR and upon CRISPR/Cas9-induced deletions of TVC-specific region. Nuclei of B7.5 lineage cells are labelled by Mesp>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). Nk4/Nkx2-5 expression was not affected in the epidermis (open arrowhead). Experiment performed in biological replicates. Scale bar, 20 μm. Fisher exact test, total numbers of individual halves scored per condition are shown in 'n='. Gene expression data for 6 hpf and ‘FGF-MAPK perturbation 10 hpf’ (Christiaen et al., 2008) and 8 to 20 hpf (Razy-Krajka et al., 2014) were previously published.

Figure 2.

Figure 2—figure supplement 1. Inhibition of FGF signaling (_Mesp>FgfrDN_−10 hpf) induces opening of ATM-specific elements.

Figure 2—figure supplement 1.

(A) Relation between expression and accessibility for DE genes in the top 99.25% quantile of fold change between expression in FgfrDN and control at 10 hpf (log2(FC) > 1.31). (B, C) A 10 kb region on chromosome 6 and 10 including Synpo1/2 (B) and KH.C10.125 (C) displaying chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). The newly identified enhancers are in boxed regions. ATAC-seq peaks were tested in vivo by reporter gene assay (see Supplementary file 4). Reporter gene assay in embryos at stage 23 (Hotta et al., 2007). GFP expression was detected specifically in ATM cells and not in TVCs. Nuclei of B7.5 lineage cells are labelled by Mesp>H2B::mCherry and Mesp>hCD4::mCherry (red). Numbers indicate observed/total of half-embryos scored. Scale bar = 30 μm.

Figure 2—figure supplement 2. General characterization of differential accessibility.

Figure 2—figure supplement 2.

(A) Spearman correlation of expression of bulk RNA-seq or microarray (for B7.5 6 vs. 10 hpf and Mesp>FgfrDN vs. Mesp>LacZ - 10 hpf) pairwise comparisons with ATAC-seq pairwise comparisons. Correlation was calculated based on log2(FC) of differentially expressed genes associated with differentially accessible peaks for each comparison. Differentially expressed genes in Mesp>FgfrDN vs. Mesp>LacZ at 10 hpf derived from Christiaen et al. (2008), Hand-r>FgfrDN vs. Hand-r>LacZ and Foxf>M RasCA vs. Foxf>LacZ at 18 hpf from Wang et al. (2019), FoxfCRISPR vs. ControlCRISPR at 10 hpf from the present study. (B) GSEA for ATAC-seq pairwise comparisons. A negative normalized enrichment score (NES) indicates elements annotated to that gene set are less accessible in the comparison. A positive NES indicates elements annotated to that gene set are more accessible in the comparison. (C) One-tailed hypergeometric test for enrichment of accessible elements overlapping a genomic feature. The leftmost column shows the expected distribution of accessible elements from the accessome.

Figure 2—figure supplement 3. Peakshift validation of sgRNA efficiency.

Figure 2—figure supplement 3.

(A) Chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). Transcript model is indicated as black bar. The newly identified TVC-specific enhancer is in boxed regions; in the zoom region, blue arrows indicate primers used to amplify the region between the target sites. In wild-type embryos, the resulting PCR product is ~1.1 kilobase pairs. The two single guide RNAs used to target Nk4/Nkx2-5 intronic element are in red. (B) Alignment of cloned PCR products amplified using the primers indicated in (A), from wild-type (wt) embryos, and from embryos electroporated with 25 µg Ef1α>nls::Cas9-gem::nls and 40 µg each of U6>sgRNA2 and U6>sgRNA4. Clone '01' contains a ~0.6 kilobase deletion between the approximate sites targeted by the two sgRNAs. (C) Quantification of indel-shifted electrophoresis chromatogram peaks (‘Peakshift’ assay; Gandhi et al., 2017) revealed sgRNA mutagenesis efficacies. (D) In situ hybridization for Nk4 (green) showing expression throughout the ventral head endoderm in embryos electroporated with _Mesp_>H2B::mCherry (red), _Mesp_>nls::Cas9-Gem::nls and _U6_>sgRNAs targeting Nk4 intron (Nk4_KhC8.2200-.2201CRISPR) or Tyrosinase (TyrosinaseCRISPR), used as control. In control embryos (panels on the left) Nk4 expression is essentially wild-type with a strong expression in ventral head endoderm (asterix and arrowhead) and TVCs. In Nk4_KhC8.2200-.2201CRISPR embryos (panels on the right) Nk4 expression is lost specifically in the B7.5 lineage and not in the other endogenous territories. Scale bars = 30 µm.

Figure 2—figure supplement 4. Candidate TVC-specific enhancers in vivo validation by reporter gene assay.

Figure 2—figure supplement 4.

(A–G) ATAC-seq peaks specifically accessible in TVC displaying chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). Transcript model is indicated as black bar. The newly identified enhancers are in boxed regions. ATAC-seq peaks were tested in vivo by reporter gene assay (see Supplementary file 4). (A’–G’’) Reporter gene assay in embryos at stage 23 (Hotta et al., 2007). GFP expression was detected specifically in TVC cells but not in ATMs. Nuclei of B7.5 lineage cells are labelled by Mesp>H2B::mCherry and Mesp>hCD4::mCherry (red). Numbers indicate observed/total of half-embryos scored. Scale bar = 30 μm.

Figure 2—figure supplement 5. Candidate TVC-specific enhancers in vivo validation by CRISPR/Cas9.

Figure 2—figure supplement 5.

(A–C) TVC-specific enhancers in Fbln, Smurf and Fgf4 loci targeted by two or three single guide RNA (sgRNAs) for CRISPR/Cas9-mediated deletions. Upper panel is showing chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). Transcript model is indicated as black bar. (D–F) Endogenous expression of Fbln (D), Smurf (E) and Fgf4 (F) visualized by in situ (green) in TyrosinaseCRISPR (right panel) and upon CRISPR/Cas9-induced deletions of TVC-specific peaks (left panel). Nuclei of B7.5 lineage cells are labelled by Mesp>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). Anterior to the left, stages indicated as 'st.'. Experiment performed in biological replicates. Scale bar = 20 μm. All the treatments were significant versus TyrosinaseCRISPR (Fisher exact test, p < 0.001); ‘n’: total numbers of embryos.

Figure 2—figure supplement 6. CRISPR validation on the TVC-specific Fgf4 enhancer.

Figure 2—figure supplement 6.

(A) Chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). Transcript model is indicated as black bar. The newly identified TVC-specific enhancer is in boxed regions; in the zoom region, blue arrows indicate primers used to amplify the region between the target sites. In wild-type embryos, the resulting PCR product is ~1.2 kilobase pairs. The two single guide RNAs used to target Fgf4 intronic element are in red. (B) Quantification of indel-shifted electrophoresis chromatogram peaks (‘Peakshift’ assay; Gandhi et al., 2017) revealed sgRNA mutagenesis efficacies. (C–D) 1% Agarose gel and alignment of cloned PCR products showing the result of Fgf4 intron-I amplification with the oligos indicated in (A) from embryos electroporated with 25 µg Ef1α>nls::Cas9-gem::nls and 40 µg each of U6>sgRNA6 and U6>sgRNA8. (C) Alignment of cloned PCR products amplified using the primers indicated in (A), from the ~1.2 kilobase band that is similar to the control. (D) Alignment of cloned PCR products amplified using the primers indicated in (A), from the ~0.5 kilobase band that corresponds to the expected deletion between the approximate sites targeted by the two sgRNAs. The intronic element that is excluded from the sites targeted by the two sgRNAs remained intact (yellow box).

Consistent with the hypothesis that FGF-MAPK-dependent, cardiopharyngeal-specific elements act as tissue-specific enhancers, they were predominantly found in intronic or intergenic regions (48% and 37%, respectively, FDR < 0.05, one-tailed hypergeometric test). Conversely, tissue-specific peaks associated with tail muscle markers were enriched in promoters, TSS and 5’UTR (one-tailed hypergeometric test, FDR < 0.05, 57%, 23% and 15%, respectively; Figure 2—figure supplement 2C). Previously characterized enhancers for TVC-specific genes Lgr4/5/6, Rhodf, Foxf, Unc5, Rgs21, Ddr, Asb2 and Gata4/5/6 showed ATAC-seq patterns consistent with cardiopharyngeal-specific accessibility (Supplementary file 3; Beh et al., 2007; Bernadskaya et al., 2019; Christiaen et al., 2008; Woznica et al., 2012). We thus leveraged differential accessibility profiles to identify novel enhancers of cardiopharyngeal gene expression. We focused on a locus containing the conserved cardiac determinant and TVC marker, Nk4/Nkx2-5 (Wang et al., 2013), and two tail-muscle specific Myosin regulatory light chain (Mrlc) genes (Kusakabe et al., 2004; Satou et al., 2001b; Sierro et al., 2006), with associated elements showing the predicted TVC- and ATM-specific accessibility patterns, respectively (Figure 2A,D). Reporter gene expression assays showed that a DNA fragment containing differentially accessible elements located in the Nk4/Nkx2-5 intron (KhC8.2200 and .2201) was sufficient to drive GFP expression specifically in cardiopharyngeal multipotent progenitors (Figure 2E). B7.5 lineage-specific CRISPR/Cas9-induced deletions of these elements reduced or eliminated Nk4/Nkx2-5 expression specifically in TVCs, thus demonstrating its role as a bona fide cardiopharyngeal enhancer (Figure 2F; Figure 2—figure supplement 3D). Extending these analyses to other loci, including Fgf4, Fzd4, Foxg-r, Fbln, Eph1, Ncaph, Hand and Smurf1/2, we identified 8 out of 15 candidate cardiopharyngeal enhancers that drove reporter expression in the multipotent progenitors (Figure 2—figure supplement 4; Supplementary file 4), and B7.5-lineage-specific CRISPR/Cas9-mediated mutagenesis targeting differentially accessible elements reduced TVC-specific expression of the neighbouring genes Fgf4, Smurf1/2 and Fbln (Figure 2—figure supplement 5; Figure 2—figure supplement 6; Supplementary file 5). Conversely, candidate ATM-specific elements activated reporter gene expression in the tail muscles, including ATM cells, but not in the cardiopharyngeal progenitors, and were located near tail muscle markers (Figure 2—figure supplement 1B,C; Supplementary file 6). Collectively, these findings indicate that genomic elements that open specifically in multipotent progenitors act as transcriptional enhancers of cardiopharyngeal gene expression and their accessibility is controlled by FGF-MAPK induction.

A Foxf-dependent code for cardiopharyngeal accessibility

Next, we harnessed chromatin accessibility patterns predictive of cardiopharyngeal enhancer activity to identify enriched sequence motifs, and thus candidate regulators of chromatin accessibility and gene expression. For this, we performed a one-tailed hypergeometric test for enrichment of known motifs. We complemented this analysis by calculating differential accessibility of motifs using chromVAR (Schep et al., 2017), which was developed to analyze sequence motifs associated with cell-type-specific accessibility (Figure 3A). Naive Mesp+ mesoderm-specific elements, which closed between 6 and 10 hpf, were enriched in motifs for Homeodomain, T-box and Ets families of transcription factors (TF), consistent with documented roles for Lhx3/4, Tbx6 and Ets homologs in B7.5 blastomeres (Figure 3A; Davidson et al., 2006; Davidson et al., 2005; Satou et al., 2001a). Candidate tail muscle-specific elements, which opened upon FgfrDN misexpression, were similar to naive _Mesp_+ mesoderm, and enriched in motifs for the basic helix-loop-helix (bHLH) family of TFs, which includes Mesp and Mrf/MyoD, a conserved muscle-specific transcription regulator that promotes tail muscle differentiation (Christiaen et al., 2008; Meedel et al., 2007; Razy-Krajka et al., 2014; Tolkin and Christiaen, 2016). By contrast, motifs for Zinc Finger, Fox/Forkhead and nuclear receptor families of TFs were enriched among candidate cardiopharyngeal-specific elements, revealing a typical mesendodermal signature for early cardiopharyngeal progenitors (Cusanovich et al., 2018).

Figure 3. Foxf is required for cardiopharyngeal-specific chromatin accessibility.

(A) Motif accessibility between libraries from chromVAR (Schep et al., 2017). Motifs were obtained and associated with Ciona transcription factors (TFs) as described in the Materials and methods. Deviations were computed for FGF signaling-dependent regions at 10 hpf and B7.5 replicates at 6 and 10 hpf. We calculated the differential accessibility of all motifs between conditions and time points. Only the most significant motif is shown for each TF. (B) Expression of transcription factors over time compared to enrichment of corresponding TF motifs in condition specific peak sets. log2(odds ratio) values (log2(OR), see Materials and methods) are shown for motifs that are significantly enriched in a peak set (one-tailed hypergeometric test, FDR < 0.05). Only TFs expressed in the B7.5 lineage are shown. (**C**) Differential expression of Foxf target genes (DE) vs. differential chromatin accessibility (DA) in _FoxfCRISPR_. ρ is the Spearman correlation of expression and accessibility for DA regions associated with DE genes. (**D**) Association between expression of Foxf target genes and accessibility of proximal regions which were both TVC-specific and closed in _FoxfCRISPR_ as in Figure 2C. (**E**) A 3.6 kb region on chromosome 14 displaying expression profiles of RNA-seq and chromatin accessibility profiles of ATAC-seq normalized tag count. Foxf core binding site (GTAAACA) is displayed as blue line. The boxed region indicates a newly identified TVC-specific enhancer in _Hand_ locus. Red arrow indicates a TVC-specific enhancer showing closed chromatin in _FoxfCRISPR_ ATAC-seq. (**F**) Enhancer-driven in vivo reporter expression (green) of tested ATAC-seq peaks. TVCs marked with _Mesp_>H2B::mCherry (red). Numbers indicate observed/total of half-embryos scored. Experiment performed in biological replicates. Scale bar, 30 μm. Gene expression data for 6 hpf and ‘FGF-MAPK perturbation 10 hpf’ (Christiaen et al., 2008), and from 8 to 20 hpf (Razy-Krajka et al., 2014) were previously published.

Figure 3.

Figure 3—figure supplement 1. Conservation of DNA-recognition motifs of Foxf proteins.

Figure 3—figure supplement 1.

(A) Conservation of DNA-recognition motifs among Homo sapiens Hnf-3/forkhead and tunicates and amphioxus Foxf proteins. The protein sequences were retrieved from NCBI (www.ncbi.nlm.nih.gov), Ensembl (www.ensembl.org) and Aniseed (www.aniseed.cnrs.fr/aniseed) databases using Foxf of Ciona robusta as query sequence (see Supplementary file 1). Alignment of conserved amino-acid residues in Foxf orthologues was obtained using Clustal Omega (MView algorithm) (McWilliam et al., 2013). (B) Amino acid sequence alignment of DNA-recognition motif between Homo sapiens Hnf-3/forkhead (AAT74628.1) and Ciona robusta Foxf (NP_001071710.1). Asterisks indicate the structurally equivalenced residues between HFN3γ and GH5 as in Clark et al. (1993). (C) Amino acid sequence alignment in HFN3γ (107-223) and GH5 based on structural alignment and secondary structural elements indicated as rectangle for ɑ-helices and arrows for β-strand (structurally equivalenced residues are indicated in upper case) adapted from Clark et al. (1993).

Figure 3—figure supplement 2. General characterization of Foxf CRISPR.

Figure 3—figure supplement 2.

(A) ATAC-seq peaks displaying Foxf and Gata4/5/6 loci with expression profiles from RNA-seq and chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). Known TVC-specific enhancers of Foxf (Beh et al., 2007) and Gata4/5/6 (Woznica et al., 2012) are in the boxed regions. (B) Differential gene expression determined by RNA-seq analysis of FoxfCRISPR vs. ControlCRISPR at 10 hpf. Colored dots represent genes that are differentially expressed (FDR < 0.05) in _FoxfCRISPR_ vs. _ControlCRISPR_ at 10 hpf. Genes not differentially expressed are shown in gray. (**C**) GSEA showing normalized enrichment score (NES) of gene sets in peaks ranked by difference in accessibility between _FoxfCRISPR_ vs. _ControlCRISPR_ at 10 hpf. Dark gray bars indicate significant enrichment. Light gray bars are not significant (FDR < 0.05). (**D**) Change in accessibility between ATAC-seq samples. Colored dots represent peaks associated with genes that are differentially accessible (FDR < 0.05 and |log2(FC)| > 0.5). (E) Association between expression and accessibility of Foxf up-regulated genes and peaks which were both closed in FoxfCRISPR and associated with TVC-specific peaks as in Figure 2C. Gene expression data for 6 hpf and ‘FGF-MAPK perturbation 10 hpf’ (Christiaen et al., 2008), and from 8 to 20 hpf (Razy-Krajka et al., 2014) were previously published.

Figure 3—figure supplement 3. Conserved binding motifs in TVC-specific Nk4/Nkx2-5 enhancer.

Figure 3—figure supplement 3.

(A) Conservation of TVC-accessible elements of Nk4/Nk2-5 locus between Ciona robusta and Ciona savignyi. The transcript model is shown in black. Conservation score was calculated using mVISTA (http://genome.lbl.gov/vista/mvista/submit.shtml). Pink peaks indicate conserved non-coding sequences (> 65% identity per 80 bp). (B) Alignment of DNA region corresponding to ‘_KhC8.2201_’ peak between Ciona savignyi and Ciona robusta (bottom sequence). Only the core nucleotides are shown for putative Forkhead and Gata binding sites. (C) Match scores for TF binding motifs present in the boxed region (~1.2 Kb) including the I-intronic element of both Ciona robusta and Ciona savignyi Nk4/Nk2-5 transcripts. Only the highest scoring match for each TF is shown.

Figure 3—figure supplement 4. Foxf loss-of-function (FoxfCRISPR) caused closing of TVC-specific enhancers.

Figure 3—figure supplement 4.

(A–E) Foxf target gene loci chromatin accessibility profiles showing ATAC-seq (normalized by total sequencing depth) in the indicated conditions. Foxf core binding sites (GTAAACA) are displayed as blue lines. The transcript model is shown in black.

Figure 3—figure supplement 5. Conserved binding motifs in TVC-specific Hand enhancer.

Figure 3—figure supplement 5.

(A) Conservation of TVC-accessible elements of Hand locus between Ciona robusta and Ciona savignyi. The transcript model is shown in black. Conservation score was calculated using mVISTA (http://genome.lbl.gov/vista/mvista/submit.shtml). Violet peaks indicate conserved sequences (> 65% identity per 80 bp). (B) Alignment of DNA region corresponding to ‘_KhC14.806_’ peak between Ciona robusta and Ciona savignyi (bottom sequence). Putative Forkhead binding sites are evidenced in green box (only the core nucleotides are shown). (C) Match scores for TF binding motifs present in the boxed region (~0.8 Kb) upstream of the coding ATG in both Ciona robusta and Ciona savignyi Hand transcripts. Only the highest scoring match is shown for each TF.

Figure 3—figure supplement 6. Conserved binding motifs in Smurf enhancer.

Figure 3—figure supplement 6.

(A) Conservation of TVC-accessible elements of Smurf locus between Ciona robusta and Ciona savignyi. The transcript model is shown in black. Conservation score was calculated using mVISTA (http://genome.lbl.gov/vista/mvista/submit.shtml). Pink peaks indicate conserved non-coding sequences (> 65% identity per 80 bp). (B) Alignment of DNA region corresponding to ‘_KhC4.606_’ peak between Ciona robusta and Ciona savignyi (bottom sequence). Only the core nucleotides are shown for putative Forkhead and GATA binding sites. (C) Match scores for TF binding motifs present in the boxed region including the II-intronic element and exon III-IV of both Ciona robusta and Ciona savignyi Smurf transcripts. Only the highest scoring match is shown for each TF.

Figure 3—figure supplement 7. Accessible elements annotated to Foxf target genes.

Figure 3—figure supplement 7.

(A) Intersection of elements associated with Foxf target genes at 10 hpf (see Material and methods) and differentially accessible elements. Peaks closed in FoxfCRISPR are peaks more accessible in ControlCRISPR vs. FoxfCRISPR at 10 hpf. A two-tailed binomial test was performed on each intersection against the null hypothesis that the intersection’s constituent sets are independent (see Statistics under Materials and methods). (B) Subset of (A) showing only elements associated with Foxf target genes. Peaks closed at 6 hpf are peaks more accessible in Mesp> LacZ 10 hpf vs. Mesp > LacZ 6 hpf. (C) ATAC-seq peaks (113 regions) associated with Foxf target genes showing differential accessibility (FDR < 0.05, see Material and methods) over time. The log2(FC)for each time point is versus the average (avg) of all control samples.

Combined with motif enrichment analyses, temporal gene expression profiles (Razy-Krajka et al., 2014) identified candidate _trans_-acting regulators of cardiopharyngeal-specific accessibility and/or activity (Figure 3B). For example, regions specifically accessible in tail muscle- or naive Mesp+ mesoderm were enriched in homeobox, Ets and T-box motifs, consistent with early expression of Mrf/MyoD, Lhx3/4, Ets1/2, and Tbx6, respectively. Similarly, the increased accessibility of motifs for K50 Paired homeodomain proteins in naive Mesp+ mesoderm indicated a possible role for Otx, which is expressed early in B7.5 blastomeres (Figure 3A,B; Hudson et al., 2003). Cardiopharyngeal-specific enrichment for Fox/Forkhead and Zinc Finger motifs pointed to several known factors, including Foxf, one of the first genes activated in multipotent progenitors upon induction by FGF-MAPK, prior to Gata4/5/6 (Beh et al., 2007; Christiaen et al., 2008; Ragkousi et al., 2011). Moreover, GATA and Forkhead proteins are founding members of a group of TFs known as pioneers, which can bind their target sites in closed chromatin and promote accessibility (Cirillo et al., 2002; Zaret and Carroll, 2011). Protein sequence alignments indicated the presence of key residues, conserved between the DNA binding domains of Foxf and the classic pioneer FOXA, which mimic linker histone H1 in its ability to displace DNA-bound nucleosomes (Figure 3—figure supplement 1B,C; Clark et al., 1993). Finally, the Foxf enhancer was accessible in naive _Mesp_+ founder cells, suggesting that it is poised for activation, unlike the intronic Gata4/5/6 enhancer (Figure 3—figure supplement 2A). Consistent with a role for Fox and GATA proteins in opening and activating the Nk4/Nkx2-5 enhancer, we found putative cognate binding sites in the newly identified element to be conserved with the closely related species, C. savignyi (Figure 3—figure supplement 3B,C). Taken together, these analyses identified a putative code for cardiopharyngeal-specific accessibility and enhancer activity, which comprise motifs for candidate DNA binding factors of the Forkhead and GATA and identified Foxf as candidate determinant of cardiopharyngeal accessibility.

To test if Foxf contributes to establishing cardiopharyngeal accessibility and gene expression profiles, we used reagents for B7.5 lineage-specific loss-of-function by CRISPR/Cas9-mediated mutagenesis (FoxfCRISPR; Gandhi et al., 2017), and performed ATAC- and RNA-seq on FACS-purified cells isolated from tailbud embryos at 10 hpf. RNA-seq confirmed that CRISPR/Cas9 mutagenesis inhibited Foxf itself and other TVC-expressed genes, including effectors of collective cell migration such as Ddr, consistent with previous microarray data (Bernadskaya et al., 2019; Christiaen et al., 2008) (Figure 3D; Figure 3—figure supplement 2B). Out of 52 differentially expressed genes (Figure 3—figure supplement 2B; Supplementary file 7), seven down-regulated genes were previously annotated as primed pan-cardiac markers, including Hand, Gata4/5/6 and Fzd4 (Wang et al., 2019). Down-regulated genes also included primed pharyngeal muscle markers, such as Rhod/f (Figure 3D; Figure 3—figure supplement 2B; Christiaen et al., 2008; Razy-Krajka et al., 2014), suggesting that Foxf promotes the onset of both the cardiac and pharyngeal muscle programs in multipotent progenitors, a feature known as multilineage transcriptional priming (Razy-Krajka et al., 2014; Wang et al., 2019).

Consistent with the effects of Foxf mutagenesis on gene expression, regions closed in FoxfCRISPR samples included known cardiopharyngeal enhancers for Gata4/5/6 and Ddr (Figure 3D; Figure 3—figure supplement 4; Bernadskaya et al., 2019; Christiaen et al., 2008; Woznica et al., 2012), newly identified enhancers for Eph1, Smurf1/2 and Fzd4, and a novel enhancer of Hand expression (Hand_ KhC14.805 -. 807; Figure 3C–F, Supplementary file 8). These differentially accessible elements contain several, evolutionary conserved, putative Fox binding sites (Figure 3—figure supplement 5; Figure 3—figure supplement 6). We identified two conserved putative Forkhead binding sites in the minimal STVC-specific enhancer from the Tbx1/10 locus (termed T12; Razy-Krajka et al., 2018), which were necessary for reporter gene expression (Figure 5—figure supplement 1D,E). Moreover, loss-of-function of Foxf (FoxfCRISPR) drastically reduced T12 enhancer activity (Figure 5—figure supplement 1). These results are consistent with the hypothesis that Foxf acts directly on the minimal Tbx1/10 enhancer to promote its activity in the second multipotent cardiopharyngeal progenitors.

Notably, 98% (40/41) of the regions with diminished accessibility following Foxf inhibition, and located near a candidate Foxf target gene, were also accessible in multipotent progenitor cells (Figure 3D; Figure 3—figure supplement 7A). Moreover, 22% (600/2,728) of the predicted multipotent progenitor-specific elements were closed upon Foxf inhibition, and gene set enrichment analysis indicated that Foxf loss-of-function generally decreased the accessibility of cardiopharyngeal-specific elements (Figure 3—figure supplement 2C; Figure 3—figure supplement 7B). Finally, 18 of 41 (44%) of the _Foxf_-dependent elements associated with candidate Foxf targets were closed in 6 hpf founder cells and appear to open specifically in the cardiopharyngeal progenitors by 10 hpf (Figure 3—figure supplement 2D; Figure 3—figure supplement 7B,C; Supplementary file 9). This dynamic is consistent with a requirement for Foxf activity following its activation in the TVCs, immediately after division of the naive _Mesp_+ progenitors. Taken together, these results indicate that, in newborn multipotent progenitors, FGF-MAPK signaling upregulates Foxf (Beh et al., 2007; Christiaen et al., 2008), which is in turn required to open a substantial fraction of cardiopharyngeal-specific elements for gene expression in multipotent progenitors, including for such essential determinants as Gata4/5/6 and Hand (Figure 4H).

Figure 4. Cardiopharyngeal lineage-specific accessibility profiles and decoupling between enhancer accessibility and activity for de novo expressed genes.

(A) Differentially expressed (DE) genes vs. differentially accessible (DA) peaks in response to FGF-MAPK perturbation in the fate-restricted cells. ρ is the Spearman correlation of expression and accessibility for DA peaks associated with DE genes. (B) Relationship between accessibility and expression of de novo pan-cardiac genes as in Figure 2C. DE genes in either condition are shown on the left. (C) Time-dependent ATAC-seq peaks associated with de novo expressed pan-cardiac genes. The accessibility of these peaks is shown for 6, 10 and 18 hpf vs. the average accessibility in the controls (LacZ) and upon FGF-MAPK perturbations at either 10 or 18 hpf. Peaks were classified as ‘Open in ASM’ (less accessible in FgfrDN vs. M-RasCA or LacZ at 18 hpf), ‘Open in Heart’ (less accessible in _M-RasCA_vs. _FgfrDN_or LacZ at 18 hpf), ‘Closed in _FoxfCRISPR_’ (less accessible in FoxfCRISPR vs. ControlCRISPR), or ‘Open in TVC’ (less accessible in FgfrDN vs. MekS216D,S220E or LacZ at 10 hpf). Only regions changing accessibility between 6 and 10 hpf, or 10 and 18 hpf are shown. (D) A 6 kb region on chromosome four displaying expression profiles of RNA-seq and chromatin accessibility profiles of ATAC-seq normalized tag count. Peak ID refers to elements tested for reporter assay in vivo. The newly identified enhancer in Lrp4/8 locus is in the boxed region. (E) Enhancer-driven in vivo reporter expression (green) of tested ‘_KhC4.137_’ peak. TVCs marked with Mesp>H2B::mCherry (red). Numbers indicate observed/total of half-embryo scored. Zoom on cardiopharyngeal cell lineage (panel on the right). (F) Endogenous expression of Lrp4/8 visualized by in situ (green) in TyrosinaseCRISPR and upon CRISPR/Cas9-induced deletion of ATAC-seq peaks. Nuclei of B7.5 lineage cells are labelled by Mesp>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). _Mesp_-driven hCD4::mCherry accumulates at the cell membrane as revealed by anti mCherry antibody (Blue). Experiment performed in biological replicates. Scale bar = 10 μm. (G) Fisher exact test; n is the total number of individual embryo halves scored per condition. (H) Summary model: patterns of chromatin accessibility dynamics and gene expression during early cardiopharyngeal fate specification.

Figure 4.

Figure 4—figure supplement 1. General characterization of FGF-MAPK perturbation.

Figure 4—figure supplement 1.

(A–D) Intersections of differentially expressed genes from bulk RNA-seq and scRNA-seq (Wang et al., 2019) comparisons in cardiopharyngeal restricted cells using Venn Diagram (A–B) and UpSet plots (C–D). FGF-MAPK activated genes at 18 hpf are defined as the intersection of genes downregulated in FgfrDN vs. LacZ at 18 hpf and genes downregulated in _FgfrDN_vs. M-RasCA at 18 hpf (A). FGF-MAPK inhibited genes at 18 hpf are defined as the intersection of genes downregulated in M-RasCA vs. LacZ at 18 hpf and genes downregulated in M-RasCA vs. FgfrDN at 18 hpf (FDR < 0.05 and |log2(FC)| > 1). (E–F) GSEA of significantly enriched gene sets shows opposite trends in FgfrDN vs. control and _M-RasCA_vs. control (LacZ). Dark gray bars indicate significant enrichment. Light gray bars are not significant (FDR < 0.05).

Figure 4—figure supplement 2. Differential accessibility in response to FGF/MAPK perturbation at 18 hpf.

Figure 4—figure supplement 2.

(A) Differentially expressed (DE) genes vs. differentially accessible (DA) peaks using bulk RNA-seq and ATAC-seq in Foxf>M RasCA vs. LacZ at 18 hpf. (B) Relationship between expression and accessibility of DE genes associated with DA peaks. The top axis shows the scale for log2(FC) of bulk RNA-seq (black dots). The bottom axis shows the log2(FC) scale for ATAC-seq (colored diamonds). Top 50 genes inhibited and activated by FGF-MAPK (based on log2(FC)) along with STVC genes. Heatmaps show gene expression over time. (C) A 15 kb region on chromosome five displaying expression profiles from RNA-seq and chromatin accessibility profiles from ATAC-seq (normalized by total sequencing depth). Transcript model is indicated as black bars, percentage of conservation score between C. savignyi and C. robusta obtained obtained from WASHU browser (Brozovic et al., 2018) (yellow peaks), accessome (light grey bars), TVC- (green bar) and ATM-specific peaks (dark gray bars). Peak accessible specifically in the heart (KhC5.1641) of Mmp21 locus is in the boxed region. (D–E) Enhancer-driven in vivo reporter expression (green) of a ~ 3 kb region upstream the coding ATG and including the DA peak (KhC5.1641). (D) dorsal view, (E) lateral view; GFP is detected in the mesenchyme (asterisk) surrounding the cardiopharyngeal progenitors (D) and in the epidermis (orange arrow) (E). B7.5 cells are marked with _Mesp_>H2B::mCherry (red). Dotted line: ventral midline. Numbers indicate observed/total of half-embryos scored. Scale bar = 25 µm.

Figure 4—figure supplement 3. Accessibility of elements annotated to de novo ASM genes.

Figure 4—figure supplement 3.

(A) Relationship between accessibility and expression of de novo ASM genes in either condition indicated. (B) Time-dependent ATAC-seq peaks associated with de novo expressed ASM genes. The accessibility of these peaks is shown for 6, 10 and 18 hpf vs. the average (avg.) accessibility in the controls and upon FGF-MAPK perturbation either at 10 and 18 hpf. Peaks are clustered based on their accessibility as ‘Open in ASM’ less accessible in Hand-r>FgfrDN vs. Foxf>M RasCA or Hand-r>LacZ at 18 hpf, ‘Open in Heart’ less accessible in Foxf>M RasCA vs. Hand-r>FgfrDN or Hand-r>LacZ at 18 hpf, ‘Closed in _FofxCRISPR_’ less accessible in FoxfCRISPR vs. ControlCRISPR, ‘Open in TVC’ less accessible in Mesp>FgfrDN vs. Mesp>MekS216D,S220E or Mesp>LacZ at 10 hpf. Only peaks changing accessibility between 6 hpf and 10 hpf or 10 hpf and 18 hpf are shown.

Figure 4—figure supplement 4. Accessibility of binding motifs over time for elements annotated to de novo-expressed genes.

Figure 4—figure supplement 4.

(A) TF binding motifs enriched in peaks associated with de novo cardiac and pharyngeal expressed genes parsed based on primed and de novo accessibility (see Material and methods). log2(OR) (see Materials and methods) are shown for motifs significantly enriched (one-tailed hypergeometric test, p < 0.05) in the indicated peak classes. Motif accessibility from chromVAR (Schep et al., 2017) is shown for all peaks associated with a de novo-expressed cardiac or ASM gene. Only the motif with the highest log2(OR) for each TF is shown. (B) Conservation of TVC-accessible peaks closed in FoxfCRISPR in Lrp4/8 locus between C. savignyi and C. robusta. The gene body is shown in black. Conservation score was obtained from WASHU browser (Brozovic et al., 2018). (C) Alignment of conserved region of ‘_KhC8.137_’ peak between C. robusta and C. savignyi. (D) C. robusta Hand locus, exon IV sequence highlighted in bold and black. Only the core nucleotides are shown for candidate Pitx binding sites.

Chromatin accessibility in late heart vs. pharyngeal muscle precursors

Besides controlling coherent chromatin opening, enhancer activity and gene expression in multipotent cardiopharyngeal progenitors, FGF-Foxf inputs also appeared to open regions associated with later de novo-expressed heart and pharyngeal muscle markers (Figure 2A,B; Supplementary file 10; Razy-Krajka et al., 2018; Razy-Krajka et al., 2014; Wang et al., 2019; Wang et al., 2013). Accessibility patterns were also better correlated between 10 and 18 hpf (Figure 1B), suggesting a decoupling between early accessibility and late heart- vs. pharyngeal muscle-specific expression in late fate-restricted precursors. To identify accessibility patterns underlying the heart vs. pharyngeal muscle fate choices, we compared bulk RNA-seq (Wang et al., 2019) and ATAC-seq datasets obtained from cardiopharyngeal lineage cells isolated from 18 hpf larvae, following the same defined perturbations of FGF-MAPK signaling (Figure 1A; Figure 4A; Figure 4—figure supplement 1A–D; Figure 4—figure supplement 2A,B; Davidson et al., 2006; Razy-Krajka et al., 2018). Among cardiac and pharyngeal muscle markers, we identified 35 FGF-MAPK-regulated genes associated with one or more elements showing consistent differential accessibility (Figure 4A,B; Figure 4—figure supplement 2A,B; Supplementary file 11). This indicated that, at least for a subset of cardiopharyngeal marker genes, FGF-MAPK-dependent changes in gene expression follow corresponding changes in chromatin accessibility in early heart and pharyngeal muscle precursors.

Gene-level inspection of differential accessibility associated with either inhibition or activation of gene expression revealed that only a fraction of associated elements was either closing or opening upon perturbation of FGF-MAPK signaling (Figure 4B; Figure 4—figure supplement 2B). For example, the first heart lineage marker Matrix metalloproteinase 21/Mmp21 (Wang et al., 2019) was associated with multiple upstream and intronic elements, but only some of these elements were differentially accessible following either gain or loss of FGF-MAPK function (Figure 4—figure supplement 2A–C), and a 3 kb fragment containing the upstream differentially accessible element sufficed to drive reporter gene expression throughout the cardiopharyngeal lineage, but not specifically in the first heart precursors (Figure 4—figure supplement 2C–E). Similarly, reporter gene expression assays showed that DNA fragments containing differentially accessible elements located ~0.5 kb upstream of the coding region of the de novo-expressed gene Tmtc2 (KhC2.3468), and upstream of KH.C1.1093_ZAN (KhC3.47, KhC3.46), were sufficient to drive GFP expression in both cardiac and pharyngeal muscle progenitors, consistent with the notion that electroporated plasmids are not ‘chromatinized’ and thus constitutively accessible (Figure 4—figure supplement 2F–G; Figure 4—figure supplement 3C–D). This suggested that, for genes like Mmp21, Tmtc2 and Zan, cell-type-specific accessibility determines cardiac vs. pharyngeal muscle-specific gene expression.

Remarkably, the vast majority (91%, 356 genes out of 391, Supplementary file 12) of differentially expressed genes were not associated with differentially accessible elements (Figure 4A,B; Figure 4—figure supplement 2A,B). Specifically, out of 30 de novo-expressed pan-cardiac genes that were also differentially expressed upon FGF-MAPK perturbation at 18 hpf, only 8 (27% ± 8%, SE) were associated with one differentially accessible element following perturbation of FGF-MAPK signaling (Figure 4B). Similarly, out of 23 de novo-expressed pharyngeal muscle genes, which were also differentially expressed upon FGF-MAPK perturbation at 18 hpf, 11 were associated with one differentially accessible element following perturbation of FGF-MAPK signaling (Figure 4—figure supplement 3A). This suggested that most differential gene expression in early heart and pharyngeal muscle precursors arise from differential _cis_-regulatory activity of elements that are otherwise accessible throughout the cardiopharyngeal mesoderm. In keeping with this hypothesis, accessible regions associated with de novo-expressed pan-cardiac and pharyngeal muscle markers tended to open between 6 and 10 hpf, in a pattern consistent with FGF- and Foxf-dependent cardiopharyngeal-specific accessibility (Figure 4C; Figure 4—figure supplement 3B). These observations suggest that _cis_-regulatory elements controlling cell-type-specific de novo gene expression open in multipotent progenitors, prior to becoming active in fate-restricted precursors. Such decoupling between enhancer accessibility and activity has been observed in other developmental contexts, including early cardiogenesis in mammals (Paige et al., 2012; Wamstad et al., 2012).

As a proof of principle, we analyzed the Lrp4/8 locus, which harbors two intronic elements (KhC4.137 and KhC4.144) that opened upon TVC induction in an FGF- and Foxf-dependent manner, prior to Lrp4/8 upregulation in cardiac progenitors (Wang et al., 2019), and were not differentially accessible at 18 hpf (Figure 4C,D). Of the two regions, only KhC4.137 was sufficient to drive GFP expression in heart precursors indicating enhancer activity, and illustrating the decoupling between early and broad accessibility and late, cell-type-specific, activity (Figure 4E). Reporter gene expression and CRISPR/Cas9-mediated mutagenesis assays followed by FISH indicated that KhC4.137 is both necessary and sufficient to activate gene expression in heart precursors (Figure 4E,F), showing that it acts as a bona fide enhancer, and demonstrating a specific case of decoupling between early and broad accessibility and late, cell-type-specific, activity.

To identify candidate regulators of late accessibility and/or activity, we parsed accessible elements associated with de novo-expressed heart and pharyngeal muscle markers into pre-accessible/primed or de novo-accessible elements and discovered sequence motifs enriched in each category (Figure 4—figure supplement 4A; Supplementary file 13). Putative binding sites for SMAD and homeodomain proteins such as Smad4 and Pitx respectively were enriched among pre-accessible elements associated with cardiac markers, and found in the primed elements regulating Lrp4/8 upregulation (Figure 4—figure supplement 4A,B), suggesting a specific role in transcriptional activation, consistent with conserved roles for Pitx2 and BMP-SMAD signaling during heart development (Figure 4—figure supplement 4C,D; Nowotschin et al., 2006; Schultheiss et al., 1997). Motifs for known regulators of cardiac development, including Meis (Desjardins and Naya, 2016; Paige et al., 2012), were over-represented among de novo-accessible elements associated with cardiac markers, suggesting roles in establishing accessibility and/or regulating enhancer activity (Figure 4—figure supplement 4A). Notably, GATA motifs were enriched in primed accessible elements associated with cardiac markers, consistent with conserved roles for GATA factors as pioneer factors, and during cardiac development (Pikkarainen et al., 2004). Among motifs enriched in accessible elements associated with de novo-expressed pharyngeal muscle markers, the presence of ETS-, bHLH, and EBF-family motifs is consistent with established roles for FGF-MAPK, Hand-r, Mrf and Ebf in pharyngeal muscle specification (Razy-Krajka et al., 2018; Razy-Krajka et al., 2014; Stolfi et al., 2014; Stolfi et al., 2010; Wang et al., 2013). Notably, the enrichment of Ebf motifs among de novo-accessible elements associated with de novo-expressed pharyngeal muscle markers is reminiscent of the ability of EBF-family factors to interact with nucleosome-bound cognate sites, suggestive of a pioneering activity in committed pharyngeal muscle precursors (Boller et al., 2016; Buenrostro et al., 2013). In summary, this analysis identified distinct combinations of established and putative _trans_-acting factors differentially controlling the accessibility and/or activity of _cis_-regulatory elements that govern heart- vs. pharyngeal-muscle-specific gene expression (Figure 4H).

Combinatorial _cis_-regulatory control of cardiopharyngeal determinants

The above analyses focused on one-to-one associations between accessible elements and neighboring genes to uncover candidate _trans_-acting inputs controlling gene expression through defined elements. However, most genes are associated with multiple accessible regions, especially developmental regulators (Figure 1—figure supplement 3D–F), presumably exposing diverse motifs for transcription factor binding. Moreover, distinct elements associated with the same neighboring gene often exhibited different accessibility dynamics. For instance, the loci of several de novo-expressed heart and pharyngeal muscle markers contained both primed-accessible and de novo-accessible elements (Figure 4C; Figure 4—figure supplement 3B; Figure 4—figure supplement 4A,B). This suggested that individual genes respond to a variety of regulatory inputs mediated through separate _cis_-regulatory elements.

To explore this possibility, we focused on Tbx1/10 and Ebf, two established determinants of cardiopharyngeal fates (Razy-Krajka et al., 2018; Razy-Krajka et al., 2014; Stolfi et al., 2014; Stolfi et al., 2010; Tolkin and Christiaen, 2016; Wang et al., 2013). Both loci contained multiple accessible regions, including elements open already in the naive _Mesp_+ mesoderm (e.g. Ebf_ KhL24.35/36), cardiopharyngeal-lineage-specific elements that open prior to gene activation but after induction of multipotent progenitors (e.g. Ebf_ KhL24.34), and elements that open de novo in fate-restricted pharyngeal muscle precursors, where the gene is activated (e.g. Ebf_ KhL24.37) (Figure 5A,B). Previous reporter gene expression assays identified the latter element, Ebf_ KhL24.37, as a weak minimal enhancer with pharyngeal muscle-specific activity (Wang et al., 2013). CRISPR/Cas9-mediated mutagenesis assays followed by FISH indicated that each one of these elements is necessary for proper activation of Ebf in pharyngeal muscle progenitors (Figure 5C–E; Figure 5—figure supplement 2B,C). Consistently, we found that targeted deletions of individual accessible elements upstream of Ebf induced pharyngeal muscle precursor migration defects. When targeting Ebf_ KhL24.37, 37 ± 5% (SE) of 28 hpf larvae showed such defects (n = 101, Figure 5—figure supplement 3). These observations indicate that each accessible _cis_-regulatory element upstream of Ebf is necessary for proper expression and subsequent pharyngeal muscle morphogenesis.

Figure 5. Combinations of _cis_-regulatory elements with distinct chromatin accessibility profiles are required for Ebf transcription in pharyngeal-muscle precursors.

(A) A 12 kb region of the scaffold L24 displaying expression profiles of RNA-seq and chromatin accessibility profiles of ATAC-seq (normalized tag count) in the Ebf locus. sgRNAs used to target ATAC-seq peaks are shown in red; intronic antisense riboprobes are shown in orange (B) Schematic representation showing sequential opening of _cis_-regulatory elements required for Ebf activation in pharyngeal muscle founder cells, and maintenance by auto-regulation in committed precursor. (C) Schematic representation of Ebf cis_-regulatory elements targeted for CRISPR/Cas9-mediated deletions. Shapes represent binding sites located in the regulatory elements and differentially accessible over time. (D) Proportions of larva halves showing the indicated Ebf transcription patterns, in indicated experimental conditions; all the treatments were significant versus Tyrosinase (Fisher exact test, p < 0.001). (**E**) Endogenous expression of _Ebf_ visualized by in situ (green) in _TyrosinaseCRISPR_ and upon CRISPR/Cas9-induced deletion of ATAC-seq peaks as indicated, at stage 25 (**E**) and 27 (**F**) based on Hotta et al. (2007). For stage 25, an anti-sense riboprobe for the full length cDNA was used, whereas for stage 27 an intronic anti-sense riboprobe targeting the first three introns of _Ebf_ transcript (orange lines) as previously used in Wang et al. (2013). Nuclei of B7.5 lineage cells are labelled by _Mesp>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). _Mesp_-driven hCD4::mCherry accumulates at the cell membrane as revealed by anti mCherry antibody (Blue). Scale bar = 10 µm.

Figure 5.

Figure 5—figure supplement 1. Combinations of _cis_-regulatory elements with distinct chromatin accessibility profiles are required for Tbx1/10 transcription in pharyngeal-muscle precursors.

Figure 5—figure supplement 1.

(A) An 11 kb region on chromosome seven displaying expression profiles of RNA-seq and chromatin accessibility profiles ATAC-seq normalized tag count in Tbx1/10 locus. (B) STVC-specific enhancer (T12) (Razy-Krajka et al., 2018) driven in vivo reporter expression (green) in ASMFs and SHPs at stage 26 (Hotta et al., 2007). Nuclei of B7.5 lineage cells are labelled by _Mesp_>H2B::mCherry (red). Scale bar, 20 μm. T12 enhancer tested alone or fused to the intronic element (T12+KhC7.914). Statistical significance of the difference in reporter expression was tested using a Fisher exact test (p < 0.001) (**C**) Motif scores in each experimentally validated peak in the _Tbx1/10_ locus. Only the highest match score is shown for each motif. (**D**) Sequence alignment of _Tbx1/10_ enhancer (_T12_) between _Ciona robusta_/_Ciona savignyi_. Conserved blocks in the orange boxes with putative Forkhead binding sites. In blue is highlighted the single guide RNA (sgRNA#2) used to target CRISPR/Cas9 system, with the PAM domain in red; the point mutations induced in two conserved putative Forkhead binding sites (Fox1 and Fox2) are in bold and red after the asterisks. (**E**) Proportion of larvae expressing both GFP and mCherry in the STVC progeny when co-electroporated wild-type and mutant _Tbx1/10_ reporters lacking the indicated putative Forkhead binding sites and _Mesp_>H2B::mCherry in comparison to the control. (F) Proportions of larvae halves showing GFP expressed in the ASMFs and SHPs in embryos electroporated with _Mesp_>Cas9 along with single guide RNAs targeting Tyrosinase (controlCRISPR) as well as Foxf (FoxfCRISPR). (G) Schematic representation of regulatory elements in Tbx1/10 locus as displayed in ATAC-seq profiles targeted for CRISPR/Cas9-mediated deletions. (H) Proportions of larvae halves showing the indicated Tbx1/10 transcription patterns, in indicated experimental conditions (Fisher exact test, p < 0.001). (**I**) Endogenous expression of _Tbx1/10_ visualized by in situ (green) in _TyrosinaseCRISPR_ (left panel) and upon CRISPR/Cas9-induced deletion of TVC-specific peaks (right panel) at stage 25 according to Hotta et al. (2007). Nuclei of B7.5 lineage cells are labelled by _Mesp_>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). Nuclei of B7.5 lineage cells are labelled by _Mesp_>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). Scale bar = 10 µm. Experiment performed in biological replicates. Total numbers of individual halves scored per condition are shown in 'n='. All the treatments were significant versus Control (Tyrosinase) (Fisher exact test, p < 0.001).

Figure 5—figure supplement 2. Ebf regulatory regions showing differentially accessibility over time contain distinct binding motifs.

Figure 5—figure supplement 2.

(A) Motif scores in each experimentally validated peak in the Ebf locus. Only the highest match score is shown for each motif. (B) Endogenous expression of Ebf visualized by in situ (green) in TyrosinaseCRISPR and upon CRISPR/Cas9-induced deletions of TVC-specific peaks at stage 24 according to Hotta et al. (2007). Nuclei of B7.5 lineage cells are labelled by _Mesp_>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). Nuclei of B7.5 lineage cells are labelled by _Mesp_>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). Scale bar = 10 µm. (C) Proportions of larvae halves showing the indicated Ebf transcription patterns, in indicated experimental conditions. Experiment performed in biological replicates. All the treatments were significant versus Control (Tyrosinase) (Fisher exact test, p < 0.001).

Figure 5—figure supplement 3. CRISPR/Cas9-mediated deletions on individual accessible elements upstream of Ebf caused phenotypic impact on pharyngeal muscle precursors morphogenesis.

Figure 5—figure supplement 3.

(A) Schematic representation of _Ebf cis_-regulatory elements targeted for CRISPR/Cas9-mediated deletions; the shapes of the distinct _cis_-regulatory elements are as in Figure 5C. (B) Proportions of larva halves showing GFP-driven STVC-specific enhancer of Tbx1/10 in indicated experimental conditions; all the treatments were significant versus Tyrosinase (Fisher exact test, p < 0.001; ‘n’ is the total number of individual halves scored per condition.). (**C**) Example of an embryo at 28 hpf showing GFP expression only in the ASM (solid arrowhead) and SHP (arrow) but not in the FHP (open arrowheads), where _Tbx1/10_ enhancer is not active (_TyrosinaseCRISPR_, first panel on the left). Targeted deletions in _KhL24.37_ peak induced ASMP cell migration defects. B7.5 lineage cells are labelled by _Mesp_>LacZ. Scale bar = 10 µm.

Figure 5—figure supplement 4. Intronic and distal enhancer accessibility in the Tbx1/10 locus tested by dCas9-KRAB.

Figure 5—figure supplement 4.

(A) ASM-specific enhancer of Ebf (Ebf-3348/–178) (Wang et al., 2013) driven in vivo reporter expression (green) in embryos at stage 27 (Hotta et al., 2007) electroporated with Mesp>dCas9 KRAB along with single guide RNAs targeting Tyrosinase (control, left panel) as well as intronic (KhC7.914) and distal (KhC7.909) of Tbx1/10 locus (right panel) as in Figure 5—figure supplement 1G. Nuclei of B7.5 lineage cells are labelled by Mesp_>H2B::mCherry (red). White asterisk indicated central nervous system (CNS). Scale bar, 50 μm. (B) Proportions of larvae halves showing GFP expressed in the ASMFs in indicated experimental conditions (Fisher exact test, p < 0.001). (**C**) Double in situ hybridization of _Ebf_ (green) and _Tbx1/10_ (red) on embryos at stage 27 electroporated with _Mesp>dCas9 KRAB along with single guide RNAs targeting Tyrosinase (control, left panel) as well as intronic (KhC7.914) and distal (KhC7.909) of Tbx1/10 locus. White asterisk indicated central nervous system (CNS). Scale bar, 20 μm. (D) Proportions of larvae halves expressing both Tbx1/10 and Ebf in the ASMFs in indicated experimental conditions. Experiment performed in biological replicates. Total numbers of individual halves scored per condition are shown in 'n='. All the treatments were significant versus Control (Tyrosinase) (Fisher exact test, p < 0.001).

Consistent with the established roles of Hand-r, Tbx1/10 and Ets-mediated FGF-MAPK signaling in activating Ebf, the primed cardiopharyngeal-specific element (KhL24.34) contained Fox and bHLH motifs, and the more distal de novo-accessible minimal enhancer (KhL24.37) also contained putative Ets and RORγ binding sites, whereas the constitutively accessible elements (KhC24.35 and .36) contained primarily CREB and T-box binding sites (Figure 5C; Figure 5—figure supplement 2A). Tbx1/10 showed a similar logic, whereby a constitutively accessible upstream element (KhC7.909) acts as an enhancer of cardiopharyngeal expression (Razy-Krajka et al., 2018), and whose activity also requires a primed cardiopharyngeal-specific intronic element (KhC7.914) (Figure 5—figure supplement 1A–F). As a complement, and to more directly test the importance of enhancer accessibility, we targeted the intronic and distal elements in the Tbx1/10 locus using dCas9::KRAB (Klann et al., 2017), which recruits deacetylases and presumably closes chromatin (Sripathy et al., 2006; Groner et al., 2010; Schultz et al., 2002; Reynolds et al., 2012; Thakore et al., 2015). We oserved loss of Tbx1/10 expression and function, as evaluated by expression of its target, Ebf (Wang et al., 2013) (Figure 5—figure supplement 4).

Of note, Ebf expression is maintained by auto-regulation (Razy-Krajka et al., 2018), which requires separate intronic elements that harbor putative Ebf binding sites and open later (Figure 5B, Figure 5—figure supplement 2A). Together with Ebf's potent myogenic and anti-cardiogenic effects (Razy-Krajka et al., 2014; Stolfi et al., 2014; Stolfi et al., 2010; Tolkin and Christiaen, 2016), this auto-regulatory logic catalyzes the pharyngeal muscle fate, stressing the importance of spatially and temporally accurate onset of expression to avoid ectopic ASM specification at the expense of cardiac identities, especially in the second heart lineage. These observations suggest that pharyngeal muscle fate specification relies on ‘combined enhancers’, characterized by a combination of _trans_-acting inputs mediated by distinct elements with variable dynamics of accessibility, to control the onset of Ebf expression in the cardiopharyngeal mesoderm (Figure 5F).

To test whether ‘combined enhancers’ drive spatially and temporally accurate expression in pharyngeal muscle progenitors, we built a reporter containing multiple copies of the minimal, but weak, Ebf enhancer (KhL24.37) (Wang et al., 2013). Two and three copies of the KhL24.37 element (2x and 3x KhL24.37) significantly increased reporter gene expression in pharyngeal muscle precursors (43 ± 3% SE for 2x KhL24.37; 58 ± 2% SE for 3x KhL24.37), compared to a single copy construct) (14 ± 3% SE for 1x KhL24.37) (Figure 6A–C), restoring reporter gene expression to levels similar to ‘full length’ upstream element encompassing all combined enhancers, with endogenous genomic spacing (Ebf-full length −3348 /- 178) (Wang et al., 2013) (80 ± 1%, SE). Remarkably, unlike the ‘full length’ combined enhancers, the 3x KhL24.37 construct induced precocious reporter gene expression in the STVCs (89 ± 3%, SE) (n = 95, Figure 6D–F) causing an ectopic GFP expression in the second heart lineage (13 ± 2%, SE) (n = 218, Figure 6B,C). To test whether spacing between accessible elements could affect transcriptional output, we built a concatemer of KhL24.37, .36, .35, and .34 elements without endogenous spacer sequences. This construct increased the proportion of embryos with ASM cells expressing the reporter to 92 ± 2% (SE, n = 130; Figure 6C), but it did not induce ectopic expression in the second heart lineage, supporting the notion that combined enhancers drive high but spatially and temporally accurate expression.

Figure 6. Multiple copies of a weak Ebf enhancer drive ectopic reporter gene expression.

Figure 6.

(A) Schematic representation of cardiopharyngeal lineage cells at Stage 27 (Hotta et al., 2007); First heart precursors (FHPs, red and open arrowheads), second heart precursors (SHPs, orange and arrows), inner ASM precursors and derivatives (iASMPs, violet and solid arrowhead), outer ASM precursors and derivatives (oASMPs, dark blue and solid arrowhead) (see (Razy-Krajka et al., 2014); black bars link sister cells. (B) Lineage tracing in individual larvae expressing single (1x) and multiple copy (3x) of _Ebf cis_-regulatory element, KhL24.37, driving H2B::mCherry (red) in cardiopharynegal progenitors at stage 27; B7.5 lineage is marked with Mesp>H2B::GFP (green). The single copy of KhL24.37 element drives H2B::mCherry reporter expression specifically in the ASMPs (upper left panel, in white); three copies of KhL24.37 (3x L24.37) drives expression in ASMPs and induces ectopic expression SHPs (lower left panel, in white). Experiment performed in biological replicate. Scale bar = 30 µm. (C) Proportions of embryos expressing H2B::mCherry in indicated cell-type progenitors by the indicated _cis_-regulaory elements. The ‘full length’ upstream region encompassing all combined enhancers with endogenous spacing (Ebf_-full length −3348 /- 178) (Wang et al., 2013) as well as the concatemer of KhL24.37, .36, .35, and .34 elements, lacking endogenous spacer sequences were used as controls. Statistical analysis using a Fisher exact test showed all comparisons with either control to be significant (p < 0.01); ‘n’ is the total number of individual halves scored per condition. (**D**) Schematic representation of cardiopharyngeal lineage cells at Stage 24 (Hotta et al., 2007); First heart precursors (FHPs, red and open arrowheads), secondary TVC (arrows). (**E**) Expression of GFP visualized by in situ hybridization on embryos at Stage 24 electroporated with single, multi copies and full-length of _Ebf cis_-regulatory element. Nuclei of B7.5 lineage cells are labelled by _Mesp>nls::LacZ and revealed with an anti beta-galactosidase antibody (red). _Mesp_-driven hCD4::mCherry accumulates at the cell membrane as revealed by anti mCherry antibody (Blue). Scale bar = 10 µm. (F) Proportions of embryos showing GFP-driven by the indicated _Ebf cis_-regulatory element (Fisher exact test, p < 0.01).

We reasoned that high but precocious activation of Ebf should suffice to trigger the autoregulatory loop, and cause ectopic pharyngeal muscle specification in the cells that normally form the second heart lineage, as observed previously (Razy-Krajka et al., 2014; Stolfi et al., 2010). To directly test this possibility, we used different combinations of the _cis_-regulatory elements upstream of Ebf to drive expression of a functional Ebf cDNA, and assayed the effects endogenous Ebf expression and cardiopharyngeal fates (Figure 7). As expected, one copy of the KhL24.37 element (1x KhL24.37) failed to cause ectopic activation of the endogenous locus, as evaluated using intronic probes, whereas using the KhL24.37 trimer (3x KhL24.37) to drive expression of the Ebf cDNA sufficed to activate the endogenous locus ectopically in ~40% of the embryos, as shown by the presence of nascent transcripts in 4 out of 6 nuclei per side, instead of 2 (Figure 7A,B). This observation is consistent with our model that the upstream elements mediate high activating inputs for the onset of Ebf expression, whereas maintenance upon commitment to a pharyngeal muscle identity relies on autoregulation (Figure 5F; Razy-Krajka et al., 2018).

Figure 7. Multimer of one weak Ebf enhancer drives ectopic pharyngeal muscle fate specification.

Figure 7.

(A) Targeted expression of an Ebf cDNA by three copies of KhL24.37 element induces expression of endogenous Ebf in four cells compared to the control where Ebf is detected only in the two ASMPs. Proportions of embryos expressing Ebf in ASMPs only or in ASMPs and SHPs in the indicated conditions. The single _KhL24.37 cis_-regulatory was used as control. (B) Expression of Ebf visualized by in situ hybridization on embryos at Stage 26 electroporated with single or three copies KhL24.37 element driving Ebf cDNA. Intron-specific probes show nascent Ebf transcripts in ASMPs. Nuclei of B7.5 lineage cells are labelled by _Mesp_>nls::LacZ and revealed with an anti beta-galactosidase antibody (Red). _Mesp_-driven hCD4::mCherry accumulates at the cell membrane as revealed by anti mCherry antibody (Blue). Scale bar = 10 µm. (C) Proportions of embryos showing GFP-driven STVC-specific enhancer of Tbx1/10 (Fisher exact test, p < 0.01). The ‘full length’ upstream element encompassing all combined enhancers with endogenous spacing (Ebf-full length −3348 / -178) as well as the ‘full length’ without endogenous spacing (_KhL24.37_, ._36_, ._35_, ._34_) and the single copy _KhL24.37_ element driving _Ebf_ were used as a positive controls. (**D**) Example of an embryo at 28 hpf showing GFP expression only in the ASMP (solid arrowhead) and SHP (arrow) but not in the FHP (open arrowheads), where _Tbx1/10_ enhancer is not active. Targeted expression of _Ebf_ cDNA by three copies of _KhL24.37_ element induces cells that normally form the second heart lineage to migrate alongside the pharyngeal muscle. Nuclei of B7.5 lineage cells are labelled by _Mesp_>H2B::mCherry. Scale bar = 15 µm. (E) Proposed role of combinatorial logics in fostering both precision and spatial and temporal accuracy of cell fate choices. The _Ebf_-full length, with or without endogenous spacers, fosters spatial/temporal accurate and precise cell fate choice, whereas the single copy of KhL24.37 (1x L24.37) gives spatiotemporally accurate reporter expression but is likely insufficient to induce a precise ASM fate. Multiple copies of _KhL24.37 cis_-regulatory element rescue a high reporter activity that reflects precise pharyngeal fate at the expense of spatial and temporal accuracy. The shapes of the distinct _cis_-regulatory elements are as in Figure 5C. Statistical analysis using a Fisher exact test (p < 0.01); ‘n’ is the total number of individual halves scored per condition.

These results suggested that inaccurate activation and maintenance of Ebf expression would cause ectopic pharyngeal muscle specification at the expense of the second heart lineage. We tested this possibility by analyzing pharyngeal muscle morphogenesis in stage 30 larvae, which is characterized by collective migration away from heart progenitors and formation of a ring of atrial siphon muscle precursors (Figure 7). We used the STVC-specific Tbx1/10 enhancer to visualize both the pharyngeal muscle precursors, which migrate and form a ring, and the second heart precursors, which remain associated with the _Tbx1/10_-negative first heart lineage (Figure 7D; Razy-Krajka et al., 2018). Remarkably, Ebf misexpression using three copies of the KhL24.37 element induced cells that normally form the second heart lineage to migrate alongside the pharyngeal muscles in 38% of larvae at 28 hpf (38 ± 3%, SE) (n = 265, Figure 7C–D). Importantly, neither one copy of KhL24.37 (1x KhL24.37), nor the full combined enhancers, with or without endogenous spacers, sufficed to cause substantial fate transformation of _Tbx1/10_+ second heart lineage into migratory pharyngeal muscle precursors (Figure 7C–D). These results indicate that driving expression of an Ebf cDNA by multimerizing a weak Ebf enhancer sufficed to cause ectopic activation of the endogenous locus and transformation of second heart lineage cells into migratory pharyngeal muscle precursors.

Taken together, these results provide evidence for a decoupling between two essential aspects of transcriptional activation, whereby multiple copies of a single regulatory element enable robust gene expression, compatible with precise fate specification, albeit at the expense of spatial and temporal accuracy; whereas individual regulatory elements appear to integrate distinct _trans_-acting inputs controlling enhancer accessibility and/or activity, thus increasing the repertoire of regulatory inputs controlling developmental gene expression (Figure 7E). In other words, while multiple elements are required for proper activation of cell fate determinants, such as Tbx1/10 and Ebf, in a manner reminiscent of super- and shadow enhancers (Lagha et al., 2012), we propose that combined enhancers foster spatially and temporally accurate cell fate decisions.

Discussion

We characterized the accessible genome of the tunicate Ciona, with a special focus on the cardiopharyngeal lineage that produces heart and pharyngeal muscles. As seen in other systems, less than 10% of the Ciona genome is accessible, and distributed across thousands of short regions, most of which are stably accessible across time and lineages, especially promoter regions. By contrast, developmentally regulated regions either closed upon induction of multipotent progenitors or opened specifically in the cardiopharyngeal lineage in response to FGF-MAPK signaling and Foxf activity. The latter elements were predominantly found in intergenic and intronic regions, and near cardiopharyngeal markers, consistent with their function as transcriptional enhancers. Similarly to other Forkhead factors (Zaret and Carroll, 2011), Ciona Foxf is required to open cardiopharyngeal elements for either immediate or later activation in multipotent or fate-restricted progenitors, respectively. Notably, Foxf homologs play deeply conserved roles in visceral muscles specification (Jakobsen et al., 2007; Scimone et al., 2018; Zaffran et al., 2001), including during heart development in mammals (Hoffmann et al., 2014). GATA motifs are also over-represented among cardiopharyngeal-specific elements, consistent with a conserved role for GATA homologs in heart development (Holtzinger and Evans, 2007; Molkentin et al., 1997; Qian and Bodmer, 2009; Reiter et al., 1999; Sorrentino et al., 2005; Zhao et al., 2008). As combinations of FOX and GATA inputs play well-established roles in early endoderm specification (Cirillo et al., 2002), we speculate that cardiopharyngeal regulatory programs were built upon an ancestral endomesodermal chromatin landscape during Olfactores evolution.

The majority of cell-type-specific markers expressed de novo are associated with ‘primed accessible’ elements, as observed in numerous systems including cardiac differentiation of embryonic stem cells (Paige et al., 2012; Wamstad et al., 2012), and consistent with the role of pioneer factors in establishing competence for subsequent activation (Zaret and Carroll, 2011). In the case of Tbx1/10 and Ebf, spatially and temporally accurate activation is essential to permit the emergence of first and second cardiac, and pharyngeal muscle lineages (Razy-Krajka et al., 2018; Razy-Krajka et al., 2014; Wang et al., 2019; Wang et al., 2013). We found that several elements, exhibiting distinct accessibility dynamics, are required for proper activation of both Tbx1/10 and Ebf. Specifically, a minimal distal enhancer proved necessary, but not sufficient for Ebf activation in newborn pharyngeal muscle precursors. By contrast, multiple copies of the same element sufficed to restore high transcriptional activity, but caused precocious activation in the _Tbx1/10_+ multipotent progenitors (aka STVCs, Figure 1A and Figure 6E), ectopic GFP expression in the second heart lineage (Figure 6), and eventually heart-to-pharyngeal muscle fate transformation when used to express Ebf itself (Figure 7). We propose that, whereas the activity of multiple elements with similar spatio-temporal transcriptional outputs permits precise and robust gene activation (Bentovim et al., 2017; Lagha et al., 2012), the modular organization of combined enhancers increases the repertoire of regulatory inputs, acting through both accessibility and activity, to control gene activation. Together with canalizing mechanisms, including positive auto-regulatory feedbacks, such multi-level combinatorial inputs achieve exquisite spatio-temporal control while permitting strong activation, thus ensuring both precise and accurate developmental fate choices (Figure 7E).

Materials and methods

Animals and electroporations

Gravid wild Ciona intestinalis type A, now called Ciona robusta (Pennati et al., 2015), were obtained from M-REP (Carlsbad, CA, USA), and kept under constant light to avoid spawning. Gametes from several animals were collected separately for in vitro cross-fertilization followed by dechorionation and electroporation as previously described (Christiaen et al., 2009a)​​, and cultured in filtered artificial seawater (FASW) in agarose-coated plastic Petri dishes at 18°C. Different quantities of plasmids were electroporated depending on the constructs: the amount of fluorescent reporter DNA (​_Mesp>nls::lacZ_​, Mesp>hCD4::mCherry_​, ​_Mesp>tagRFP_​, ​_MyoD905>eGFP and ​_Hand-r>tagBFP_​) and NLS::lacZ was typically 50 ​μg, but only 15 µg for Mesp>H2B::mCherry. For perturbation constructs (Mesp>FgfrDN, Mesp>MekS216D,S220E, Foxf>M RasCA, Hand-r>FgfrDN), 70 µg were usually electroporated, except for Mesp>nls::Cas9-Gem::nls (30 µg) and pairs of U6>sgRNA plasmids (25 µg each).

Molecular cloning

Putative enhancers were amplified from Ciona robusta genome using primers containing specific sequence tails (Supplementary file 4) for subcloning into a vector upstream of a basal promoter from Zfpm (aka Friend of GATA/Fog [Rothbächer et al., 2007]) driving expression of green fluorescent protein (GFP) fused to cytoplasmic unc76 (Stolfi et al., 2010).

The ‘full length’ Ebf enhancer with non-endogenous spacing in the genome (KhL24.37, .36, .35, .34) was generated by synthesizing DNA fragment (Twist Bioscience) and subcloning into the full-length reporter plasmids. nls::dCas9-KRAB::nls was derived from ‘pLV-dCas9-KRAB-PGK-HygR’ (Klann et al., 2017) (Addgene plasmid: #83890) and inserted downstream of the promoter Mesp (Davidson et al., 2005).

CRISPR/Cas9-mediated mutagenesis of ATAC-seq peaks

Two to four single guide RNAs (sgRNA) with Doench scores (http://crispor.tefor.net; Haeussler et al., 2016) higher than 60 were designed to induce deletions in selected accessible elements using CRISPR/Cas9 in the B7.5 lineage as described (Gandhi et al., 2017). sgRNAs targeting non-overlapping sequences per gene are listed in Supplementary file 5. The efficiency of sgRNAs was evaluated using the peakshift method as described (Gandhi et al., 2017)​​ (Figure 2—figure supplement 3B–C; Figure 2—figure supplement 6C). CRISPR/Cas9-mediated deletions were also evaluated by PCR-amplification directly from embryo lysates following electroporated with Eef1a>nls::Cas9-Gem::nls (Figure 2—figure supplement 6C). sgRNAs were expressed using the ​Ciona robusta ​U6 promoter (Stolfi et al., 2014) (Figure 2—figure supplement 3D; Figure 2—figure supplement 6C–D). For each peak, three or four guide RNAs were used in combination with 25 ​μg of each expression plasmid. 25 ​μg of ​Mesp>nls::Cas9-Gem::nls and Mesp>nls::dCas9-KRAB::nls plasmids were co-electroporated with guide RNA expression plasmids for B7.5 lineage-specific CRISPR/Cas9-mediated mutagenesis. Two guide RNAs were used to mutagenize Tyrosinase, which is not expressed in the cardiopharyngeal lineage and thus used to control the specificity of the CRISPR/Cas9 system (Wang et al., 2019).

Fluorescent ​in situ Hybridization-Immunohistochemistry (FISH-IHC) in ​Ciona embryos

FISH-IHC were performed as previously described (Christiaen et al., 2009b; Razy-Krajka et al., 2018)​. Embryos were harvested and fixed at desired developmental stages for 2 hr in 4% MEM-PFA and stored in 75% ethanol at −20°C. Antisense RNA probes were synthesized as described (Racioppi et al., 2014). For Ebf FISH-IHC, an anti-sense riboprobe targeting the full length cDNA was used at Stage 24 (Figure 5—figure supplement 2B) and 25 (Figure 5E, upper panel) where an anti-sense riboprobe targeting the first three intronic elements in the Ebf transcript as previously used in Wang et al. (2013) was used for stage 27 (Figure 5E, lower panel; Figure 6E). The template for Fbln anti-sense riboprobe was PCR-amplified with the following oligos: Fbln_pb_fw (5’ TTGCGCTAAGTCATGACAGC 3’), Fbln_pb_rev (5’CATTTGCCGATTCAGCTATGT3’). In vitro antisense RNA synthesis was performed using T7 RNA Polymerase (Roche, Cat. No. 10881767001) and DIG RNA Labeling Mix (Roche, Cat. No. 11277073910). Anti-Digoxigenin-POD Fab fragment (Roche, IN) was first used to detect the hybridized probes, then the signal was revealed using Tyramide Signal Amplification (TSA) with Fluorescein TSA Plus Evaluation Kits (Perkin Elmer, MA). Anti–β-galactosidase monoclonal mouse antibody (Promega) was co-incubated with anti-mCherry polyclonal rabbit antibody (Bio Vision, Cat. No. 5993–100) for immunodetection of ​Mesp>nls::lacZ and ​Mesp>hCD4::mCherry products respectively. Goat anti-mouse secondary antibodies coupled with AlexaFluor-555 and AlexaFluor-633 were used to detect β-galactosidase-bound mouse antibodies and mCherry-bound rabbit antibodies after the TSA reaction. FISH samples were mounted in ProLong Gold Antifade Mountant (ThermoFisher Scientific, Waltham, MA. Catalog number P36930).

Imaging

Images were acquired with an inverted Leica TCS SP8 X confocal microscope, using HC PL APO 63×/1.30 objective. Z-stacks were acquired with 1 μm z-steps. Maximum projections were processed with maximum projection tools from the LEICA software LAS-AF.

Cell dissociation and FACS

Sample dissociation and FACS were performed as previously described (Christiaen et al., 2016; Wang et al., 2018)​. Embryos and larvae were harvested at 6, 10, 15, 18 and 20 hpf in 5 ml borosilicate glass tubes (Fisher Scientific, Waltham, MA. Cat.No. 14-961-26) and washed with 2 ml calcium- and magnesium-free artificial seawater (CMF-ASW: 449 mM NaCl, 33 mM Na​2S​O4, 9 mM KCl, 2.15 mM NaHCO​3, ​ 10 mM Tris-Cl pH 8.2, 2.5 mM EGTA). Embryos and larvae were dissociated in 2 ml 0.2% trypsin (w/v, Sigma, T- 4799) in CMF-ASW by pipetting with glass Pasteur pipettes. The dissociation was stopped by adding 2 ml filtered ice cold 0.05% BSA CMF-ASW. Dissociated cells were passed through 40 μm cell-strainer and collected in 5 ml polystyrene round-bottom tube (Corning Life Sciences, Oneonta, New York. REF 352235). Cells were collected by centrifugation at 800 g for 3 min at 4°C, followed by two washes with ice cold 0.05% BSA CMF-ASW. Cell suspensions were filtered again through a 40 μm cell-strainer and stored on ice. Cell suspensions were used for sorting within 1 hr. B7.5 lineage cells were labeled by ​Mesp>tagRFP reporter. The B-line mesenchyme cells were counter-selected using ​MyoD905>GFP as described​ (Christiaen et al., 2016; Wang et al., 2018)​. The TVC-specific Hand-r>tagBFP reporter was used in a 3-color FACS scheme for positive co-selection of TVC-derived cells (Wang et al., 2018), in order to minimize the effects of mosaicism. Dissociated cells were loaded in a BD FACS AriaT​M cell sorter. 488 nm laser, FITC filter was used for GFP; 407 nm laser, 561 nm laser, DsRed filter was used for tagRFP and Pacific BlueTM filter was used for tagBFP. The nozzle size was 100 μm. eGFP+ mesenchyme cells were collected for downstream ATAC-seq analysis, whereas tagRFP+, tagBFP+ and eGFP- cardiopharyngeal lineage cells were collected for both ATAC- and RNA-seq analyses.

Preparation and sequencing of ATAC-seq library

ATAC-seq was performed with minor modifications to the published protocol (Buenrostro et al., 2015). 4,000 cells obtained by FACS were centrifuged at 800 x g for 4 min at 4°C for each sample. Cells were resuspended in 5 µL of cold lysis buffer (10 mM Tris-HCl [pH 7.4], 10 nM NaCl, 3 nM MgCl2, 0,1% v/v Igepal CA-360 [Sigma-Aldrich] and incubated on ice for 5 min. After centrifugation of the cells at 500 x g for 10 min at 4°C, the supernatant was discarded. The tagmentation reaction was performed at 37°C for 30 min using Nextera Sample Preparation kit (Illumina) with the addition of a tagmentation-stop step by the addition of EDTA to a final concentration of 50 nM and incubation at 50°C for 30 min (Hockman et al., 2019). After tagmentation, transposed DNA fragments were amplified using the following PCR condition: 1 cycle of 72°C for 5 min and 98°C for 30 s, followed by 12 to 14 cycles of 98°C for 10 s, 63°C for 30 s and 72°C for 1 min. Amplified libraries were purified using PCR purification MinElute kit (Qiagen), and the quality of of the purified library was assessed by using 2100 Bioanalyzer (Agilent Technologies) and a High Sensitivity DNA analysis Kit (Agilent Technologies) to confirm a period pattern in the size of the amplified DNA. After size selection by using Agencourt Ampure XP beads (Beckman Coulter, Cat#A633880) with a bead-to-sample ratio of 1.8:1, the libraries were sequenced as paired-end 50 bp by using the HiSeq 2500 platform (Illumina) according to the manufacturer’s instructions. An input control was also generated by using 10 ng of C. robusta genomic DNA.

ATAC-seq data analysis

Alignment of ATAC-seq reads

Raw reads were preprocessed by FastQC (version 0.11.2, http://www.bioinformatics.babraham.ac.uk/projects/fastqc) and adaptors were trimmed using Trim Galore (version 0.4.4, http://www.bioinformatics.babraham.ac.uk/projects/trim_galore) before being aligned to Ciona robusta genome (joined scaffold (KH), http://ghost.zool.kyoto-u.ac.jp/datas/JoinedScaffold.zip) using Bowtie2 (version 2.3.2, Langmead and Salzberg, 2012) with the parameters --no-discordant -p 20 -X 1000. Reads with mapping quality score > 30 were kept for downstream analysis using SAMtools (version 1.2, Li et al., 2009). Mitochondrial reads were removed using NGSutils 0.5.9 (Breese and Liu, 2013). At least three replicates for each sample were merged for peak calling by MACS2 (version 2.7.9) (Zhang et al., 2008) (--nomodel --nolambda --gsize 1.2e8). To correct for nonspecific sequencing biases, we subtracted sequenced gDNA from these libraries (Toenhake et al., 2018).We defined the accessome by concatenating all narrow peaks from MACS2 using bedtools merge (2.26.0) (Quinlan and Hall, 2010). Peaks were filtered for length > 50 bp, to obtain a final number of 56,090 peaks.

Accessome annotation

We annotated each peak to all transcripts from the C. robusta genome (http://ghost.zool.kyoto-u.ac.jp/datas/KH.KHGene.2013.gff3.zip) within 10 kb using GenomicRanges (version 1.28.6, Lawrence et al., 2013). We also identified whether each peak overlapped with the genomic feature categories ‘TSS, 5’ UTR, intron, exon, 3’ UTR, TTS, or intergenic region’ (Figure 1—figure supplement 1C). The KH2013 gene IDs are identical to the KH2012 gene IDs after merging all transcripts of each gene. Only the genomic feaatures differ between the two annotations. We used a TSS-seq dataset of high-density TSS regions, or TSS clusters (TSCs) (Yokomori et al., 2016) from Ciona larvae as guidelines for core promoter length. We used the mean TSC length plus two standard deviations to define a window of 107 bp upstream and downstream of the start of the 5’ UTR of any transcript as the TSS. We defined a window of 1 kb upstream from the TSS to be the promoter, which we further divided into two 500 bp windows. We defined the TTS to be a window 200 bp upstream and downstream from the end of the 3’ UTR of any transcript. We defined regions not covered by any feature to be intergenic. We performed a two-tailed binomial test for enrichment of accessible elements in each of these features, as well as the previously published TSCs. The fraction of the genome covered by each feature was used as the null probability. We tested for differential expression of promoters by extending the TSS 1 kb upstream, then aligning bulk RNA-seq reads to this putative promoter region.

Differential accessibility analysis

We evaluated the significant change in chromatin accessibility between different samples using DESeq2 (Love et al., 2014) with the parameters (method = ‘LRT’, alpha = 0.05, cooksCutoff = FALSE) for differential accessibility analysis. Libraries with fewer than 500,000 reads were removed. We performed a likelihood ratio test for four models:

  1. the 10 hpf model, tested for difference in accessibility between any of the conditions at 10 hpf (Mesp>FgfrDN, Mesp>MekS216D,S220E, U6>Foxf.1, U6>Foxf.4, listed in Supplementary file 5) and controls (Mesp>nls::LacZ and U6>Ngn.p1, U6>Ngn.p2 used for FoxfCRISPR, listed in Supplementary file 5).
  2. the 15–20 hpf model, tested for difference in accessibility between Hand-r>FgfrDN and control Hand-r>LacZ at any of the later time points (15, 18 or 20 hpf), Foxf>M RasCA and Foxf>nls::LacZ for the 18 hpf samples, versus a null model where accessibility is only dependent on time.
  3. the time model, tested for any change in accessibility between control B7.5 samples (Mesp>nls::LacZ and Hand-r>nls::LacZ) at 6, 10, 15, 18 and 20 hpf versus a null model dependent where accessibility is not changing.
  4. the control vs. mesenchyme model, tested for any difference in accessibility between B7.5 control (Mesp>nls::LacZ and Hand-r>nls::LacZ) and mesenchymal (MyoD905>GFP) samples versus a null model dependent at 10, 15, 18 and 20 hpf.

Cell-type-specific accessibility

We defined cell-type-specific accessibility based on the peak sets:

We tested for enrichment of peaks overlapping with each category of genomic feature in each of these peak sets (Figure 2—figure supplement 2C).

Definition of gene sets

We used cardiopharyngeal lineage cell-type-specific gene sets (cardiopharyngeal markers), including primed and de novo-expressed Cardiac and ASM/pharyngeal muscle as well as secondary TVC (STVC), first heart precursor (FHP), second heart precursor (SHP) markers defined by scRNA-seq (Wang et al., 2019). Anterior tail muscles (ATM) marker genes were derived from publicly available expression data (Brozovic et al., 2018) (Supplementary file 6). To these we added the gene sets either activated or inhibited by MAPK at 10 or 18 hpf as defined by microarray (Christiaen et al., 2008) and Bulk RNA-seq analysis (Wang et al., 2019).

Gene set enrichment analysis (GSEA)

We performed GSEA on elements associated with these gene sets using fgsea 1.2.1 (Sergushichev, 2016) with parameters: minSize = 5, maxSize = Inf, nperm = 10000. Because the gene sets used for GSEA need not be exclusive, we can use ATAC-seq peaks as ‘genes’ and define a gene set as all peaks annotated to genes in the set. We measured enrichment for peaks associated with the gene sets in open peaks (indicated by a positive enrichment score, ES) or closed peaks (indicated by a negative ES) in each comparison.

For GSEA of the ATAC-seq data, peaks were ranked by log2(FC) for a specific pairwise comparison, and all peaks annotated to a gene in a given gene set (as described above) were considered members of that gene set. The normalized enrichment score (NES) is a measure of whether members of the gene set are enriched at the top or bottom of a ranked list (Figure 1D). NES is calculated by keeping a running tally of whether or not the peak at each position of the list is a member of the gene set and comparing the maximum value to that of a null distribution calculated from a random permutation of peaks. For example, In Figure 1D,a positive NES indicates that peaks associated with a gene set are more open at 6 hpf than at 10 hpf, while a negative NES indicates that peaks associated with a gene set are more accessible at 10 hpf than at 6 hpf.

Ciona robusta motif inference

We obtained inferred C. intestinalis (former name of C. robusta) transcription factor (TF) binding motifs from CIS-BP (Weirauch et al., 2014). The TF binding motifs were inferred based on similarity of DNA-binding domain sequences between C. robusta TFs and those of TFs with known binding motifs. A strong correlation has been observed between conserved DNA-binding domains and bound motifs (Weirauch et al., 2014).

Motif selection

To the inferred motifs we added SELEX-seq motifs for C. robusta from Nitta et al. (2019). We then searched the C. robusta genome for orthologs to known motifs from the HOMER database (Heinz et al., 2010). All motifs with orthologous C. robusta TFs were added to the set of candidate TF motifs. This set of candidate TF motifs consisted of 1,745 motifs associated with 294 C. robusta TFs. Of these, 1,228 motifs could be unambiguously associated with one of 237 TFs. To obtain motifs for the remaining 57 TFs, to minimize redundancy, we remove from consideration all associations to TFs with an assigned motif, then assign all unambiguous motifs to a TF. This iteration continues until no more unambiguous motifs can be assigned. Motifs for the remaining TFs are assigned by repeating the process for motifs associated with only two TFs, then three TFs, until each TF has an assigned motif. Our final set of motifs had 1,487 members, with the median TF associated with two motifs.

For Figure 1—figure supplement 2E, we used previously published eukaryotic core promoter motifs (Haberle and Stark, 2018). For Figure 1—figure supplement 2F, correlations were performed on the inferred C. robusta DNA-binding motifs.

Motif analyses of ATAC-seq

We compiled the C. robusta genome (http://ghost.zool.kyoto-u.ac.jp/datas/JoinedScaffold.zip) as a BSGenome (version 1.44.2) package using ‘forgeBSgenomeDataPkg’. We searched for all motifs present in accessible elements using motifmatchr (version 0.99.8) and TFBSTools (version 1.14.2, Tan and Lenhard, 2016). Background nucleotide frequencies were computed for the whole accessome. Matches were returned if they passed a cutoff of p < 5*10−5. This cutoff was used to set a score threshold using the standard method of comparing the probability of observing the sequence given the PWM versus the probability of observing the sequence given the background nucleotide sequence of the accessible elements (Staden, 1989). We tested for motif enrichment in a peak set versus the accessome using a one-tailed hypergeometric test. The odds ratio represents the probability that an outcome (a peak containing the indicated motif) will occur if the peak is contained in the indicated peak set, compared to the probability of the outcome occurring in an element randomly selected from the accessome. For Figure 4—figure supplement 4A, the test indicates the probability that an element will contain a motif given when the element is accessible.

Primed or de novo accessible elements associated with de novo cardiac expressed genes were considered cardiac accessible if they had log2(FC) < 0 in Foxf>M RasCA vs. Foxf>nls::LacZ or log2(FC) > 0 in Hand-r>FgfrDN vs. Hand-r>nls::LacZ. Primed or de novo elements associated with de novo pharyngeal muscle (ASM) genes were considered ASM accessible if they had log2(FC) > 0 in Foxf>M RasCA vs. Foxf>nls::LacZ or log2(FC) < 0 in Hand-r>FgfrDN vs. Hand-r>nls::LacZ. The peak sets tested were primed cardiac accessible peaks associated with de novo cardiac expressed genes, de novo cardiac accessible peaks associated with de novo cardiac expressed genes, primed ASM accessible peaks associated with de novo ASM expressed genes, and de novo ASM accessible peaks associated with de novo ASM expressed genes. To ensure that the elements’ accessibility is not reflecting expression of primed genes, elements were removed if they were also associated with a gene expressed earlier in the cardiopharyngeal lineage. Elements associated with de novo ASM genes were removed if they were also associated with a cardiac or FGF-MAPK inhibited at 18 hpf gene. Elements associated with de novo cardiac genes were removed if they were also associated to an ASM or FGF-MAPK activated at 18 hpf gene.

Accessibility of motif sites in differentially accessible (FDR < 0.05) peaks in the 10 hpf model was calculated using chromVAR (version 0.99.3, Schep et al., 2017). Peaks were resized to 200 bp. Deviations were computed between all 10 hpf conditions as well as the 6 hpf controls. Replicates were grouped by both conditions and time. We tested for significantly differentially accessible motifs (FDR < 0.01) using the ‘differentialDeviations’ function from chromVAR.

chromVAR algorithm

The accessibility of a motif is calculated by summing the fragment counts for all accessible elements containing the motif. An expected accessibility of the motif is calculated from the total counts for a motif from all replicates normalized to the library size of each replicate. For each biological replicate, the raw accessibility deviation of each motif is defined as the difference of the observed counts and the expected counts divided by the expected counts. A normally distributed background accessibility deviation for each motif is calculated by iteratively sampling accessible elements from the same replicate. The deviation Z-score for each motif is given by the difference of the raw accessibility deviation and the mean of the background accessibility deviation divided by the standard deviation of the background accessibility deviation. For a more thorough description of the algorithm underlying chromVAR, see Schep et al. (2017).

Bulk RNA-seq library preparation, Sequencing and Analyses

We used bulk RNA-seq performed following defined perturbations of FGF-MAPK signaling and FACS (Wang et al., 2019). To profile transcriptomes of FACS-purified cells from FoxfCRISPR and control samples, 1,000 cells were directly sorted in 100 ​μl​ lysis buffer from the RNAqueous-Micro Total RNA Isolation Kit (Ambion). For each condition, samples were obtained in two biological replicates. The total RNA extraction was performed following the manufacturer’s instruction. The quality and quantity of total RNA was checked using Agilent RNA Screen Tape (Agilent) using 4200 TapeStation system. RNA samples with RNA Integrity Number (RIN) > 9 were kept for downstream cDNA synthesis. 250–2000 pg of total RNA was loaded as a template for cDNA synthesis using the SMART-Seq v4 Ultra Low Input RNA Kit (Clontech) with template switching technology. RNA-Seq Libraries were prepared and barcoded using Ovation Ultralow System V2 1–16 (NuGen). Up to six barcoded samples were pooled in one lane of the flow cell and sequenced by Illumina NextSeq 500 (MidOutput run). Paired-end 75 bp length reads were obtained from all the bulk RNA-seq libraries. Bulk RNA-seq libraries were aligned using STAR 2.5.3a (Zhang et al., 2018) with the parameters ‘--runThreadN 20 --bamRemoveDuplicatesType UniqueIdentical’. Counts were obtained using htseq-count (0.6.1p1) (Anders et al., 2015). Differential expression was calculated using edgeR for the 10 hpf conditions and DESeq2 for the conditions at 18 hpf.

Cell-type-specific expression gene sets

We tested for differential expression in the bulk RNA-seq data for the pairwise comparisons FoxfCRISPR vs. controlCRISPR, Foxf>M RasCA vs. control (Foxf>nls::LacZ), Hand-r>FgfrDN vs. control (Hand-r>nls::LacZ), and Foxf>M RasCA vs. Hand-r>FgfrDN using DESeq2 1.16.1 (Love et al., 2014). We defined genes downregulated in FoxfCRISPR vs. controlCRISPR (log2(FC) < −0.75 and FDR < 0.05) as Foxf targets. We defined ‘MAPK-inhibited genes at 18 hpf’ as the intersection of genes downregulated in Foxf>M RasCA vs. control and genes downregulated in Foxf>M RasCA vs. Hand-r>FgfrDN (log2(FC) < −1 and FDR < 0.05) as. We defined the intersect of genes downregulated in Hand-r>FgfrDN vs. control (log2(FC) < −1 and FDR < 0.05) and genes upregulated in Foxf>M RasCA vs. Hand-r>FgfrDN (log2(FC) > 1 and FDR < 0.05) as ‘MAPK activated genes’ (Figure 4—figure supplement 1A–D). For Figure 4B and Figure 4—figure supplement 3A, we defined the union of differentially expressed genes in either Foxf>M RasCA vs. control or Hand-r>FgfrDN vs. control as ‘differentially expressed upon FGF-MAPK perturbation at 18 hpf.’

We integrated microarray data to obtain gene sets for MAPK perturbation at 10 hpf (Christiaen et al., 2008). We defined genes downregulated in Mesp>FgfrDN vs. control (log2(FC) < −1 and p_ < 0.05) as ‘MAPK activated genes at 10 hpf.’ We defined genes upregulated in _Mesp>FgfrDN vs. control (Mesp>nls::LacZ) (log2(FC) > 1 and p < 0.05) as ‘MAPK inhibited genes at 10 hpf.’

From this same data set we obtained genes downregulated in control 6 hpf vs. 10 hpf (log2(FC) < −1 and _p_ < 0.05) and upregulated in _LacZ_ 6 hpf vs. 10 hpf (log2(FC) > 1 and p < 0.05).

Statistics

We used R (version 3.6.1) (Rizzo, 2016) to perform all statistical analysis. An alpha level of 0.05 was adapted for statistical significance throughout the analyses. We corrected for multiple hypothesis testing using false discovery rate (FDR) where indicated. All other _p_-values are unadjusted. For statistics presented in FISH panels, a Fisher’s exact test was run for each condition vs. the control. Embryos were classified as wild type, reduced expression, or no expression. A _p_-value less than 0.01 was considered significant (Figure 2F; Figure 2—figure supplement 5; Figure 3F; Figure 4G; Figure 5—figure supplement 1B,E). The two-tailed binomial tests in Figure 1—figure supplements 1C, 2D and 3G assume accessibility to not be feature-specific and dependent only on the fraction of the genome in base pairs covered by a feature. If promoter sites comprise 10% of the genome, we expect 10% of accessible elements will overlap with promoter sites by chance. For Figure 3—figure supplement 7A and B, we performed a two-tailed binomial test on the joint probability of each intersection against the null hypothesis that the probability of each set is independent. The null probability was calculated as the product of the probabilities of all constituent sets of an intersection. This is similar to the method used to calculate the deviation score in the original UpSet plot (Lex et al., 2014).

Software

All computations were performed on an x86_64-centos-linux-gnu platform. In addition to software specified elsewhere in the Materials and methods section, we created a SQLite database using RSQLite 2.1.1 (Müller et al., 2018), ComplexHeatmap (2.0.0) (Gu et al., 2016), circlize (0.4.6) (Gu et al., 2014), latticeExtra (0.6–28) (Sarkar and Andrews, 2016) and Integrative Genomics Viewer (Robinson et al., 2011).

Code availability

All analyses were done with publicly available software. All scripts for data analysis are available at https://github.com/ChristiaenLab/ATAC-seq. (Wiechecki and Racioppi, 2019; https://github.com/elifesciences-publications/ATAC-seq).

Data availability

All sequencing data were deposited on GEO with accession GSE126691.

Acknowledgements

We are grateful to Karen Lam, Sana Badri, Emily Miraldi and Richard Bonneau for help processing ATAC-seq data in the early phase of this study. We thank Dayanne M Castro for help and discussion on ATAC-seq data analysis. We thank Emily Huang for TVC-specific enhancers characterization, Tina Jiang for Tbx1/10 enhancer validation, Elizabeth Zelid and Carly Vaccaro for help with reporter assays. We thank Alberto Stolfi for assistance on cloning and constructive input and discussions. We are grateful to Wei Wang for sharing single cell and bulk RNA-seq data as well as Florian Razy-Krajka for sharing MekS216D,S220E prior to publication and to Pui Leng Ip for support with the FACS. We thank Esteban Mazzoni and members of the Christiaen lab for discussions, in particular Yelena Bernadskaya and Keaton Schuster for comments on the manuscript and figures. We thank Tatjana Sauka-Spengler for advice with ATAC-seq protocol. CR has been supported by a long-term fellowship ALTF 1608–2014 from EMBO. This work was funded by awards R01 HL108643 from NIH/NHLBI, R01 HD096770 from NIH/NICHD, and 15CVD01 from the Leducq Foundation to LC.

Funding Statement

The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication.

Contributor Information

Claudia Racioppi, Email: cr1636@nyu.edu.

Lionel Christiaen, Email: lc121@nyu.edu.

Brad Davidson, Swarthmore College, United States.

Didier Y Stainier, Max Planck Institute for Heart and Lung Research, Germany.

Funding Information

This paper was supported by the following grants:

Additional information

Competing interests

No competing interests declared.

Author contributions

Conceptualization, Data curation, Formal analysis, Funding acquisition, Validation, Investigation, Visualization, Writing—original draft, Writing—review and editing.

Data curation, Software, Formal analysis, Validation, Investigation, Visualization, Writing—original draft, Writing—review and editing.

Conceptualization, Resources, Supervision, Funding acquisition, Methodology, Writing—original draft, Project administration, Writing—review and editing.

Additional files

Supplementary file 1. Description of ATAC-seq dataset, with sample name, conditions, type of cells, number of mapped reads for all the biological replicates.

Supplementary file 2. Differentially FGF-MAPK-regulated genes (microarray) associated with differentially accessible peaks (ATAC-seq) in the pairwise comparison _Mesp_>FgfrDN vs. _Mesp_>LacZ – 10 hpf; genomic coordinates of each peak are reported.

Supplementary file 3. Known enhancers for TVC-specific genes characterized in previous studies; genomic coordinates of each peak are reported.

Supplementary file 4. Primer sequences used for cloning putative regulatory regions for functional enhancer assay.

Supplementary file 5. Primer sequences of single guide RNAs used to induce CRISPR/Cas9-mediated deletions in ATAC-seq peak as well as in Foxf (FoxfCRISPR) and to induce point mutations in Tbx1/10 STVC-specific enhancer; genomic coordinates and peak IDs are reported for each peak validated.

Supplementary file 6. ATM gene selected based on expression data deposited in ANISEED (https://www.aniseed.cnrs.fr).

Supplementary file 7. Differentially expressed genes from RNA-seq on FACS-purified cells following CRISPR/Cas9-induced loss of function of Foxf.

Supplementary file 8. Differentially accessible peaks from ATAC-seq on FACS-purified cells following CRISPR/Cas9-induced loss of function of Foxf (FoxfCRISPR vs. ControlCRISPR at 10 hpf).

Sheet1: closed in FoxfCRISPR log2(FC) < −0.5 and FDR < 0.05); sheet2: open in _FoxfCRISPR_ log2(FC) > 0.5 and FDR < 0.05.

Supplementary file 9. Differentially accessible peaks from ATAC-seq, closed in the founder cells (LacZ 6 hpf < LacZ 10 hpf) and closed in FoxfCRISPR vs. _ControlCRISPR_at 10 hpf.

Supplementary file 10. Differentially accessible peaks closed in FgfrDN at 10 hpf and associated with de novo cardiac (sheet 1) and pharyngeal muscle (sheet 2) expressed genes.

Supplementary file 11. Differentially FGF-MAPK-regulated genes associated with differentially accessible peaks in MAPK activated −18 hpf (sheet1) and MAPK inhibited −18 hpf (sheet 2); genomic coordinates of each peak are reported.

Supplementary file 12. FGF-MAPK-regulated genes associated with ATAC-seq peaks that are non-differentially accessible in MAPK activated −18 hpf and MAPK inhibited −18 hpf.

Supplementary file 13. Accessible elements associated with de novo expressed heart and pharyngeal muscle markers into pre-accessible/primed or de novo accessible elements; genomic coordinates of each peak are reported.

Supplementary file 14. FASTA sequences including the Foxf proteins (with their accession numbers) used for Figure 3—figure supplement 2.

Transparent reporting form

Data availability

All sequencing data were deposited on GEO with accession GSE126691.

The following dataset was generated:

Racioppi C, Wiechecki K, Christiaen L. 2019. Genome-wide maps of chromatin accessibility governing the transition from founder cells to distinct fate-restricted cardiopharyngeal precursors. NCBI Gene Expression Omnibus. GSE126691

The following previously published datasets were used:

Wang W, Niu X, Stuart T, Jullian E, Mauck WM, Kelly RG, Satija R, Christiaen L. 2019. A single cell transcriptional roadmap for cardiopharyngeal fate diversification. NCBI Gene Expression Omnibus. GSE99846

Christiaen L, Davidson B, Kawashima T, Powell W, Nolla H, Vranizan K, Levine M. 2008. Transcription profiling of B7.5 lineage cells in the tunicate Ciona intestinalis to study cardiac cell migration. ArrayExpress. E-MEXP-1478

Razy-Krajka F, Lam K, Wang W, Stolfi A, Joly M, Bonneau R, Christiaen L. 2014. Collier/OLF/EBF-dependent transcriptional dynamics control muscle specification from cardiopharyngeal progenitors. NCBI Gene Expression Omnibus. GSE54746

References

  1. Anders S, Pyl PT, Huber W. HTSeq--a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31:166–169. doi: 10.1093/bioinformatics/btu638. [DOI] [PMC free article] [PubMed] [Google Scholar]
  2. Barolo S. Shadow enhancers: frequently asked questions about distributed cis-regulatory information and enhancer redundancy. BioEssays. 2012;34:135–141. doi: 10.1002/bies.201100121. [DOI] [PMC free article] [PubMed] [Google Scholar]
  3. Beh J, Shi W, Levine M, Davidson B, Christiaen L. FoxF is essential for FGF-induced migration of heart progenitor cells in the ascidian Ciona intestinalis. Development. 2007;134:3297–3305. doi: 10.1242/dev.010140. [DOI] [PubMed] [Google Scholar]
  4. Bentovim L, Harden TT, DePace AH. Transcriptional precision and accuracy in development: from measurements to models and mechanisms. Development. 2017;144:3855–3866. doi: 10.1242/dev.146563. [DOI] [PMC free article] [PubMed] [Google Scholar]
  5. Bernadskaya YY, Brahmbhatt S, Gline SE, Wang W, Christiaen L. Discoidin-domain receptor coordinates cell-matrix adhesion and collective polarity in migratory cardiopharyngeal progenitors. Nature Communications. 2019;10:57. doi: 10.1038/s41467-018-07976-3. [DOI] [PMC free article] [PubMed] [Google Scholar]
  6. Boller S, Ramamoorthy S, Akbas D, Nechanitzky R, Burger L, Murr R, Schübeler D, Grosschedl R. Pioneering activity of the C-Terminal domain of EBF1 shapes the chromatin landscape for B cell programming. Immunity. 2016;44:527–541. doi: 10.1016/j.immuni.2016.02.021. [DOI] [PubMed] [Google Scholar]
  7. Breese MR, Liu Y. NGSUtils: a software suite for analyzing and manipulating next-generation sequencing datasets. Bioinformatics. 2013;29:494–496. doi: 10.1093/bioinformatics/bts731. [DOI] [PMC free article] [PubMed] [Google Scholar]
  8. Brozovic M, Dantec C, Dardaillon J, Dauga D, Faure E, Gineste M, Louis A, Naville M, Nitta KR, Piette J, Reeves W, Scornavacca C, Simion P, Vincentelli R, Bellec M, Aicha SB, Fagotto M, Guéroult-Bellone M, Haeussler M, Jacox E, Lowe EK, Mendez M, Roberge A, Stolfi A, Yokomori R, Brown CT, Cambillau C, Christiaen L, Delsuc F, Douzery E, Dumollard R, Kusakabe T, Nakai K, Nishida H, Satou Y, Swalla B, Veeman M, Volff JN, Lemaire P. ANISEED 2017: extending the integrated ascidian database to the exploration and evolutionary comparison of genome-scale datasets. Nucleic Acids Research. 2018;46:D718–D725. doi: 10.1093/nar/gkx1108. [DOI] [PMC free article] [PubMed] [Google Scholar]
  9. Buenrostro JD, Giresi PG, Zaba LC, Chang HY, Greenleaf WJ. Transposition of native chromatin for fast and sensitive epigenomic profiling of open chromatin, DNA-binding proteins and nucleosome position. Nature Methods. 2013;10:1213–1218. doi: 10.1038/nmeth.2688. [DOI] [PMC free article] [PubMed] [Google Scholar]
  10. Buenrostro JD, Wu B, Chang HY, Greenleaf WJ. ATAC-seq: a method for assaying chromatin accessibility Genome-Wide. Current Protocols in Molecular Biology. 2015;109: 21.29.1– 21.29.9. doi: 10.1002/0471142727.mb2129s109. [DOI] [PMC free article] [PubMed] [Google Scholar]
  11. Christiaen L, Davidson B, Kawashima T, Powell W, Nolla H, Vranizan K, Levine M. The transcription/migration interface in heart precursors of Ciona intestinalis. Science. 2008;320:1349–1352. doi: 10.1126/science.1158170. [DOI] [PubMed] [Google Scholar]
  12. Christiaen L, Wagner E, Shi W, Levine M. Isolation of sea squirt (Ciona) gametes, fertilization, dechorionation, and development. Cold Spring Harbor Protocols. 2009a;2009:pdb.prot5344. doi: 10.1101/pdb.prot5344. [DOI] [PubMed] [Google Scholar]
  13. Christiaen L, Wagner E, Shi W, Levine M. Whole-mount in situ hybridization on sea squirt (Ciona intestinalis) embryos. Cold Spring Harbor Protocols. 2009b;2009:pdb.prot5348. doi: 10.1101/pdb.prot5348. [DOI] [PubMed] [Google Scholar]
  14. Christiaen L, Wagner E, Shi W, Levine M. Isolation of individual cells and tissues from electroporated sea squirt (Ciona) Embryos by Fluorescence-Activated cell sorting (FACS) Cold Spring Harbor Protocols. 2016;2009:pdb.prot5349. doi: 10.1101/pdb.prot5349. [DOI] [PubMed] [Google Scholar]
  15. Cirillo LA, Lin FR, Cuesta I, Friedman D, Jarnik M, Zaret KS. Opening of compacted chromatin by early developmental transcription factors HNF3 (FoxA) and GATA-4. Molecular Cell. 2002;9:279–289. doi: 10.1016/S1097-2765(02)00459-8. [DOI] [PubMed] [Google Scholar]
  16. Clark KL, Halay ED, Lai E, Burley SK. Co-crystal structure of the HNF-3/fork head DNA-recognition motif resembles histone H5. Nature. 1993;364:412–420. doi: 10.1038/364412a0. [DOI] [PubMed] [Google Scholar]
  17. Cooley J, Whitaker S, Sweeney S, Fraser S, Davidson B. Cytoskeletal polarity mediates localized induction of the heart progenitor lineage. Nature Cell Biology. 2011;13:952–957. doi: 10.1038/ncb2291. [DOI] [PMC free article] [PubMed] [Google Scholar]
  18. Cusanovich DA, Reddington JP, Garfield DA, Daza RM, Aghamirzaie D, Marco-Ferreres R, Pliner HA, Christiansen L, Qiu X, Steemers FJ, Trapnell C, Shendure J, Furlong EEM. The cis-regulatory dynamics of embryonic development at single-cell resolution. Nature. 2018;555:538–542. doi: 10.1038/nature25981. [DOI] [PMC free article] [PubMed] [Google Scholar]
  19. Daugherty AC, Yeo RW, Buenrostro JD, Greenleaf WJ, Kundaje A, Brunet A. Chromatin accessibility dynamics reveal novel functional enhancers in C. elegans. Genome Research. 2017;27:2096–2107. doi: 10.1101/gr.226233.117. [DOI] [PMC free article] [PubMed] [Google Scholar]
  20. Davidson B, Shi W, Levine M. Uncoupling heart cell specification and migration in the simple chordate Ciona intestinalis. Development. 2005;132:4811–4818. doi: 10.1242/dev.02051. [DOI] [PubMed] [Google Scholar]
  21. Davidson B, Shi W, Beh J, Christiaen L, Levine M. FGF signaling delineates the cardiac progenitor field in the simple chordate, Ciona intestinalis. Genes & Development. 2006;20:2728–2738. doi: 10.1101/gad.1467706. [DOI] [PMC free article] [PubMed] [Google Scholar]
  22. Davidson EH. The regulatory genome: gene regulatory networks in development and evolution. Elsevier; 2010. [Google Scholar]
  23. Desjardins CA, Naya FJ. The function of the MEF2 family of transcription factors in cardiac development, cardiogenomics, and direct reprogramming. Journal of Cardiovascular Development and Disease. 2016;3:26. doi: 10.3390/jcdd3030026. [DOI] [PMC free article] [PubMed] [Google Scholar]
  24. Diogo R, Kelly RG, Christiaen L, Levine M, Ziermann JM, Molnar JL, Noden DM, Tzahor E. A new heart for a new head in vertebrate cardiopharyngeal evolution. Nature. 2015;520:466–473. doi: 10.1038/nature14435. [DOI] [PMC free article] [PubMed] [Google Scholar]
  25. Frankel N, Davis GK, Vargas D, Wang S, Payre F, Stern DL. Phenotypic robustness conferred by apparently redundant transcriptional enhancers. Nature. 2010;466:490–493. doi: 10.1038/nature09158. [DOI] [PMC free article] [PubMed] [Google Scholar]
  26. Frankel N. Multiple layers of complexity in cis-regulatory regions of developmental genes. Developmental Dynamics : An Official Publication of the American Association of Anatomists. 2012;241:1857–1866. doi: 10.1002/dvdy.23871. [DOI] [PubMed] [Google Scholar]
  27. Gandhi S, Haeussler M, Razy-Krajka F, Christiaen L, Stolfi A. Evaluation and rational design of guide RNAs for efficient CRISPR/Cas9-mediated mutagenesis in Ciona. Developmental Biology. 2017;425:8–20. doi: 10.1016/j.ydbio.2017.03.003. [DOI] [PMC free article] [PubMed] [Google Scholar]
  28. Ganot P, Kallesoe T, Reinhardt R, Chourrout D, Thompson EM. Spliced-leader RNA trans splicing in a chordate, Oikopleura dioica, with a compact genome. Molecular and Cellular Biology. 2004;24:7795–7805. doi: 10.1128/MCB.24.17.7795-7805.2004. [DOI] [PMC free article] [PubMed] [Google Scholar]
  29. Gaulton KJ, Nammo T, Pasquali L, Simon JM, Giresi PG, Fogarty MP, Panhuis TM, Mieczkowski P, Secchi A, Bosco D, Berney T, Montanya E, Mohlke KL, Lieb JD, Ferrer J. A map of open chromatin in human pancreatic islets. Nature Genetics. 2010;42:255–259. doi: 10.1038/ng.530. [DOI] [PMC free article] [PubMed] [Google Scholar]
  30. Groner AC, Meylan S, Ciuffi A, Zangger N, Ambrosini G, Dénervaud N, Bucher P, Trono D. KRAB-zinc finger proteins and KAP1 can mediate long-range transcriptional repression through heterochromatin spreading. PLOS Genetics. 2010;6:e1000869. doi: 10.1371/journal.pgen.1000869. [DOI] [PMC free article] [PubMed] [Google Scholar]
  31. Gu Z, Gu L, Eils R, Schlesner M, Brors B. Circlize implements and enhances circular visualization in R. Bioinformatics. 2014;30:2811–2812. doi: 10.1093/bioinformatics/btu393. [DOI] [PubMed] [Google Scholar]
  32. Gu Z, Eils R, Schlesner M. Complex heatmaps reveal patterns and correlations in multidimensional genomic data. Bioinformatics. 2016;32:2847–2849. doi: 10.1093/bioinformatics/btw313. [DOI] [PubMed] [Google Scholar]
  33. Haberle V, Stark A. Eukaryotic core promoters and the functional basis of transcription initiation. Nature Reviews Molecular Cell Biology. 2018;19:621–637. doi: 10.1038/s41580-018-0028-8. [DOI] [PMC free article] [PubMed] [Google Scholar]
  34. Haeussler M, Schönig K, Eckert H, Eschstruth A, Mianné J, Renaud JB, Schneider-Maunoury S, Shkumatava A, Teboul L, Kent J, Joly JS, Concordet JP. Evaluation of off-target and on-target scoring algorithms and integration into the guide RNA selection tool CRISPOR. Genome Biology. 2016;17:148. doi: 10.1186/s13059-016-1012-2. [DOI] [PMC free article] [PubMed] [Google Scholar]
  35. Hastings KE. SL trans-splicing: easy come or easy go? Trends in Genetics. 2005;21:240–247. doi: 10.1016/j.tig.2005.02.005. [DOI] [PubMed] [Google Scholar]
  36. He A, Gu F, Hu Y, Ma Q, Ye LY, Akiyama JA, Visel A, Pennacchio LA, Pu WT. Dynamic GATA4 enhancers shape the chromatin landscape central to heart development and disease. Nature Communications. 2014;5:4907. doi: 10.1038/ncomms5907. [DOI] [PMC free article] [PubMed] [Google Scholar]
  37. Heinz S, Benner C, Spann N, Bertolino E, Lin YC, Laslo P, Cheng JX, Murre C, Singh H, Glass CK. Simple combinations of lineage-determining transcription factors prime cis-regulatory elements required for macrophage and B cell identities. Molecular Cell. 2010;38:576–589. doi: 10.1016/j.molcel.2010.05.004. [DOI] [PMC free article] [PubMed] [Google Scholar]
  38. Hnisz D, Abraham BJ, Lee TI, Lau A, Saint-André V, Sigova AA, Hoke HA, Young RA. Super-enhancers in the control of cell identity and disease. Cell. 2013;155:934–947. doi: 10.1016/j.cell.2013.09.053. [DOI] [PMC free article] [PubMed] [Google Scholar]
  39. Hockman D, Chong-Morrison V, Green SA, Gavriouchkina D, Candido-Ferreira I, Ling ITC, Williams RM, Amemiya CT, Smith JJ, Bronner ME, Sauka-Spengler T. A genome-wide assessment of the ancestral neural crest gene regulatory network. Nature Communications. 2019;10:4689. doi: 10.1038/s41467-019-12687-4. [DOI] [PMC free article] [PubMed] [Google Scholar]
  40. Hoffmann AD, Yang XH, Burnicka-Turek O, Bosman JD, Ren X, Steimle JD, Vokes SA, McMahon AP, Kalinichenko VV, Moskowitz IP. Foxf genes integrate tbx5 and hedgehog pathways in the second heart field for cardiac septation. PLOS Genetics. 2014;10:e1004604. doi: 10.1371/journal.pgen.1004604. [DOI] [PMC free article] [PubMed] [Google Scholar]
  41. Holtzinger A, Evans T. Gata5 and Gata6 are functionally redundant in zebrafish for specification of cardiomyocytes. Developmental Biology. 2007;312:613–622. doi: 10.1016/j.ydbio.2007.09.018. [DOI] [PMC free article] [PubMed] [Google Scholar]
  42. Hong JW, Hendrix DA, Levine MS. Shadow enhancers as a source of evolutionary novelty. Science. 2008;321:1314. doi: 10.1126/science.1160631. [DOI] [PMC free article] [PubMed] [Google Scholar]
  43. Hotta K, Mitsuhara K, Takahashi H, Inaba K, Oka K, Gojobori T, Ikeo K. A web-based interactive developmental table for the ascidian Ciona intestinalis , including 3D real-image embryo reconstructions: I. from fertilized egg to hatching larva. Developmental Dynamics. 2007;236:1790–1805. doi: 10.1002/dvdy.21188. [DOI] [PubMed] [Google Scholar]
  44. Hudson C, Darras S, Caillol D, Yasuo H, Lemaire P. A conserved role for the MEK signalling pathway in neural tissue specification and posteriorisation in the invertebrate chordate, the ascidian Ciona intestinalis. Development. 2003;130:147–159. doi: 10.1242/dev.00200. [DOI] [PubMed] [Google Scholar]
  45. Jakobsen JS, Braun M, Astorga J, Gustafson EH, Sandmann T, Karzynski M, Carlsson P, Furlong EEM. Temporal ChIP-on-chip reveals biniou as a universal regulator of the visceral muscle transcriptional network. Genes & Development. 2007;21:2448–2460. doi: 10.1101/gad.437607. [DOI] [PMC free article] [PubMed] [Google Scholar]
  46. Jänes J, Dong Y, Schoof M, Serizay J, Appert A, Cerrato C, Woodbury C, Chen R, Gemma C, Huang N, Kissiov D, Stempor P, Steward A, Zeiser E, Sauer S, Ahringer J. Chromatin accessibility dynamics across C. elegans development and ageing. eLife. 2018;7:e37344. doi: 10.7554/eLife.37344. [DOI] [PMC free article] [PubMed] [Google Scholar]
  47. Kaplan N, Razy-Krajka F, Christiaen L. Regulation and evolution of cardiopharyngeal cell identity and behavior: insights from simple chordates. Current Opinion in Genetics & Development. 2015;32:119–128. doi: 10.1016/j.gde.2015.02.008. [DOI] [PMC free article] [PubMed] [Google Scholar]
  48. Khan A, Mathelier A, Zhang X. Super-enhancers are transcriptionally more active and cell type-specific than stretch enhancers. Epigenetics. 2018;13:910–922. doi: 10.1080/15592294.2018.1514231. [DOI] [PMC free article] [PubMed] [Google Scholar]
  49. Klann TS, Black JB, Chellappan M, Safi A, Song L, Hilton IB, Crawford GE, Reddy TE, Gersbach CA. CRISPR-Cas9 epigenome editing enables high-throughput screening for functional regulatory elements in the human genome. Nature Biotechnology. 2017;35:561–568. doi: 10.1038/nbt.3853. [DOI] [PMC free article] [PubMed] [Google Scholar]
  50. Klemm SL, Shipony Z, Greenleaf WJ. Chromatin accessibility and the regulatory epigenome. Nature Reviews Genetics. 2019;20:207–220. doi: 10.1038/s41576-018-0089-8. [DOI] [PubMed] [Google Scholar]
  51. Kusakabe T, Yoshida R, Ikeda Y, Tsuda M. Computational discovery of DNA motifs associated with cell type-specific gene expression in Ciona. Developmental Biology. 2004;276:563–580. doi: 10.1016/j.ydbio.2004.09.037. [DOI] [PubMed] [Google Scholar]
  52. Kvon EZ, Kazmar T, Stampfel G, Yáñez-Cuna JO, Pagani M, Schernhuber K, Dickson BJ, Stark A. Genome-scale functional characterization of Drosophila developmental enhancers in vivo. Nature. 2014;512:91–95. doi: 10.1038/nature13395. [DOI] [PubMed] [Google Scholar]
  53. Lagha M, Bothma JP, Levine M. Mechanisms of transcriptional precision in animal development. Trends in Genetics. 2012;28:409–416. doi: 10.1016/j.tig.2012.03.006. [DOI] [PMC free article] [PubMed] [Google Scholar]
  54. Lam DD, de Souza FS, Nasif S, Yamashita M, López-Leal R, Otero-Corchon V, Meece K, Sampath H, Mercer AJ, Wardlaw SL, Rubinstein M, Low MJ. Partially redundant enhancers cooperatively maintain mammalian pomc expression above a critical functional threshold. PLOS Genetics. 2015;11:e1004935. doi: 10.1371/journal.pgen.1004935. [DOI] [PMC free article] [PubMed] [Google Scholar]
  55. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nature Methods. 2012;9:357–359. doi: 10.1038/nmeth.1923. [DOI] [PMC free article] [PubMed] [Google Scholar]
  56. Lawrence M, Huber W, Pagès H, Aboyoun P, Carlson M, Gentleman R, Morgan MT, Carey VJ. Software for computing and annotating genomic ranges. PLOS Computational Biology. 2013;9:e1003118. doi: 10.1371/journal.pcbi.1003118. [DOI] [PMC free article] [PubMed] [Google Scholar]
  57. Lescroart F, Kelly RG, Le Garrec J-F, Nicolas J-F, Meilhac SM, Buckingham M. Clonal analysis reveals common lineage relationships between head muscles and second heart field derivatives in the mouse embryo. Development. 2010;137:3269–3279. doi: 10.1242/dev.050674. [DOI] [PubMed] [Google Scholar]
  58. Lescroart F, Chabab S, Lin X, Rulands S, Paulissen C, Rodolosse A, Auer H, Achouri Y, Dubois C, Bondue A, Simons BD, Blanpain C. Early lineage restriction in temporally distinct populations of Mesp1 progenitors during mammalian heart development. Nature Cell Biology. 2014;16:829–840. doi: 10.1038/ncb3024. [DOI] [PMC free article] [PubMed] [Google Scholar]
  59. Lescroart F, Wang X, Lin X, Swedlund B, Gargouri S, Sànchez-Dànes A, Moignard V, Dubois C, Paulissen C, Kinston S, Göttgens B, Blanpain C. Defining the earliest step of cardiovascular lineage segregation by single-cell RNA-seq. Science. 2018;359:1177–1181. doi: 10.1126/science.aao4174. [DOI] [PMC free article] [PubMed] [Google Scholar]
  60. Lex A, Gehlenborg N, Strobelt H, Vuillemot R, Pfister H. UpSet: visualization of intersecting sets. IEEE Transactions on Visualization and Computer Graphics. 2014;20:1983–1992. doi: 10.1109/TVCG.2014.2346248. [DOI] [PMC free article] [PubMed] [Google Scholar]
  61. Li L, Zhu Q, He X, Sinha S, Halfon MS. Large-scale analysis of transcriptional cis-regulatory modules reveals both common features and distinct subclasses. Genome Biology. 2007;8:R101. doi: 10.1186/gb-2007-8-6-r101. [DOI] [PMC free article] [PubMed] [Google Scholar]
  62. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, Marth G, Abecasis G, Durbin R, 1000 Genome Project Data Processing Subgroup The sequence alignment/Map format and SAMtools. Bioinformatics. 2009;25:2078–2079. doi: 10.1093/bioinformatics/btp352. [DOI] [PMC free article] [PubMed] [Google Scholar]
  63. Long HK, Prescott SL, Wysocka J. Ever-Changing landscapes: transcriptional enhancers in development and evolution. Cell. 2016;167:1170–1187. doi: 10.1016/j.cell.2016.09.018. [DOI] [PMC free article] [PubMed] [Google Scholar]
  64. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biology. 2014;15:550. doi: 10.1186/s13059-014-0550-8. [DOI] [PMC free article] [PubMed] [Google Scholar]
  65. Madgwick A, Magri MS, Dantec C, Gailly D, Fiuza U-M, Guignard L, Hettinger S, Gomez-Skarmeta JL, Lemaire P. Evolution of embryonic cis-regulatory landscapes between divergent Phallusia and Ciona ascidians. Developmental Biology. 2019;448:71–87. doi: 10.1016/j.ydbio.2019.01.003. [DOI] [PubMed] [Google Scholar]
  66. Mavrich TN, Ioshikhes IP, Venters BJ, Jiang C, Tomsho LP, Qi J, Schuster SC, Albert I, Pugh BF. A barrier nucleosome model for statistical positioning of nucleosomes throughout the yeast genome. Genome Research. 2008a;18:1073–1083. doi: 10.1101/gr.078261.108. [DOI] [PMC free article] [PubMed] [Google Scholar]
  67. Mavrich TN, Jiang C, Ioshikhes IP, Li X, Venters BJ, Zanton SJ, Tomsho LP, Qi J, Glaser RL, Schuster SC, Gilmour DS, Albert I, Pugh BF. Nucleosome organization in the Drosophila genome. Nature. 2008b;453:358–362. doi: 10.1038/nature06929. [DOI] [PMC free article] [PubMed] [Google Scholar]
  68. Mayran A, Khetchoumian K, Hariri F, Pastinen T, Gauthier Y, Balsalobre A, Drouin J. Pioneer factor Pax7 deploys a stable enhancer repertoire for specification of cell fate. Nature Genetics. 2018;50:259–269. doi: 10.1038/s41588-017-0035-2. [DOI] [PubMed] [Google Scholar]
  69. McWilliam H, Li W, Uludag M, Squizzato S, Park YM, Buso N, Cowley AP, Lopez R. Analysis tool web services from the EMBL-EBI. Nucleic Acids Research. 2013;41:W597–W600. doi: 10.1093/nar/gkt376. [DOI] [PMC free article] [PubMed] [Google Scholar]
  70. Meedel TH, Chang P, Yasuo H. Muscle development in Ciona intestinalis requires the b-HLH myogenic regulatory factor gene Ci-MRF. Developmental Biology. 2007;302:333–344. doi: 10.1016/j.ydbio.2006.09.043. [DOI] [PMC free article] [PubMed] [Google Scholar]
  71. Molkentin JD, Lin Q, Duncan SA, Olson EN. Requirement of the transcription factor GATA4 for heart tube formation and ventral morphogenesis. Genes & Development. 1997;11:1061–1072. doi: 10.1101/gad.11.8.1061. [DOI] [PubMed] [Google Scholar]
  72. Müller K, Wickham H, James DA, Falcon S. RSQLite: “SQLite” Interface for R. 2.1.12018 https://cran.r-project.org/web/packages/RSQLite/vignettes/RSQLite.html
  73. Nitta KR, Vincentelli R, Jacox E, Cimino A, Ohtsuka Y, Sobral D, Satou Y, Cambillau C, Lemaire P. High-Throughput protein production combined with high- Throughput SELEX identifies an extensive atlas of Ciona robusta transcription factor DNA-Binding specificities. Methods in Molecular Biology. 2019;2025:487–517. doi: 10.1007/978-1-4939-9624-7_23. [DOI] [PubMed] [Google Scholar]
  74. Nowotschin S, Liao J, Gage PJ, Epstein JA, Campione M, Morrow BE. Tbx1 affects asymmetric cardiac morphogenesis by regulating Pitx2 in the secondary heart field. Development. 2006;133:1565–1573. doi: 10.1242/dev.02309. [DOI] [PubMed] [Google Scholar]
  75. Paige SL, Thomas S, Stoick-Cooper CL, Wang H, Maves L, Sandstrom R, Pabon L, Reinecke H, Pratt G, Keller G, Moon RT, Stamatoyannopoulos J, Murry CE. A temporal chromatin signature in human embryonic stem cells identifies regulators of cardiac development. Cell. 2012;151:221–232. doi: 10.1016/j.cell.2012.08.027. [DOI] [PMC free article] [PubMed] [Google Scholar]
  76. Pennati R, Ficetola GF, Brunetti R, Caicci F, Gasparini F, Griggio F, Sato A, Stach T, Kaul-Strehlow S, Gissi C, Manni L. Morphological differences between larvae of the Ciona intestinalis species complex: hints for a valid taxonomic definition of distinct species. PLOS ONE. 2015;10:e0122879. doi: 10.1371/journal.pone.0122879. [DOI] [PMC free article] [PubMed] [Google Scholar]
  77. Perry MW, Boettiger AN, Bothma JP, Levine M. Shadow enhancers foster robustness of Drosophila gastrulation. Current Biology. 2010;20:1562–1567. doi: 10.1016/j.cub.2010.07.043. [DOI] [PMC free article] [PubMed] [Google Scholar]
  78. Perry MW, Boettiger AN, Levine M. Multiple enhancers ensure precision of gap gene-expression patterns in the Drosophila embryo. PNAS. 2011;108:13570–13575. doi: 10.1073/pnas.1109873108. [DOI] [PMC free article] [PubMed] [Google Scholar]
  79. Pikkarainen S, Tokola H, Kerkelä R, Ruskoaho H. GATA transcription factors in the developing and adult heart. Cardiovascular Research. 2004;63:196–207. doi: 10.1016/j.cardiores.2004.03.025. [DOI] [PubMed] [Google Scholar]
  80. Pott S, Lieb JD. What are super-enhancers? Nature Genetics. 2015;47:8–12. doi: 10.1038/ng.3167. [DOI] [PubMed] [Google Scholar]
  81. Qian L, Bodmer R. Partial loss of GATA factor Pannier impairs adult heart function in Drosophila. Human Molecular Genetics. 2009;18:3153–3163. doi: 10.1093/hmg/ddp254. [DOI] [PMC free article] [PubMed] [Google Scholar]
  82. Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26:841–842. doi: 10.1093/bioinformatics/btq033. [DOI] [PMC free article] [PubMed] [Google Scholar]
  83. Racioppi C, Kamal AK, Razy-Krajka F, Gambardella G, Zanetti L, di Bernardo D, Sanges R, Christiaen LA, Ristoratore F. Fibroblast growth factor signalling controls nervous system patterning and pigment cell formation in Ciona intestinalis. Nature Communications. 2014;5:4830. doi: 10.1038/ncomms5830. [DOI] [PMC free article] [PubMed] [Google Scholar]
  84. Ragkousi K, Beh J, Sweeney S, Starobinska E, Davidson B. A single GATA factor plays discrete, lineage specific roles in ascidian heart development. Developmental Biology. 2011;352:154–163. doi: 10.1016/j.ydbio.2011.01.007. [DOI] [PMC free article] [PubMed] [Google Scholar]
  85. Razy-Krajka F, Lam K, Wang W, Stolfi A, Joly M, Bonneau R, Christiaen L. Collier/OLF/EBF-dependent transcriptional dynamics control pharyngeal muscle specification from primed cardiopharyngeal progenitors. Developmental Cell. 2014;29:263–276. doi: 10.1016/j.devcel.2014.04.001. [DOI] [PMC free article] [PubMed] [Google Scholar]
  86. Razy-Krajka F, Gravez B, Kaplan N, Racioppi C, Wang W, Christiaen L. An FGF-driven feed-forward circuit patterns the cardiopharyngeal mesoderm in space and time. eLife. 2018;7:e29656. doi: 10.7554/eLife.29656. [DOI] [PMC free article] [PubMed] [Google Scholar]
  87. Reiter JF, Alexander J, Rodaway A, Yelon D, Patient R, Holder N, Stainier DY. Gata5 is required for the development of the heart and endoderm in zebrafish. Genes & Development. 1999;13:2983–2995. doi: 10.1101/gad.13.22.2983. [DOI] [PMC free article] [PubMed] [Google Scholar]
  88. Reynolds N, Salmon-Divon M, Dvinge H, Hynes-Allen A, Balasooriya G, Leaford D, Behrens A, Bertone P, Hendrich B. NuRD-mediated deacetylation of H3K27 facilitates recruitment of polycomb repressive complex 2 to direct gene repression. The EMBO Journal. 2012;31:593–605. doi: 10.1038/emboj.2011.431. [DOI] [PMC free article] [PubMed] [Google Scholar]
  89. Rizzo ML. Statistical Computing with R. CRC Press; 2016. [Google Scholar]
  90. Robinson JT, Thorvaldsdóttir H, Winckler W, Guttman M, Lander ES, Getz G, Mesirov JP. Integrative genomics viewer. Nature Biotechnology. 2011;29:24–26. doi: 10.1038/nbt.1754. [DOI] [PMC free article] [PubMed] [Google Scholar]
  91. Rosa-Garrido M, Karbassi E, Monte E, Vondriska TM. Regulation of chromatin structure in the cardiovascular system. Circulation Journal. 2013;77:1389–1398. doi: 10.1253/circj.CJ-13-0176. [DOI] [PMC free article] [PubMed] [Google Scholar]
  92. Rothbächer U, Bertrand V, Lamy C, Lemaire P. A combinatorial code of maternal GATA, ets and beta-catenin-TCF transcription factors specifies and patterns the early ascidian ectoderm. Development. 2007;134:4023–4032. doi: 10.1242/dev.010850. [DOI] [PubMed] [Google Scholar]
  93. Sarkar D, Andrews F. latticeExtra: Extra Graphical Utilities Based on Lattice. 2016 https://cran.r-project.org/web/packages/latticeExtra/index.html
  94. Satou Y, Imai KS, Satoh N. Early embryonic expression of a LIM-homeobox gene Cs-lhx3 is downstream of beta-catenin and responsible for the endoderm differentiation in Ciona savignyi embryos. Development. 2001a;128:3559–3570. doi: 10.1242/dev.128.18.3559. [DOI] [PubMed] [Google Scholar]
  95. Satou Y, Takatori N, Yamada L, Mochizuki Y, Hamaguchi M, Ishikawa H, Chiba S, Imai K, Kano S, Murakami SD, Nakayama A, Nishino A, Sasakura Y, Satoh G, Shimotori T, Shin-I T, Shoguchi E, Suzuki MM, Takada N, Utsumi N, Yoshida N, Saiga H, Kohara Y, Satoh N. Gene expression profiles in Ciona intestinalis tailbud embryos. Development. 2001b;128:2893–2904. doi: 10.1242/dev.128.15.2893. [DOI] [PubMed] [Google Scholar]
  96. Satou Y, Hamaguchi M, Takeuchi K, Hastings KE, Satoh N. Genomic overview of mRNA 5'-leader trans-splicing in the ascidian Ciona intestinalis. Nucleic Acids Research. 2006;34:3378–3388. doi: 10.1093/nar/gkl418. [DOI] [PMC free article] [PubMed] [Google Scholar]
  97. Schep AN, Wu B, Buenrostro JD, Greenleaf WJ. chromVAR: inferring transcription-factor-associated accessibility from single-cell epigenomic data. Nature Methods. 2017;14:975–978. doi: 10.1038/nmeth.4401. [DOI] [PMC free article] [PubMed] [Google Scholar]
  98. Schultheiss TM, Burch JB, Lassar AB. A role for bone morphogenetic proteins in the induction of cardiac myogenesis. Genes & Development. 1997;11:451–462. doi: 10.1101/gad.11.4.451. [DOI] [PubMed] [Google Scholar]
  99. Schultz DC, Ayyanathan K, Negorev D, Maul GG, Rauscher FJ. SETDB1: a novel KAP-1-associated histone H3, lysine 9-specific methyltransferase that contributes to HP1-mediated silencing of euchromatic genes by KRAB zinc-finger proteins. Genes & Development. 2002;16:919–932. doi: 10.1101/gad.973302. [DOI] [PMC free article] [PubMed] [Google Scholar]
  100. Scimone ML, Wurtzel O, Malecek K, Fincher CT, Oderberg IM, Kravarik KM, Reddien PW. foxF-1 controls specification of Non-body wall muscle and phagocytic cells in planarians. Current Biology. 2018;28:3787–3801. doi: 10.1016/j.cub.2018.10.030. [DOI] [PMC free article] [PubMed] [Google Scholar]
  101. Sergushichev A. An algorithm for fast preranked gene set enrichment analysis using cumulative statistic calculation. bioRxiv. 2016 doi: 10.1101/060012. [DOI]
  102. Sierro N, Kusakabe T, Park KJ, Yamashita R, Kinoshita K, Nakai K. DBTGR: a database of tunicate promoters and their regulatory elements. Nucleic Acids Research. 2006;34:D552–D555. doi: 10.1093/nar/gkj064. [DOI] [PMC free article] [PubMed] [Google Scholar]
  103. Sorrentino RP, Gajewski KM, Schulz RA. GATA factors in Drosophila heart and blood cell development. Seminars in Cell & Developmental Biology. 2005;16:107–116. doi: 10.1016/j.semcdb.2004.10.005. [DOI] [PubMed] [Google Scholar]
  104. Sripathy SP, Stevens J, Schultz DC. The KAP1 corepressor functions to coordinate the assembly of de novo HP1-demarcated microenvironments of heterochromatin required for KRAB zinc finger protein-mediated transcriptional repression. Molecular and Cellular Biology. 2006;26:8623–8638. doi: 10.1128/MCB.00487-06. [DOI] [PMC free article] [PubMed] [Google Scholar]
  105. Staden R. Methods for calculating the probabilities of finding patterns in sequences. Bioinformatics. 1989;5:89–96. doi: 10.1093/bioinformatics/5.2.89. [DOI] [PubMed] [Google Scholar]
  106. Stolfi A, Gainous TB, Young JJ, Mori A, Levine M, Christiaen L. Early chordate origins of the vertebrate second heart field. Science. 2010;329:565–568. doi: 10.1126/science.1190181. [DOI] [PMC free article] [PubMed] [Google Scholar]
  107. Stolfi A, Gandhi S, Salek F, Christiaen L. Tissue-specific genome editing in Ciona embryos by CRISPR/Cas9. Development. 2014;141:4115–4120. doi: 10.1242/dev.114488. [DOI] [PMC free article] [PubMed] [Google Scholar]
  108. Tan G, Lenhard B. TFBSTools: an R/bioconductor package for transcription factor binding site analysis. Bioinformatics. 2016;32:1555–1556. doi: 10.1093/bioinformatics/btw024. [DOI] [PMC free article] [PubMed] [Google Scholar]
  109. Thakore PI, D'Ippolito AM, Song L, Safi A, Shivakumar NK, Kabadi AM, Reddy TE, Crawford GE, Gersbach CA. Highly specific epigenome editing by CRISPR-Cas9 repressors for silencing of distal regulatory elements. Nature Methods. 2015;12:1143–1149. doi: 10.1038/nmeth.3630. [DOI] [PMC free article] [PubMed] [Google Scholar]
  110. Toenhake CG, Fraschka SA, Vijayabaskar MS, Westhead DR, van Heeringen SJ, Bártfai R. Chromatin Accessibility-Based characterization of the gene regulatory network underlying plasmodium falciparum Blood-Stage development. Cell Host & Microbe. 2018;23:557–569. doi: 10.1016/j.chom.2018.03.007. [DOI] [PMC free article] [PubMed] [Google Scholar]
  111. Tolkin T, Christiaen L. Rewiring of an ancestral Tbx1/10-Ebf-Mrf network for pharyngeal muscle specification in distinct embryonic lineages. Development. 2016;143:3852–3862. doi: 10.1242/dev.136267. [DOI] [PMC free article] [PubMed] [Google Scholar]
  112. Vandenberghe AE, Meedel TH, Hastings KE. mRNA 5'-leader trans-splicing in the chordates. Genes & Development. 2001;15:294–303. doi: 10.1101/gad.865401. [DOI] [PMC free article] [PubMed] [Google Scholar]
  113. Wamstad JA, Alexander JM, Truty RM, Shrikumar A, Li F, Eilertson KE, Ding H, Wylie JN, Pico AR, Capra JA, Erwin G, Kattman SJ, Keller GM, Srivastava D, Levine SS, Pollard KS, Holloway AK, Boyer LA, Bruneau BG. Dynamic and coordinated epigenetic regulation of developmental transitions in the cardiac lineage. Cell. 2012;151:206–220. doi: 10.1016/j.cell.2012.07.035. [DOI] [PMC free article] [PubMed] [Google Scholar]
  114. Wang W, Razy-Krajka F, Siu E, Ketcham A, Christiaen L. NK4 antagonizes Tbx1/10 to promote cardiac versus pharyngeal muscle fate in the ascidian second heart field. PLOS Biology. 2013;11:e1001725. doi: 10.1371/journal.pbio.1001725. [DOI] [PMC free article] [PubMed] [Google Scholar]
  115. Wang W, Racioppi C, Gravez B, Christiaen L. Purification of fluorescent labeled cells from dissociated Ciona embryos. Advances in Experimental Medicine and Biology. 2018;1029:101–107. doi: 10.1007/978-981-10-7545-2_9. [DOI] [PMC free article] [PubMed] [Google Scholar]
  116. Wang W, Niu X, Stuart T, Jullian E, Mauck WM, Kelly RG, Satija R, Christiaen L. A single-cell transcriptional roadmap for cardiopharyngeal fate diversification. Nature Cell Biology. 2019;21:674–686. doi: 10.1038/s41556-019-0336-z. [DOI] [PMC free article] [PubMed] [Google Scholar]
  117. Weirauch MT, Yang A, Albu M, Cote AG, Montenegro-Montero A, Drewe P, Najafabadi HS, Lambert SA, Mann I, Cook K, Zheng H, Goity A, van Bakel H, Lozano JC, Galli M, Lewsey MG, Huang E, Mukherjee T, Chen X, Reece-Hoyes JS, Govindarajan S, Shaulsky G, Walhout AJM, Bouget FY, Ratsch G, Larrondo LF, Ecker JR, Hughes TR. Determination and inference of eukaryotic transcription factor sequence specificity. Cell. 2014;158:1431–1443. doi: 10.1016/j.cell.2014.08.009. [DOI] [PMC free article] [PubMed] [Google Scholar]
  118. Whyte WA, Orlando DA, Hnisz D, Abraham BJ, Lin CY, Kagey MH, Rahl PB, Lee TI, Young RA. Master transcription factors and mediator establish super-enhancers at key cell identity genes. Cell. 2013;153:307–319. doi: 10.1016/j.cell.2013.03.035. [DOI] [PMC free article] [PubMed] [Google Scholar]
  119. Wiechecki KA, Racioppi C. Github; 2019. https://github.com/ChristiaenLab/ATAC-seq [Google Scholar]
  120. Woznica A, Haeussler M, Starobinska E, Jemmett J, Li Y, Mount D, Davidson B. Initial deployment of the cardiogenic gene regulatory network in the basal chordate, Ciona intestinalis. Developmental Biology. 2012;368:127–139. doi: 10.1016/j.ydbio.2012.05.002. [DOI] [PMC free article] [PubMed] [Google Scholar]
  121. Yokomori R, Shimai K, Nishitsuji K, Suzuki Y, Kusakabe TG, Nakai K. Genome-wide identification and characterization of transcription start sites and promoters in the tunicate Ciona intestinalis. Genome Research. 2016;26:140–150. doi: 10.1101/gr.184648.114. [DOI] [PMC free article] [PubMed] [Google Scholar]
  122. Zaffran S, Küchler A, Lee HH, Frasch M. Biniou (FoxF), a central component in a regulatory network controlling visceral mesoderm development and midgut morphogenesis in Drosophila. Genes & Development. 2001;15:2900–2915. doi: 10.1101/gad.917101. [DOI] [PMC free article] [PubMed] [Google Scholar]
  123. Zaidi S, Choi M, Wakimoto H, Ma L, Jiang J, Overton JD, Romano-Adesman A, Bjornson RD, Breitbart RE, Brown KK, Carriero NJ, Cheung YH, Deanfield J, DePalma S, Fakhro KA, Glessner J, Hakonarson H, Italia MJ, Kaltman JR, Kaski J, Kim R, Kline JK, Lee T, Leipzig J, Lopez A, Mane SM, Mitchell LE, Newburger JW, Parfenov M, Pe'er I, Porter G, Roberts AE, Sachidanandam R, Sanders SJ, Seiden HS, State MW, Subramanian S, Tikhonova IR, Wang W, Warburton D, White PS, Williams IA, Zhao H, Seidman JG, Brueckner M, Chung WK, Gelb BD, Goldmuntz E, Seidman CE, Lifton RP. De novo mutations in histone-modifying genes in congenital heart disease. Nature. 2013;498:220–223. doi: 10.1038/nature12141. [DOI] [PMC free article] [PubMed] [Google Scholar]
  124. Zaret KS, Carroll JS. Pioneer transcription factors: establishing competence for gene expression. Genes & Development. 2011;25:2227–2241. doi: 10.1101/gad.176826.111. [DOI] [PMC free article] [PubMed] [Google Scholar]
  125. Zeitlinger J, Zinzen RP, Stark A, Kellis M, Zhang H, Young RA, Levine M. Whole-genome ChIP-chip analysis of dorsal, twist, and snail suggests integration of diverse patterning processes in the Drosophila embryo. Genes & Development. 2007;21:385–390. doi: 10.1101/gad.1509607. [DOI] [PMC free article] [PubMed] [Google Scholar]
  126. Zhang Y, Liu T, Meyer CA, Eeckhoute J, Johnson DS, Bernstein BE, Nusbaum C, Myers RM, Brown M, Li W, Liu XS. Model-based analysis of ChIP-Seq (MACS) Genome Biology. 2008;9:R137. doi: 10.1186/gb-2008-9-9-r137. [DOI] [PMC free article] [PubMed] [Google Scholar]
  127. Zhang P, Hung L-H, Lloyd W, Yeung KY. Hot-starting software containers for STAR aligner. GigaScience. 2018;7:giy092. doi: 10.1093/gigascience/giy092. [DOI] [PMC free article] [PubMed] [Google Scholar]
  128. Zhao R, Watt AJ, Battle MA, Li J, Bondow BJ, Duncan SA. Loss of both GATA4 and GATA6 blocks cardiac myocyte differentiation and results in acardia in mice. Developmental Biology. 2008;317:614–619. doi: 10.1016/j.ydbio.2008.03.013. [DOI] [PMC free article] [PubMed] [Google Scholar]

In the interests of transparency, eLife publishes the most substantive revision requests and the accompanying author responses.

Acceptance summary:

This article provides fundamental insights regarding regulatory mechanisms that dictate precise gene expression and thus ensure proper cell fate specification. By leveraging an exceptionally well-characterized heart and pharyngeal muscle lineage gene network in an invertebrate chordate model, this study provides profound insights into principles governing reliable and accurate network deployment. The authors employ in-depth, comprehensive chromatin accessibility profiling and extensive functional manipulations to generate substantive support for a novel "combined enhancer" model. According to this model, precise regulation of crucial transcription factors is mediated by heterologous combinations of regulatory elements with distinct inputs and chromatin accessibility profiles. Strikingly, these heterologous elements function synergistically as a regulatory unit to generate a single, precise spatio-temporal expression pattern. The authors convincingly demonstrate that removal of individual elements within these combinatorial units or deployment of homologous combinations disrupt expression and can lead to improper specification. Should future studies indicate that such combinatorial regulatory elements are broadly deployed, this would have wide-ranging implications regarding gene network architecture as well as the potential impact of non-coding mutations in the context of evolution and disease.

Decision letter after peer review:

Thank you for submitting your article "Combinatorial chromatin dynamics foster accurate cardiopharyngeal fate choices" for consideration by eLife. Your article has been reviewed by three peer reviewers, and the evaluation has been overseen by a Reviewing Editor and Didier Stainier as the Senior Editor. The reviewers have opted to remain anonymous.

The reviewers have discussed the reviews with one another and the Reviewing Editor has drafted this decision to help you prepare a revised submission.

There was a strong consensus among the reviewers that this paper contains findings that are of broad general interest and are appropriate for publication in eLife once revisions have been made to more clearly align their data with their conclusions. In particular, the reviewers were impressed by the broad scope of their ATAC-seq analysis, the innovative application of this technique to Ciona, some of the more in-depth findings regarding specific regulatory elements and some of the insights provided by this analysis in regards to more general principles of regulatory logic.

However, the reviewers had a number of critical concerns including a lack of clarity in some key definitions and figures, some questions regarding interpretation of specific assays and experiments and a lack of sufficient support for some of their major conclusions including the proposed combined enhancer model.

The was also a consensus that these concerns can be addressed in two ways. The authors can either modify their conclusions to better match their current data or they can conduct further experiments to more robustly support their current conclusions.

Reviewer #1:

In this study, Racioppi, Wiechecki and Christiaen have deployed an impressive series of genome-wide analyses complemented with CRISPR/Cas9-mediated genome editing. The study addresses chromatin accessibility profiles during early development of cardiopharyngeal lineages in Ciona. In this lineage, three rounds of asymmetric cell divisions, which are controlled by differential activation of FGF-MAPK signals in a reiterative manner, lead to generation of four cell types, first heart precursor (FHP), second heart precursor (SHP), atrial siphon muscle precursor (ASMP) and anterior tail muscle (ATM). The Christiaen laboratory is one of the major players in deciphering the underlying transcriptional and signalling control of the cardiopharyngeal lineage segregation. The current study reveals another layer of such regulatory mechanisms, in particular, chromatin accessibility, by deploying ATAC-seq on isolated cardiopharyngeal lineage cells. The study contains a huge amount of bioinformatic analyses, which, to be honest, are very hard for me to digest. Nonetheless, I really appreciate the last part of their study highlighted in Figures 5, Figure 5—figure supplements 1 and 2, which focuses on chromatin accessibility profiles of two key genes of cardiopharyngeal fate choices, Ebf and Tbx1/10. It reveals convincingly that these genes are controlled by multiple _cis_-regulatory elements with distinct temporal chromatin accessibility profiles and distinct TF-binding motif compositions. This finding indeed represents the main message of the current study.

I have a few comments (questions):

  1. For Figure 5A and Figure 5—figure supplement 2A, it would be informative to include the ATAC-seq dataset of FoxF(CRISPR)-10 hpf.

  2. I assume that changes in chromatin accessibility between FGFR(DN)-10 hpf and FoxF(CRISPR)-10 hpf would exhibit a certain degree of correlation since FoxF is transcriptionally activated by FGF signals in TVCs. It would be nice if the authors could conduct a correlation test. There seems to be such a correlation for the peak ID KhL20.169 within Gata4/5/6 locus (Figure 3—figure supplement 2A and Figure 3—figure supplement 4D).

  3. Why do the authors call MyoD905>GFP cells as mesenchymal cells? Aren't they muscle cells (Christiaen et al., 2008)?

  4. The authors mention "cardiopharyngeal markers" that include TVC progenitor markers, primed cardiac markers, primed ASM markers, de novo cardiac markers, de novo ASM markers and ATM markers. The definition for these different markers is supposed to be described in supplementary text. I imagine that the authors are meaning that in the "Gene Set Enrichment Analysis (GSEA)" section but it doesn't help me to understand the definition. A clearer definition, together with gene lists for the different categories of cardiopharyngeal markers as a supplementary table, would be helpful.

Reviewer #2:

This manuscript by Racioppi et al. performs the first ATAC-seq experiments in the early heart cell lineages of Ciona in order to find enhancers that are accessible in sublineages of the heart. This paper creates a very valuable resource for the heart and ciona communities. I don't feel that the authors make the most of their data and the value of this work is skipped over to focus on the combined enhancer. Overall I find there are too many interpretations/extrapolations without sufficient evidence and these detract from the excellent quality and substantial value of this dataset for understanding heart development. I would recommend that the authors focus on their data and what it tells us about the regulatory networks and developmental programs of the heart rather than the concept of combined enhancers.

Combined enhancers:

A major claim of the paper, that "combined enhancers" foster spatially and temporally accurate fate choices, by increasing the repertoire of regulatory inputs that control gene expression, through either accessibility and/or activity" is only supported by a study of only one locus shown in Figure 6, and I have issues with how the experiment is designed. The authors do an experiment in Figure 6 where they multimerize a weak enhancer (1x, 2x, 3x) and find that when 3 copies of the weak enhancer are put next to each other in a certain way, they find ectopic expression by reporter assay. They then say that this, coupled with evidence that many weak enhancers in the genomic context do not cause ectopic expression, show that these "combined enhancers" foster spatially and temporally accurate fate choices. However, this is comparing apples (3 identical enhancers with non-endogenous spacing together on a plasmid reporter) to oranges (3 unique endogenous enhancers with endogenous spacing in the genome). Since spacing between enhancers and genomic context can affect transcriptional output, I would believe this more if they put the 3 unique weak enhancers from the genomic context into a synthetic construct with spacings like the 3x multimerized enhancer and showed that these constructs did not yield ectopic expression. Even if this experiment were to show that the different inputs led to specific expression, I don't understand how this type of enhancer is different from enhancers using different inputs turning on at different times and space during a developmental program.

Definition of Shadow enhancer: The authors state that "shadow enhancer promotes robust transcription through the actions of multiple elements mediating similar regulatory inputs". Based on the literature I wasn't able to find a consensus that shadow enhancers have to use similar inputs, indeed it is very hard to pin down exactly what the inputs for many shadow enhancers are. Perry, 2011/Hong, 2008.

Reviewer #3:

In this manuscript, Racioppi et al. investigate the temporal dynamics of chromatin accessibility during the different steps of cardiopharyngeal progenitor specification in larvae of the chordate Ciona. They sample different cell types at different time points post fertilization, ranging from Mesp1+ founder cells after 6 hours to first and second heart progenitors, as well as atrial syphon precursors and anterior tail muscle cells after 18 hours. They show that global accessibility decreases with time, and that most enhancers are opened early during the process, in multipotent progenitors, even though their associated genes are mostly activated later. Using CRISPR-mediated deletions and reporter assays, they demonstrate that the enhancers they identify are necessary and sufficient to activate gene expression with temporal and spatial specificity. They then identify putative regulators of stage- and cell-type specific enhancer accessibility and confirm Foxf as a driver of enhancer accessibility in the cardiac lineage. Then, they study the locus of Ebf, where multiple enhancers are necessary for its expression in atrial muscle precursors. Finally, they study combinatorial effects of Ebf enhancers using a reporter assay and find that multiple copies of the same enhancer increase transcription efficiency, at the cost of precocious activation, whereas one copy of a combined enhancer construct recapitulates gene expression with greater temporal fidelity.

In this well-conducted study, the authors present solid data to back their conclusions, which are interesting for the field of development. Some points could be addressed to strengthen their claims:

- Does their reporter assay recapitulate enhancer accessibility or activity? This is unclear because the loss of the genomic context may cause an enhancer to be more easily activated. What is the timing of activation of the reporter gene in regard to the timing of chromatin accessibility and gene transcription activation?

- For the enhancers they define as being dependent on Foxf, does the knockout of Foxf impede their activation in a reporter assay? Or are enhancers containing mutated Foxf sites inactive (endogenously or in a reporter essay)?

- The authors could use histone marks to assess if the decoupling between accessibility and gene activation is also seen between gene activation and accumulation of histone marks typical of active genetic elements, such as histone 3 lysine 27 acetylation, if it is technically possible for them with their low number of cells.

- Figure 4A-C and the corresponding text are unclear. It is difficult to understand how the later comparison of FGF-MAPK perturbation relates to the rest of the paper.

Reviewer #1:

[…] I have a few comments (questions):

  1. For Figure 5A and Figure 5—figure supplement 2A, it would be informative to include the ATAC-seq dataset of FoxF(CRISPR)-10 hpf.

We have included FoxfCRISPR ATAC-seq data at 10 hpf in Figure 5B, and in Figure 5—figure supplement 2A.

  1. I assume that changes in chromatin accessibility between FGFR(DN)-10 hpf and FoxF(CRISPR)-10 hpf would exhibit a certain degree of correlation since FoxF is transcriptionally activated by FGF signals in TVCs. It would be nice if the authors could conduct a correlation test. There seems to be such a correlation for the peak ID KhL20.169 within Gata4/5/6 locus (Figure 3—figure supplement 2A and Figure 3—figure supplement 4D).

This is correct. We observed several peaks showing similar change in chromatin accessibility following FGF perturbation (FgfrDN) or tissue-specific Foxf loss-of-function (FoxfCRISPR), like peaks KhC14.806,.807 and.808 associated to Hand in Figure 3E, peaks KhC4.137 and.144 in the Lrp4/8 locus in Figure 4D, and peaks KhC1.479 in the Eph1 locus, KhC9.3113 in the Ddr1/2 locus, KhC4.606 in the Smurf1/2 locus, KhC6.1052 and KhC6.1046 in the Fzd4 locus shown in Figure 3—figure supplement 4.

Consistent with this observation, we found that chromatin accessibility changes induced by Mesp>FgfrDN and FoxfCRISPR compared to controlsat 10 hpf are correlated, with a Spearman correlation coefficient (𝜌) of 0.49. We added this additional information in Figure 3—figure supplement 2F as a scatter plot to compare the correlation of differentially accessible regions in “Mesp>FgfrDN and “FoxfCRISPR at 10 hpf” ATAC-seq samples.

  1. Why do the authors call MyoD905>GFP cells as mesenchymal cells? Aren't they muscle cells (Christiaen et al., 2008)?

The _cis_-regulatory DNA from the myogenic differentiation factor MyoD (termed Myod905) was used to express green fluorescent protein in the mesenchymal cells (Christiaen at al., 2008). This is because the enhancer does not contain all control elements of Myod/Mrf and “leaks” strongly in the mesenchyme. This enhancer is active in the tail muscle cells as well. Nevertheless, mesenchymal cells are smaller and more numerous than tail muscle cells, they are thus more efficiently isolated by FACS than tail muscle cells. There is a possibility that muscle cells are sorted along with mesenchymal cells but they will be underrepresented in the FACSorted cell population. In Christiaen et al., 2008, the Myod905>YFPVenus+ population expressed higher levels of mesenchyme markers.

In Christiaen et al., 2008 (Supplementary Figure 1), we described a similar trend for the two sister cells, TVC and ATM. TVC cells are smaller than ATMs and they are more efficiently sorted by FACS rather than their sister cells ATMs (Christiaen et al., 2008, Supporting material section 3.2).

In order to be more accurate in defining ATAC-seq samples isolated from mesenchymal cells at different time points, we replaced “Mesenchyme” ATAC-seq samples with “MyoD905>GFP”.

4) The authors mention "cardiopharyngeal markers" that include TVC progenitor markers, primed cardiac markers, primed ASM markers, de novo cardiac markers, de novo ASM markers and ATM markers. The definition for these different markers is supposed to be described in Supplementary text. I imagine that the authors are meaning that in the "Gene Set Enrichment Analysis (GSEA)" section but it doesn't help me to understand the definition. A clearer definition, together with gene lists for the different categories of cardiopharyngeal markers as a supplementary table, would be helpful.

In the Materials and methods under the section “Definition of gene sets”, we clarified the definition of gene sets used for "Gene Set Enrichment Analysis (GSEA)" by listing the cell-type-specific gene sets including primed and de novo-expressed Cardiac and ASM/pharyngeal muscle, as defined from scRNA-seq by Wang et al., 2019, as well as anterior tail muscle (ATM) gene sets.

More specifically, we replaced:

“Gene Set Enrichment Analysis (GSEA)

We define gene sets as described in the supplementary text. […] We measured enrichment for peaks associated to the gene categories defined by both bulk RNA-seq and scRNA-seq (Wang et al., 2019) in open peaks (indicated by a positive enrichment score, ES) or closed peaks (indicated by a negative ES) in each comparison.”

With:

“Definition of gene sets

We used cardiopharyngeal lineage cell-type-specific gene sets (cardiopharyngeal markers), including primed and de novo-expressed Cardiac and ASM/pharyngeal muscle as well as secondary TVC (STVC), first heart precursor (FHP), second heart precursor (SHP) markers defined by scRNA-seq (Wang et al., 2019). […] To these we added the gene sets either activated or inhibited by MAPK at 10 or 18 hpf as defined by microarray (Christiaen et al., 2008) and Bulk RNA-seq analysis (Wang et al., 2019).”

“Gene Set Enrichment Analysis (GSEA)

We performed GSEA on elements associated to these gene sets using fgsea 1.2.1 (Sergushichev, 2016) with parameters: minSize=5, maxSize=Inf, nperm=10000. Because the gene sets used for GSEA need not be exclusive, we can use ATAC-seq peaks as “genes” and define a gene set as all peaks annotated to genes in the set. We measured enrichment for peaks associated to the gene sets in open peaks (indicated by a positive enrichment score, ES) or closed peaks (indicated by a negative ES) in each comparison.”

Furthermore, we added “We performed bulk RNA-seq following defined perturbations of FGF-MAPK signaling and FACS (Wang et al., 2019).” as the first sentence under Bulk RNA-seq Library Preparation, Sequencing and Analyses.

We deleted:

“To run GSEA, we used gene sets including primed and de novo-expressed Cardiac and ASM/pharyngeal muscle, as well as secondary TVC (STVC), first heart precursor (FHP), second heart precursor (SHP) markers defined by scRNA-seq(Wang et al., 2019). Anterior tail muscles (ATM) marker genes were derived from publicly available expression data (Brozovic et al., 2018) (Tables S46). The gene sets either activated or inhibited by MAPK at 10 or 18 hpf were derived from bulk RNA-seq following defined perturbations of FGF-MAPK signaling and FACS (Wang et al., 2019).”

Finally, we replaced “Supplementary text” with “Materials and methods”.

Reviewer #2:

This manuscript by Racioppi et al. performs the first ATAC-seq experiments in the early heart cell lineages of Ciona in order to find enhancers that are accessible in sublineages of the heart. This paper creates a very valuable resource for the heart and ciona communities. I don't feel that the authors make the most of their data and the value of this work is skipped over to focus on the combined enhancer.

We appreciate this comment and indeed tried to not underplay the findings that pertain to the general description of accessibility patterns, and the roles of FGF-MAPK signaling and Foxf in establishing these patterns, as well as defining novel enhancers. Figures 1 to 4 are focused on this, and hopefully do it justice. Arguably, this dataset contains more information to be mined and followed-up on.

Overall I find there are too many interpretations/extrapolations without sufficient evidence and these detract from the excellent quality and substantial value of this dataset for understanding heart development. I would recommend that the authors focus on their data and what it tells us about the regulatory networks and developmental programs of the heart rather than the concept of combined enhancers.

Our initial rationale was that the general concepts of decoupling early accessibility from late enhancer activity, of a forkhead factor acting as a possible pioneer had become classic concepts in the past few years. We also wanted to circle back to the regulation of actual gene expression, especially heart vs. pharyngeal muscle-specific gene expression, and do justice to the novel observations that genes like Ebf and Tbx1/10 were flanked by elements that showed different dynamics of accessibility. This seemed to offer an opportunity for more novel and yet broadly relevant findings, which we report as “combined enhancers”.

We added new experimental data (Figures 6 and 7), which we think strengthen our conclusions: indeed, we could experimentally decouple the level of expression of Ebf from its spatio-temporal accuracy and thus cause ectopic activation of the endogenous locus and fate transformation of the second heart lineage into pharyngeal muscles. These new results provide more direct experimental support for the claim that “combined enhancers foster accurate fate choices” to paraphrase the title.

We considered the option of splitting the paper, but practical constraints and the rationale explained above led us to keep the paper largely as is, while strengthening conclusions.

Combined enhancers:

A major claim of the paper, that "combined enhancers" foster spatially and temporally accurate fate choices, by increasing the repertoire of regulatory inputs that control gene expression, through either accessibility and/or activity" is only supported by a study of only one locus shown in Figure 6.

We reported a similar pattern of variable dynamics of accessibility for enhancer elements for both Tbx1/10 and Ebf, but could only perform the more extensive analysis on Ebf indeed. This part of the paper stands as a proof of a rather novel concept, and it is conceivable that future studies will be needed to probe its generality. We suspect that it will matter particularly for developmental control genes like Ebf, which mediate commitment and rely on autoregulation for long-term maintenance.

I have issues with how the experiment is designed. The authors do an experiment in Figure 6 where they multimerize a weak enhancer (1x, 2x, 3x) and find that when 3 copies of the weak enhancer are put next to each other in a certain way, they find ectopic expression by reporter assay. They then say that this, coupled with evidence that many weak enhancers in the genomic context do not cause ectopic expression, show that these "combined enhancers" foster spatially and temporally accurate fate choices. However, this is comparing apples (3 identical enhancers with non-endogenous spacing together on a plasmid reporter) to oranges (3 unique endogenous enhancers with endogenous spacing in the genome). Since spacing between enhancers and genomic context can affect transcriptional output, I would believe this more if they put the 3 unique weak enhancers from the genomic context into a synthetic construct with spacings like the 3x multimerized enhancer and showed that these constructs did not yield ectopic expression.

We extensively addressed this point experimentally and basically showed that the concatemer of distinct elements without endogenous spacer sequences behave like the full length upstream sequence, in that it drives accurate expression in the pharyngeal muscles, whereas the multimer of one distal weak enhancer causes ectopic activation, and ultimately fate transformation.

We added the following paragraph in the Results section: “To test whether spacing between accessible elements could affect transcriptional output, we built a concatemer of KhL24.37,.36,.35, and.34 elements without endogenous spacer sequences. This construct increased the proportion of embryos with ASM cells expressing the reporter to 92% ± 2% (SE, n=130; Figure 6C), but it did not induce ectopic expression in the second heart lineage”, and: “neither one copy of KhL24.37 (1x KhL24.37), nor the full combined enhancers, with or without endogenous spacers, sufficed to cause substantial fate transformation of _Tbx1/10_+ second heart lineage into migratory pharyngeal muscle precursors (Figure 7C-D)”.

Even if this experiment were to show that the different inputs led to specific expression, I don't understand how this type of enhancer is different from enhancers using different inputs turning on at different times and space during a developmental program.

In this case, the key difference is that these elements contribute to the onset of Ebf expression for ~ 1 hour in the pharyngeal muscle precursors (as we showed in Razy-Krajka et al., 2018). In other words, their actions are at the exact same time and place in development. This is thus different from classic modular _cis_-regulatory system paradigms such as Drosophila Even-skipped, and many other examples.

Definition of Shadow enhancer: The authors state that " shadow enhancer promotes robust transcription through the actions of multiple elements mediating similar regulatory inputs" Based on the literature I wasn't able to find a consensus that shadow enhancers have to use similar inputs, indeed it is very hard to pin down exactly what the inputs for many shadow enhancers are. Perry, 2011/Hong, 2008.

Three main features characterized shadow enhancers: (1) they each drive a pattern of transcription resembling that of a previously identified “primary” enhancer that is more proximal to the promoter being regulated; (2) they bind the same TFs as the primary enhancer, suggesting a similar regulatory logic; and (3) they are located either within an intron of, or on the far side of a neighboring gene (Barolo et al., 2012).

Based on the literature, Mike Levine’s lab was the first to coin the term “Shadow Enhancer” for “remote secondary enhancers mapping far from the target gene and mediating activities overlapping the primary enhancer”. This paper is cited in the Introduction (Hong et al., 2008).

Moreover, we cited the “Zeitlinger et al., 2007” paper where Levine’s group described the identification of shadow enhancers by chromatin immunoprecipitation (ChIP)-chip assays.

The corresponding author was still a postdoc in the Levine lab when these studies were conducted and is thus familiar with the issue. Basically the concept of shadow enhancer focuses primarily on the transcriptional output, and the fact that shadow enhancers integrate the same inputs is generally implied, and was an important part of their discovery.

More recently, studies are beginning to uncover that shadow enhancers integrate distinct inputs. This is essentially what we are saying for combined enhancers, except that we show differences in the dynamics of chromatin accessibility, which is particularly innovative here.

Reviewer #3:

[…] In this well-conducted study, the authors present solid data to back their conclusions, which are interesting for the field of development. Some points could be addressed to strengthen their claims:

- Does their reporter assay recapitulate enhancer accessibility or activity? This is unclear because the loss of the genomic context may cause an enhancer to be more easily activated.

For several elements analyzed in this study, we observed that reporter gene assays recapitulate enhancer activity and not accessibility. We added additional information to support this point on Figure 4—figure supplement 2.

In the third paragraph of the subsection “Chromatin accessibility in late heart vs. pharyngeal muscle precursors”, we described that very few de novo-expressed pan-cardiac genes (30) were also differentially expressed upon FGF-MAPK perturbation at 18 hpf and that only 8 out of these 30 genes were associated with differentially accessible element following perturbation of FGF-MAPK signaling (Figure 4A and B). In this group of 8 elements, we identified the KhC5.1641 peak that is hosted in the 5’UTR of the de novo expressed cardiac marker, Mmp21, gene locus (Figure 4—figure supplement 2C).

We experimentally tested a ~3kb region upstream the coding ATG of Mmp21 encompassing the KhC5.1641 peak by gene reporter assay (Figure 4—figure supplement 2C). Even if Mmp21 is expressed only in the first heart precursors (Wang et al., 2009), the ~3kb region was able to activate GFP in both cardiac (first and second heart precursors) and pharyngeal muscle cell precursors in embryos at stage 25 embryos (18hpf) (Figure 4—figure supplement 2D and E).

Similarly, we identified an element (peak ID: KhC2.3468) upstream the coding region of the de novo expressed pharyngeal muscle marker, Tmct2. This peak was exclusively accessible in the pharyngeal cell upon perturbation of FGF-MAPK signaling (Figure 4—figure supplement 2F). However, by reporter assay, it was able to activate GFP expression in the whole B7.5 cell lineage, including ATMs and cardiac and pharyngeal cell muscles (Figure 4—figure supplement 2G).

We additionally add another example regarding two accessible elements upstream the coding region of the primed-expressed cardiac marker, KH.C1.1093_ZAN that were able to activate GFP in both cardiac and pharyngeal cell muscles (Figure 4—figure supplement 2H, I).

We added the following paragraph in the Results section:

“Similarly, reporter gene expression assays showed that DNA fragments containing differentially accessible elements located ~0.5 kb upstream of the coding region of the de novo-expressed gene Tmtc2 (KhC2.3468), and upstream of KH.C1.1093_ZAN (KhC3.47, KhC3.46), were sufficient to drive GFP expression in both cardiac and pharyngeal muscle progenitors, consistent with the notion that electroporated plasmids are not “chromatinized” and thus constitutively accessible (Figure 4—figure supplement 2F-G; Figure 4—figure supplement 3C-D). This suggested that, for genes like Mmp21, Tmtc2 and Zan, cell-type-specific accessibility determines cardiac vs. pharyngeal muscle-specific gene expression.”

What is the timing of activation of the reporter gene in regard to the timing of chromatin accessibility and gene transcription activation?

We observed different timing of the reporter gene activation in regard to the timing of chromatin accessibility and gene transcription. We often observed a decoupling between chromatin accessibility and chromatin activity for de novo expressed genes. For example, the previously characterized STVC-specific enhancer located upstream the starting coding of Tbx1/10 (Razy-Krajka et al., 2018) was already accessible in the naive Mesp+ mesoderm at 6 hpf preceding the gene transcription which starts at 15hpf (confirmed by in situ and RNA-seq data) (Figure 5—figure supplement 2A). We describe this in the Results section: “Tbx1/10 showed a similar logic, whereby a constitutively accessible upstream element (KhC7.909) acts as an enhancer of cardiopharyngeal expression (Razy-Krajka et al., 2018)”. Similarly, the TVC-specific enhancer of Foxf (Beh et al., 2007) is accessible in founder cells at 6 hpf anticipating the gene expression which starts immediately after the founder cells divide in two sister cells, TVC and ATM, at 8 hpf (Figure 3—figure supplement 2A). We describe this in the Results section: “Finally, the Foxf enhancer was accessible in naive Mesp+ founder cells, suggesting that it is poised for activation, unlike the intronic Gata4/5/6 enhancer (Figure 3—figure supplement 2A).”

The pharyngeal muscle determinant Ebf gene locus is surrounded by multiple regions showing temporal dynamics in chromatin accessibility during cardiopharyngeal development (Figure 5A). However, the minimal enhancer driving Ebf in the pharyngeal muscles (KhL24.37 peak) opens only at 15hpf, when the gene turns on in the ASMPs. We describe this in the Results section: “Both loci contained multiple accessible regions, including elements open already in the naive _Mesp_+ mesoderm (e.g. Ebf_KhL24.35/36), cardiopharyngeal-lineage-specific elements that open prior to gene activation but after induction of multipotent progenitors (e.g. Ebf_KhL24.34), and elements that open de novo in fate-restricted pharyngeal muscle precursors, where the gene is activated (e.g. Ebf_KhL24.37) (Figure 5A, B). Previous reporter gene expression assays identified the latter element, Ebf_KhL24.37, as a weak minimal enhancer with pharyngeal muscle-specific activity (Wang et al., 2013).”.

We also observed some elements opening simultaneously with the gene transcription activation. For example, TVC-specific enhancers controlling two primed-expressed genes as GATA4/5/6 (Woznica et al., 2012) (Figure 3—figure supplement 4D) as well as Hand (Figure 3E) are accessible at 10 hpf when the genes are transcribed in the multipotent progenitors.

- For the enhancers they define as being dependent on Foxf, does the knockout of Foxf impede their activation in a reporter assay? Or are enhancers containing mutated Foxf sites inactive (endogenously or in a reporter essay)?

We added new data (Figure 5—figure supplement 2) to support that loss-of-function of Foxf in the cardiopharyngeal lineage reduce the enhancer elements that we have characterized in this work.

We identified two conserved Forkhead binding sites in the minimal STVC-specific enhancer controlling Tbx1/10 activation (termed T12) and induced point mutation in this binding site (“mFoxf.1:TAAACA in TAACCA” and “mFoxf.2:GTTTA in GTGTA”; we listed the oligos used for point mutation in Supplementary file 5).

We tested the enhancer activity by reporter assay using T12 encompassing the two point mutations in the Forkhead binding sites and T12 wt. We found that the enhancer activity is drastically reduced in ~92% of larvae expressing T12-mFoxf.2 and ~88% of larvae expressing T12-mFoxf.1 compared to the wt T12 element (Figure 5—figure supplement 2D, E).

Moreover, in embryos where we induced loss-of-function of Foxf (FoxfCRISPR), we observed T12 enhancer activity drastically reduced compared to the CRISPR control embryos. ~84% of embryos were not expressing GFP in STVC progeny in embryos at stage 26 compared to the control embryos where ~28% of embryos were expressing GFP in both ASMPs and SHPs, accordingly to the normal activity of T12 (Figure 5——figure supplement 2F).

We added the following paragraph in the Results section: “We identified two conserved putative Forkhead binding sites in the minimal STVC-specific enhancer from the Tbx1/10 locus (termed T12, (Razy-Krajka et al., 2018), which were necessary for reporter gene expression (Figure 5—figure supplement 2D, E). Moreover, loss-of-function of Foxf (FoxfCRISPR) drastically reduced T12 enhancer activity (Figure 5—figure supplement 2). These results are consistent with the hypothesis that Foxf acts directly on the minimal Tbx1/10 enhancer to promote its activity in the second multipotent cardiopharyngeal progenitors.”

- The authors could use histone marks to assess if the decoupling between accessibility and gene activation is also seen between gene activation and accumulation of histone marks typical of active genetic elements, such as histone 3 lysine 27 acetylation, if it is technically possible for them with their low number of cells.

This would be fantastic, and our lab has been testing the recently developed Cut&Tag approach, but this is not established yet, and would extend far beyond the scope of this paper, since we carefully avoided making unsubstantiated claims about histone modifications.

- Figure 4A-C and the corresponding text are unclear. It is difficult to understand how the later comparison of FGF-MAPK perturbation relates to the rest of the paper.

In a nutshell, we observed very few differences in accessibility between fate-restricted heart and pharyngeal muscle precursors, which contrasts with differences in gene expression, and is consistent with the idea that pan-cardiopharyngeal patterns of accessibility are established early, upon induction of progenitors, and that differential activity of otherwise accessible elements is a key driver of differential gene expression, with the noted exceptions.

We added more information regarding the pan-cardiac and pharyngeal muscle genes that are differentially expressed upon FGF signaling perturbation.

We added the following paragraph in the Results section: “Similarly, out of 23 de novo-expressed pharyngeal muscle genes that were also differentially expressed upon FGF-MAPK perturbation at 18 hpf, 11 were associated with one differentially accessible element following perturbation of FGF-MAPK signaling (Figure 4—figure supplement 3A).”

These sub-set of genes are surrounded by few chromatin accessible regions showing differential accessibility in response to FGF signaling molecular perturbations. Even if few peaks were significantly changing accessibility, we identified DNA fragments located ~0.5 Kb upstream the coding region of the de novo-expressed Tmtc2 (KhC2.3468) and upstream of the primed-expressed KH.C1.1093_ZAN (KhC3.47, KhC3.46) were sufficient to drive GFP expression in both cardiac and pharyngeal muscle progenitors (Figure 4—figure supplement 2F-G; Figure 4—figure supplement 3C-D). These data along with the validation of Lrp4/8 associated peaks in Figure 4D-G, reinforced the concept of a decoupling between enhancer activity and accessibility.

Data Citations

  1. Racioppi C, Wiechecki K, Christiaen L. 2019. Genome-wide maps of chromatin accessibility governing the transition from founder cells to distinct fate-restricted cardiopharyngeal precursors. NCBI Gene Expression Omnibus. GSE126691
  2. Wang W, Niu X, Stuart T, Jullian E, Mauck WM, Kelly RG, Satija R, Christiaen L. 2019. A single cell transcriptional roadmap for cardiopharyngeal fate diversification. NCBI Gene Expression Omnibus. GSE99846 [DOI] [PMC free article] [PubMed]
  3. Christiaen L, Davidson B, Kawashima T, Powell W, Nolla H, Vranizan K, Levine M. 2008. Transcription profiling of B7.5 lineage cells in the tunicate Ciona intestinalis to study cardiac cell migration. ArrayExpress. E-MEXP-1478
  4. Razy-Krajka F, Lam K, Wang W, Stolfi A, Joly M, Bonneau R, Christiaen L. 2014. Collier/OLF/EBF-dependent transcriptional dynamics control muscle specification from cardiopharyngeal progenitors. NCBI Gene Expression Omnibus. GSE54746 [DOI] [PMC free article] [PubMed]

Supplementary Materials

Supplementary file 1. Description of ATAC-seq dataset, with sample name, conditions, type of cells, number of mapped reads for all the biological replicates.

Supplementary file 2. Differentially FGF-MAPK-regulated genes (microarray) associated with differentially accessible peaks (ATAC-seq) in the pairwise comparison _Mesp_>FgfrDN vs. _Mesp_>LacZ – 10 hpf; genomic coordinates of each peak are reported.

Supplementary file 3. Known enhancers for TVC-specific genes characterized in previous studies; genomic coordinates of each peak are reported.

Supplementary file 4. Primer sequences used for cloning putative regulatory regions for functional enhancer assay.

Supplementary file 5. Primer sequences of single guide RNAs used to induce CRISPR/Cas9-mediated deletions in ATAC-seq peak as well as in Foxf (FoxfCRISPR) and to induce point mutations in Tbx1/10 STVC-specific enhancer; genomic coordinates and peak IDs are reported for each peak validated.

Supplementary file 6. ATM gene selected based on expression data deposited in ANISEED (https://www.aniseed.cnrs.fr).

Supplementary file 7. Differentially expressed genes from RNA-seq on FACS-purified cells following CRISPR/Cas9-induced loss of function of Foxf.

Supplementary file 8. Differentially accessible peaks from ATAC-seq on FACS-purified cells following CRISPR/Cas9-induced loss of function of Foxf (FoxfCRISPR vs. ControlCRISPR at 10 hpf).

Sheet1: closed in FoxfCRISPR log2(FC) < −0.5 and FDR < 0.05); sheet2: open in _FoxfCRISPR_ log2(FC) > 0.5 and FDR < 0.05.

Supplementary file 9. Differentially accessible peaks from ATAC-seq, closed in the founder cells (LacZ 6 hpf < LacZ 10 hpf) and closed in FoxfCRISPR vs. _ControlCRISPR_at 10 hpf.

Supplementary file 10. Differentially accessible peaks closed in FgfrDN at 10 hpf and associated with de novo cardiac (sheet 1) and pharyngeal muscle (sheet 2) expressed genes.

Supplementary file 11. Differentially FGF-MAPK-regulated genes associated with differentially accessible peaks in MAPK activated −18 hpf (sheet1) and MAPK inhibited −18 hpf (sheet 2); genomic coordinates of each peak are reported.

Supplementary file 12. FGF-MAPK-regulated genes associated with ATAC-seq peaks that are non-differentially accessible in MAPK activated −18 hpf and MAPK inhibited −18 hpf.

Supplementary file 13. Accessible elements associated with de novo expressed heart and pharyngeal muscle markers into pre-accessible/primed or de novo accessible elements; genomic coordinates of each peak are reported.

Supplementary file 14. FASTA sequences including the Foxf proteins (with their accession numbers) used for Figure 3—figure supplement 2.

Transparent reporting form

Data Availability Statement

All sequencing data were deposited on GEO with accession GSE126691.

All sequencing data were deposited on GEO with accession GSE126691.

The following dataset was generated:

Racioppi C, Wiechecki K, Christiaen L. 2019. Genome-wide maps of chromatin accessibility governing the transition from founder cells to distinct fate-restricted cardiopharyngeal precursors. NCBI Gene Expression Omnibus. GSE126691

The following previously published datasets were used:

Wang W, Niu X, Stuart T, Jullian E, Mauck WM, Kelly RG, Satija R, Christiaen L. 2019. A single cell transcriptional roadmap for cardiopharyngeal fate diversification. NCBI Gene Expression Omnibus. GSE99846

Christiaen L, Davidson B, Kawashima T, Powell W, Nolla H, Vranizan K, Levine M. 2008. Transcription profiling of B7.5 lineage cells in the tunicate Ciona intestinalis to study cardiac cell migration. ArrayExpress. E-MEXP-1478

Razy-Krajka F, Lam K, Wang W, Stolfi A, Joly M, Bonneau R, Christiaen L. 2014. Collier/OLF/EBF-dependent transcriptional dynamics control muscle specification from cardiopharyngeal progenitors. NCBI Gene Expression Omnibus. GSE54746