Zebrafish are vertebrate animals, and all experiments and animal husbandry were conducted in accordance with the Dana-Farber Cancer Institute’s Institutional Animal Care and Use Committee (IACUC)-approved animal research protocol #02-107. The stocks of wild-type and mutant lines were maintained in accordance with a previously reported protocol (10).

One-cell fertilized zebrafish embryos were injected with various amounts of mRNAs encoding zinc finger nucleases (ZFNs) harboring DD/RR or EL/KK variant FokI nuclease domains, as previously described (7). Site-specific ZFN function was verified by PCR on genomic DNA (gDNA) using tet2-specific primers and then sequencing. The resulting tet2^m/m zebrafish line, which has been described and characterized previously (7–9), was used in this study.

For adult stages, fish were individually genotyped using one forward (5′ ATCTCCAAGGTCTTGCAACCTA 3′) and one reverse (5′ ATACAAGCCCTCATCCACTGAT 3′) primer, as previously shown (7). The forward and reverse primers were located in the sixth and eighth exons of tet2, respectively. The PCR products, 668 bp in length, were sequenced with the forward primer.

The cell suspension preparation and flow cytometry were conducted as previously described (7). In short, 4-month- and 15-month-old wild-type tet2^wt/wt and tet2^m/m fish were anesthetized. The anesthesia of the fish was carried out prior to kidney removal to harvest the kidney marrow. The fish were anesthetized by immersion in a small fish tank containing 200 mg/L of tricaine methanesulfonate dissolved in fish water. The working solution (200 mg/L of tricaine methanesulfonate dissolved in fish water) was made fresh by diluting a stock solution of tricaine methanesulfonate (4 mg/mL) in deionized water buffered with Tris-HCL to pH 7.0, which was stored in a freezer at −20°C. Tricaine methanesulfonate is an ester-type local anesthetic agent that acts systemically when absorbed through the gills, and thereafter is distributed throughout the body via the blood and acts on the peripheral and central nervous systems. After the fish were completely anesthetized, the kidney containing the whole kidney marrow was surgically removed. We then euthanized the fish by transferring them to a small fish tank immersed in an ice bath with a higher concentration of tricaine methanesulfonate (of 500 mg/L), which was prepared fresh by diluting the stock solution as previously described. The fish were kept in this euthanizing dose of tricaine methanesulfonate for at least 30 minutes on ice and then frozen in plastic bags at −20°C in a freezer designated as the morgue, prior to disposal by our Animal Resource Facility.

The removed kidneys were dissected and placed in ice-cold 0.9 × phosphate-buffered saline (PBS) containing 5% fetal calf serum (FCS). Single-cell suspensions were generated by aspiration, followed by mild teasing of the kidney on a 40-μm nylon mesh filter with a pipette tip. Progenitor cells were sorted using a JF Aria II Cell Sorter (BD Biosciences) with high forward scatter/side scatter (FSC/SSC), as previously described.

The enhanced reduced representation of bisulfite sequencing (ERRBS) was conducted as previously described (11, 12). Briefly, high-molecular-weight DNA was isolated from the progenitor cells using Qiagen’s AllPrep™ extraction kit, per the manufacturer’s instructions. The DNA was then digested with the MspI restriction enzyme. Fragments were purified and subjected to end repair, A-tailing, and ligation of methylated Illumina adapters. This was followed by bisulfite conversion and PCR amplification. The libraries were sequenced on an Illumina 2500 using a 50-bp single-end read approach, per the manufacturer’s instructions. The sequencing data were aligned with Zv9 as the reference genome using our in-house processing pipeline (11). The CpG sites with at least 10× coverage were used for downstream analysis. The sequencing statistics are included in Supplementary Table 1A .

For each CpG, the interquartile range (IQR) across all samples was calculated using the “iqr” function in R (version 3.3.0; The R Foundation for Statistical Computing, Vienna, Austria). CpGs with an IQR greater than the 90th percentile were annotated as high-variance CpGs.

The RNA was isolated from the progenitor cells using Qiagen’s AllPrep extraction kit, per the manufacturer’s instructions. RNA-seq libraries were prepared using TruSeq™ RNA-Seq by poly(A) enrichment (Illumina) and sequenced on a HiSeq^® 2000 system (Illumina) using a 50-bp single-end approach, per the manufacturer’s instructions. Sequencing statistics are included in Supplementary Table 1B .

