Notice:
If you happen to see a question you know the answer to, please do chime in and help your fellow community members. We encourage our fourm members to be more involved, jump in and help out your fellow researchers with their questions. GATK forum is a community forum and helping each other with using GATK tools and research is the cornerstone of our success as a genomics research community.We appreciate your help!

Test-drive the GATK tools and Best Practices pipelines on Terra


Check out this blog post to learn how you can get started with GATK and try out the pipelines in preconfigured workspaces (with a user-friendly interface!) without having to install anything.

about combine GVCFs

after filtering the vcf files with ISHet Filter, I am using the command CombineGVCFs as follows:
.././gatk CombineGVCFs -R reference.fasta --variant filtered1.g.vcf --variant filtered2.g.vcf --variant filtered3.snps.g.vcf -O combined.vcf

The input vcf files have thousands of SNPs but after combining these files, there is no SNP in the OUTPUT file. why is so? Kindly guide me .

I am using vGATK4.1.0.0

Answers

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    @navneetsekhon

    Can you please confirm that the ref files used with combineGVCFs match the ref file used to generate the vcf files?

  • navneetsekhonnavneetsekhon Member
    Yes I have used the same reference file for both the steps.
  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin
    edited April 4

    @navneetsekhon

    Thats odd. CombineGVCFs is merely combining the individual gvcfs into one mutli-sample gvcf so it should not filter any variants.

    Can you try to run CombineGVCF without the filtering step. The purpose of gvcf files is to contain every genomic location information and not just the variants. Hence should not be filtered. See this doc for more info: https://gatkforums.broadinstitute.org/gatk/discussion/11004/gvcf-genomic-variant-call-format

  • navneetsekhonnavneetsekhon Member
    Problem is still the same. I combined without filtering but there was nothing under header in the output file
  • navneetsekhonnavneetsekhon Member
    kindly guide me how to remove indels using GATK4.1.0.0 pipeline.
  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @navneetsekhon

    I checked with the dev team and this error has been seen with outdated vcf index files.
    Try to delete the index files and remake them and then Rerun CombineGVCFs.

    tabix -p vcf test.g.vcf

  • navneetsekhonnavneetsekhon Member
    I am trying to do base quality score recalibration (BQSR) using following 2 commands.

    ./gatk T PrintReads -R reference.fasta -I dedup_reads.bam -BQSR recalibration_report.grp -o output.bam

    ./gatk BaseRecalibrator -R reference.fasta -I dedup_reads.bam -knownSites latest_dbsnp.vcf - Q recal_data.table

    But these commands are not working. KIndly guide me. I am using v4.1.0.0 version and my data is GWAS data

    Please update me with the right command
  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @navneetsekhon

    You can find the right commands here in this tool docs for PrintReads and BaseRecalibrator.

  • navneetsekhonnavneetsekhon Member
    Kindly tell me how to combineGVCFs. I am still not able to use this command. I am writing here all the commands which I used for my samples for your reference:
    Alignment – Map to Reference

    bwa mem -M -R '@RG\tID:sample_name\tLB:sample_name\tPL:ILLUMINA\tPM:HISEQ\tSM:sample_name' reference.fasta sample_name_1_clean.fq.gz sample_name_2_clean.fq.gz > aligned_reads.sam

    #Sort SAM file by coordinate, convert to BAM

    java -jar picard.jar SortSam INPUT=aligned_reads.sam OUTPUT=sorted_reads.bam SORT_ORDER=coordinate

    mark Duplicates

    java -jar picard.jar MarkDuplicates INPUT=sorted_reads.bam OUTPUT=dedup_reads.bam METRICS_FILE=metrics.txt

    #Build BAM Index

    java -jar picard.jar BuildBamIndex INPUT=dedup_reads.bam

    # Call Variants

    ./gatk HaplotypeCaller -R ref -I dedup_reads.bam -o raw_variants.vcf

    #Extract SNPs and Indels

    ./gatk SelectVariants -R reference.fasta -V raw_variants.vcf -O raw_snps.vcf

    #Filter SNP

    ./gatk VariantFiltration -R reference.fasta -V raw_snps.vcf --genotype-filter-expression "isHet ==1" --genotype-filter-name "isHetFilter" -O filtered_snps.vcf

    ./gatk SelectVariants -V filtered_snps.vcf --set-filtered-gt-to-nocall -O trioGGVCF_VF_SV.vcf

    #Consolidate GVCFs

    ./gatk CombineGVCFs -R Bnigreference.fasta --variant trioGGVCF_VF_SV.vcf --variant trioGGVCF_2VF_SV.vcf -O cohort.g.vcf


    I am not able to get combines GVCF files. Kindly help!!!!
  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    @navneetsekhon

    To produce gvcf from haplotypecaller you need to use -ERC GVCF option. See tool docs for more info on that: https://software.broadinstitute.org/gatk/documentation/tooldocs/current/org_broadinstitute_hellbender_tools_walkers_haplotypecaller_HaplotypeCaller.php
    Once you generate the gvcf then you can use CombineGVCFs to combine them.
    Please read the tools docs for information on how to use these tools.

  • navneetsekhonnavneetsekhon Member

    Thank you dear Bhanu. Your suggestions solved my problem. Thanks alot for always replying promptly.

    I combined GVCF files and then used GenotypeGVCFs from the output file of CombineGVCFs file. The output vcf file from GenotypeGVCFs contain indels.

    1. Kindly tell me how to remove those indels using GATK pipeline??
    2. Please let me know that how to filter this vcf file based on genotypic quality, minor allele frequency and allele count ?
  • navneetsekhonnavneetsekhon Member

    Please reply. I am waiting for your answer.

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    HI @navneetsekhon

    You can remove indels and filter variants using SelectVariants. Please read this tool doc for information on its usage: https://software.broadinstitute.org/gatk/documentation/tooldocs/current/org_broadinstitute_hellbender_tools_walkers_variantutils_SelectVariants.php

  • navneetsekhonnavneetsekhon Member

    Hi,

    I am using the following command for filtration of the file (output of genotypeGVCFs file)

    ../gatk VariantFiltration -V genotyped_snps.vcf --FilterExpression "QD < 2.0 || FS > 60.0 || MQ < 40.0 || -12.5 || ReadPosRankSum < -8.0" --filterName "my_snp_filter" -O hardfilter_SNPs.vcf

    When I run this command the following, following error id displayed:

    A USER ERROR has occurred: FilterExpression is not a recognized option.

    Please suggest the possible wayout !

  • navneetsekhonnavneetsekhon Member

    Hi,

    I am using the following command for filtration of the file (output of genotypeGVCFs file)

    ../gatk VariantFiltration -V genotyped_snps.vcf --FilterExpression "QD < 2.0 || FS > 60.0 || MQ < 40.0 || -12.5 || ReadPosRankSum < -8.0" --filterName "my_snp_filter" -O hardfilter_SNPs.vcf

    When I run this command the following, error is displayed:

    ******A USER ERROR has occurred: FilterExpression is not a recognized option.******

    Please suggest the possible wayout !

  • navneetsekhonnavneetsekhon Member

    Hi Bhanu

    the following command finally worked:

    ../gatk VariantFiltration -V genotyped_snps.vcf --Filter-expression "QD < 2.0 || FS > 60.0 || MQ < 40.0 || -12.5 || ReadPosRankSum < -8.0" --filter-name "my_snp_filter" -O hardfilter_SNPs.vcf

    But, the no.of SNPs after hard filtering remains same.
    It seems that hard filtering only marks the SNPs that are passed through filters. It is right?
    Kindly tell how to remove the failed SNPs?

  • TwesidaveTwesidave Member

    Hello,

    I would like to run the CombineGVCFs command...is there a way I can supply the input gVCF files in one go rather than writing the -V flag 124 times? Thank you.

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @navneetsekhon

    With VariantFiltration records are hard-filtered by changing the value in the FILTER field to something other than PASS. Filtered records will be preserved in the output unless their removal is requested in the command line. You can use --excludeFiltered argument in SelectVariants to exclude sites that have been marked as filtered (i.e. have anything other than . or PASS in the FILTER field) will be excluded from the output.

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @Twesidave

    You can provide a list of gvcfs files to CombineGVCFs in a .list file in the -V argument instead of writing the -V flag 124 times.

Sign In or Register to comment.