STARR foundation award (MQZ), NIH grant ES017166 (MQZ), NSF grant DMS0906593 (HX).

1. Preparing Input Files for BCP Analysis
<ol>
 <li>Align the short reads produced from sequencing runs (ChIP and input libraries) to the appropriate reference genome using the preferred short read alignment software. The mapped locations should be converted to the 6 column browser extensible data (BED) format13 (UCSC genome browser, <a href="http://genome.ucsc.edu/" target="_blank">http://genome.ucsc.edu/</a>), a tab-delimited line per mapped read indicating the mapped chromosome, start position (0-based), end position (half-open), read name, score (optional), and strand.</li>
</ol>
2a. Diffuse Read Profiles: Preprocessing ChIP Read Densities for Detection of Enriched Islands in Diffuse Data
<ol>
 <li>Extend the ChIP and input mapped locations to a predetermined fragment length, i.e. the fragment size targeted during enzyme digestion or sonication of the DNA, usually around 200 bp. Fragment counts are then aggregated in adjacent bins. By default, bin size is set to the estimated fragment length of 200 bp.</li>
 <li>Any possible change-points in a set of bins with identical read counts will most likely fall at the outer most boundaries. Accordingly, it is improbable that a change point will occur at an internal boundary between two bins with the same read counts. So, group adjacent bins, with identical reads per bin, into a single block, i.e. bedGraph format13.</li>
</ol>
2b. Punctate Read Profiles: Preprocessing ChIP and Input BED Files for Detection of Peaks in Punctate Data
<ol>
 <li>Aggregate overlapping reads for plus and minus strand ChIP reads separately. The strand specific read densities should form a bimodal profile of plus and minus peaks. Choose plus/minus pairs of the most enriched peaks and use the distance between their summits as an estimate for the library fragment length.</li>
 <li>Shift the ChIP and input reads half the fragment length to the center and recalculate the read density of the shifted and merged plus and minus strand reads. This methodology for estimating the fragment length was adopted from Zhang, et al.3. Positions with identical merge counts should be grouped into blocks, similar to step 2a.2.</li>
</ol>
3. Estimate the Posterior Mean Read Density of Each Block using our BCMIX Approximation
<ol>
 <li>The read density of each block is modeled as a Poisson distribution, Pois(&theta;t), with a mean parameter following a mixture of Gamma distributions, &Gamma;(&alpha;,&beta;), and a prior probability of a change point occurring at any block boundary of p. Conditioning Pois(&theta;t) on G(&alpha;,&beta;) effectively renders the model an infinite state HMM. Estimate the hyper-parameters, &alpha;, &beta;, and p, using maximum posterior likelihood.</li>
 <li>Explicitly calculate the Bayes estimates for each block, &theta;t, as E(&theta;t|&gamma;Z). Replace the more traditional but time consuming forward and backward filters often used in HMMs, with the more computationally efficiently Bounded Complexity Mixture approximation to estimate posterior means, &theta;c. The resulting posterior means will be "smoothed" into an approximate piecewise constant profile so blocks with identical, &theta;c, should be further blocked together with updated boundary coordinates.</li>
</ol>
4a. Diffuse Read Profiles: Post-process Posterior Means into Segments of Diffuse Enrichment
<ol>
 <li>Use the number of input reads per each new&theta;c block as the background rate, Pois(&lambda;a) and determine enrichment using a simple hypothesis test based on whether the ChIP posterior mean, &theta;c, exceeds some threshold &delta;. The 90th-quantile is the default d and is appropriate in most cases.</li>
 <li>Merge adjacent &theta;c blocks that exceed the enrichment into a single region and report merge coordinates in simple BED format. Alternatively, one can report the &theta;c for each block in bedGraph format to preserve the high-resolution details of the read density estimates.</li>
</ol>
4b. Punctate Read Profiles: Post-process Posterior Means into Peak Candidates
<ol>
 <li>Define the background rate, Pois(&lambda;a), as the average of all read counts (&gamma;2) and identify all blocks which exceed the threshold, d. Since punctate peaks are expected to be more substantially enriched, the default &delta; is set to the 99th-quantile of Pois(&lambda;a).</li>
 <li>Set the block with the maximal &theta;c as the candidate peak summit and adjoin flanking blocks that share a similar read density (&plusmn;1 read count to allow for slight variation). This adjoined region is defined as a candidate binding site.</li>
 <li>Calculate&lambda;2 as the average read counts in the ChIP candidate binding site and hypothesis test this versus input background were the null hypothesis, H0, is that &lambda;1 &ge;&lambda;2 and reject H0 based on a p-value threshold. Output candidate peaks in BED format.</li>
</ol>

<ol>
	<li>Park, P. J. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=ChIP-seq:+advantages+and+challenges+of+a+maturing+technology.">ChIP-seq: advantages and challenges of a maturing technology.</a> Nat. Rev. Genet. 10, 669-680 (2009).</li><li>Barski, A., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=High-resolution+profiling+of+histone+methylations+in+the+human+genome.">High-resolution profiling of histone methylations in the human genome.</a> Cell. 129, 823-837 (2007).</li><li>Zhang, Y., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=Model-based+Analysis+of+ChIP-Seq+(MACS).">Model-based Analysis of ChIP-Seq (MACS).</a> Genome Biol. 9, R137 (2008).</li><li>Zang, C., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=A+clustering+approach+for+identification+of+enriched+domains+from+histone+modification+ChIP-Seq+data.">A clustering approach for identification of enriched domains from histone modification ChIP-Seq data.</a> Bioinformatics. 25, 1952-1958 (2009).</li><li>Jothi, R., Cuddapah, S., Barski, A., Cui, K., Zhao, K. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=Genome-wide+identification+of+in+vivo+protein-DNA+binding+sites+from+ChIP-Seq+data.">Genome-wide identification of in vivo protein-DNA binding sites from ChIP-Seq data.</a> Nucleic Acids Res. 36, 5221-5231 (2008).</li><li>Qin, Z. S., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=HPeak:+an+HMM-based+algorithm+for+defining+read-enriched+regions+in+ChIP-Seq+data.">HPeak: an HMM-based algorithm for defining read-enriched regions in ChIP-Seq data.</a> BMC Bioinformatics. 11, 369 (2010).</li><li>Song, Q., Smith, A. D. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=Identifying+dispersed+epigenomic+domains+from+ChIP-Seq+data.">Identifying dispersed epigenomic domains from ChIP-Seq data.</a> Bioinformatics. 27, 870-871 (2011).</li><li>Spyrou, C., Stark, R., Lynch, A. G., Tavar&#233;, S. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=BayesPeak:+Bayesian+analysis+of+ChIP-seq+data.">BayesPeak: Bayesian analysis of ChIP-seq data.</a> BMC Bioinformatics. 10, 299 (2009).</li><li>Lai, T., Xing, H. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=A+simple+Bayesian+approach+to+multiple+change-points.">A simple Bayesian approach to multiple change-points.</a> Statistica Sinica. , (2011).</li><li>Robertson, G., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=Genome-wide+profiles+of+STAT1+DNA+association+using+chromatin+immunoprecipitation+and+massively+parallel+sequencing.">Genome-wide profiles of STAT1 DNA association using chromatin immunoprecipitation and massively parallel sequencing.</a> Nat. Methods. 4, 651-657 (2007).</li><li>Stitzel, M. L., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=Global+epigenomic+analysis+of+primary+human+pancreatic+islets+provides+insights+into+type+2+diabetes+susceptibility+loci.">Global epigenomic analysis of primary human pancreatic islets provides insights into type 2 diabetes susceptibility loci.</a> Cell Metab. 12, 443-455 (2010).</li><li>Bernstein, B. E., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=The+NIH+Roadmap+Epigenomics+Mapping+Consortium.">The NIH Roadmap Epigenomics Mapping Consortium.</a> Nat. Biotechnol. 28, 1045-1048 (2010).</li><li>Karolchik, D., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=The+UCSC+Table+Browser+data+retrieval+tool.">The UCSC Table Browser data retrieval tool.</a> Nucleic Acids Res. 32, 493-496 (2004).</li><li>Matys, V., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=TRANSFAC:+transcriptional+regulation,+from+patterns+to+profiles.">TRANSFAC: transcriptional regulation, from patterns to profiles.</a> Nucleic Acids Res. 31, 374-378 (2003).</li><li>Portales-Casamar, E., et al. <a target="_blank" href="http://www.ncbi.nlm.nih.gov/entrez/query.fcgi?db=PubMed&cmd=Search&doptcmdl=Citation&defaultField=Title+Word&term=JASPAR+2010:+the+greatly+expanded+open-access+database+of+transcription+factor+binding+profiles.">JASPAR 2010: the greatly expanded open-access database of transcription factor binding profiles.</a> Nucleic Acids Res. 38, D105-D110 (2010).</li></ol>

No conflicts of interest declared.

We set out to develop a model for analyzing ChIPseq data that could identify both punctate and diffuse data structures equally well. Until now, regions of enrichment, particularly diffuse regions, which reflect the presupposed expectation of large island size, have been difficult to identify. To address these problems, we utilized the most recent advances in HMM technology, which possess many advantages over existing heuristic models and less innovative HMMs.
Our model makes use of a Bayesian...

<table border="1">
 <tbody>
 <tr>
 <td>Name of the reagent</td>
 <td>Company</td>
 <td>Catalogue number</td>
 <td>Comments (optional)</td>
 </tr>
 <tr>
 <td>Linux-based workstation</td>
 <td></td>
 <td></td>
 <td></td>
 </tr>
</tbody>
</table>

a novel bayesian change-point algorithm for genome-wide analysis of diverse chipseq data types

ChIPseq is a widely used technique for investigating protein-DNA interactions. Read density profiles are generated by using next-sequencing of protein-bound DNA and aligning the short reads to a reference genome. Enriched regions are revealed as peaks, which often differ dramatically in shape, depending on the target protein1. For example, transcription factors often bind in a site- and sequence-specific manner and tend to produce punctate peaks, while histone modifications are more pervasive and are characterized by broad, diffuse islands of enrichment2. Reliably identifying these regions was the focus of our work.
Algorithms for analyzing ChIPseq data have employed various methodologies, from heuristics3-5 to more rigorous statistical models, e.g. Hidden Markov Models (HMMs)6-8. We sought a solution that minimized the necessity for difficult-to-define, ad hoc parameters that often compromise resolution and lessen the intuitive usability of the tool. With respect to HMM-based methods, we aimed to curtail parameter estimation procedures and simple, finite state classifications that are often utilized.
Additionally, conventional ChIPseq data analysis involves categorization of the expected read density profiles as either punctate or diffuse followed by subsequent application of the appropriate tool. We further aimed to replace the need for these two distinct models with a single, more versatile model, which can capably address the entire spectrum of data types.
To meet these objectives, we first constructed a statistical framework that naturally modeled ChIPseq data structures using a cutting edge advance in HMMs9, which utilizes only explicit formulas-an innovation crucial to its performance advantages. More sophisticated then heuristic models, our HMM accommodates infinite hidden states through a Bayesian model. We applied it to identifying reasonable change points in read density, which further define segments of enrichment. Our analysis revealed how our Bayesian Change Point (BCP) algorithm had a reduced computational complexity-evidenced by an abridged run time and memory footprint. The BCP algorithm was successfully applied to both punctate peak and diffuse island identification with robust accuracy and limited user-defined parameters. This illustrated both its versatility and ease of use. Consequently, we believe it can be implemented readily across broad ranges of data types and end users in a manner that is easily compared and contrasted, making it a great tool for ChIPseq data analysis that can aid in collaboration and corroboration between research groups. Here, we demonstrate the application of BCP to existing transcription factor10,11 and epigenetic data12 to illustrate its usefulness.

BCP excels at identifying regions of broad enrichment in histone modification data. As a point of reference, we previously compared our results to those of SICER3, an existing tool which has demonstrated strong performance. To best illustrate BCP's advantages, we examined a histone modification that had been well studied to establish a foundation for assessing success rates. With this in mind, we then analyzed H3K36me3, since it has been shown to associate strongly with actively transcribed gene bodies (<stron...

Watch this Scientific Journal Video about A Novel Bayesian Change-point Algorithm for Genome-wide Analysis of Diverse ChIPseq Data Types at JoVE.com

A Novel Bayesian Change-point Algorithm for Genome-wide Analysis of Diverse ChIPseq Data Types

Our Bayesian Change Point (BCP) algorithm builds on state-of-the-art advances in modeling change-points via Hidden Markov Models and applies them to chromatin immunoprecipitation sequencing (ChIPseq) data analysis. BCP performs well in both broad and punctate data types, but excels in accurately identifying robust, reproducible islands of diffuse histone enrichment.

ChIPseq is a widely used technique for investigating protein-DNA interactions. Read density profiles are generated by using next-sequencing of ...

a-novel-bayesian-change-point-algorithm-for-genome-wide-analysis

Research

JoVE Journal

Biology

11.2K Views.  Stony Brook University. Our Bayesian Change Point (BCP) algorithm builds on state-of-the-art advances in modeling change-points via Hidden Markov Models and applies them to chromatin immunoprecipitation sequencing (ChIPseq) data analysis. BCP performs well in both broad and punctate data types, but excels in accurately identifying robust, reproducible islands of diffuse histone enrichment.