Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

No Locus contains enough reads to extract read length #33

Closed
seunghun23 opened this issue Aug 29, 2018 · 9 comments
Closed

No Locus contains enough reads to extract read length #33

seunghun23 opened this issue Aug 29, 2018 · 9 comments
Assignees

Comments

@seunghun23
Copy link

Hello,
Thank you for the nice tool!

I just installed GangSTR and tried to run it on a BAM file, but I got following error message
[GangSTR-1.4] ERROR: No Locus contains enough reads to extract read length. (Possible mismatch in chromosome names)

Thanks for the help in advance!

@nmmsv
Copy link
Collaborator

nmmsv commented Sep 2, 2018

Hello,
Can you please check if the ref name (chromosome names) match between alignment (bam) file, the regions (bed) file, and the reference genome (fasta/fa) that you're using? (i.e., if the convention is chr1:100-200 or just simply 1:100-200)
Also, is the bam file whole-genome or is it just limited to a region? How large is that region?
Thanks for reporting this issue!
Nima

@nmmsv nmmsv self-assigned this Sep 2, 2018
@eudesbarbosa
Copy link

Hi,

I had the same error: "[GangSTR-1.4] ERROR: No Locus contains enough reads to extract read length. (Possible mismatch in chromosome names)"
Everything seems to be working fine when I use the whole genome, but when I extract a region and run it again, I receive the error message.

BAM Region: 4:3076000-3078000
RepeatId:* HTT
TargetRegion: 4:3076604-3076660
OS: CentOS Linux 7 (Core)

Is there an argument I should pass?

Thanks for the help!

Eudes

@nmmsv
Copy link
Collaborator

nmmsv commented Nov 19, 2018

Hi,
In order to find the read length, GangSTR tries to find 10 consecutive reads with the same length in a 20k window around the repeat. Based on what you describe, this shouldn't be an issue in your run.
But you can go around this problem by manually setting the read length directly from command line (use --readlength ).
Please let me know if that didn't work.
Best,
Nima

@eudesbarbosa
Copy link

Hi,

Thanks for the quick reply. Even manually defining it (--readlength) I stumble on the same error =/
For the whole genome, it runs regardless of the parameter.

Process:
$ samtools view file.bam | head -n 10000 | awk '{print length($10)}' | sort -n | uniq -c #10000 150
$ samtools view -b file.bam "4:3076000-3078000" > region.bam
$ samtools view region.bam | awk '{print length($10)}' | sort -n | uniq -c #414 150
$ GangSTR --bam region.bam --ref genome.fa --regions HTT_hg19.bed --out test --readlength 150

Thanks for the help so far,
Cheers.
Eudes

@nmmsv
Copy link
Collaborator

nmmsv commented Nov 21, 2018

Hi Eudes,
Can you please check if the bed file that you're using HTT_hg19.bed has chromosome names formatted similar to your bam file? It seems like your bam file uses single numbers as chromosome IDs (4 as opposed to chr4).
If the chromosome IDs in bed and bam match, can you please run GangSTR in --verbose mode and send the log? That will help me narrow down the source of the problem.
Thanks,
Nima

@eudesbarbosa
Copy link

Hi, Nima.

The Fasta file I'm using has single numbers/character as identifiers, so it was easier to just modify the repeat definitions.

I attached the log for the whole bam (whole_bam.log) and restricted region (region_bam.log)

$ GangSTR --bam whole.bam --ref genome.fa --regions HTT_hg19.bed --out test_whole --verbose
$ GangSTR --bam region.bam --ref genome.fa --regions HTT_hg19.bed --out test _region --readlength 150 --verbose

Let me know if you need something else,
Cheers.
Eudes

@nmmsv
Copy link
Collaborator

nmmsv commented Nov 27, 2018

Thank you for sending the logs.
It seems to me that the error is in the computation of coverage and mean and standard deviation of insert size. This is expected as the number of reads in the smaller bam file is probably not enough to pass the threshold for computing these parameters. The threshold is set to 300 in the released version, and it is not exposed to the user. While it is possible to lower this threshold, I think a better solution is to supply coverage and mean and standard deviation of insert size as parameters in input. You can use the values computed in the whole genome log.
Parameters to set:
--insertmean <float>: Fragment length mean. Default: -1
--insertsdev <float>: Fragment length standard deviation. Default: -1
--coverage <float>: Average coverage. must be set for exome/targeted data. Default: -1

@eudesbarbosa
Copy link

Hi, Nima.

It worked! Thanks a lot for the help =)

--
Eudes

@nmmsv
Copy link
Collaborator

nmmsv commented Nov 27, 2018

Perfect!
No problem.

@nmmsv nmmsv closed this as completed Nov 27, 2018
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants