add --skip-pls flag to gemini load.
this will not add the gt_phred_ll_het/homref/homalt BLOB columns
to the database.
MODIFY
docs/content/history.rst
docs/content/history.rst
@@ -19,6 +19,8 @@ Release History
#. Allow specifying a --where clause for `gene_wise` to filter considered variants. (Thanks Uma for suggestion).
#. Fix bug in propagating vep-extras to variant_impacts table (Thanks Sergey for reporting).
#. Update clinvar annotation.
+#. Add new `gemini load` flag `--skip-pls` which can give size reductions for the database and speed Improvements
+ for querying and loading. If you never use `gt_phred_ll_homalt/het/homref` then you can load with this flag.
0.19.1
======

MODIFY
docs/content/preprocessing.rst
docs/content/preprocessing.rst
@@ -95,6 +95,9 @@ variants. By default GEMINI will leave these variants in the database during
loading but they can be filtered out during the loading step by passing the
``--passonly`` flag to load.
+You can create a smaller, faster database if you dont need the genotype likelihoods format
+each sample by passing the `--skip-pls` flag.
+
================================
Using multiple CPUs for loading

MODIFY
gemini/database.py
gemini/database.py
@@ -99,10 +99,18 @@ def get_path(path):
return path
-def create_tables(path, effect_fields=None):
+def create_tables(path, effect_fields=None, pls=True):
"""
Create our master DB tables
"""
+ if pls:
+ pls = """\
+ gt_phred_ll_homref blob,
+ gt_phred_ll_het blob,
+ gt_phred_ll_homalt blob,
+"""
+ else:
+ pls = ""
if effect_fields:
effect_string = "".join(e + (" float,\n" if e.endswith("_num") else " TEXT,\n") for e in effect_fields)
else:
@@ -129,9 +137,7 @@ def create_tables(path, effect_fields=None):
gt_alt_depths blob,
gt_quals blob,
gt_copy_numbers blob,
- gt_phred_ll_homref blob,
- gt_phred_ll_het blob,
- gt_phred_ll_homalt blob,
+ %s,
call_rate float,
max_aaf_all float,
in_dbsnp bool,
@@ -273,7 +279,7 @@ def create_tables(path, effect_fields=None):
gnomad_num_het int,
gnomad_num_hom_alt int,
gnomad_num_chroms int,
- %s""" % effect_string.rstrip(","),
+ %s""" % (pls, effect_string.rstrip(",")),
variant_impacts="""
variant_id integer,
@@ -371,7 +377,7 @@ def create_tables(path, effect_fields=None):
for table in db:
db[table] = db[table].strip().strip(",").split(",\n")
- db[table] = [x.strip().split() for x in db[table]]
+ db[table] = [x.strip().split() for x in db[table] if x.strip()]
cols = [sql.Column(c[0], lookup[c[1].lower()]) for c in db[table]]
if table != "variant_impacts":

MODIFY
gemini/gemini_load.py
gemini/gemini_load.py
#!/usr/bin/env python
-from __future__ import division
+from __future__ import division, print_function
# native Python imports
import os.path
@@ -107,21 +107,24 @@ def merge_chunks(chunks, db, kwargs):
cleanup_temp_db_files(chunks)
return db
-def get_merge_chunks_cmd(chunks, db, tempdir=None, vcf=None, anno_type=None):
+def get_merge_chunks_cmd(chunks, db, tempdir=None, vcf=None, anno_type=None, skip_pls=False):
chunk_names = ""
for chunk in chunks:
chunk_names += " --chunkdb " + chunk
tempdir_string, vcf_string, annotype_string = "", "", ""
+ skip_pls_string = ""
if tempdir is not None:
tempdir_string = " --tempdir " + tempdir
if vcf is not None:
vcf_string = " --vcf " + vcf
if anno_type is not None:
annotype_string = " -t " + anno_type
+ if skip_pls:
+ skip_pls_string = "--skip-pls"
return ("gemini merge_chunks {chunk_names} {tempdir_string} "
- "{vcf_string} {annotype_string} --db {db}").format(**locals())
+ "{vcf_string} {annotype_string} {skip_pls_string} --db {db}").format(**locals())
def finalize_merged_db(tmp_db, db):
@@ -160,7 +163,7 @@ def merge_chunks_multicore(chunks, args):
tmp_dbs = get_temp_dbs(len(sub_merges), os.path.dirname(sub_merges[0][0]))
for sub_merge, tmp_db in zip(sub_merges, tmp_dbs):
cmd = get_merge_chunks_cmd(sub_merge, tmp_db, tempdir=args.tempdir, vcf=args.vcf,
- anno_type=args.anno_type)
+ anno_type=args.anno_type, skip_pls=args.skip_pls)
procs.append(subprocess.Popen(cmd, shell=True, stderr=sys.stderr))
wait_until_finished(procs)
cleanup_temp_db_files(chunks)
@@ -232,6 +235,8 @@ def load_chunks_multicore(grabix_file, args):
if args.passonly is True:
passonly = "--passonly"
+ skip_pls = "--skip-pls" if args.skip_pls else ""
+
skip_info_string = ""
if args.skip_info_string is True:
skip_info_string = "--skip-info-string"
@@ -310,6 +315,8 @@ def load_chunks_ipython(grabix_file, args, view):
if args.skip_info_string is True:
skip_info_string = "--skip-info-string"
+ skip_pls = "--skip-pls" if args.skip_pls else ""
+
vcf, _ = os.path.splitext(os.path.basename(grabix_file))
chunk_steps = get_chunk_steps(grabix_file, args)
total_chunks = len(chunk_steps)
@@ -320,6 +327,7 @@ def load_chunks_ipython(grabix_file, args, view):
"chunk_dir": chunk_dir,
"vcf": vcf,
"grabix_file": grabix_file,
+ "skip_pls": skip_pls,
"no_genotypes": no_genotypes,
"no_load_genotypes": no_load_genotypes,
"skip_gerp_bp": skip_gerp_bp,
@@ -358,7 +366,7 @@ def cleanup_temp_db_files(chunk_dbs):
def gemini_pipe_load_cmd():
grabix_cmd = "grabix grab {grabix_file} {start} {stop}"
gemini_load_cmd = ("gemini load_chunk -v - {anno_type} {ped_file}"
- " {no_load_genotypes} {no_genotypes}"
+ " {no_load_genotypes} {no_genotypes} {skip_pls}"
" {skip_gerp_bp} {skip_gene_tables} {skip_cadd}"
" {passonly} {skip_info_string} {test_mode} {tempdir}"
" -o {start} {chunk_dir}{vcf}.chunk{chunk_num}.db")

MODIFY
gemini/gemini_load_chunk.py
gemini/gemini_load_chunk.py
@@ -354,7 +354,7 @@ class GeminiLoader(object):
db_path = self.args.db if not hasattr(self.args, 'tmp_db') else self.args.tmp_db
if os.path.exists(db_path):
os.remove(db_path)
- self.c, self.metadata = database.create_tables(db_path, effect_fields or [])
+ self.c, self.metadata = database.create_tables(db_path, effect_fields or [], not self.args.skip_pls)
session = self.c
if session.bind.name == "sqlite":
self.c.execute('PRAGMA synchronous=OFF')
@@ -527,6 +527,7 @@ class GeminiLoader(object):
gt_quals = var.gt_quals
#gt_copy_numbers = np.array(var.gt_copy_numbers, np.float32) # 1.0 2.0 2.1 -1
gt_copy_numbers = None
+ if not self.args.skip_pls:
gt_phred_ll_homref = var.gt_phred_ll_homref
gt_phred_ll_het = var.gt_phred_ll_het
gt_phred_ll_homalt = var.gt_phred_ll_homalt
@@ -561,9 +562,6 @@ class GeminiLoader(object):
gt_phases=pack_blob(gt_phases), gt_depths=pack_blob(gt_depths),
gt_ref_depths=pack_blob(gt_ref_depths), gt_alt_depths=pack_blob(gt_alt_depths),
gt_quals=pack_blob(gt_quals), gt_copy_numbers=pack_blob(gt_copy_numbers),
- gt_phred_ll_homref=pack_blob(gt_phred_ll_homref),
- gt_phred_ll_het=pack_blob(gt_phred_ll_het),
- gt_phred_ll_homalt=pack_blob(gt_phred_ll_homalt),
call_rate=call_rate, in_dbsnp=bool(in_dbsnp),
rs_ids=rs_ids,
@@ -695,6 +693,10 @@ class GeminiLoader(object):
)
+ if not self.args.skip_pls:
+ variant['gt_phred_ll_homref'] = pack_blob(gt_phred_ll_homref)
+ variant['gt_phred_ll_het'] = pack_blob(gt_phred_ll_het)
+ variant['gt_phred_ll_homalt'] = pack_blob(gt_phred_ll_homalt)
variant['max_aaf_all'] = max(-1,
variant['aaf_esp_ea'],

MODIFY
gemini/gemini_main.py
gemini/gemini_main.py
@@ -193,6 +193,10 @@ def main():
action='store_true',
help='Load in test mode (faster)',
default=False)
+ parser_load.add_argument('--skip-pls',
+ action='store_true',
+ help='dont create columns for phred-scaled genotype likelihoods',
+ default=False)
def load_fn(parser, args):
from gemini import gemini_load
if args.vcf != "-":
@@ -291,6 +295,10 @@ def main():
action='store_true',
help='Load in test mode (faster)',
default=False)
+ parser_loadchunk.add_argument('--skip-pls',
+ action='store_true',
+ help='dont create columns for phred-scaled genotype likelihoods',
+ default=False)
parser_loadchunk.add_argument('--tempdir', dest='tempdir',
default=tempfile.gettempdir(),
help='Local (non-NFS) temp directory to use for working around SQLite locking issues '
@@ -328,6 +336,10 @@ def main():
action='store_true',
help='Create all database indexes. If multiple merges are used to create a database, only the last merge '
'should create the indexes.')
+ parser_mergechunks.add_argument('--skip-pls',
+ action='store_true',
+ help='dont create columns for phred-scaled genotype likelihoods',
+ default=False)
def mergechunk_fn(parser, args):
from gemini import gemini_merge_chunks

MODIFY
gemini/gemini_merge_chunks.py
gemini/gemini_merge_chunks.py
@@ -173,7 +173,9 @@ def merge_db_chunks(args):
if os.path.exists(args.db):
os.remove(args.db)
- gemini_db.create_tables(args.db, gemini_load_chunk.get_extra_effects_fields(args) if args.vcf else [])
+ gemini_db.create_tables(args.db,
+ gemini_load_chunk.get_extra_effects_fields(args) if args.vcf else
+ [], not args.skip_pls)
main_conn = sqlite3.connect(args.db)
main_conn.isolation_level = None

MODIFY
requirements.txt
requirements.txt
@@ -4,7 +4,7 @@ geneimpacts>=0.1.3
cython>=0.22.1
sqlalchemy>=1
pysam>=0.6
-cyvcf2>=0.4.2
+cyvcf2>=0.7.2
PyYAML>=3.10
pybedtools>=0.6.2
jinja2>=2.7.1

MODIFY
test/data-setup.sh
test/data-setup.sh
set -eo pipefail
gemini load -v test.auto_rec.vcf tt.db --no-genotypes && rm tt.db
+gemini load --cores 2 --skip-gene-tables --test-mode --skip-pls -v test.x-linked.vcf --skip-gerp-bp --skip-cadd -t VEP -p test.x-linked.ped test.x_linked.db
gemini load --skip-gene-tables --test-mode -v test.x-linked.vcf --skip-gerp-bp --skip-cadd -t VEP -p test.x-linked.ped test.x_linked.db
gemini load --skip-gene-tables --test-mode -t VEP -v test-vep-extra.vcf test.vep.extra.db &
gemini load --skip-gene-tables --test-mode --skip-gerp --skip-cadd -v test.PLs.vcf test.PLs.db &

add --skip-pls flag to gemini load.

this will not add the gt_phred_ll_het/homref/homalt BLOB columns
to the database.