Skip to main content
Bioinformatics

Return to Answer

updated to account for random record sampling
Source Link
gringer ♦
gringer ♦
  • 15.1k
  • 5
  • 24
  • 83

Properly parsing FASTQ is a nightmare. It's not possible to start at a random position inTo emphasise the file and guarantee that a seek forwardissue (or backwards) to an '@' symbol will arrive atas outlined in the start of a FASTQ header. That's a problem if log(1/n) time complexityupdate), consider an input file with one record that is desired5 million bases long, because the problem meansand one records that it's not possible to read through the entire file to get a sampleis 100 bases long. You want an equal probability of sequencesselecting either of these two records. Any random seek methods will overwhelmingly pick out the long record.

IfI expect that indexing the allowed FASTQ format were more strict, i.e. 4 lines perlocations of record: [header, sequence, header, quality], then it would be possible to do this: seek to a random location in the file, seek to find starts is really the next line that contains only a + (preceded by a DNA sequence to exclude + in the quality string), seek backwards two linesworkable option here, and there'sparticularly if multi-line records are possible. Create an index file containing the start locations of a FASTQeach record (as identically sized integers, e.

If the file were indexedg. 64-bit), then it would be possible to do the shuffling onsample from the index and use thatfile (which is an identical length for each record) to select reads fromfetch the file..start location. butI'd envisage that this file would also require a complexity of over nonly contain the start locations; any additional metadata (because at least indexing is neededincluding sequence name) would require seeking in the original file.

Once the indexing is done, the file can be compressed with bgzip, with specific offsets retrieved using the samtools view-b has theand -s <float> option for subsampling aligned reads, but that's in BAM formatoptions. However, rather than FASTQ format (andI expect that compression would still go through the whole file before arriving at the last read)not be particularly efficient if multiple random records were desired.

Properly parsing FASTQ is a nightmare. It's not possible to start at a random position in the file and guarantee that a seek forward (or backwards) to an '@' symbol will arrive at the start of a FASTQ header. That's a problem if log(n) time complexity is desired, because the problem means that it's not possible to read through the entire file to get a sample of sequences.

If the allowed FASTQ format were more strict, i.e. 4 lines per record: [header, sequence, header, quality], then it would be possible to do this: seek to a random location in the file, seek to find the next line that contains only a + (preceded by a DNA sequence to exclude + in the quality string), seek backwards two lines, and there's the start of a FASTQ record.

If the file were indexed, then it would be possible to do the shuffling on the index and use that to select reads from the file... but that would also require a complexity of over n (because at least indexing is needed).

samtools view has the -s <float> option for subsampling aligned reads, but that's in BAM format, rather than FASTQ format (and would still go through the whole file before arriving at the last read).

To emphasise the issue (as outlined in the 1/n update), consider an input file with one record that is 5 million bases long, and one records that is 100 bases long. You want an equal probability of selecting either of these two records. Any random seek methods will overwhelmingly pick out the long record.

I expect that indexing the locations of record starts is really the only workable option here, particularly if multi-line records are possible. Create an index file containing the start locations of each record (as identically sized integers, e.g. 64-bit), then sample from the index file (which is an identical length for each record) to fetch the start location. I'd envisage that this file would only contain the start locations; any additional metadata (including sequence name) would require seeking in the original file.

Once the indexing is done, the file can be compressed with bgzip, with specific offsets retrieved using the -b and -s options. However, I expect that compression would not be particularly efficient if multiple random records were desired.

Source Link
gringer ♦
gringer ♦
  • 15.1k
  • 5
  • 24
  • 83

Properly parsing FASTQ is a nightmare. It's not possible to start at a random position in the file and guarantee that a seek forward (or backwards) to an '@' symbol will arrive at the start of a FASTQ header. That's a problem if log(n) time complexity is desired, because the problem means that it's not possible to read through the entire file to get a sample of sequences.

If the allowed FASTQ format were more strict, i.e. 4 lines per record: [header, sequence, header, quality], then it would be possible to do this: seek to a random location in the file, seek to find the next line that contains only a + (preceded by a DNA sequence to exclude + in the quality string), seek backwards two lines, and there's the start of a FASTQ record.

If the file were indexed, then it would be possible to do the shuffling on the index and use that to select reads from the file... but that would also require a complexity of over n (because at least indexing is needed).

samtools view has the -s <float> option for subsampling aligned reads, but that's in BAM format, rather than FASTQ format (and would still go through the whole file before arriving at the last read).