HOMER

annotatePeaks.pl <peak/BED file> <genome> [options] > <output file>

i.e. annotatePeaks.pl ERpeaks.txt hg18 > outputfile.txt

The great thing about annotatePeaks.pl is that it combines many features into a single location.  It is the primary program to investigate how sequencing reads, sequence motifs, and other information and annotations interact.

Three primary options are available to specify types of data that can be processed by annotatePeaks.pl:


By default, these data types are processed relative to each peak/region provided in the primary input file.  There are a bunch of options that help fine tune how each type of data is considered by the program covered below.

However, annotatePeaks.pl can take the same input data and do other things, such as make histograms and heatmaps, allowing you to explore the data in a different way.

annotatePeaks.pl accepts HOMER peak files or BED files:

HOMER peak files should have at minimum 5 columns (separated by TABs, additional columns will be ignored):

BED files should have at minimum 6 columns (separated by TABs, additional columns will be ignored)

In theory, HOMER will accept BED files with only 4 columns (+/- in the 4th column), and files without unique IDs, but this is NOT recommended.  For one, if you don't have unique IDs for your regions, it's hard to go back and figure out which region contains which peak.

Mac Users: If using a EXCEL to prepare input files, make sure to save files as a "Text (Windows)" if running MacOS - saving as "Tab delimited text" in Mac produces problems for the software.  Otherwise, you can run the script "changeNewLine.pl <filename>" to convert the Mac-formatted text file to a Windows/Dos/Unix formatted text file.

If errors occur, it is likely that the file is not in the correct format, or the first column is not actually populated with unique identifiers.

TSS Mode

Pretty much everything explained in the following sections depends heavily on the "-size <#>" parameter.  A couple of quick notes:


For example, if you peaks are actually transcription start sites, you might want to specify "-size -500,100" to perform the analysis upstream -500 bp to +100 bp downstream.  If your peaks/regions are actually "transcript" regions, specifying "-size given" will count reads along the entire transcript.  If it doesn't make sense, watch Delta Force I and II back to back.  That should numb the brain enough to get it.

annotatePeaks.pl is useful program for cross-referencing data from multiple experiments.  In order to count the number of tags from different sequencing experiments, you must first create tag directories for each of these experiments.  Once created, tag counts from these directories in the vicinity of your peaks can be added by specifying "-d <tag directory 1> <tag directory 2> ...".  You can specify as many tag directories as you like.  Tag totals for each directory will be placed in new columns starting on column 18.  For example:


HOMER automatically normalizes each directory by the total number of mapped tags such that each directory contains 10 million tags.  This total can be changed by specifying "-norm <#>" or by specifying "-noadj", which will skip this normalization step.

The other important parameter when counting tags is to specify the size of the region you would like to count tags in with "-size <#>".  For example, "-size 1000" will count tags in the 1kb region centered on each peak, while "-size 50" will count tags in the 50 bp region centered on the peak (default is 200).  The number of tags is not normalized by the size of the region.

One last thing to keep in mind is that in order to fairly count tags, HOMER will automatically center tags based on their estimated ChIP-fragment lengths.  This is can be overridden by specifying a fixed ChIP-fragment length using "-len <#>" or "-fragLength <#>".  This is important to consider when trying to count RNA tags, or things such as 5' RNA CAGE/TSS-Seq, where you may want to specify "-len 0" so that HOMER doesn't try to move the tags before counting them.

HOMER can also quantify signal in bedGraph and WIG files using "-bedGraph <bedgraph file1> ..." and "-wig <wiggle file1> ...", respectively.

X-Y scatter plots are a great way to present information, and by counting tag densities from different tag directories, you can visualize the relative levels of different sequencing experiments.  Because of the large range of values, it can be difficult to appreciate the relationship between data sets without log transforming the data (or sqrt to stay Poisson friendly).  Also, due to the digital nature of tag counting, it can be hard to properly assess the data from a X-Y scatter plot since may of the data points will have the same values and overlap.  To assist with these issues, you can specify "-log" or "-sqrt" to transform the data.  These functions will actually report "log(value+1+rand)" and "sqrt(value+rand)", respectively, where rand is a random "fraction of a tag" that adds jitter to your data so that data points with low tag counts will not have exactly the same value. For example, lets look at the distribution of H3K4me1 and H3K4me3 near Oct4 peaks in mouse embryonic stem cells:

Figuring out which peaks have instances of motifs found with findMotifsGenome.pl is very easy.  Simply use "-m <motif file1> <motif file 2>..." with annotatePeaks.pl.  (Motif files can be concatenated into a single file for ease of use)  This will search for each of these motifs near each peak in your peak file.  Use "-size <#>" to specify the size of the region around the peak center you wish to search.  Found instances of each motif will be reported in additional columns of the output file.  For example:

Each instance of the motif is specified in the following format (separated by commas):

There are also a bunch of motif specific options for specialized analysis:

This feature may seem slightly out of place, but since annotatePeaks.pl is the workhorse of HOMER, you can add "-mbed <filename>" in conjunction with "-m <motif file 1> [motif file 2] ..." to produce a BED file describing motif positions near peaks that can be loaded as a custom track in the UCSC genome browser.  In the example below, you would load motif.bed as a custom track:

