[Bioc-sig-seq] WIG from coverage

Seidel, Christopher CWS at stowers.org
Wed Sep 30 01:28:24 CEST 2009


In an earlier thread: BED to WIG format conversion (Thu 9/17/2009 4:50 PM) Michael Lawrence mentioned:

Ignoring efficiency, it's pretty straight-forward. In the development version of IRanges, one can convert an Rle (or RleList for multiple chromosomes) as output to coverage() to a RangedData. In the development version of rtracklayer, this is done implicitly:

> export(cov, "coverage.wig")

as a way to create a wig file from a coverage() result. When I try this, my wig file looks like the following:

track name="R Track" type=wiggle_0
.       0       2       168
.       2       4       175
.       4       5       177
.       5       7       178
.       7       8       182
.       8       12      185
.       12      14      188
.       14      15      191

This isn't compatible with any of the UCSC wig formats that I'm aware of (is it?). Shouldn't rtracklayer be putting out a variableStep format for this? Am I doing something wrong?

By the way, here is my work flow for going from solexa export to coverage:

# read in some IP data
aln <- readAligned(spath, pattern="s_1_export.txt", type="SolexaExport", filter=filt)
# convert to GenomeData object
s1 <- as(aln,"GenomeData")

# read in some WCE data
aln <- readAligned(spath, pattern="s_3_export.txt", type="SolexaExport", filter=filt)
# convert to GenomeData object
s3 <- as(aln,"GenomeData")

# create GenomeDataList from ip and wce
ypd <- GenomeDataList(list(ip=s1, wce=s3))

# extend reads
ip1 <- extendReads(ypd$ip$chrI, seqLen=150)
wce1 <- extendReads(ypd$wce$chrI, seqLen=150)

# find coverage for chromosome 1
ip1.cov <- coverage(ip1, width=230208)
wce1.cov <- coverage(wce1, width=230208)

# make a track
export(ip1.cov, "ip1_coverage.wig")

> sessionInfo()
R version 2.10.0 Under development (unstable) (2009-09-03 r49555)

[1] C

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] rtracklayer_1.5.13 RCurl_1.2-0        bitops_1.0-4.1     chipseq_0.1.25
[5] ShortRead_1.3.34   lattice_0.17-25    BSgenome_1.13.11   Biostrings_2.13.39
[9] IRanges_1.3.75

loaded via a namespace (and not attached):
[1] Biobase_2.5.6 XML_2.6-0     grid_2.10.0   hwriter_1.1


Chris Seidel, Ph.D.

More information about the Bioc-sig-sequencing mailing list