FASTQ files were checked for quality using FASTQC (version 0.10.1). The reads were then aligned using STAR aligner (version 2.3.0) (13) in single-end mode. After sorting and indexing, the aligned files were used with BEDTools (version 2.16.1) (14) and the Ensembl zebrafish transcriptome (Zv9) to generate read counts for genes for all samples. We calculated the normalized expression values using DESeq2 normalization. The normalized expression values were averaged across the two replicates in each of the wild-type and mutant categories. Standard deviation was used to calculate the bounds for error bars.

Principal component analysis (PCA) was conducted using the “prcomp” function in R. The first two PCs were plotted using ggplot2 in R.

The change in DNA CpG methylation between any given groups for a given region was calculated using methylKit (15). Regions with a minimum coverage of 10 reads and a maximum coverage greater than the 90th percentile of read counts were excluded from the analysis. Regions with an absolute methylation difference of at least 15% at a statistical significance of less than 0.05 were annotated as differentially methylated regions (DMRs).

Differential gene expression analysis of any two comparisons was conducted using DESeq2 (version 1.14.1) (16) package in R (version 3.3.0). The batch correction variable accounting for variation between replicates was added in the design equation. For each comparison, genes with greater than five reads across the samples were retained for analysis. Genes with an absolute value of log-fold change greater than 1 and a q-value less than 0.05 were annotated as differentially expressed genes (DEGs).

To identify significantly enriched pathways for a given set of genes (either DEGs or DMRs), we converted the zebrafish gene names to their human orthologs using BioMart from the Ensembl database (v70). Overlap of the human orthologs with known pathways was assessed using “enrichR” in the R package (17–19). Significantly enriched pathways were identified at a p-value less than 0.05.

We identified a higher density of 75%–100% DNA CpG methylation in 4-month-old tet2^m/m fish than in 4-month-old tet2^wt/wt fish (Wilcoxon rank-sum test p < 0.05; Figure 2B ), which was consistent with the role of tet2 in catalyzing the removal of CpG methylation (20). To examine the differences in DNA methylation at the level of specific gene promoters upon tet2 loss, we compared the overall methylation levels of promoter regions in 4-month-old tet2^m/m fish with those in 4-month-old tet2 ^wt/wt fish, as visualized in the volcano plot shown in Figure 2C . Of the 4,096 promoter regions analyzed, 98 were significantly hypermethylated in 4-month-old tet2^m/m fish, and 20 were hypomethylated [ Figure 2C and Supplementary Table 2A ; methylation difference of at least 15% and statistical significance (q-value) of less than 0.05].

We also examined the genes and pathways that were affected by differential CpG methylation in tet2-mutant hematopoietic progenitor cells. We conducted gene set enrichment analysis (GSEA) to identify transcription factor gene targets with promoters that were significantly hypermethylated or hypomethylated in the HSPCs of 4-month-old tet2^m/m fish compared with those of tet2 ^wt/wt fish, based on the CHEA Transcription Factor Target Dataset (21, 22). Of the 98 hypermethylated promoters in 4-month-old tet2^m/m fish, 28 were predicted to be targets of SUZ12, a component of the polycomb repressive complex 2 (PRC2), which mediates H3K27 trimethylation (H3K27me3) (23). Another 21 gene promoters with increased methylation were predicted to be targets of the metal-response-element-binding transcription factor 2 (MTF2) (p-value < 0.05). MTF2 is a metal-response-element-binding transcription factor that recruits PRC2 to specific genes (24). Furthermore, 19 of the hypermethylated genes are common targets of SUZ12 and MTF2, including the transcription factors NKX-2, DBMX1, TSC22D3, and NR2F2. Therefore, the pathways affected by increased methylation due to tet2 loss involve both the recruitment and functional activity of PRC2 ( Figure 2D , Supplementary Table 2B ). Hypomethylated gene promoters in 4-month-old tet2^m/m fish were the targets of the master transcription factors MYB, FOXM1, SCL (or TAL1), and KLF5, which are involved in the control of cell proliferation ( Figure 2E , Supplementary Table 2C ) (25, 26). Each of these master transcription factors is upregulated in cancers of different types.

We also assessed differences in the methylation levels of CpG islands (CGIs) between 4-month-old tet2^m/m fish and 4-month-old tet2^wt/wt fish. Of the 4,350 CpG islands assessed, 97 were hypermethylated and 15 were hypomethylated in 4-month-old tet2^m/m fish, similar to the pattern of promoter methylation ( Figure 2F , Supplementary Table 2D ). Notably, nine of the CpG islands were within gene promoter regions and showed a consistent pattern of differential methylation as the promoter region in 4-month-old tet2^m/m fish ( Supplementary Table 2D ).

