Chipseq peak analysis using ChIPseeker

Question

[Mike]

Hi ,

I am using ChIPseeker tool for my Chipseq data analysis, But my heatmap of ChiP binding to TSS regions is almost blank “attached image”,

I used macs2 output file.

I used following commands:

txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene


peak <- readPeakFile("chip_vs_input_summits.bed")
covplot(peak, weightCol="V5")



##Profile of ChIP peaks binding to TSS regions
promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
tagMatrix <- getTagMatrix(peak, windows=promoter)
tagHeatmap(tagMatrix, xlim=c(-3000, 3000), color="red")

[Mike]

Sorry in above case I used peak summit file, now I tried with *peaks.narrowPeak macs2 output files,

awk 'BEGIN {OFS="\t"} {print $1,$2,$3,$4,$7,$9}'  chip_vs_input_peaks.narrowPeak >chip_vs_input.bed

when I used this file it gives some warnings:

 promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
tagMatrix <- getTagMatrix(peak, windows=promoter)
Warning message:
In .Seqinfo.mergexy(x, y) :
  Each of the 2 combined objects has sequence levels not in the other:
  - in 'x': chr20, chr21, chr22, chrY, chrM, chr1_gl000191_random, chr1_gl000192_random, chr4_ctg9_hap1, chr4_gl000193_random, chr4_gl000194_random, chr6_apd_hap1, chr6_cox_hap2, chr6_dbb_hap3, chr6_mann_hap4, chr6_mcf_hap5, chr6_qbl_hap6, chr6_ssto_hap7, chr7_gl000195_random, chr8_gl000196_random, chr8_gl000197_random, chr9_gl000198_random, chr9_gl000199_random, chr9_gl000200_random, chr9_gl000201_random, chr11_gl000202_random, chr17_ctg5_hap1, chr17_gl000203_random, chr17_gl000204_random, chr17_gl000205_random, chr17_gl000206_random, chr18_gl000207_random, chr19_gl000208_random, chr19_gl000209_random, chr21_gl000210_random, chrUn_gl000211, chrUn_gl000212, chrUn_gl000213, chrUn_gl000214, chrUn_gl000215, chrUn_gl000216, chrUn_gl000217, chrUn_gl000218, chrUn_gl000219, chrUn_gl000220, chrUn_gl000221, chrUn_gl000222, chrUn_gl000223, chrUn_gl000224, chrUn_gl000225, chrUn_gl000226, chrUn_gl000227, chrUn_gl000228, chrUn_gl0 [... truncated]

[Guangchuang Yu]

maybe send me your bed file and source code for me to reproduce the issue.

[Mike]

Thanks,

this is my macs2 command:

macs2 callpeak -t chip.bam -c input.bam -f BAM -B -g hs -n rep1

and input bed file:

http://www.filedropper.com/showdownload.php/sample1peaks

or

http://s000.tinyupload.com/index.php?file_id=66463745037323321210

Answer 1

[Guangchuang Yu]

Warning message:
In .Seqinfo.mergexy(x, y) :
  Each of the 2 combined objects has sequence levels not in the other:

The warning message is OK, as you can check:

> seqlevels(txdb)
 [1] "chr1"                  "chr2"                  "chr3"
 [4] "chr4"                  "chr5"                  "chr6"
 [7] "chr7"                  "chr8"                  "chr9"
