F BS-seq and AZD0865 web oxBS-seq data, and thus we analyzed BS-seq and oxBSseq data separately.Binomial test with conversion efficiencyWe used the binomial test with the conversion efficiencies (BSEff = 0.99) as described in the supplement of [32] to quantify the presence of 5mC and 5hmC for each CpG. Since [32] does not provide a way to handle replicates, we combined the replicate-specific p values using Fisher’s method. We used this strategy to analyze both wild-type and knockout conditions separately. The obtained p values therefore provide a proxy for the amount of 5mC and 5hmC; low p values correspond to high amounts of cytosine modifications. Using a p value threshold we can decide the presence of 5mC and/or 5hmC in both conditions and call a difference in methylation modification levels, which we defined by using the minimum of the two p values. Finally, by sliding the p value threshold from 0 to 1 we can then generate the ROC graph and the AUC score as illustrated in Figure S8d in Additional file 4.Simulation of dataEuropean Bioinformatics Institute ArrayExpress Archive (E-MTAB-1042). Bismark v0.7.12 [63] was used to align the BS and oxBS reads against the mm9 reference genome. The alignment was done using the single-end Bowtie 2 [64] backend with the following parameters: -N 1 20. The “bismark_methylation_extractor” script distributed with the Bismark aligner was used to extract the number of unconverted and converted read-outs for each cytosine with the following parameters: utoff 5 edGraph ounts. The PCR primers given in [32] were aligned against the mm9 reference genome and the locations of the CCGG sites within the loci were extracted. The methylation levels of the second cytosine PubMed ID:https://www.ncbi.nlm.nih.gov/pubmed/26740125 within the CCGG sites were estimated using Lux ( = (0.8, 0.8, 0.8)). The Booth et al. estimates and glucMS-qPCR measurements were taken from [32].Integrative analysis of BS-seq, TAB-seq, and fCAB-seq dataThe counts of unconverted read-outs out of N read-outs from BS-seq and oxBS-seq experiments are assumed to be binomially distributed random variables with the derived emission probabilities. The experimental parameters and methylation levels are varied as indicated. Downsampling was done by sampling data from binomial distributions defined by the parameters estimated from the complete data. That is, for a given cytosine and BS-seq experiment we calculated the fraction of unconverted read-outs, NBS,C/NBS. This value was used as the success probability parameter, i.e., the probability of observing “C”. Using the defined binomial distribution, we sampled a number of “C” read-outs out of N read-outs. The same procedure was used for oxBS-seq but in that case we calculated the fractions NoxBS,C/NoxBS.Kernel density estimation in the open two-dimensional simplexA kernel density estimator was applied to data prior to ternary plotting. To deal with compositional data correctly we utilized a published method based on the use of the isometric log-ratio normal kernel (iln) [75].Comparison with glucMC-qPCR dataThe raw BS-seq and oxBS-seq data sets were downloaded from the European Molecular Biology Laboratory-First, we derived the statistical model for the simultaneous and integrative analysis of BS-seq, TAB-seq, and fCAB-seq data. The derivation of BS-seq/TAB-seq/fCABseq model followed the same principle as the aforementioned derivation of the BS-seq/oxBS-seq model. Briefly, for a given cytosine, we used a Dirichlet random variable of order four to model proportion.