We next compared the DNA CpG methylation patterns in the gene promoter regions of 15-month-old tet2^m/m fish to those of 4-month-old tet2^m/m fish ( Figure 3A ). Of the 3,873 promoters analyzed, we found that 55 were significantly hypermethylated and 74 were significantly hypomethylated ( Figure 3A , Supplementary Table 3A ) at 15 months. We conducted GSEA to identify the transcription factor gene targets with promoters that had acquired significant hypermethylation or hypomethylation at 15 months, based on the CHEA Transcription Factor Target Dataset (18, 22). Hypermethylated gene promoters were identified to be targets of CEBPB, RING1B, SOX2, SETDB1, and TCF4 ( Figure 3B , Supplementary Table 3B ). In addition, hypomethylated gene promoters were identified to be targets of EZH2, CDX2, SOX9, RING1B, and SUZ12 ( Figure 3C , Supplementary Table 3C ). In addition to promoters, we also identified changes in DNA methylation in CGIs ( Figure 3D , Supplementary Figure 3D ). Of the 4,336 CGIs analyzed, we identified 76 that were significantly hypermethylated and 68 that were significantly hypomethylated in the HSPCs of 15-month-old tet2^m/m fish compared with 4-month-old tet2^m/m fish.

We then assessed if the changes in the DNA CpG methylation acquired by HSPCs in 15-month-old tet2^m/m fish were a continuation of the differences identified between 4-month-old tet2^m/m and 4-month-old tet2^wt/wt fish ( Supplementary Figure 1 ). We found that 89 of the gene promoters that were unchanged (or not significantly differentially methylated) in the 4-month-old fish were newly differentially methylated in the 15-month-old tet2^m/m fish compared with the 4-month-old tet2^m/m fish ( Supplementary Figure 1 ). Interestingly, the new DMRs arising in 15-month-old tet2^m/m fish were approximately evenly divided between significant hypermethylation and hypomethylation. Among the 44 newly hypermethylated gene promoters, 15 were predicted to be either CEBPB or RING1B targets, whereas the targets of the 45 newly hypomethylated gene promoters included CDX2, YY1, CEBPB, ERA, and CREB1 ( Supplementary Figure 1 and Supplementary Table 4 ).

Furthermore, 66 of the gene promoters that were hypermethylated in the HSPCs of 4-month-old tet2^m/m fish compared with those of 4-month-old tet2^wt/wt fish were not differentially methylated in the 15-month-old tet2^m/m compared with the 4-month-old tet2^m/m fish.

Twenty-eight promoters were hypermethylated in the HSPCs of 4-month-old tet2^m/m compared with those of 4-month-old tet2^wt/wt fish; these then lost methylation in 15-month-old tet2^m/m fish, whereas seven hypomethylated promoters in the HSPCs of 4-month-old tet2^m/m fish, compared with 4-month-old tet2^wt/wt fish, gained methylation in the HSPCs of 15-month-old tet2^m/m fish. Interestingly, the promoters of TSC22D3 and NKX2-5, which were hypermethylated in 4-month-old tet2^m/m fish, become hypomethylated in 15-month-old tet2^m/m fish.

We next sought to determine if the DNA CpG methylation changes identified in HSPCs between 4- and 15-month-old tet2^m/m fish were distinct from those observed during the normal aging process. We compared the DNA CpG methylation patterns in the gene promoter regions of HSPCs of 15-month-old tet2^wt/wt fish to those of the 4-month-old tet2^wt/wt fish ( Figure 4A , Supplementary Table 5A ). Of the 4,145 promoter regions analyzed, 99 promoters were hypermethylated and 44 were hypomethylated in the 15-month-old tet2^wt/wt fish ( Figure 4A , and bottom row Figure 4B ). Interestingly, most of these changes affect different genes than those that differed between the 4-month-old tet2^m/m fish and 4-month-old tet2^wt/wt fish (middle row, Figure 4B ). Of the promoters identified as different between the 15-month-old tet2^m/m fish and 4-month-old tet2^m/m fish (top row, Figure 4B ), most were different from the promoters that were differentially methylated during normal aging (bottom row, Figure 4B ). Thus, of the small percentage of promoters with altered levels DNA methylation during normal aging, most of these alterations appear to depend on tet2 because they are lost in tet2^m/m HSPCs ( Supplementary Table 5B ).

