12
votes
Difference between BWA-backtrack and BWA-MEM
TL;DR:
BWA-backtrack is based on backtracking.
This approach is appropriate only when the dissimilarity between the reads and the reference is low,
or when you want to find all best hits or enumerate ...
- 1,899
11
votes
Random access on a FASTQ file
Arbitrary record access in constant time
To get a random record in constant time, it is sufficient to get an arbitrary record in constant time.
I have two solutions here: One with ...
- 2,236
9
votes
Why does cutadapt remove low quality bases from the ends of reads only?
If you check the read QC statistics of an Illumina run in e.g. fastQC, you will see that at the end of the read the quality decreases. This is because of exhaustion of chemicals at the end of the run. ...
- 3,571
9
votes
Difference between BWA-backtrack and BWA-MEM
To quote the Introduction to BWA on sourceforge:
BWA is a software package for mapping low-divergent sequences against a large reference genome, such as the ...
- 663
9
votes
Accepted
Extracting all reads from bam file which match read IDs in another file
It is still slow but grep has a -f option to take in a file
samtools view inbam.bam | grep -f read_names.txt > read_locs.txt
- 2,574
8
votes
Filtering step for read counts data
If you have your counts in a data.frame called counts, something like this might work:
...
- 3,571
7
votes
Accepted
Random access on a FASTQ file
As wkretzsch suggested this was worthy of an actual answer, I feel the obvious solution is missing here; index the FASTQ.
Index it
As much as I typically hesitate ...
- 772
7
votes
Accepted
What is local realignment and what is the problem it solves?
SNPs are likely to be created and InDels are likely to be missed. Suppose you have a read, ACTGACTGACTGTAC and you align it to a reference sequence ...
- 19.6k
7
votes
Accepted
How to remove all BAM read groups from all reads (not just the header)?
You might have to manually strip those auxiliary tags off:
...
- 19.6k
7
votes
Accepted
Why does cutadapt remove low quality bases from the ends of reads only?
I am commenting on this part:
The algorithm only removes low quality bases from the end until it reaches a good quality base. If there is a bad quality base beyond that, it is not trimmed.
...
- 6,485
6
votes
Counting repeated kmers sequences that match at least x % of reads sequence
The following javascript does what you want. You need node.js to run it. It should be easy to translate the code to Python. I have not carefully tested it. Use with caution.
EDIT (response to new ...
- 6,485
6
votes
Accepted
How to check whether all BAM read contain defined read groups?
According to the man page, running samtools stats --split RG <file1.bam> should produce summary statistics separated by read group. If it doesn't produce a ...
- 13.8k
5
votes
Random access on a FASTQ file
You could shuffle the FASTQ once and then read sequences off the top of the file as you need them:
...
- 2,236
5
votes
Accepted
What is the right way of calculating a Phred score by hand?
Not sure if my explanation is any good but let's try...
First, let's convince ourselves that converting phred to probability is the right thing to do (i.e. your Way 2).
Imagine a read of length 10000 ...
- 689
4
votes
Accepted
Why are there missing calls in a VCF file from exome sequencing?
Missing variant calls due to lack of coverage shouldn't happen in the targeted capture region and I'd think most of these would come from off-target regions where some samples had reads mapped. I'd ...
- 393
4
votes
How to check whether all BAM read contain defined read groups?
using samjdk and invoking the function getReadGroup()
getReadGroup() returns The SAMReadGroupRecord from the SAMFileHeader for this SAMRecord, or null if 1) this record has no RG tag, or 2) the ...
- 1,501
4
votes
Accepted
Select top 100 genes ranked by variance in read counts
I assume by "reverse sort of variance" you mean "highest variance". Assuming you made a matrix out of that (set the row names to the first column and then remove it) and called it ...
- 19.6k
4
votes
What is and how to detect a dephased read
Phasing in Illumina data means that the read you get is shifted one base before or after what the sequence really should be. They must be expecting some constant 6 bases at the end of each read, and ...
- 1,882
4
votes
Extract reads from bam files by their @RG
Typically I use samtools for operations like this. Specifically I use samtools view with either -r or ...
- 2,574
3
votes
While opening a Bam file with the SeqAn library I get 'seqan::FileOpenError'
Remove the getAbsolutePath call — your path is already absolute. getAbsolutePath garbles it, as you can see in the error message....
- 4,845
3
votes
Accepted
Scaling by linear regression against the number of reads
I don't know if this question has been solved already, but what they try to do is equalize the depth of sequencing for each cell. Therefore, they scale for the total number of reads. If you regress ...
- 200
3
votes
Filtering step for read counts data
While this answers explains how to do it I want to address when and why and which thresholds to do it.
Filtering the genes with low counts is usually done because the counts are not reliable it would ...
- 4,693
3
votes
Why are there missing calls in a VCF file from exome sequencing?
I've just been generating data like this, so can tell you about why/how missing calls are created in my dataset. There are two main reasons:
1. Sequencing failure
When reads don't map across the ...
- 13.8k
3
votes
Random access on a FASTQ file
One possibility is to:
reformat the data such that each record is a single line containing the read description, bases, and quality scores
pad out each record to a maximum length in each field such ...
- 31
3
votes
Reads mapped to exonic, intronic and intergenic regions
As @DevonRyan mentioned, it's very likely that those samples were degraded, which is good justification for excluding them from subsequent analysis.
- 5,060
3
votes
If a gene is expressed at a level of 1/1200 compared to the average gene, how is probability 50:50 that we have a read mapped to it?
I actually disagree with that.. I guess I write this down as a discussion.
The expected value for an average gene is 1200. For this gene at 1/1200 expression, you expect 1 read.
However, because of ...
- 1,678
3
votes
If fastp output is not a good measure of FASTQ correctness, what is?
On a regular basis, using FastQC is quite enough as to assess the quality of FastQ-formatted data. It gives you plenty of details as can be seen in the screenshots you shared.
As others have commented ...
- 827
2
votes
Random access on a FASTQ file
I wrote a tool called sample that you can use to do random sampling without reading the entire file into memory.
It can be used where GNU ...
- 3,125
2
votes
Random access on a FASTQ file
One of the most thorough treatments of this question (or a similar question: grabbing a random subset of reads) was given by Jared Simpson in a blog post a few years ago. http://simpsonlab.github.io/...
- 5,060
2
votes
Genome assembly from error-prone reads
For short reads, the typical and the most widely used solution is to correct away sequencing errors before assembly. You can correct errors with k-mer spectrum, a trie or multi-alignment. There are ...
- 6,485
Only top scored, non community-wiki answers of a minimum length are eligible