1
$\begingroup$
library(ChIPseeker)
library(TxDb.Hsapiens.UCSC.hg38.knownGene)
txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene
library(clusterProfiler)
library(org.Hs.eg.db)

peaks <- readPeakFile("MACS_folder/pol2_chip/Pol_2pausing/pol2_diff_peak.bed")

print(peaks)


#peakAnno <- annotatePeak("merge/blast.bed", tssRegion=c(-3000, 3000),
#                        TxDb=txdb, annoDb="org.Hs.eg.db")

peakAnno <- annotatePeak(peaks, tssRegion=c(-1000, 1000),
                         TxDb=txdb, annoDb="org.Hs.eg.db")

So this is basically what im doing to annotate the peaks.

Now the issue comes with annotation for example I see TLR4 annotated to this

118082896   118083202   307 *   shp1_R1_peak_679    Exon (ENST00000413759.1/ENST00000413759.1, exon 1 of 1) 9   117708197   117708979   783 1   7099    ENST00000490685.1   374699  ENSG00000136869 TLR4

Now The peak is 679 as you can see from the IGV plot as well as annotated file. But when i see in IGV browser TLR4 is far away from the peak where i get the annotation.

So now is it an issue with the code im using or due to the inherent nature of the data.

Any help or suggestion would be really appreciated.

TLR4

$\endgroup$
1
$\begingroup$

I had a similar case just recently.

The tssRegion argument ChIPseeker is used to define what is considered a promoter, which is frequently the question at hand. However, ChIPseeker will not perform any kind of filtering on the annotated peak file.

In the results there is a column annotation that contains information what kind of genomic region the peak is in. If it is within your tssRegion window it will be labeled as a promoter.
In your case the annotation is 1st Exon = Exon (ENST00000413759.1/ENST00000413759.1, exon 1 of 1) 9. This tells you that the peak is not part of the promoter region you defined and located within the first exon of the indicated transcript. The peak is really far away from TLR4 (374699 bp) but I think all other sequences nearby are pseudogenes and this is why TLR4 is selected as the closest gene.

You can get an overview over the annotations using plotAnnoPie(peakAnno)

tl;dr: you need to filter the results of the ChIPseeker annotation to get the genomic regions you want.

$\endgroup$
3
  • $\begingroup$ Glad this is a reproducible problem and now its more confusing for some peak the genes assigned are perfectly fine matches to the coordinate given but for others its would be in the intro or some exon but would show as annotation being promoter. " ChIPseeker annotation to get the genomic regions you want." how shall i do that ? simple distanetoTSS filter? $\endgroup$
    – kcm
    Dec 2 '20 at 11:12
  • 1
    $\begingroup$ you can either directly use the distance to TSS to or define the promoter region via tssRegion and then filter all peaks where the annotation contains "Promoter" i.e. grepl("Promoter", annotated.df$annotation) $\endgroup$
    – PPK
    Dec 2 '20 at 11:56
  • $\begingroup$ i did annotate again with homer it seems pretty accurate as compared to chipseeker ran homer with default parameters it works fine ..need to find why the same can;t be seen with chipseeker $\endgroup$
    – kcm
    Dec 3 '20 at 9:23

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.