SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
How to call SNP using known SNP information nkwuji Bioinformatics 1 07-23-2013 09:16 PM
Ti/Tv ratio for known and novel SNP call sets deviates seqpig100 Bioinformatics 1 06-27-2011 09:32 AM
Using dbsnps or use samtools to call from 1000 genomes southan Bioinformatics 4 06-27-2011 02:00 AM
SOLiD SNP CALL Problem pr0t3us Bioinformatics 2 03-17-2011 07:56 AM
Any way to check if a SNP is present in 1000 genomes yuanzhi Bioinformatics 3 10-09-2010 08:23 AM

Reply
 
Thread Tools
Old 03-22-2011, 02:13 AM   #1
zhanglu295
Junior Member
 
Location: HongKong

Join Date: Mar 2011
Posts: 8
Unhappy 1000 genome SNP call

I want to call SNPs from the 1000 genome project original .bam files

Such as the data we can download from:

ftp://ftp.1000genomes.ebi.ac.uk/vol1...096/alignment/

I used samtools to call SNPs from the bam file named "HG00096.mapped.ILLUMINA.bwa.GBR.low_coverage.20101123.bam" in this folder.

but the pileup file is empty, the snp call procedure should be correct, is the bam file incompleted or other potential problems?

Please give me some advise.

Thank you very much

Eric
zhanglu295 is offline   Reply With Quote
Old 03-22-2011, 02:37 AM   #2
ulz_peter
Senior Member
 
Location: Graz, Austria

Join Date: Feb 2010
Posts: 219
Default

Could you give the exact command of you SNP calling step?
ulz_peter is offline   Reply With Quote
Old 03-22-2011, 03:08 AM   #3
ulz_peter
Senior Member
 
Location: Graz, Austria

Join Date: Feb 2010
Posts: 219
Default

Just a quick guess:
Samtools homepage recommends snp calling with the following command:

1. samtools mpileup -uf ref.fa aln1.bam aln2.bam | bcftools view -bvcg - > var.raw.bcf
2. bcftools view var.raw.bcf | vcfutils.pl varFilter -D100 > var.flt.vcf

One possible cause for your problem could lie in the reference sequence. 1000genomes use reference sequences called 1,2,3,4,5,6,etc for the chromosomes
If you use chr1,chr2,chr3 (UCSC-style) reference sequence, that could possibly lead to an empty pileup file.
ulz_peter is offline   Reply With Quote
Old 03-22-2011, 06:53 AM   #4
zhanglu295
Junior Member
 
Location: HongKong

Join Date: Mar 2011
Posts: 8
Default

Thank you very much for your help.

Because I just want to call SNP for each individual separately, so the command used is:

samtools pileup -vcf ref.fa aln.bam > raw.pileup

I think it should be the issue of reference sequence, which I used for snp call is UCSC hg19. Maybe I can test whether this problem still exist after changing the reference sequence.
zhanglu295 is offline   Reply With Quote
Old 03-22-2011, 08:11 AM   #5
colindaven
Senior Member
 
Location: Germany

Join Date: Oct 2008
Posts: 415
Default

I think you're right, the error is in the name of the reference sequence. Check to see the chromosome names in the alignment and reference are _exactly_ the same. We have had this problem frequently.
colindaven is offline   Reply With Quote
Old 03-30-2011, 09:57 AM   #6
heekim1
Junior Member
 
Location: Houston

Join Date: Dec 2010
Posts: 1
Default mpileup generates only chr1 in a bcf file

Quote:
Originally Posted by colindaven View Post
I think you're right, the error is in the name of the reference sequence. Check to see the chromosome names in the alignment and reference are _exactly_ the same. We have had this problem frequently.
Hi,
I had a similar situation.
I have a bam file containing read maps across entire chromosomes.
I did mpileup,

samtools mpileup -uf Homo_sapiens_assembly18.fasta s_4.merged.sorted.rmdup.bam | bcftools view -cvbg - > s_4.merged.sorted.rmdup.bam.raw.bcf

by using bcftools view s_4.merged.sorted.rmdup.bam.raw.bcf > s_4.merged.sorted.rmdup.bam.raw.vcf,

I realized that the .vcf file has SNP information on only chr1.
there is no chr2, chr3, ... or chrX

Do you have any idea about this problem?

As you suggested, i checked the chromosome name in both reference sequence and bam file. the name is identical.

Here is a part of my .fai derived from the reference sequence:
chrM 16571 6 50 51
chr1 247249719 16915 50 51
chr2 242951149 252211635 50 51
chr3 199501827 500021813 50 51
chr4 191273063 703513683 50 51
chr5 180857866 898612214 50 51
:
:

Here is a part of my .bam file:
HWI-EAS276_0022_FC70B81AAXX:4:91:4774:3269#0/1 0 chr1 12060 255 34M * 0 0 CTGGAGTGGAGTTTTCCTGTGGAGAGGAGCCATG BB=B=DD;DBBBBBDEEDD?ABABEDFEFGGG@G XA:i:0 MD:Z:34 NM:i:0

:
:
HWI-EAS276_0022_FC70B81AAXX:4:29:16345:8488#0/1 16 chr2 34145 255 34M * 0 0 TCATAGTTCTGCTAGACTTCTCTGAGGTGAGCTA @IGDIGIHDIIIIHIIIIGIIIIIIIIHGIIIII XA:i:0 MD:Z:34 NM:i:0
HW
:
:
HWI-EAS276_0022_FC70B81AAXX:4:63:1614:2851#0/1 0 chr3 90279 255 34M * 0 0 TTTTATAAGGGGCTTTTCCCCCTTTGCTCAGCAC IIIIHIDIIIHHIIIIIIIIIIDIIIIHBHHHGH XA:i:0 MD:Z:34 NM:i:0


Thank you

Hee
heekim1 is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 12:13 AM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO