steinmann / peakzilla Goto Github PK

PEAKZILLA
---------
Peakzilla identifies sites of enrichment and transcription factor binding sites from transcription factor ChIP-seq and ChIP-exo experiments at hight accuracy and resolution. It is designed to perform equally well for data from any species. All necessary parameters are estimated from the data. Peakzilla is suitable for both single and paired end data from any sequencing platform.

Note that peakzilla is not suited for the identification of broad regions of enrichment (e.g. ChIP-seq for histone marks), we recommand using MACS instead: Zhang et al. Model-based Analysis of ChIP-Seq (MACS). Genome Biol (2008) 9(9):R137

USAGE
-----
Options:
  -h, --help            show this help message and exit
  -m N_MODEL_PEAKS, --model_peaks=N_MODEL_PEAKS
                        number of most highly enriched regions used to
                        estimate peak size: default = 200
  -c ENRICHMENT_CUTOFF, --enrichment_cutoff=ENRICHMENT_CUTOFF
                        minimum cutoff for fold enrichment: default = 2
  -s SCORE_CUTOFF, --score_cutoff=SCORE_CUTOFF
                        minimum cutoff for peak score: default = 1
  -e, --gaussian        use empirical model estimate instead of gaussian
  -p, --bedpe           input is paired end and in BEDPE format
  -l LOG, --log=LOG     directory/filename to store log file to: default =
                        log.txt
  -n, --negative        write negative peaks to negative_peaks.tsv


DEPENDENCIES
------------
* Runs on OSX and Linux based operating systems
* Requires Python 2.5 or greater
* For 2x better perormance use PyPy instead of CPython


INPUT FORMAT
------------
Peakzilla accepts BED formated alignments as input.

For converstion to BED format and working with BED files and alignments in
general I highly reccommend:

* bowtie (http://bowtie-bio.sourceforge.net/)
* SAMtools (http://samtools.sourceforge.net/)
* bedtools (http://code.google.com/p/bedtools/)


WORKFLOW EXAMPLE
----------------
# use bowtie to map uniquely mappable reads to the genome
bowtie -p4 -m1 --sam genome_index input.fastq input.sam
bowtie -p4 -m1 --sam genome_index chip.fastq chip.sam

# convert to BAM format
samtools view -bS input.sam > input.bam
samtools view -bS chip.sam > chip.bam

# convert to BED format
bamToBed -i input.bam > input.bed
bamToBed -i chip.bam > chip.bed

# run peakzilla
python peakzilla.py chip.bed input.bed > chip_peaks.tsv

# Comparison of 2 datasets
#    Determine significant peaks with a score threshold of 10
python peakzilla.py -s 10 chip1.bed input1.bed > chip1_s10_peaks.tsv
#    Determine enriched regions with a score threshold of 2
python peakzilla.py -s 2 chip2.bed input2.bed > chip2_s2_peaks.tsv
#    Overlap significant peaks from chip1 with enriched regions from chip2
intersectBed -a chip1_s10_peaks.tsv -b chip2_s2_peaks.tsv > intersect_peaks.tsv

For example datasets as well as an example of a computational pipeline for the comparative analysis of ChIP-seq datasets, please refer to our publication: Bardet AF et al. A computational pipeline for comparative ChIP-seq analyses. Nature Protocols (2011) 7(1):45-61 (http://www.starklab.org/data/bardet_natprotoc_2011/)

OPTIONS
-------
One of peakzilla's design goals is to learn all the necessary information
from the data. The usage of the options should therefore not be required.

OUTPUT FORMAT
-------------
* Results are printed as a table of tab delimited values to stdout
* Logs are appended to logs.txt in the current directory or a custom directory/filename specified by the -l option
* Enriched regions in the control sample are written to negative_peaks.tsv or a custom directory/filename specified by the -n option
* Columns represent Chromosome / Start / End / Name / Summit / Score / ChIP / Control / FoldEnrichment / DistributionScore / FDR (%)