Expression problems in bcftools view

[pishotta]

Commands like:
"bcftools view -i 'ref="A"' p3.vcf.gz" work fine, yielding lines like

>> chr1    883619  .       A       .       999     .       DP=44;MQSB=1;MQ0F=0...
>> chr1    883622  .       A       .       999     .       DP=45;MQSB=1;MQ0F=0...
>> chr1    883625  .       A       G       141     .       DP=45;VDB=6.66888e-18...
>> chr1    883627  .       A       .       999     .       DP=46;MQSB=1;MQ0F=0;AN=2;DP4=1,45,0,0...
>> chr1    883629  .       A       .       999     .       DP=46;MQSB=1;MQ0F=0;AN=2;DP4=1,44,0,0...

"bcftools view -i 'alt="G"' p3.vcf.gz" yields every line with either a "G" or a "." in 
the ALT column.  Not the behavior I expect/want.

>> chr1    908306  .       G       .       999     .      DP=143;MQSB=0.20004...
>> chr1    908309  .       G       .       27.283  .      DP=139;VDB=0.478589;SGB=-0.651104...
>> chr1    908310  .       C       G       7.91518 .      DP=112...
>> chr1    908311  .       C       .       999     .      DP=139;MQSB=0.774446...
>> chr1    908314  .       C       .       999     .      DP=137;MQSB=0.701554...

"bcftools view -i 'alt="G" && ref="A"' p3.vcf.gz" prints the header then dumps core.  Ouch!
"bcftools view -v snps -i 'alt="G" && ref="A"' p3.vcf.gz" seem to work OK.
The "." in the ALT column would seem to be the problem.

Have tried many syntactic variations to no avail.  Using Version: 1.2-1-gf0a7dcd, and 
htslib 1.2-1-gb0742b0.  Peculiar.

-- Fred P.

[pd3]

Hi Fred, can you run make clean in both htslib and bcftools and recompile? If the problem persist, can you please send me a small test case? I am not able to reproduce the problem with my data.

[pishotta]

At the moment, I can't recompile with the latest release.
[root@ubccan-hpmysql ~/bcftools]# date
Tue Feb  3 08:00:32 PST 2015
[root@ubccan-hpmysql ~]# rm -r bcftools htslib samtools
[root@ubccan-hpmysql ~]# git clone --branch=develop git://github.com/samtools/htslib.git
Initialized empty Git repository in /root/htslib/.git/
remote: Counting objects: 5654, done.
remote: Compressing objects: 100% (19/19), done.
remote: Total 5654 (delta 9), reused 4 (delta 1)
Receiving objects: 100% (5654/5654), 4.82 MiB | 1.09 MiB/s, done.
Resolving deltas: 100% (3591/3591), done.
[root@ubccan-hpmysql ~]# git clone --branch=develop git://github.com/samtools/bcftools.git
Initialized empty Git repository in /root/bcftools/.git/
remote: Counting objects: 5951, done.
remote: Compressing objects: 100% (52/52), done.
remote: Total 5951 (delta 24), reused 9 (delta 5)
Receiving objects: 100% (5951/5951), 5.54 MiB | 3.18 MiB/s, done.
Resolving deltas: 100% (3856/3856), done.
[root@ubccan-hpmysql ~]# git clone --branch=develop git://github.com/samtools/samtools.git
Initialized empty Git repository in /root/samtools/.git/
remote: Counting objects: 6871, done.
remote: Compressing objects: 100% (29/29), done.
remote: Total 6871 (delta 10), reused 0 (delta 0)
Receiving objects: 100% (6871/6871), 8.53 MiB | 262 KiB/s, done.
Resolving deltas: 100% (4017/4017), done.
[root@ubccan-hpmysql ~]# cd bcftools/
[root@ubccan-hpmysql ~/bcftools]# make
cd ../htslib && make lib-static
make[1]: Entering directory `/root/htslib'
gcc -g -Wall -O2 -I. -c -o kfunc.o kfunc.c
gcc -g -Wall -O2 -I. -c -o knetfile.o knetfile.c
gcc -g -Wall -O2 -I. -c -o kstring.o kstring.c
gcc -g -Wall -O2 -I. -c -o bgzf.o bgzf.c
gcc -g -Wall -O2 -I. -c -o faidx.o faidx.c
In file included from faidx.c:33:
htslib/faidx.h:94: error: wrong number of arguments specified for ‘__deprecated__’ attribute
make[1]: *** [faidx.o] Error 1
make[1]: Leaving directory `/root/htslib'
make: *** [../htslib/libhts.a] Error 2

