SAMTOOLS problem

[lapti]

Dear all,

I have two different problems related to samtools.

1. When I try to create a pileup using pileup option in samtools and sorted bam file, I get segmentation fault. I'm using following command.

samtools pileup -cv -f all_chr.fa bamfile.sorted.bam

2. when I use pysam to find different variations at particular coordinates of the human genome I get following errors.

Traceback (most recent call last):
File "twobitReader.py", line 5, in <module>
samfile = pysam.Samfile("Brain_DNA.sorted.bam","rb")
File "csamtools.pyx", line 266, in csamtools.Samfile.__cinit__ (pysam/csamtools.c:2759)
File "csamtools.pyx", line 376, in csamtools.Samfile._open (pysam/csamtools.c:3871)
ValueError: invalid BAM binary header (is this a BAM file?)

This problem is only related to files which are more than 10 GBs in size. I have tested with smaller files and it works.

--
regards,
Anmol

[michalkovac]

Dear All,

Since this thread has been created for the problems regarding SAMtools I'd like to take advantage of it and add mine, which just popped up today.

As of 15th November, the manual entry for SAMtools m/pileup command has been updated by replacing following lines

samtools pileup -vcf hsRef.fa - > var-X.raw.txt
samtools.pl varFilter -D100 var-X.raw.txt > var-X.flt.txt
awk '($3=="*"&&$6>=50)||($3!="*"&&$6>=20)' var-X.flt.txt > var-X.final.txt

with

samtools mpileup -ugf ref.fa aln.sorted.bam | bcftools view -bvcg - > var.raw.bcf
bcftools view var.raw.bcf | vcfutils.pl varFilter -D 100 > var.flt.vcf

Mine problem with the latter pipeline lays in not getting the whole source parsed as it always stops at certain position (file size?) regardless what source file I use.

This is what I get:

awk '{print $1}' normal.var.final.vcf | sort | uniq -c
46343 chr1
30048 chr10
1 chr10_random
30111 chr11
6 chr11_random
27205 chr12
19945 chr13
114 chr13_random

and without using mpileup-bcf command:

awk '{print $1}' ./Finals/2007.02142B.normal.var.final.ST.plp | sort | uniq -c
22053 chr1
14314 chr10
7 chr10_random
13908 chr11
5 chr11_random
12635 chr12
8127 chr13
96 chr13_random
7786 chr14
8239 chr15
8 chr15_random
10247 chr16
1 chr16_random
9566 chr17
308 chr17_random
6534 chr18
5 chr18_random
8030 chr19
86 chr19_random
16 chr1_random
21137 chr2
7047 chr20
5105 chr21
292 chr21_random
4800 chr22
81 chr22_random
3 chr2_random
16190 chr3
3 chr3_random
16118 chr4
124 chr4_random
13458 chr5
1 chr5_h2_hap1
13203 chr6
45 chr6_cox_hap1
56 chr6_qbl_hap2
29 chr6_random
13484 chr7
293 chr7_random
11937 chr8
23 chr8_random
11056 chr9
506 chr9_random
3556 chrX
2 chrX_random
3110 chrY

Could anyone hint why the whole file has not been parsed? Or are there any ways to transform plp format into vcf?

Many thanks for help.
Michal

p.s. segmentation fault pops up after running samtools mpileup -ugf ref.fa aln.sorted.bam (and not being piped into bcftools) and the output size exceeds 10gb.