YuLab-SMU / ChIPseeker

:dart: ChIP peak Annotation, Comparison and Visualization
https://onlinelibrary.wiley.com/share/author/GYJGUBYCTRMYJFN2JFZZ?target=10.1002/cpz1.585
220 stars 74 forks source link

annotatePeak issue with genomicAnnotationPriority order #142

Closed Liufan-YIN closed 3 years ago

Liufan-YIN commented 3 years ago

Hi, Dear Prof. Yu. I am sure that I have used R version 4.0 and the lastest verion of ChIPseeker from Bioconductor. I have a problem when I use annotatePeak function about genomicAnnotationPriority.

When it is set as Promoter is the first order: tssRegion=c(-2000,0) genomicAnnotationPriority = c("Promoter","5UTR","3UTR","Exon","Intron","Downstream","Intergenic") The result show the peak belongs to Promoter(1-2kb) and distanceToTSS is -1037

suppressMessages(library("ChIPseeker")) 
suppressMessages(library("TxDb.Athaliana.BioMart.plantsmart28"));
txdb=TxDb.Athaliana.BioMart.plantsmart28;
saple1.loc <- readPeakFile("upload_file.bed")
[upload_file.bed.zip](https://github.com/YuLab-SMU/ChIPseeker/files/6035140/upload_file.bed.zip)

saple1.loc.PeakAnno <- annotatePeak(saple1.loc,TxDb =txdb,tssRegion=c(-2000,0),assignGenomicAnnotation = TRUE,genomicAnnotationPriority = c("Promoter","5UTR","3UTR","Exon","Intron","Downstream","Intergenic"))
saple1.loc.PeakAnno.DF.H2AZ <- as.data.frame(saple1.loc.PeakAnno) 
subset(saple1.loc.PeakAnno.DF.H2AZ,saple1.loc.PeakAnno.DF.H2AZ[,2]==20321001) 
#     seqnames    start      end width strand  V4  V5           V6       V7
#3539        1 20321001 20321999   999      * 406 117 1.153546e-63 2.625169
#               V8       annotation geneChr geneStart  geneEnd geneLength
#3539 6.046514e-63 Promoter (1-2kb)       1  20323036 20328009       4974
#     geneStrand    geneId transcriptId distanceToTSS
#3539          1 AT1G54440  AT1G54440.2         -1037

However when it is set as Promoter is the second order after 5UTR : tssRegion=c(-2000,0) genomicAnnotationPriority = c("5UTR","Promoter","3UTR","Exon","Intron","Downstream","Intergenic") The result show distanceToTSS of the peak is still -1037 and the peak should still belongs to Promoter(1-2kb), But the result showed the peak belongs to Distal Intergenic.

saple1.loc.PeakAnno <- annotatePeak(saple1.loc,TxDb =txdb,tssRegion=c(-2000,0),assignGenomicAnnotation = TRUE,genomicAnnotationPriority = c("5UTR","Promoter","3UTR","Exon","Intron","Downstream","Intergenic"))
saple1.loc.PeakAnno.DF.H2AZ <- as.data.frame(saple1.loc.PeakAnno) 
subset(saple1.loc.PeakAnno.DF.H2AZ,saple1.loc.PeakAnno.DF.H2AZ[,2]==20321001) 
#     seqnames    start      end width strand  V4  V5           V6       V7
#3539        1 20321001 20321999   999      * 406 117 1.153546e-63 2.625169
#               V8        annotation geneChr geneStart  geneEnd geneLength 
#3539 6.046514e-63 Distal Intergenic       1  20323036 20328009       4974 
#     geneStrand    geneId transcriptId distanceToTSS                      
#3539          1 AT1G54440  AT1G54440.2         -1037 

the track of this peak on igv is showed:

trac

GuangchuangYu commented 3 years ago

@MingLi-929 can you help testing this issue?

GuangchuangYu commented 3 years ago

fixed in v >=1.26.1