[jkbonfield]

It looks like we accidentally ended up with a requirement for gcc 4.5 or newer. Until it is fixed, either use a newer gcc or modify the offending line in htslib/htslib/hts_defs.h.

[jmarshall]

Note btw that that's the current develop rather than "the latest release".

If you checkout the latest htslib release (cd htslib; git checkout 1.2.1) you'll avoid this problem too.

[mp15]

Oops, I didn't realise the attribute's message syntax was a late addition. I'll add a pull to fix it for legacy compilers.

Sent from my iPhone

On 3 Feb 2015, at 16:46, John Marshall notifications@github.com wrote:

Note btw that that's the current develop rather than "the latest release".

If you checkout the latest htslib release (cd htslib; git checkout 1.2.1) you'll avoid this problem too.

—
Reply to this email directly or view it on GitHub.

The Wellcome Trust Sanger Institute is operated by Genome Research
Limited, a charity registered in England with number 1021457 and a
company registered in England with number 2742969, whose registered
office is 215 Euston Road, London, NW1 2BE.

[pishotta]

OK, compile problems resolved, thanks, but with the attached test file:

--> bcftools view -i 'ref="A"' test.vcf.gz <--< works OK, but
--> bcftools view -i 'alt="A"' test.vcf.gz <--< dumps core.

It's a RHEL6 linux system.
[root@ubccan-hpmysql /vcf_repo]# uname -a
Linux ubccan-hpmysql 2.6.32-504.8.1.el6.x86_64 #1 SMP Fri Dec 19 12:09:25 EST 2014 x86_64 x86_64 x86_64 GNU/Linux

[pishotta]

Update...
I recompiled and re-ran the above commands on a different system under RHEL7...
> [root@ubccan-s-tardis bcftools]# uname -a
> Linux ubccan-s-tardis.can.ubc.ca 3.10.0-123.13.2.el7.x86_64 #1 SMP Thu Dec 18 14:09:13 UTC 2014 x86_64 x86_64 x86_64 GNU/Linux
...with the same results.  Core dump trying to check the ALT field.  I then replaced all the
"." characters in the ALT column of the vcf file with an "A", and the command completed
successfully.  Continuing to experiment, I note that the problem does not occur if there
is a "." in the REF column, only when there is one in the ALT column.  I grant that most
vcf files would not contain a "." in the ALT column if they only store variants, but for some
downstream statistical purposes I am also storing homozygous reference calls in them,
which does put a "." in the ALT column and which is a valid (if uncommon) use.  They were, 
in fact, generated by "bcftools call" with option settings to do so.
Perhaps the fact that a "." is a regex wildcard character is contributing to the problem?
-- Fred P.

[pishotta]

A test file that can be used to demonstrate this behavior can be found here:

http://ubccan-public.can.ubc.ca/download/test.vcf.gz

[pd3]

Thank you for the test case, now fixed.

[pishotta]

That fixed most things, thanks. No segfaults, but

bcftools view -i 'alt="."' test.vcf.gz

still fails to return rows where there is a “.” in the ALT column of the vcf file.
Thought it might be an “escape” problem on the command line (bash), but
various versions of putting back-slashes before the dot, didn’t help.

Test file is still at http://ubccan-public.can.ubc.ca/download/test.vcf.gz

-- Fred P.

From: pd3 [mailto:notifications@github.com]
Sent: Thursday, February 05, 2015 12:41 AM
To: samtools/bcftools
Cc: Pishotta, Frederick
Subject: Re: [bcftools] Expression problems in bcftools view (#201)

Thank you for the test case, now fixed.

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/201#issuecomment-73012149.

[jmarshall]

In file included from faidx.c:33:
htslib/faidx.h:94: error: wrong number of arguments specified for ‘__deprecated__’ attribute

This borkage has now been fixed on htslib's develop by samtools/htslib@a233385.