In order to find the nearest peak from another set of peaks, use "-p <peak file 1> [peak file 2] ...".  This will add columns to the output spreadsheet that will specify the nearest peak ID and the distance to that peak.  If all you want is the distance (so you can sort this column), add the option "-pdist" to the command.  Otherwise, if you prefer to count the number of peaks in the peak file found within the indicated regions (i.e. with "-size <#>"), add "-pcount".

HOMER can be used to make histograms that document sequencing library and motif densities relative to specific positions in the genome.  This can be done near peaks, subsets of peaks, or near promoters, exon junctions or anywhere else you find interesting.  To make histograms, use the annotatePeaks.pl program but add the parameters "-hist <#>" to produce a tab delimited text file that can then be visualized using EXCEL or other data visualization software.

Basic usage:

annotatePeaks.pl <peak file> <genome> -size <#> -hist <#> -d <tag directory 1> [tag directory2] ... -m <motif 1> <motif 2> ... >  <output matrix file>

i.e. annotatePeaks.pl ERpeaks.txt hg18 -size 6000 -hist 25 -d MCF7-H3K4me1/ MCF7-H3K4me2/ MCF7-H3K4me3/ > outputfile.txt

Running this command is very similar to creating annotated peak files - in fact, most of the data can be used to make both types of files - hence the reason for combining this functionality in the same command.  Be default, HOMER normalizes the output histogram such that the resulting units are per bp per peak, on top of the standard total mapped tag normalization of 10 million tags.

For each tag directory or motif, HOMER will output 3 columns in the histogram.  In the case of tag directories, the first column will indicate ChIP-Fragment Coverage, which is calculated by extending tags by their estimated ChIP-fragment length, and is analogous to the profiles made for the UCSC Genome Browser.  The 2nd and 3rd columns report the density of 5' and 3' aligned tags, and are independent of fragment length.  For example, lets look at H3K4me2 distribution near Androgen Receptor (AR) peaks before and after 16 hours of treatment with testosterone (dht):

annotatePeaks.pl ARpeaks.txt hg18 -size 4000 -hist 10 -d H3K4me2-control/ H3K4me2-dht-16h/ > outputfile.txt

Opening outputfile.txt with EXCEL, we see:



Graphing columns B and E while using column A for the x-coordinates, we get the following:


However, if we graph only the 5' and 3' tags that come from the H3K4me2-dht-16h directory (columns F and G):


Here we can see how the 5' and 3' reads from the H3K4me2 marked nucleosomes are distributed near the AR sites.

Making histograms out of motif occurrences is very similar to sequencing tag distributions.  Run the annotatePeaks.pl program with "-hist <#>" and "-m <motif file>" (you can also find motif densities and tag densities at the same time):


Graphing outputfile.txt with EXCEL:

One cool analysis strategy is to center peaks on a specific motif.  For example, by centering peak for the Androgen Receptor Transcription Factor on the ARE motif (GNACANNNTGTNC), you can map the spacial relationship between the motif and other sequence features and sequencing reads.

To center peaks on a motif, run annotatePeaks.pl with the following options:

annotatePeaks.pl <peak file> <genome> -size <#> -center <motif file> > newpeakfile.txt

i.e. annotatePeaks.pl ARpeaks.txt hg18r -size 200 -center are.motif > areCenteredPeaks.txt

Now the idea is to use the new peak file to perform analysis:

HOMER is not capable of generating actual Heatmaps per se, but it will generate the data matrix (similar to a gene expression matrix) than can then be visualized using standard gene expression heatmap tools.  For example, I will generate a heatmap data matrix file using HOMER, and then open it with Cluster 3.0 (Micheal Eisen/de Hoon) to cluster it and/or visualize it with Java Tree View (by Alok J. Saldanha).  In reality, you can use any clustering and/or heatmap visualization software (i.e. R).

Basic usage (add "-ghist" when making a histogram):

annotatePeaks.pl <peak file> <genome> -size <#> -hist <#> -ghist -d <tag directory 1> [tag directory2] ... > <output matrix file>

i.e. annotatePeaks.pl ARpeaks.txt hg18 -size 6000 -hist 25 -ghist -d H3K4me2-control/ H3K4me2-dht-16h/ > outputfile.txt

Running this command is very similar to making histograms with annotatePeaks.pl.  In fact, a heatmap isn't really all that different from a histogram - basically, instead of averaging all of the data from each peak, we keep data from each peak separate and visualize it all together in a heatmap.  The key difference when making a heat map or a histogram is that you must add "-ghist" when making a heatmap.

The resulting file is a tab-delimited text file where the first row is a header file and the remaining rows represent each peak from the input peak file.  The first set of columns will describe the read densities from the first tag directory as a function of distance from the center of the peak, with bin sizes corresponding to the parameter used with "-hist <#>".  After the first block of columns, a second block will start over with read densities from the 2nd tag directory, and so on.  If you would like to cluster this file to help organize the patterns, make sure you only cluster the "genes" (i.e. rows).

-strand <+|-|both> (only look for motif etc. on + or - strand, default both)
-fragLength <#> (Fragment length, default=auto, might want to set to 0 for RNA)
-pc <#> (maximum number of tags to count per bp, default=0 [no maximum])
-cons (Retrieve conservation information for peaks/sites - creates new column for this information)
-CpG (Calculate CpG/GC content)
-norevopp (do not search for motifs on the opposite strand [works with -center too])
-norm <#> (normalize tags to this tag count, default=1e7, 0=average tag count in all directories)