This is a test version of Biostars. For the public version, visit https://www.biostars.org.
No annotation peaks when running homer

Hi there, I'm trying to use homer to annotate my peaks, what I'm doing is

perl annotatePeaks.pl peaks_cntr1_peaks.narrowPeak hg38 > annotated_peaks_ctr1

and I've got this message

Peak file =peaks_cntr1_peaks.narrowPeak
    Genome = hg38
    Organism = human
sh: 1: bed2pos.pl: not found
sh: 1: checkPeakFile.pl: not found
sh: 1: cleanUpPeakFile.pl: not found
    Reading Positions...
    -----------------------
    Finding Closest TSS...
sh: 1: annotateRelativePosition.pl: not found
sh: 1: assignGenomeAnnotation: not found
readline() on closed filehandle IN at annotatePeaks.pl line 3381.
rm: cannot remove '0.395676903267489.tmp': No such file or directory
    NOTE: If this part takes more than 2 minutes, there is a good chance
        your machine ran out of memory: consider hitting ctrl+C and rerunning
        the command with "-noann"
    To capture annotation stats in a file, use "-annStats <filename>" next time
sh: 1: assignGenomeAnnotation: not found
readline() on closed filehandle IN at annotatePeaks.pl line 3418.
    Counting Tags in Peaks from each directory...
    Organism: human
    Loading Gene Informaiton...
    Outputing Annotation File...
    Done annotating peaks file

The script generates the file but when I open it, I only have the header with no further information

PeakID (cmd=peaks_cntr1_peaks.narrowPeak hg38)  Chr Start   End Strand  Peak Score  Focus Ratio/Region Size Annotation  Detailed Annotation Distance to TSS Nearest PromoterID  Entrez ID   Nearest Unigene Nearest Refseq  Nearest Ensembl Gene Name   Gene Alias  Gene Description    Gene Type

Any idea?

PS1. This is the script for macs2 :

`macs2 callpeak -t CNTR_1.bam -c Input_1.bam  -f BAM -g hg -s 36 -m 8 30 -p 0.00001 -n peaks_cntr1 --outdir peaks`

PS2. This is the flagstats result:

40076034 + 0 in total (QC-passed reads + QC-failed reads)
0 + 0 secondary
0 + 0 supplementary
0 + 0 duplicates
25088016 + 0 mapped (62.60% : N/A)
0 + 0 paired in sequencing
0 + 0 read1
0 + 0 read2
0 + 0 properly paired (N/A : N/A)
0 + 0 with itself and mate mapped
0 + 0 singletons (N/A : N/A)
0 + 0 with mate mapped to a different chr
0 + 0 with mate mapped to a different chr (mapQ>=5)
homer peaks tss
• 2,124 views • link
written by Lila M Lila M

1 answer

Try raising your p-value to make sure that they aren't all just being filtered out by that threshold (which seems pretty low).

Edit: I may have misunderstood. Do you have peaks below that threshold? Serendipitously, I'm also running homer (and chipSeeker) today for my dataset. Upon a closer look, it seems like you're missing perl scripts, but I don't see any reference to those in the HOMER page for annotatePeaks.pl (http://homer.ucsd.edu/homer/ngs/annotation.html).

• 0 views • link
written by geneticatt • geneticatt

Hi, I have around 20,000 peaks. I see your point about missing a perl scrip, are those messages appears eventualy? Thanks Best

Edit Problem solved, I've forgot to edit my ~/.bash_profile

• 0 views • link
written by Lila M Lila M

Log in to answer this question.

More posts like this