[10] "chr10"                 "chr11"                 "chr12"
[13] "chr13"                 "chr14"                 "chr15"
[16] "chr16"                 "chr17"                 "chr18"
[19] "chr19"                 "chr20"                 "chr21"
[22] "chr22"                 "chrX"                  "chrY"
[25] "chrM"                  "chr1_gl000191_random"  "chr1_gl000192_random"
[28] "chr4_ctg9_hap1"        "chr4_gl000193_random"  "chr4_gl000194_random"
[31] "chr6_apd_hap1"         "chr6_cox_hap2"         "chr6_dbb_hap3"
[34] "chr6_mann_hap4"        "chr6_mcf_hap5"         "chr6_qbl_hap6"
[37] "chr6_ssto_hap7"        "chr7_gl000195_random"  "chr8_gl000196_random"
[40] "chr8_gl000197_random"  "chr9_gl000198_random"  "chr9_gl000199_random"
[43] "chr9_gl000200_random"  "chr9_gl000201_random"  "chr11_gl000202_random"
[46] "chr17_ctg5_hap1"       "chr17_gl000203_random" "chr17_gl000204_random"
[49] "chr17_gl000205_random" "chr17_gl000206_random" "chr18_gl000207_random"
[52] "chr19_gl000208_random" "chr19_gl000209_random" "chr21_gl000210_random"
[55] "chrUn_gl000211"        "chrUn_gl000212"        "chrUn_gl000213"
[58] "chrUn_gl000214"        "chrUn_gl000215"        "chrUn_gl000216"
[61] "chrUn_gl000217"        "chrUn_gl000218"        "chrUn_gl000219"
[64] "chrUn_gl000220"        "chrUn_gl000221"        "chrUn_gl000222"
[67] "chrUn_gl000223"        "chrUn_gl000224"        "chrUn_gl000225"
[70] "chrUn_gl000226"        "chrUn_gl000227"        "chrUn_gl000228"
[73] "chrUn_gl000229"        "chrUn_gl000230"        "chrUn_gl000231"
[76] "chrUn_gl000232"        "chrUn_gl000233"        "chrUn_gl000234"
[79] "chrUn_gl000235"        "chrUn_gl000236"        "chrUn_gl000237"
[82] "chrUn_gl000238"        "chrUn_gl000239"        "chrUn_gl000240"
[85] "chrUn_gl000241"        "chrUn_gl000242"        "chrUn_gl000243"
[88] "chrUn_gl000244"        "chrUn_gl000245"        "chrUn_gl000246"
[91] "chrUn_gl000247"        "chrUn_gl000248"        "chrUn_gl000249"
> seqlevels(peak)
 [1] "chr1"         "chr1_random"  "chr10"        "chr11"        "chr12"
 [6] "chr13"        "chr13_random" "chr14"        "chr15"        "chr16"
[11] "chr17"        "chr18"        "chr19"        "chr2"         "chr3"
[16] "chr4"         "chr5"         "chr6"         "chr7"         "chr8"
[21] "chr8_random"  "chr9"         "chrUn_random" "chrX"         "chrY"
[26] "chrY_random"

I can produce the tagHeatmap using the sample bed file you provided.

[Mike]

Thanks you so much Guangchuang Yu,

what are chr1/13/8_random and chrUn_random, why some chr (chr20, ch21, chr22 ) are missing in my data.

Thanks.

[Mike]

sorry one more problem,

promoter <- getPromoters(TxDb=txdb, upstream=3000, downstream=3000)
tagMatrixList <- lapply(files, getTagMatrix, windows=promoter)
plotAvgProf(tagMatrixList, xlim=c(-3000, 3000))

I get this error

Error in plotAvgProf.internal(tagMatrix, xlim = xlim, xlab = xlab, ylab = ylab,  : 
  tagMatrix should be a named list...

how to give list of bed files in plotAvgProf

but below command is works:

 tagHeatmap(tagMatrixList, xlim=c(-3000, 3000), color=NULL)

And when I use :

peakAnnoList <- lapply(files, annotatePeak, TxDb=txdb, tssRegion=c(-3000, 3000), verbose=FALSE)
 plotAnnoBar(peakAnnoList)

error:

Error in eval(expr, envir, enclos) : object '.id' not found

[Guangchuang Yu]

just named your files.

names(files) <- basename(files)

or

names(files) <- paste0('file_', seq_along(files))

something like this.

[Guangchuang Yu]

I don't know either.