samtools mpileup error
0
0
Entering edit mode
4.3 years ago
teresayahu • 0

I am using samtools function mpileup to get reads in sites. The command I used is:

samtools mpileup -l hg19.position -f hg19.fa -Q 20 -q 20 -I file.bam

in the command, hg19.position is position-based site file.

When running these is an error printed on screen: ...

[E::faidx_fetch_seq] The sequence "1" not found
[E::faidx_fetch_seq] The sequence "1" not found
[E::faidx_fetch_seq] The sequence "1" not found
[E::faidx_fetch_seq] The sequence "1" not found
[E::faidx_fetch_seq] The sequence "1" not found
[E::faidx_fetch_seq] The sequence "1" not found

...

May I ask the reason of this error and how could I fix it? Thanks! Ya

Assembly software error next-gen sequence • 4.5k views
ADD COMMENT
1
Entering edit mode

Have you indexed the hg19.fa file using the same version of samtools as per your command above and using the samtools faidx command?

ADD REPLY
0
Entering edit mode
  • I'm getting the same error as the OP "[E::faidx_fetch_seq] The sequence "1" not found"

I tried the index command you suggested,

samtools faidx hg19.fa

After pressing Enter there's a few second pause, followed by a new line:

ben@DESKTOP-EMRLM4V:~$ samtools faidx hg19.fa

ben@DESKTOP-EMRLM4V:~$

... Does this mean the fa file is now indexed?

If so, why do I keep getting the same error, "[E::faidx_fetch_seq] The sequence "1" not found", when executing my initial command:

bcftools mpileup -Ou -f hg19.fa I19139.hg19.bam | bcftools call -vmO z -o Nubians.vcf.gz

Thanks for any help,

ADD REPLY
0
Entering edit mode

For added context, here's a view of some lines from the BAM file using the samtools view H command:

@HD     VN:1.5  GO:none SO:coordinate
@SQ     SN:1    LN:249250621
@SQ     SN:2    LN:243199373
@SQ     SN:3    LN:198022430
@SQ     SN:4    LN:191154276
@SQ     SN:5    LN:180915260
@SQ     SN:6    LN:171115067
@SQ     SN:7    LN:159138663
@SQ     SN:8    LN:146364022
@SQ     SN:9    LN:141213431
@SQ     SN:10   LN:135534747
@SQ     SN:11   LN:135006516
@SQ     SN:12   LN:133851895
@SQ     SN:13   LN:115169878
@SQ     SN:14   LN:107349540
@SQ     SN:15   LN:102531392
@SQ     SN:16   LN:90354753
@SQ     SN:17   LN:81195210
@SQ     SN:18   LN:78077248
@SQ     SN:19   LN:59128983
@SQ     SN:20   LN:63025520
@SQ     SN:21   LN:48129895
@SQ     SN:22   LN:51304566
@SQ     SN:X    LN:155270560
@SQ     SN:Y    LN:59373566
@SQ     SN:MT   LN:16569
ADD REPLY
0
Entering edit mode

I would say you're bam or your reference using the 'chr1' notation, not just '1'

ADD REPLY
0
Entering edit mode

Thanks! both the bam and reference files are using the '1'. Is there a way that I can change bam notation from '1' to 'chr1'?

ADD REPLY
0
Entering edit mode

Please use ADD COMMENT/ADD REPLY when responding to existing posts to keep threads logically organized.

ADD REPLY
0
Entering edit mode

what about hg19.position ??

ADD REPLY

Login before adding your answer.

Traffic: 2905 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6