Integration of ATAC- and RNA-sequencing identifies chromatin and transcriptomic signatures in classical and non-classical zebrafish osteoblasts and indicates mechanisms of entpd5a regulation
- Institute of Cardiovascular Organogenesis and Regeneration, Faculty of Medicine, Universität Münster, Germany
- Flow Cytometry Unit, Max Planck Institute for Molecular Biomedicine, Germany
- Department of Cancer and Genomic Sciences, School of Clinical Sciences, College of Medicine and Health, University of Birmingham, United Kingdom
Figures
Figure 1
entpd5a but not osterix is a marker for a non-classical osteoblast population in zebrafish.
(A) Schematic of the embryonic notochord structure (adapted from Grotmol et al., 2005). The ‘notochord sheath’ comprises the three outer layers, specifically the basal lamina (grey), the medial layer of fibrillar collagen (light blue), and the outer layer of loosely arranged matrix (lighter blue). Alizarin red staining of mineralised bone reveals non-classical osteoblasts along the trunk that are negative for osterix (B, Bi) and positive for entpd5a (D) at 6 dpf, as well as classical osteoblasts of the head and mature vertebrae (arrowheads) that are positive for osterix (C, Ci) (21 dpf). (E) entpd5a and col9a2 label bone-producing and intersegmental NSCs, respectively. cl, cleithrum; CNS, central nervous system; da, dorsal arch; no, notochord; NSCs, notochord sheath cells; op, operculum; VCs, vacuolated cells. All lateral views, anterior to the left. (B, Bi, D, E) Scale bars: 100 μm. (C, Ci) Scale bars: 200 μm.
© 2005, John Wiley and Sons. Panel A was modfied with permission from Grotmol et al., 2005. Further reproduction of this panel/illustration/figure would need permission from the copyright holder.
Figure 2
The complex expression pattern of entpd5a is dynamically regulated during zebrafish development.
(A) Schematic of the BAC driving expression of Gal4FF in the entpd5a:Gal4FF transgenic line. The BAC contains the entire open reading frame of entpd5a and part of the open reading frames of adjacent genes, coq6 and aldh6a1. (B, C, D) GFP expression is detected in osteoblasts (arrow) and (partially) in cartilage (asterisk; D, inset) making up the head skeleton. (C) Strong GFP expression is seen in the notochord and the cleithrum, but also in a subset of CNS neurons. (E–J) Using the entpd5a:Kaede photoconversion line, we first detect entpd5a expression at the 15 somite stage (E). Following the same embryo, active expression of the gene continues until prior to 24 hpf (F–H). Between 24 hpf and until notochord segmentation takes place, entpd5a is only actively expressed in the ventral-most notochord sheath cells (NSCs) of the notochord, while expression in the remaining NSCs and vacuolated cells is turned off (I, J). (I–J, insets) Cross section of the notochord at 48 hpf and at 72 hpf in the position indicated by the white vertical line. (B) Ventral view; (C–J) lateral views, anterior towards the left. (B–J) Scale bars: 100 μm. (D, I, J) Insets scale bars: 50 μm.
Figure 3 with 1 supplement
Quality assessment of assay for transposase-accessible chromatin (ATAC)-seq of skeletogenic cells.
(A, B) Trunk and head, respectively, of transgenic fish at 15 dpf, used for collection of FAC-sorted cells. Dashed lines indicate the site where the cut was made to separate head from trunk tissue. GFP+ cells indicate mineralising cells (A, notochord sheath cells [NSCs]; B, cranial osteoblasts) while mCherry+ cells indicate (A) intersegmental NSCs and (B) head cartilage. (C) Principal component analysis (PCA), with red circle indicating mCherry+ cell samples and green circle GFP+ cell samples. (D) Heatmap indicating correlation amongst individual replicates. (E, F) Volcano plots showing for (E) head and (F) trunk the identified significantly accessible chromatin regions. Positive log fold change indicates regions open in mCherry+ cells, negative log fold change indicates regions open in GFP+ cells. (A, B) Lateral views, anterior to the left. Scale bars: 100 μm.
Figure 3—figure supplement 1
Gating strategy for FAC-sorting cells for assay for transposase-accessible chromatin (ATAC)-seq and RNA-seq.
Cells from (A) heads and (C) trunks of GFP-; mCherry- embryos were used as gating controls. Cells from (B) heads and (D) trunks of GFP+; mCherry+ siblings were then sorted.
Figure 4 with 2 supplements
Validation of potential enhancer regions in the entpd5a promoter indicates distinct regulation of classical and non-classical osteoblasts.
(A) Chromatin accessibility profiles of representative samples of each cell type. Highlighted are six differentially accessible regions between entpd5a+ and entpd5a- cells. Chromatin regions 4–6 lie proximal to the entpd5a open reading frame (ORF), regions 2–3 are positioned within introns 1 and 3, respectively, whereas region 1 is found in the 3’ end of the coq6 gene. The three construct schematics, (entpd5a(1.1):GFP, entpd5a(2.2):GFP, entpd5a(2.2-introns):GFP) indicate the sequences cloned for generation of the corresponding transgenic lines. (B) View of the UCSC browser showing a total of 10 regions of open chromatin within the 40 kb surrounding the ORF of entpd5a, including regions 1–6 shown in (A) and regions 7–10 located further upstream. ORFs of entpd5a and downstream coq6 are shown in purple. The TSS of entpd5a is marked as a green bar in the 5’UTR. The ZF_carp_phastcons track indicates conservation amongst zebrafish, goldfish, and two species of carp. Peaks identified by macs2 as significant compared to background are shown as boxes for each sample (red tracks for cartilage and intersegmental cells and green tracks for osteoblasts and segments). Boxes darker in colour indicate peaks with higher score values. Peaks of interest 1–10 are highlighted and it is indicated using asterisks whether a peak was tested for enhancer function by inserting the sequence in a construct driving GFP and/or by deletion in the BAC construct. (C) Embryo from a stable transgenic line with region 4 placed upstream of GFP, imaged at 6 dpf. GFP is expressed in osteoblasts (marked using entpd5a:pkRed), but not in the notochord. Insets demonstrate GFP and pkRed overlap in the operculum (op), but not in pkRed+ segmented notochord sheath cells (NSCs). (D, F) Lateral and (E) ventral head view of stable line where GFP is driven both by region 4 and by regions 2 and 3. Cranial osteoblasts are expressing GFP (D, E), and segmentation is clearly observed in the notochord (D, F). (G) Quantification of fluorescence in segments 1–6 vs. intersegmental regions as shown in (F). p-Value = 0.003 using an unpaired two-sample t-test (n = 6 segments, n = 5 intersegmental regions). cl, cleithrum; no, notochord; op, operculum. (C, D, F) Lateral views; (E) ventral view. Anterior to the left in all images. (C–F) Scale bars: 100 μm. Operculum (C) inset scale bar: 20 μm.
Figure 4—figure supplement 1
No additional enhancers are present in the sequence within peak number 5.
(A) Chromatin accessibility profiles of representative samples of each cell type. In squares are the peak region covered by the entpd5a(2.2):GFP construct and the one covered by the entpd5a(5.7):GFP construct. (B, C) GFP driven in a stable line by the 5.7 kb upstream of the start codon of entpd5a, is only expressed in cranial osteoblasts (B), but not in the notochord (C) or in cartilage of the head skeleton (B, inset). (B) Ventral view; (C) lateral view. Anterior towards the left. (D, E) Scale bars: 100 μm. (D, inset) Scale bar: 20 μm.
Figure 4—figure supplement 2
GFP expression in the notochord under the control of the introns is first detected in notochord progenitor cells.
(A, B) GFP and pkRed expression (under the control of the complete BAC) appear to completely overlap in notochord progenitor cells at 24 hpf. Scale bar: 100 μm.
Figure 5
Deletion of open chromatin further upstream or downstream of the entpd5a(2.2-introns) domain does not affect the entpd5a expression pattern.
(A) Chromatin accessibility profiles along the length of the entpd5a BAC, starting from the ATG of entpd5a. Representative tracks for each cell population are shown. Peaks of interest are labelled 4–10, following the labelling in Figure 4. The region highlighted in red indicates the 21.4 kb deleted in the BAC construct entpd5a(Δ21):pkRed. (B–E) Embryo stably expressing the entpd5a(Δ21):pkRed construct shows normal segmentation pattern at 6 dpf (B–C, asterisks), as well as normal osteoblast and chondrocyte (D, E, insets) expression of pkRed. (F–H) Two representative entpd5a:pkRed+ embryos injected with entpd5a(2.2-introns):GFP show mosaic GFP expression along notochord sheath cells (NSCs) associated with pkRed+ segments at 5 dpf (F, G). In contrast, representative embryos injected with entpd5a(2.2-coq6):GFP only show background GFP at the same stage (H, I). (J) Schematic of the entpd5a(2.2-coq6):GFP construct, integrating peak number 1 downstream of the GFP. (K) Proportion of embryos injected with GFP constructs with observable GFP expression in NSCs. Bar heights show the pooled percentage (entpd5a(2.2-introns):GFP: 28/97 embryos from 2 independent injections; entpd5a(2.2-coq6):GFP: 1/70 embryos from 3 independent injections). Error bars show the mean absolute deviation of the percentages across independent injections. Groups were compared by Fisher’s exact test on pooled embryo counts (p = 6.4 × 10⁻⁷). c, chondrocyte; cl, cleithrum; op, operculum; VC, vacuolated cell. (B, C, F–I) Lateral views; (D, E) ventral views. Heads positioned towards the left. (B–I) Scale bars 100 μm. (D, E) Insets’ scale bar: 20 μm.
Figure 6
Sequential deletions of the entpd5a(2.2):GFP construct highlight a 422 bp region containing active entpd5a enhancers.
(A) Schematic of the initial construct (2.2 kb cloned upstream of GFP) and the subsequent deleted constructs. In the final construct, the purple star depicts the –31 bp deletion and the light blue star the –37 bp deletion. The osteoblast assay for transposase-accessible chromatin (ATAC) peak position is indicated. The 422 bp of interest is highlighted in red. (B–I) Representative entpd5a:pkRed embryos at 3 dpf injected with the respective deletion constructs, or as a control with the complete construct (B). The cleithra (cl) for each embryo are depicted without transmitted light. (J) Graph indicating, across all experimental repeats, the percentage of successfully injected embryos (n) for each injected construct in which GFP+ cells were observed in the cleithrum. Successfully injected embryos were identified based on the presence of GFP+ cells in background tissues. Lateral views, anterior to the left. (B–I) Scale bar: 100 μm.
Figure 7 with 1 supplement
RNA-seq analysis of skeletogenic populations of the head and trunk.
(A, B) Heatmaps for each tissue, indicating the top differentially regulated genes. The colours on the main heatmap indicate the Z-score value, while the second heatmap indicates the corresponding log fold change, and the third the average expression of each gene. (C) Venn diagram indicating the overlap of total differentially expressed genes (DEGs) in all four cell populations. Osteogenic and cartilage tissues share 872 and 765 upregulated genes, respectively, while the non-contrasted osteoblasts vs. intersegmental and cartilage vs. segments share 43 and 7 genes, respectively. (D) Principal component analysis (PCA) plot showing the clustering of cell populations by tissue (head vs. trunk) and by entpd5a expression (positive and negative cell populations in blue and orange, respectively).
Figure 7—figure supplement 1
Chondrocyte markers maintain an open chromatin configuration in osteoblasts and are actively expressed in entpd5a+ cells in levels lower than that found in entpd5a- cells.
(A–D) Views of the col9a2 (A), sox9a (B), acana (C), and ccn6 (D) promoters in representative assay for transposase-accessible chromatin (ATAC) samples for each of our cell populations, with the scale bars indicating the maximum normalised ATAC signal (read counts) in each sample. (E) Normalised col9a2, sox9a, acana, and ccn6 transcript counts in the three entpd5a+ (green spots) and the three entpd5a- (red circles) samples sequenced in the head and in the trunk. (F) Corresponding log2 fold change values for each gene when entpd5a+ vs. entpd5a- cell populations are compared within each tissue.
Figure 8
Integration of assay for transposase-accessible chromatin (ATAC)-seq and RNA-seq data highlights differentially expressed genes (DEGs) with accessible promoter regions in distinct cell populations.
(A, B) Bar graphs indicate the proportion of ATAC peaks found within 5 kb, 10 kb, and 20 kb of the TSS in each cell population. In (A) peaks were searched in all of the Zv11 genes and in (B) only in DEGs. (C–E) Venn diagrams indicating for each cell population the overlap between genes in Zv11 associated with ATAC peaks within 20 kb of the TSS, as well as the genes differentially regulated.
Figure 9
Distinct predicted binding sites could facilitate transcription factor binding of the entpd5a promoter, aiding cell type-specific regulation of the gene’s expression in classical and non-classical osteoblasts.
Our analysis suggests that, with the exception of Runx2, distinct transcriptional regulators are functioning in (A) classical vs. (B) non-classical osteoblasts. In silico analyses predict the presence and distribution of binding sites (squares, coloured according to transcription factor specificity) of the corresponding regulators in the three non-coding regions of entpd5a, which are shown to have regulatory function. The width of each square represents the length of sequence of the respective binding site. The red box indicates the location of the 37 bp deleted using entpd5a(r37):GFP.
Tables
Table 1
Revised list of tissues with detected entpd5a in the developing zebrafish embryo.
| Tissue | 15 s – 24 hpf | 24–72 hpf | >4 dpf |
|---|---|---|---|
| Head | No expression | Osteoblasts | Osteoblasts |
| Chondrocytes | Chondrocytes | ||
| Neurons | |||
| Trunk | NPCs | Ventral NSCs | Segmented NSCs (from anterior to posterior) |
| Neurons | Ventral NSCs (posterior) |
Table 2
Candidate transcription factors (TFs) regulating bone and cartilage development.
| Cell population | Enriched within | Homer prediction | HOMER discovery method | p-Value | % of Targets sequences with motif | % of Background sequences with motif | Zebrafish gene |
|---|---|---|---|---|---|---|---|
| Osteoblasts | 0–5 kb | DLX1 | Known motif | 1.00E-03 | 44.29% | 24.22% | dlx1a |
| 0–5 kb | DLX2 | Known motif | 1.00E-02 | 42.86% | 27.51% | dlx2a | |
| 0–5 kb | Dlx3 | Known motif | 1.00E-03 | 25.71% | 11.67% | dlx3b | |
| 0–5 kb | DLX5 | Known motif | 1.00E-05 | 34.29% | 13.49% | dlx5a | |
| 5–20 kb | HAND2 | De novo discovery (0.790) | 1.00E-08 | 11.49% | 0.76% | hand2 | |
| 5–20 kb | HLF | De novo discovery (0.623) | 1.00E-09 | 32.18% | 8.65% | hlfb | |
| 0–5 kb | IsI1 | Known motif | 1.00E-02 | 44.29% | 29.56% | isl1 | |
| 5–20 kb | MEIS2 | De novo discovery (0.750) | 1.00E-07 | 10.34% | 0.81% | meis2a | |
| 0–5 kb | NFAT | Known motif | 1.00E-06 | 37.14% | 13.07% | nfatc4 | |
| 0–5 kb | PRDM1 | De novo discovery (0.732) | 1.00E-02 | 8.57% | 2.05% | prdm1a | |
| 5–20 kb | RBPJ | Known motif | 1.00E-02 | 8.05% | 2.68% | rbpjb | |
| 5–20 kb | RUNX1 | Known motif | 1.00E-02 | 25.29% | 14.84% | runx1t1 | |
| 5–20 kb | Sox14 | De novo discovery (0.743) | 1.00E-07 | 13.79% | 1.59% | sox14 | |
| 5–20 kb | Twist2 | Known motif | 1.00E-03 | 43.68% | 26.61% | twist2 | |
| Segments | 0–5 kb | Hic1 | De novo discovery (0.747+0.693) | 1.00E-03 | 33.33% | 0.00% | hic1 |
| 0–5 kb | HOXB9 | De novo discovery (0.623) | 1.00E-04 | 33.33% | 0.00% | hoxb9a | |
| 0–5 kb | Hoxc8 | De novo discovery (0.741) | 1.00E-04 | 33.33% | 0.00% | hoxc8a | |
| 5–20 kb | Hoxc9 | De novo discovery (0.731) | 1.00E-03 | 33.33% | 0.00% | hoxc9a | |
| 0–5 kb | Hoxd10 | Known motif | 1.00E-02 | 100.00% | 15.50% | hoxd10a | |
| 0–5 kb | KIf7 | De novo discovery (0.689) | 1.00E-02 | 33.33% | 0.17% | klf7a | |
| 0–5 kb | klf7b | ||||||
| 5–20 kb | NFAT | De novo discovery (0.737) | 1.00E-04 | 33.33% | 0.00% | nfatc4 | |
| 5–20 kb | NFIC | De novo discovery (0.662) | 1.00E-04 | 66.67% | 0.28% | nfic | |
| 0–5 kb | NKX2-8 | De novo discovery (0.669) | 1.00E-03 | 33.33% | 0.01% | nkx2.2a | |
| 0–5 kb | Nr5a2 | De novo discovery (0.591) | 1.00E-04 | 33.33% | 0.00% | nr5a2 | |
| 0–5 kb | Tb ×20 | De novo discovery (0.703) | 1.00E-04 | 33.33% | 0.00% | tbx20 | |
| 0–5 kb | Tcf7 | De novo discovery (0.694) | 1.00E-04 | 33.33% | 0.00% | tcf7 | |
| 0–5 kb | TEAD | Known motif | 1.00E-02 | 66.67% | 5.44% | tead3b | |
| Intersegmental | 0–5 kb | Foxa2 | Known motif | 1.00E-03 | 47.83% | 15.82% | foxa2 |
| 0–5 kb | Foxa3 | Known motif | 1.00E-04 | 34.78% | 6.46% | foxa3 | |
| 5–20 kb | MYB | De novo discovery (0.640) | 1.00E-05 | 34.48% | 5.74% | mybl1 | |
| 0–5 kb | NFAT | De novo discovery (0.660) | 1.00E-03 | 4.35% | 0.00% | nfat5b | |
| 0–5 kb | NFAT | De novo discovery (0.660) | 1.00E-03 | 4.35% | 0.00% | nfatc2a | |
| 5–20 kb | Sox8 | De novo discovery (0.816) | 1.00E-07 | 31.03% | 2.82% | sox8a | |
| 0–5 kb | Sox9 | Known motif | 1.00E-03 | 39.13% | 9.44% | sox9a | |
| 0–5 kb | TBX3 | De novo discovery (0.644) | 1.00E-02 | 4.35% | 0.01% | tbx3b | |
| Cartilage | 0–5 kb | Foxa3 | Known motif | 1.00E-10 | 22.22% | 6.16% | foxa3 |
| 0–5 kb | FOXP1 | Known motif | 1.00E-05 | 22.84% | 10.14% | foxp1a | |
| 0–5 kb | Sox8 | De novo discovery (0.755) | 1.00E-13 | 17.90% | 3.21% | sox8a | |
| 0–5 kb | Sox9 | Known motif | 1.00E-04 | 19.75% | 9.60% | sox9a | |
| 5–20 kb | Tbx20 | De novo discovery (0.698) | 1.00E-11 | 26.49% | 8.84% | tbx20 |
Table 3
Volume of TDE1 tagment DNA enzyme, total reaction volume, and incubation time of transposition reaction according to cell numbers.
| Cell number | 5–6.9k | 7–9.9k | 10–12.9k | 13–20.9k | 21–39.9k | 40–60k |
|---|---|---|---|---|---|---|
| TDE1 enzyme | 1 μl | 1.5 μl | 1.8 μl | 2.3 μl | 3 μl | 3 μl |
| Total reaction volume | 20 μl | 20 μl | 50 μl | 50 μl | 50 μl | 50 μl |
| Incubation time | 30 min | 30 min | 30 min | 30 min | 30 min | 45 min |
Additional files
-
Supplementary file 1
Quality control measurements for assay for transposase-accessible chromatin (ATAC)-seq samples.
Fragment length distribution distance (FLDD) is defined as the distance of the experiment’s FLD to a reference distribution. Negative and positive distances are associated with under- and over-transposed samples, respectively. TSS enrichment score is calculated based on the normalised fragment coverage ±2 kb around the TSS.
- https://cdn.elifesciences.org/articles/100230/elife-100230-supp1-v1.docx
-
Supplementary file 2
Primers used for cloning.
- https://cdn.elifesciences.org/articles/100230/elife-100230-supp2-v1.docx
-
MDAR checklist
- https://cdn.elifesciences.org/articles/100230/elife-100230-mdarchecklist1-v1.pdf
Download links
A two-part list of links to download the article, or parts of the article, in various formats.
Downloads (link to download the article as PDF)
Open citations (links to open the citations from this article in various online reference manager services)
Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)
Integration of ATAC- and RNA-sequencing identifies chromatin and transcriptomic signatures in classical and non-classical zebrafish osteoblasts and indicates mechanisms of entpd5a regulation
eLife 13:RP100230.
https://doi.org/10.7554/eLife.100230.3
