SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
IGV - accessing BAM file from different server adrian Bioinformatics 2 07-16-2013 09:09 AM
extremely large vcf file generated using sam tools janaahan13 RNA Sequencing 0 12-18-2012 08:43 PM
extremely large vcf file generated by GATK slowsmile Bioinformatics 2 12-09-2012 08:23 PM
Accessing .vcf.gz files on a Windows platform bnfoguy Bioinformatics 15 08-13-2012 02:32 AM
BAMseek large file viewer - now with VCF and FASTQ support BAMseek Bioinformatics 1 07-25-2011 06:57 PM

Reply
 
Thread Tools
Old 04-29-2019, 11:34 AM   #1
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default Accessing Genomic data in a large vcf.gz file

Hello all, I am very new to bioinformatics; I'm a high school intern. I have a .vcf.gz file that is 332 GB. It contains the whole genomes of around 800 dogs and all I need to do is check for a mutation for each of them. But before I can do that, I need to figure out how to crack the file open. I have heard of ways to possibly access them without decompressing the .vcf.gz because someone told me its possible that the uncompressed file could be up to three terabytes.

If anyone has any suggestions on how to proceed, I would be eternally grateful as I am nearing the end of my internship and need to figure this out.
XeroxHero69 is offline   Reply With Quote
Old 04-29-2019, 03:37 PM   #2
GenoMax
Senior Member
 
Location: East Coast USA

Join Date: Feb 2008
Posts: 6,953
Default

I think "bcftools view" command should give you what you need. Check the help page here.
GenoMax is offline   Reply With Quote
Old 04-30-2019, 03:50 AM   #3
JTH-Genome
Junior Member
 
Location: USA, NC

Join Date: Apr 2019
Posts: 6
Thumbs up Canine genomes

I liked GenoMax's suggestion, so go with that if you haven't already. Also be sure to use long read (PacBio) based assembly for your analysis against a reference genome.
If you are interested in canine genomes, here are some articles that may be of interest:
https://www.nature.com/articles/s41598-018-29190-3
https://www.youtube.com/watch?v=Ojsv...ature=youtu.be
https://www.pacb.com/blog/german-shepherd-genome/
https://www.pacb.com/blog/dog-meet-d...anine-genomes/
#SeqLife
JTH-Genome is offline   Reply With Quote
Old 04-30-2019, 09:09 AM   #4
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default

Quote:
Originally Posted by GenoMax View Post
I think "bcftools view" command should give you what you need. Check the help page here.
Thanks so much, I'll try this out
XeroxHero69 is offline   Reply With Quote
Old 04-30-2019, 11:04 AM   #5
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default

Quote:
Originally Posted by GenoMax View Post
I think "bcftools view" command should give you what you need. Check the help page here.
So I have bcftools installed. I want to check for a single-replacement mutation on a specific gene in each of the genomes in the file. Do you know of a way to do this using bcftools? I am new to bioinformatics and this is going way over my head. Thanks
XeroxHero69 is offline   Reply With Quote
Old 04-30-2019, 07:53 PM   #6
questor2010
Junior Member
 
Location: Seattle, WA

Join Date: Mar 2010
Posts: 5
Default

Quote:
Originally Posted by GenoMax View Post
I think "bcftools view" command should give you what you need. Check the help page here.
Note that for very large vcf files, using -t/-T instead of -r/-R is going to be much faster. And use the --threads option, that will help as well.
questor2010 is offline   Reply With Quote
Old 05-01-2019, 09:16 AM   #7
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default

Quote:
Originally Posted by questor2010 View Post
Note that for very large vcf files, using -t/-T instead of -r/-R is going to be much faster. And use the --threads option, that will help as well.
Thanks for your reply. Could you give me an example of the format of what a command that would use bcftools view to find a base at a position would look like?
XeroxHero69 is offline   Reply With Quote
Old 05-01-2019, 09:35 AM   #8
questor2010
Junior Member
 
Location: Seattle, WA

Join Date: Mar 2010
Posts: 5
Default

Quote:
Originally Posted by XeroxHero69 View Post
Thanks for your reply. Could you give me an example of the format of what a command that would use bcftools view to find a base at a position would look like?
Sure - here's an example:

bcftools view -f PASS --threads 8 -T target.bed -o gnomad.genomes.r2.1.sites.target.vcf.gz -O z gnomad.genomes.r2.1.sites.vcf.bgz

In this case, I'm using a bed file, instead of a single region. It is pulling out the FILTER=PASS variants that intersect the bed file into a new compressed vcf file. The source vcf file in this case is 465GB. If you have a single variant, you could use -t 1:11022. It might be best to specify a short range (1:11015-11030) if you're looking at indels - variant callers represent indels in different ways and you want to be sure you properly intersect them.
questor2010 is offline   Reply With Quote
Old 05-01-2019, 10:22 AM   #9
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default