We have previously shown that tet2-mutant fish develop trilineage dysplasia by 15 months of age, which progresses over the next 7 months to anemia, thus fulfilling the diagnostic criteria of MDS (27). The HSPCs of 15-month-old zebrafish are thus equivalent to the tet2 mutant “pre-MDS” HSPCs in humans. A previous study assessed for the effect of tet2 mutations in CD34⁺ HSPCs from MDS patients with age-matched controls (27). To determine if the DNA CpG methylation changes observed in the tet2 fish model recapitulate the observations in humans, we compared DNA CpG methylation in gene promoter regions of 15-month-old tet2^m/m to the 15-month-old tet2^wt/wt fish ( Figure 5A ). Of the 3,929 promoters analyzed, we identified 49 hypermethylated and 62 hypomethylated in the 15-month-old old tet2^m/m fish ( Figure 5A , Supplementary Table 6 ). The small number of differentially methylated promoters in HSPCs of 15-month-old tet2^m/m fish overall compared with that in human MDS patient cells presumably reflects the less-complete sequencing coverage of the zebrafish genome and the fact that 15-month-old fish did not yet have MDS, which they would not develop for another 7 months (i.e., at 22 months of age). Nevertheless, 19.8% (n = 22) of the significantly differentially methylated gene promoters in 15-month-old zebrafish HSPCs overlapped with published differentially methylated promoters of human orthologs in MDS patients compared with those in age-matched control individuals ( Figure 5B ). These genes included NKX2-5, which is predicted to regulate aberrant gene expression in human MDS (28) and CLK2, a proposed oncogene in breast cancer (29); the promoters of both of these genes were hypomethylated in 15-month-old tet2^m/m fish compared with 15-month-old tet2^wt/wt fish.

We next determined the gene expression pathways that were either upregulated or downregulated in tet2^m/m HSPCs compared with tet2^wt/wt HSPCs in 4-month-old fish. We first identified differentially expressed genes (DEGs) in the 4-month-old tet2^m/m and 4-month-old tet2^wt/wt fish ( Figure 6A , Supplementary Table 7A ). We identified 105 upregulated genes and 119 downregulated genes (absolute log2-fold change > 1 and q-value < 0.05; Figure 6A ). Pathway enrichment of the human orthologs of these zebrafish genes identified several dysregulated pathways within the hallmark gene sets from the Molecular Signature Database (MSigDB) (30). Upregulated genes were enriched in interferon alpha and gamma pathway genes and genes controlling MTORC1 signaling ( Figure 6B , Supplementary Table 7B ). Interferon-induced protein 44-like gene (ifi44l) was included in genesets for both interferon alpha and gamma pathways. Interestingly, each of the four orthologs of ifi44l in the zebrafish was upregulated by 5- to 19-fold in tet2-mutant zebrafish HSPCs compared with wild-type zebrafish HSPCs, as illustrated in Figure 7 . The other upregulated genes in tet2-mutant zebrafish HSPCs in the interferon alpha and gamma pathways are shown in Supplementary Figure 2 . In addition, the relative expression of representative upregulated genes in the mTORC1 signaling pathway is plotted in Supplementary Figure 3 .

The genes in pathways that were downregulated in tet2-mutant HSPCs include pqpig, slc4a1, and cpox, which are associated with heme metabolism ( Figure 8 , Supplementary Table 7C ). The downregulation of genes involved in heme metabolism in tet2-mutant zebrafish HSPCs is in line with the finding that the cell type found to be significantly decreased in the peripheral blood was the erythrocyte. Furthermore, by 24 months of age, most of the tet2 homozygous mutant fish tested were found to have anemia (7). Downregulated genes were also enriched for the E2F targets and G2M checkpoint ( Supplementary Figure 4 ) and interleukin 2 (IL2)/signal transducer and activator of transcription 5 (STAT5) signaling ( Supplementary Figure 5 ). It will be important in future studies to examine protein expression by Western blotting to verify that the levels of differences of the expressed proteins correspond to the differences of the RNA expression values.

Finally, we assessed the possible correlation between the gene promoters subject to differential methylation and the pattern of gene expression using GSEA (30). We found a negative association between hypermethylated gene promoters and gene expression changes observed in 4-month-old tet2^m/m fish compared with 4-month-old tet2^wt/wt fish ( Figure 9 ).