htslib in pysam.fetch on S3 Bucket

[StephanHolgerD]

Hi, I want to report a potentially problematic behaviour using pysam.fetch on AWS S3 bucket infrastructure. Using the following pseudo code on a Bam file in a S3 Bucket will create requests without a defined end range.

Code

with pysam.AlignmentFile(bamfile_S3,filepath_index=baifile_S3) as f:
      for r in f.fetch(chrom,start,end):

Request

This kind of 'open' request results in high egress costs because aws logs the whole file after the start byte as delivered, even if you stop reading the data at the end of your fetch coordinates.

Compared to the requests from IGV on S3 data (low egress costs, only the exact byte range is logged)

Request

Initially I reported this here:

pysam-developers/pysam#1215

[daviesrob]

That looks unfortunate. We'll investigate and see if we can make these requests less open-ended. It may need a bit of rework to how our http requests work though so I can't be sure how long it will take.

[StephanHolgerD]

Thx, initially I used pysam.fetch which created the problematic open end requests.
After some debugging I switched to pysam.view (more or less a wrapper around the samtools view cmd), this created clean range requests.

[StephanHolgerD]

I checked the requests from Samtools view, they are also open and will create inflated egress costs.

Hi @StephanHolgerD , did pysam.view end up producing clean range requests? I wonder if you would be able to share some code for how you implemented it as I am looking for a similar functionality and haven't found an easy solution!

[StephanHolgerD]

@apena23 hey, sorry for the late reply. I did not get clean range requesta using the htslib. Currently working on an implementation in Python

[StephanHolgerD]

https://github.com/StephanHolgerD/BGZblocksfrombam

[mfz]

Has there been any progress on this issue? We are running into similar problems after having transferred our cram files from NFS to S3 storage.

The problem is that S3 actually seems to retrieve the object from the start to the end of the file, which makes it rather slow.

Some background:

We are doing pooled calling on 500 000 cram files in windows of 20 kbp. This means that every worker needs to retrieve 20 kbp slices from 500 000 cram files. Each slices is around 200 kB. As this is much smaller than the index file, we create index files that cover small regions only, and copy those to the worker nodes. We are using htslib from C++.

Using nfs storage we can load around 10 slices/thread/second. For S3 storage this goes down to around 1 slice/thread/second. The reason seems to be that htslib only specifies the start of the range, such that S3 (at least internally) retrieves the object from start of the range to the end.

The calling has already be performed on the cram files using NFS storage and currently we have the files on both kind of storage. The current test is done to see if it is possible to store the files (around 13 PB) on S3 storage only which is significantly cheaper.

[whitwham]

I have something that seems to work but it needs tidying and testing to get it production ready. Hopefully it will not be too much longer.

[mfz]

@whitwham , have you had the chance to look at this?

[whitwham]

@mfz I should have something to test soon. Are you willing and able to test a pull request when I make it?

[mfz]

@whitwham, yes, of course

[wiedenhoeft]

Am I right in assuming that this would affect any range query via tabix as well?

[whitwham]

@mfz There is a draft PR #1958. It is not quite done but there should be enough to test, if you are still willing.