[Bioc-sig-seq] Questions about coverage visualization in wig file output
Chen-Yi Chen
ChenYiChen at lbl.gov
Wed Mar 24 00:46:41 CET 2010
Hi all,
Recently we found the greatness of this bioconductor ht-seq pacakage and
we would like to do some analysis in this package. We have aligned our
ChIP and input data by Bowtie aligner and try to output wig file for the
coverage, so we are able to use UCSC genome browser to visualize them. By
searching the documentation and work flow online, I found this great
thread in bioc-sig-sequencing archive:
https://stat.ethz.ch/pipermail/bioc-sig-sequencing/2009-September/000632.html
We basically followed the work flow in the above thread, with some
slightly modifications:
ChIPread = readAligned("./", pattern="0_Rad_H3K4_ChIP.out", type="Bowtie")
ChIP = as(ChIPread, "GenomeData")
inputread = readAligned('./', pattern='0_Rad_input.out', type='Bowtie')
input = as(inputread, "GenomeData")
gdlist = GenomeDataList(list(input=input, ChIP=ChIP))
input1 = extendReads(gdlist$input$chr1, seqLen=200)
ChIP1 = extendReads(gdlist$ChIP$chr1, seqLen=200)
library(BSgenome.Hsapiens.UCSC.hg18)
human.chromlens = seqlengths(Hsapiens)
input1.cov = coverage(input1, width=human.chromlens['chr1'])
ChIP1.cov = coverage(ChIP1, width=human.chromlens['chr1'])
inputrange1 = IRanges(start=start(input1.cov), end=end(input1.cov))
score = runValue(input1.cov)
chr1 = rep("chr1", length(inputrange1))#
export(GenomicData(inputrange1, score, chrom=chr1), "input_test.wig")
We successfully got a wig file output, however, when we try to visualize
it on UCSC genome browser, the peaks won't display correctly.
I took a glance on the output wig file, it has the format like the following:
track name="R Track" type=wiggle_0
variableStep chrom=chr1 span=1
251 7
395 7
14041 3
78204 1
...
variableStep chrom=chr1 span=2
12 2
212 9
346 7
...
variableStep chrom=chr1 span=3
343 8
404 7
...
...
fixedStep chrom=chr1 start=117727703 step=129030611 span=1990
0
0
variableStep chrom=chr1 span=1991
...
Apparently the output wig file format we got is different than what's
mentioned in the thread, and I believed those different variable steps and
fixed steps are confusing the UCSC genome browser.
Is there a way to output the wig file in the format without all those
steps, for example,
track name="R Track" type=wiggle_0
chrI 0 2 169
chrI 2 4 176
chrI 4 5 178
chrI 5 7 179
...
which mentioned in the thread?
Thanks for your help!
-Charlie-
More information about the Bioc-sig-sequencing
mailing list