Quote:
Originally Posted by questor2010 View Post
Sure - here's an example:

bcftools view -f PASS --threads 8 -T target.bed -o gnomad.genomes.r2.1.sites.target.vcf.gz -O z gnomad.genomes.r2.1.sites.vcf.bgz

In this case, I'm using a bed file, instead of a single region. It is pulling out the FILTER=PASS variants that intersect the bed file into a new compressed vcf file. The source vcf file in this case is 465GB. If you have a single variant, you could use -t 1:11022. It might be best to specify a short range (1:11015-11030) if you're looking at indels - variant callers represent indels in different ways and you want to be sure you properly intersect them.
Thanks so much. I tried this command:

bcftools view -f PASS --threads 8 -r chr9:55252802-55252802 -o 722g.990.SNP.INDEL.chrAll.vcf.gz -O z 722g.990.SNP.INDEL.chrAll.vcf.gz

and it returned:
[W::hts_idx_load2] The index file is older than the data file: 722g.990.SNP.INDEL.chrAll.vcf.gz.tbi
[W::hts_idx_load2] The index file is older than the data file: 722g.990.SNP.INDEL.chrAll.vcf.gz.tbi
[W::bgzf_read_block] EOF marker is absent. The input is probably truncated

Do you think you could point out what I did wrong?
Edit: I think i see my error, i made it output the data to the same file and now I think i ruined that data file... its only 16.1 kb now

Last edited by XeroxHero69; 05-01-2019 at 10:38 AM.
XeroxHero69 is offline   Reply With Quote
Old 05-01-2019, 10:51 AM   #10
questor2010
Junior Member
 
Location: Seattle, WA

Join Date: Mar 2010
Posts: 5
Default

Quote:
Originally Posted by XeroxHero69 View Post
Thanks so much. I tried this command:

bcftools view -f PASS --threads 8 -r chr9:55252802-55252802 -o 722g.990.SNP.INDEL.chrAll.vcf.gz -O z 722g.990.SNP.INDEL.chrAll.vcf.gz

and it returned:
[W::hts_idx_load2] The index file is older than the data file: 722g.990.SNP.INDEL.chrAll.vcf.gz.tbi
[W::hts_idx_load2] The index file is older than the data file: 722g.990.SNP.INDEL.chrAll.vcf.gz.tbi
[W::bgzf_read_block] EOF marker is absent. The input is probably truncated

Do you think you could point out what I did wrong?
Edit: I think i see my error, i made it output the data to the same file and now I think i ruined that data file... its only 16.1 kb now
Ouch. Yes - that's the problem. I hope you have another copy readily available. The index file warnings are often caused by file transfers. They transfer first for large files and then have an older timestamp. You can use "touch" to fix that, or reindex, but that takes a long time with large vcf files.
questor2010 is offline   Reply With Quote
Old 05-06-2019, 08:58 AM   #11
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default

Quote:
Originally Posted by questor2010 View Post
Ouch. Yes - that's the problem. I hope you have another copy readily available. The index file warnings are often caused by file transfers. They transfer first for large files and then have an older timestamp. You can use "touch" to fix that, or reindex, but that takes a long time with large vcf files.
So to make sure I do it right this time, what do I do to make the .vcf.bgz file?
XeroxHero69 is offline   Reply With Quote
Old 05-06-2019, 10:39 AM   #12
questor2010
Junior Member
 
Location: Seattle, WA

Join Date: Mar 2010
Posts: 5
Default

The -o switch specifies the output file name, the -O switch specifies the format:

-O, --output-type <b|u|z|v> b: compressed BCF, u: uncompressed BCF, z: compressed VCF, v: uncompressed VCF

So -O z implies that the name should be basename.vcf.gz (or basename.vcf.bgz).

If you have a single variant you just want to examine, you could use -v and output to a standard vcf file (text).
questor2010 is offline   Reply With Quote
Old 05-06-2019, 10:44 AM   #13
XeroxHero69
Member
 
Location: USA

Join Date: Apr 2019
Posts: 11
Default

Quote:
Originally Posted by questor2010 View Post
The -o switch specifies the output file name, the -O switch specifies the format:

-O, --output-type <b|u|z|v> b: compressed BCF, u: uncompressed BCF, z: compressed VCF, v: uncompressed VCF

So -O z implies that the name should be basename.vcf.gz (or basename.vcf.bgz).

If you have a single variant you just want to examine, you could use -v and output to a standard vcf file (text).
So do I need to create a new file for the output prior to running this command? how would the command look if I use -v without overwriting my file again. Sorry if these are simple questions. Thanks!
XeroxHero69 is offline   Reply With Quote
Reply

Tags
large files, vcf, vcf file

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 06:00 AM.


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