Biostar Beta. Not for public use.
Question: Calling SNPs using bcftools mpileup
0
Entering edit mode

Hi dears all, I am going to do SNP calling using bcftools mpileup for WGS data on fish. Before SNP calling I did this preprocessing as follow:

samtools sort -n F0-F06-10.bam F0-F06-10ns.bam

Add ms and MC tags for markdup to use later

samtools fixmate -m F0-F06-10ns.bam F0-F06-10-fixmate.bam

Markdup needs position order

samtools sort -o F0-F06-10PS.bam F0-F06-10-fixmate.bam

Finally mark duplicates

samtools markdup F0-F06-10PS.bam F0-F06-10-markdup.bam

Then using samtools faidx the fasta reference file got indexed and bcftools mpileup was used for SNP calling. For the last step (bcftools mpileup) I got involved with an error.

samtools faidx /home2/omid/alignment/ref/GCF_001858045.2_O_niloticus_UMD_NMBU_genomic.fa

bcftools mpileup -Ou -f GCF_001858045.2_O_niloticus_UMD_NMBU_genomic.fa.fai F0-F06-10-markdup.bam

error: [E::fai_build_core] Format error, unexpected "N" at line 1

I will be so grateful if any one can help me to solve this problem.

Cheers

ADD COMMENTlink 12 months ago jaafari.omid • 40 • updated 12 months ago finswimmer 11k
Entering edit mode
0

Hello,

the filenames used in the bcftools command doesn't match the filenames you present before. So please show us the exact command you've used.

Thanks!

fin swimmer

ADD REPLYlink 12 months ago
finswimmer
11k
Entering edit mode
0

Hello and Thanks for the comment. I have updated the commands to the exact ones that I used.

cheers,

ADD REPLYlink 12 months ago
jaafari.omid
• 40
Entering edit mode
0

You are right, When I changed the reference file to the original reference fasta file I did not see that error but the terminal got unnormal after running the command. I have attached a screen shot of my terminal. can you please let me know what can be the problem?

Cheers, Omidenter image description here

ADD REPLYlink 12 months ago
jaafari.omid
• 40
Entering edit mode
1

Hello jaafari.omid ,

please use the ADD REPLYbutton below a post you like to reply to.

In your bcftools you define that the output type is uncompressed bcf. But you neither define an output file (-o output) nor redirect it to a file (>output).

Also notice that the results of mpileup are not your final variants. This step just collect metrics about each covered site ,which are then used by bcftools call to find sites which are investigated for variants. So if you goal is finding variants, you have to do something like this:

$ bcftools mpileup -Ou -f ref.fa aln.bam | bcftools call -Ov -mv > output.vcf

fin swimmer

ADD REPLYlink 12 months ago
finswimmer
11k
Entering edit mode
0

I am really grateful for your great help. Now it worked fine. Cheers,

ADD REPLYlink 12 months ago
jaafari.omid
• 40
1
Entering edit mode

I'm quite sure that the reference file you are using in your bcftools command is the index file of your reference and not the reference file itself.

The index file must be located in the same location like the reference file. The file you specify with the -f parameter must be the reference genome file.

fin swimmer

ADD COMMENTlink 12 months ago finswimmer 11k

Login before adding your answer.

Similar Posts
Loading Similar Posts
Powered by the version 2.0