Attention:
The frontline support team will be unavailable to answer questions on April 15th and 17th 2019. We will be back soon after. Thank you for your patience and we apologize for any inconvenience!

HaplotypeCaller advanced option --force-active missing

Hello.
From reading the docs for version 4.0.11.0 (and working with GATK version 4.0.12.0)

h ttps://software.broadinstitute.org/gatk/documentation/tooldocs/current/org_broadinstitute_hellbender_tools_walkers_haplotypecaller_HaplotypeCaller.php#--disable-optimizations

I see from the option --disable-optimizations documentation :
"Setting the --force-active and --dont-trim-active-regions flags may also be necessary."

However the option --force-active seems no longer available. Because I get this error message :

"A USER ERROR has occurred: force-active is not a recognized option"

Is this option no longer relevant and is the documentation out-of-date ? or should it still be present ?
(Older versions of the documetations also mention --force-active or -forceActive).

Thank you for clarification.

Best Answer

Answers

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @rwk,

    I am looking into updating the latest version of gatk in the tool docs. I will get back to you shortly.

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin
    edited January 9

    Hi @rwk

    The latest version of 4.0.12.0 has been updated in the tool docs.
    --force-active and --dont-trim-active-regions flags are still present.

    Please send us the exact command you are using and the entire error log if the problem still persists.

  • bhanuGandhambhanuGandham Cambridge MAMember, Administrator, Broadie, Moderator admin

    Hi @rwk
    We are closing this issue for now as we have not heard from you in 2 business days. Please feel free to ask more questions if this issue has not been resolved.

  • rwkrwk Member
    Sorry for the delay,

    The command I try to launch is :

    ```
    $GATK HaplotypeCaller -R ref/ref.fasta -I bams/mother.bam -O sandbox/motherHC.vcf -bamout sandbox/motherHCdebug.bam --disable-optimizations --dont-trim-active-regions --force-active -L 20:10,000,000-10,200,000 -ip 100
    ```

    I get the following error :
    ```
    ***********************************************************************

    A USER ERROR has occurred: force-active is not a recognized option

    ***********************************************************************
    ```

    Without the --force-active option the runs finishes.

    I am running GATK 4.0.12.0

    The datasets are the test dataset from the GATK Variant Discovery Workshop.
  • dbeckerdbecker MunichMember ✭✭

    Hi,

    I have the same problem in GATK 4.0.4.0. --force-active is not mentioned in the documentation, but since you refer to this option as 'still' being available, I believe it should be there. It is only mentioned in the description of --disable-optimizations.

    Best,
    Daniel

  • AdelaideRAdelaideR Unconfirmed, Member, Broadie, Moderator admin

    @dbecker and @rwk I see that this option is not available in the most recent version of Haplotypecaller when I look at the command line, so the --disable-optimizations instructions probably need to be updated in the documentation.

    I came to this conclusion by doing a quick check of options available in the latest version of gatk.

    The simple method to do this is to run the program of interest and add "--help"

    On my local version (which I have in a docker), I ran:

    docker run broadinstitute/gatk gatk HaplotypeCaller --help
    

    And I don't see it.

    So, I will check with the developers and @bhanuGandham, who originally answered this question to determine whether the documentation or the program options need to be updated.

  • dbeckerdbecker MunichMember ✭✭

    Thanks for the clarification @AdelaideR

    As I understand it, --force-active is the only option to get the full bamout, isn't it? Is there another option?

    I have some variants in PKD1 that are called with HC but not with Sanger and some are called with Sanger but not with HC. These are based on the Pseudogene we have here, but I have to find out why exactly this happens and if we can change this by changing some parameters in the lab. Therefore I thought it might be helpful to see the whole region as HC defines it.
    Some variants that are not called by HC are visible in the bamout, which is strange. But others are in regions that are not even defined as active. I need to see how HC rearanges the reads there to understand the missing call.
    By the way... If Haplotype Caller defines active regions based on mismatches, can it be that variants that are only in the gvcf because they have been called in another sample have ./. genotype in this one because HC never looked at them? I have some Variants that do not appear in my sample but in others but every sample should have coverage. Still the have ./. genotype instead of 0/0.

    Best,
    Daniel

  • AdelaideRAdelaideR Unconfirmed, Member, Broadie, Moderator admin

    Hello @dbecker

    If I understand your question correctly, you want the bamout that contains all of the variants, not just the ones in the active regions, and all variants even if there is low coverage in some samples.

    Have you tried the following options to get the full output?

    false   If specified, we will not trim down the active region from the full region (active + extension) to just the active interval for genotyping
    
    --disable-optimizations
    false   Don't skip calculations in ActiveRegions with no variants
    
    
    --active-probability-threshold
    0.002   Minimum probability for a locus to be considered active.
    
    --enable-all-annotations
    false   Use all possible annotations (not for the faint of heart)
    
    

    There are a lot of options, found here

    As for the second question "./." means not enough coverage to make a call. For the position where it should be 0/0, how many samples are aligning to that particular location? Lowering the probability threshold may increase the number that are visible above the threshold.

    Looking at the IGV of the bamouts may be helpful to determine how much information is required to meet the threshold for "0/0" versus "./."

  • dbeckerdbecker MunichMember ✭✭

    Hi @AdelaideR,

    I want a bamout that conatins all reads in the region specified. I should not matter if parts of this region are active or not or if the have variants in them. I specifically want to check regions that don't have variants in them because I have Sanger variants for those.

    • The original bamout without any extra options contains just a few very specific parts of my specified region.
    • With --disable-optimizations these parts grow in size (horizontally) and a few new parts have reads even though they don't contain variants.
    • With --dont-trim-active-regions (additionally) these parts grow further (horizontally and in coverage).
    • Using --active-probability-threshold 0.000002 shows no further effect.
    • I tried --enable-all-annotations. This is not avialabe in 4.0.4.0 which we are using and validating for our diagnostics and if I try it in 4.0.11.0 which is also installed on our server, I get:
    A USER ERROR has occurred: Allele-specific annotations are not yet supported in the VCF mode
    

    So none of these options enables me to close all the gaps between those 'parts' of my specified region and I can't find out why some existing variants haven't been called.

    Best,
    Daniel

  • dbeckerdbecker MunichMember ✭✭

    @AdelaideR,

    Thanks! I will wait for the next version I think. We can't really update all the time since we do diagnostics and have to do a thorough validation of our pipeline with every change we make. Therefore we will update only if it is really necessary. But this question was not for our production process so i can use the next version to understand HC better.

    Best,
    Daniel

  • rwkrwk Member
    @AdelaideR
    Thank you for your response.
    I was suspecting this, seeing that the option was not available. The documentation is a bit out of sync with the version but it is good to know the option will be ported to the next version release.

    Best regards.
Sign In or Register to comment.