allow specify gt_type in gt-filter.
closes
MODIFY
docs/content/history.rst
docs/content/history.rst
@@ -6,6 +6,11 @@ Release History
===============
#. Use vcfanno for faster, more generalized variant annotation and database creation.
+0.20.1
+======
+#. Allow `--gt-filter` to specify the gt_type, e.g. `(gt_quals).(=HET).(<20).(all)`
+
+
0.20.0
======
#. Include gnomad exome allele frequencies for all populations and include these in calc of max_aaf_all.

MODIFY
docs/content/querying.rst
docs/content/querying.rst
@@ -510,6 +510,18 @@ One can also base the wildcard on multiple criteria (in this case, brown hair an
--gt-filter "(gt_types).(hair_color=='brown' and phenotype==2).(!= HET).(all)" \
extended_ped.db
+There is a special case when we want to select samples based on their genotypes. For example,
+to show only sites where HET samples have a genotype quality above 20.
+
+.. code-block:: bash
+
+ $ gemini query -q "SELECT chrom, start, end, ref, alt, gene FROM variants \
+ where gene == 'SCNN1D' limit 5" \
+ --gt-filter "(gt_quals).(=HET).(<20).(all)" test/test.query.db
+
+In this case, the 2nd value must be one of `=HET`, `=HOM_REF` or `=HOM_ALT`
+
+
Lastly, wildcards can, of course, be combined with non-wildcard criteria:
.. code-block:: bash

MODIFY
gemini/GeminiQuery.py
gemini/GeminiQuery.py
#!/usr/bin/env python
-from __future__ import absolute_import
+from __future__ import absolute_import, print_function
import os
import sys
@@ -734,7 +734,8 @@ class GeminiQuery(object):
if self.gt_filter:
try:
if 'False' == self.gt_filter: continue
- unpacked = {'sample_info': self.sample_info}
+ unpacked = {'sample_info': self.sample_info, 'HET': HET,
+ 'HOM_REF': HOM_REF, 'HOM_ALT': HOM_ALT}
for col in self.gt_cols:
if col in self.gt_filter:
unpacked[col] = row[col]
@@ -1042,20 +1043,28 @@ class GeminiQuery(object):
if token.count('.') != 3 or \
token.count('(') != 4 or \
token.count(')') != 4:
+ toks = token.split(".")
raise ValueError("Wildcard filter should consist of 4 elements. Exiting.")
-
(column, wildcard, wildcard_rule, wildcard_op) = token.split('.')
# remove the syntactic parentheses
column = column.strip('(').strip(')').strip()
+ # allow commands like (gt_quals).(=HET).(>=20).(all)
wildcard = wildcard.strip('(').strip(')').strip()
+ is_gt_expression = "=HET" == wildcard or "=HOM_REF" == wildcard or "=HOM_ALT" == wildcard
+
wildcard_rule = wildcard_rule.strip('(').strip(')').strip()
wildcard_op = wildcard_op.strip('(').strip(')').strip()
+ extra_filter = None
# collect and save all of the samples that meet the wildcard criteria
# for each clause.
# these will be used in the list comprehension for the eval expression
# constructed below.
+ if is_gt_expression:
+ self.sample_info[token_idx] = self._get_matching_sample_ids("*")
+ extra_filter = " if gt_types[sample[0]] == %s" % wildcard[1:]
+ else:
self.sample_info[token_idx] = self._get_matching_sample_ids(wildcard)
# Replace HET, etc. with 1, et.c to avoid eval() issues.
@@ -1065,16 +1074,19 @@ class GeminiQuery(object):
if wildcard_op in ["all", "any"]:
if self.variant_id_getter:
joiner = " and " if wildcard_op == "all" else " or "
- rule = joiner.join("%s[%s]%s" % (column, s[0], wildcard_rule) for s in self.sample_info[token_idx])
+ rule = joiner.join("{column}[{sample}]{rule}".format(column=column, sample=s[0],
+ rule=wildcard_rule) for s in self.sample_info[token_idx])
rule = "(" + rule + ")"
else:
- rule = wildcard_op + "(" + column + '[sample[0]]' + wildcard_rule + " for sample in sample_info[" + str(token_idx) + "])"
+ rule = wildcard_op + "(" + column + '[sample[0]]' + wildcard_rule + " for sample in sample_info[" + str(token_idx) + "]%s)"
+ rule %= (extra_filter or "")
elif wildcard_op == "none":
if self.variant_id_getter:
rule = " or ".join("%s[%s]%s" % (column, s[0], wildcard_rule) for s in self.sample_info[token_idx])
rule = "~ ((" + rule + "))"
else:
- rule = "not any(" + column + '[sample[0]]' + wildcard_rule + " for sample in sample_info[" + str(token_idx) + "])"
+ rule = "not any(" + column + '[sample[0]]' + wildcard_rule + " for sample in sample_info[" + str(token_idx) + "]%s)"
+ rule %= (extra_filter or "")
elif "count" in wildcard_op:
# break "count>=2" into ['', '>=2']
tokens = wildcard_op.split('count')
@@ -1084,7 +1096,8 @@ class GeminiQuery(object):
rule = "%s|count|%s" % (rule, count_comp.strip())
seen_count = True
else:
- rule = "sum(" + column + '[sample[0]]' + wildcard_rule + " for sample in sample_info[" + str(token_idx) + "])" + count_comp
+ rule = "sum(" + column + '[sample[0]]' + wildcard_rule + " for sample in sample_info[" + str(token_idx) + "]%s)" + count_comp
+ rule %= (extra_filter or "")
else:
raise ValueError("Unsupported wildcard operation: (%s). Exiting." % wildcard_op)

MODIFY
test/test-query.sh
test/test-query.sh
@@ -258,6 +258,7 @@ gemini query -q "select chrom, start, end, ref, alt, gene, gts.1094PC0018, gts.1
where gene == 'SCNN1D' limit 5" \
--gt-filter "gt_types.1094PC0018 == HET or gt_types.1094PC0019 == HOM_REF" test.query.db \
> obs
+
check obs exp
rm obs exp
@@ -786,3 +787,12 @@ gemini query -q "select max_aaf_all,aaf_esp_ea, aaf_esp_aa, aaf_1kg_amr, aaf_1kg
check obs exp
rm obs exp
+
+echo " query.t43..."
+echo "chr1 1219381 1219382 C G SCNN1D C/C C/C" > exp
+gemini query -q "select chrom, start, end, ref, alt, gene, gts.1094PC0018, gts.1094PC0019 \
+ from variants \
+ where gene == 'SCNN1D' limit 5" \
+ --gt-filter "(gt_quals).(=HET).(>20).(all)" test.query.db > obs
+check obs exp
+rm obs exp

allow specify gt_type in gt-filter.

closes #839