Remote bigwigs?

[jmw86069]

Have you considered supporting remote bigwig/bigbed files?

I poked around a bit, thinking CURL meant remote file access, haha. Ah well.

[jayhesselberth]

Should work now, can you test after pak::pak("rnabioco/cpp11bigwig")? Thanks

[jmw86069]

Wow this is amazing! Or so I thought at first. Caveat: This may be user error on my side, although I did try for a bit, haha.

I can get the test data to load properly:

testurl <- paste0("https://raw.githubusercontent.com/rnabioco/cpp11bigwig/main/inst/extdata/test.bw")
bwc1 <- cpp11bigwig::read_bigwig(testurl,
   chrom="1", start=2, end=250)

# and whole bigwig works as well, granted it's a small test file
bwc <- cpp11bigwig::read_bigwig(testurl)

During compilation I did observe 'libcurl found: remote file access ENABLED and confirm it is version 0.1.3.9000. Your test data worked, so that's good.

I cannot get other data to load. It,immediately crashes my R session inside Positron with error
'R 4.5.2 exited unexpectedly: exit code -1073741819'

bwurl2 <- "https://manticore.niehs.nih.gov/ucscview/archerHubs/kinyamuChipseq/hg38/GR.CST_DEX_union.hg38.bw"
bwc2a <- cpp11bigwig::read_bigwig(bwurl2, chrom="chr1", start=100000, end=100000 + 100)

So that's weird, Positron runs some form of ark as an R interpreter - which could be a factor, so I tried gasp Rterm.exe (Windows) as well, just out of curiosity. Error message there is slightly different:
Error: Failed to open file: 'https://manticore.niehs.nih.gov/ucscview/archerHubs/kinyamuChipseq/hg38/GR.CST_DEX_union.hg38.bw'

For completeness, I used curl commandline and reviewed HEAD for your test file and my test file, looking for server differences, maybe certificate or something obtuse. If there is something important different, I didn't see it.

'Failed to open file' seems like the kind of error that could appear when something else has gone wrong?
I have not tried to use libBigWig as positive control. I'm not familiar with it, so I'd be hoping there exists some kind of command-line wrapper.

I ran bigWigInfo and it showed the bigwig version 4 (didn't know there existed such a thing).

version: 4
isCompressed: yes
isSwapped: 0
primaryDataSize: 435,828,861
primaryIndexSize: 3,664,656
zoomLevels: 10
chromCount: 24
basesCovered: 3,088,269,832
mean: 1.495330
min: 0.000000
max: 454.500000
std: 1.032049

By comparison, here is bigWigInfo for test.bw:

version: 4
isCompressed: yes
isSwapped: 0
primaryDataSize: 87
primaryIndexSize: 6,212
zoomLevels: 1
chromCount: 2
basesCovered: 154
mean: 1.766883
min: 0.100000
max: 2.000000
std: 0.358110

Again, I very much appreciate your time!
Aside: I noticed your comments elsewhere about tibble and data sizes, if you had ideas about other return types (numeric or Rle of some kind) I'd be open to that. The Rle or numeric return would be expected to have the exact length (when expanded to numeric vector) of the start-to-end range.

My use case would be to query coverage for exons of a given gene of interest. Around 10 to 50 exon chunks, each around 100 to maybe 20000 width? The Rle type shaves a bit of space compared to numeric. I don't mind iterating to query each exon, since there aren't all that many. Hopefully it's pretty quick! Haha.

My driving motivation is to find suitable replacement for rtracklayer because it has been pulling the entire bigwig file just to get the seqinfo(). Some bug where it apparently cannot determine the index of the file, so it pulls the whole file then interprets the content. Apparently a known bug they attribute to the Kent code they're including in the rtracklayer package? Idk. That bug also reduces the utility of having a range-indexed bigwig file, especially since it pulls the entire file multiple times before even issuing a ranged query.

[jayhesselberth]

Should be fixed now. I also added an Rle return type if you want to kick the tires and report back.

[jmw86069]

Confirmed that it works! Thank you so much!

Now, haha, might be my last question:
I am currently looping through multiple ranges, which works fine.
The input to read_bigwig() currently only accepts one range -- would there be any benefit in accepting multiple ranges?
I totally understand this may be out of scope, and perhaps little to no benefit in speed, I don't know.

I notice in the function docs it says if given multiple chromosomes, when as='Rle', it will return RleList.
What about if given multiple ranges on the same chromosome for example?

On my side, I will optimize by reducing contiguous ranges into one large range, then split into segments after receiving the result. I think most of my queries are probably going to be this way, so I have control over that.

Thanks so much for your time!