Question

rtracklayer issue: Tracks don't show up

0

Entering edit mode

chrarnold • 0

@chrarnold-7603

Last seen 2.2 years ago

Germany

Hi there,

I am new to rtracklayer and I tried to reproduce the examples provided in the documentation. However, there seems to be an issue and I cannot reproduce it. If anyone could tell me what is wrong, it is appreciated...

This is my code, in a nutshell:

session = browserSession()
session$targets = targetTrack # targetTrack is GRanges object

The track is not found in my session:

grep ('targets', trackNames (session), ignore.case=T, v=T)

TS miRNA sites
"targetScanS"

Throws an error : Error in resolveTrackIndex(x, i) : Unknown track(s): targets
view <- browserView(session, targetTrack[1] , pack="targets")

sessionInfo()

R version 3.1.3 (2015-03-09)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 14.04.1 LTS

locale:
[1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C               LC_TIME=de_DE.UTF-8        LC_COLLATE=en_US.UTF-8     LC_MONETARY=de_DE.UTF-8
[6] LC_MESSAGES=en_US.UTF-8    LC_PAPER=de_DE.UTF-8       LC_NAME=C                  LC_ADDRESS=C               LC_TELEPHONE=C
[11] LC_MEASUREMENT=de_DE.UTF-8 LC_IDENTIFICATION=C

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

other attached packages:
[1] BiocInstaller_1.16.2      DESeq2_1.6.3              RcppArmadillo_0.4.650.1.1 Rcpp_0.11.5               Rsamtools_1.18.3
[6] Biostrings_2.34.1         XVector_0.6.0             BiocParallel_1.0.3        rtracklayer_1.26.3        GenomicRanges_1.18.4
[11] GenomeInfoDb_1.2.5        IRanges_2.0.1             S4Vectors_0.4.0           BiocGenerics_0.12.1

loaded via a namespace (and not attached):
[1] acepack_1.3-3.3         annotate_1.44.0         AnnotationDbi_1.28.2    base64enc_0.1-2         BatchJobs_1.6           BBmisc_1.9
[7] Biobase_2.26.0          bitops_1.0-6            brew_1.0-6              checkmate_1.5.2         cluster_2.0.1           codetools_0.2-11
[13] colorspace_1.2-6        DBI_0.3.1               digest_0.6.8            fail_1.2                foreach_1.4.2           foreign_0.8-63
[19] Formula_1.2-1           genefilter_1.48.1       geneplotter_1.44.0      GenomicAlignments_1.2.2 ggplot2_1.0.1           grid_3.1.3
[25] gtable_0.1.2            Hmisc_3.15-0            iterators_1.0.7         lattice_0.20-31         latticeExtra_0.6-26     locfit_1.5-9.1
[31] MASS_7.3-40             munsell_0.4.2           nnet_7.3-9              plyr_1.8.1              proto_0.3-10            RColorBrewer_1.1-2
[37] RCurl_1.95-4.5          reshape2_1.4.1          rpart_4.1-9             RSQLite_1.0.0           scales_0.2.4            sendmailR_1.2-1
[43] splines_3.1.3           stringr_0.6.2           survival_2.38-1         tools_3.1.3             XML_3.98-1.1            xtable_1.7-4
[49] zlibbioc_1.12.0

Thanks, Christian

rtracklayer bug session • 2.2k views

ADD COMMENT • link updated 9.6 years ago by Michael Lawrence ★ 11k • written 9.6 years ago by chrarnold • 0

0

Entering edit mode

I heard offline that it worked for the demo track, but not for your data. Please make your data available to me somehow, so that I can reproduce. Thanks!

ADD REPLY • link 9.6 years ago Michael Lawrence ★ 11k

0

Entering edit mode

A proxy? Interesting point, But upon testing, this does not seem to be the case. For example, using the cpneTrack data from the vignette, I can successfully load something to the session. However, with GRanges or RangedData that I constructed by myself (or even using the import function from the package that reads from a regular BED file) it simply does not work... Any other suggestions? I could send you example code that illustrates the problem

What puzzles me is the lack of any warning or error messages, it is just completely mysterious why the track is not there...

ADD REPLY • link 9.6 years ago chrarnold • 0

0

Entering edit mode

Please send the code/data.

ADD REPLY • link 9.6 years ago Michael Lawrence ★ 11k

0

Entering edit mode

Hi Michael, ok so I tested a bit further and it is just weird... It doesn't make any sense what I describe, but this is what I observe right now:

I have a data frame in BED format (100 rows only, see here: https://dl.dropboxusercontent.com/u/4001321/targets.bed)

session=NULL
session = browserSession()
restoredTrack <- import("targets.bed", asRangedData = FALSE)
track(session, "targetsRestored") =  restoredTrack
head(session[["targetsRestored"]],20)

Error in head(session[["targetsRestored"]], 20) :
error in evaluating the argument 'x' in selecting a method for function 'head': Error in track(ucscTableQuery(object, track = name, ...)) :
error in evaluating the argument 'object' in selecting a method for function 'track': Error in normArgTrack(track, trackids) : Unknown track: targetsRestored

On a related issue.

track(session, "targetsRestored") =  restoredTrack[1:100]

doesn't throw an error message, but

track(session, "targetsRestored") =  restoredTrack[1:10]

does:
Error in if (!all(unlist(lapply(split(object, strand), isStrandDisjoint)))) { :
missing value where TRUE/FALSE needed

Strangely, if I read in the file with just 10 entries, everything works (you can recapitulate that by deleting lines 11-100 in the file). However, if I add the 11th line, I get the error message as above.

A complete mystery to me... Again, with the cpneTrack from the vignette, everything is fine...

Any ideas?

THANKS!

ADD REPLY • link 9.6 years ago chrarnold • 0

score 0 · Answer 1 · 2015-04-16

0

Entering edit mode

Michael Lawrence ★ 11k

@michael-lawrence-3846

Last seen 3.0 years ago

United States

Well this is pretty surprising and embarrassing but it turns out that rtracklayer outputs invalid BED data when the strand is "*". Although not documented in the BED "spec", those "*" should become ".". I have fixed that in devel, or for now you can just set the strand to "+" before uploading. I also fixed the last problem you mentioned.

ADD COMMENT • link 9.6 years ago Michael Lawrence ★ 11k

0

Entering edit mode

Thanks Michael for looking into this.

The import function needs to have the strand in the sixth and last column, as expected. In order to use the export function, strand values must be "+", "-", or "*" and not ".", otherwise, I get a

Error in .local(x, ...) : strand values must be in '+' '-' '*'

message. Indeed, if I change the "*" strand to "+", it works!

However, "*" can then not be the only issue because when you do

export(cpneTrack, "cpneTrack.bed")

, you also have "*" as strand in the last column, and

session=NULL
session = browserSession()
export(cpneTrack, "cpneTrack.bed")
restoredTrack2 <- import("cpneTrack.bed", asRangedData = FALSE)
track(session, "cpneRestored") =  restoredTrack2
head(session[["cpneRestored"]])

works fine for me.

So, the questions I have are:

1) I don't understand why it works with https://dl.dropboxusercontent.com/u/4001321/cpneTrack.bed and not with https://dl.dropboxusercontent.com/u/4001321/targets.bed, the files seem to have an identical structure...

2) With https://dl.dropboxusercontent.com/u/4001321/targets.bed, all scores are 0 after adding it to the session (they are not zero after using import, so it must be the track function or the UCSC that converts them to 0). If I set all scores to 1, it works. If I divide all scores by their maximum so that all values are between 0 and 1, it does not work. However, it again works for cpneTrack, and I fail to understand why.

There must ba another requirement but what? According to http://genome.ucsc.edu/FAQ/FAQformat#format1, scores must be between 0 and 1000, which is true for my data (Range: 0.1418651-141.9687500). I think it would be very useful to check this and print an error or a warning for the user so he knows if the scores have been converted to 0 unless, of course, this is done silently by the UCSC.

Christian

ADD REPLY • link 9.6 years ago chrarnold • 0