The current GATK version is 3.7-0
Examples: Monday, today, last week, Mar 26, 3/26/04

Howdy, Stranger!

It looks like you're new here. If you want to get involved, click one of these buttons!

Get notifications!


You can opt in to receive email notifications, for example when your questions get answered or when there are new announcements, by following the instructions given here.

Got a problem?


1. Search using the upper-right search box, e.g. using the error message.
2. Try the latest version of tools.
3. Include tool and Java versions.
4. Tell us whether you are following GATK Best Practices.
5. Include relevant details, e.g. platform, DNA- or RNA-Seq, WES (+capture kit) or WGS (PCR-free or PCR+), paired- or single-end, read length, expected average coverage, somatic data, etc.
6. For tool errors, include the error stacktrace as well as the exact command.
7. For format issues, include the result of running ValidateSamFile for BAMs or ValidateVariants for VCFs.
8. For weird results, include an illustrative example, e.g. attach IGV screenshots according to Article#5484.
9. For a seeming variant that is uncalled, include results of following Article#1235.

Did we ask for a bug report?


Then follow instructions in Article#1894.

Formatting tip!


Wrap blocks of code, error messages and BAM/VCF snippets--especially content with hashes (#)--with lines with three backticks ( ``` ) each to make a code block as demonstrated here.

Jump to another community
Picard 2.10.2 is now available. As of 2.10.0, Picard supports NovaSeq CBCL data. Download and read release notes at https://github.com/broadinstitute/picard/releases.
**GATK4-BETA.2** is here. That's TWO, as in the second beta release. Be sure to read about the known issues before test driving. See Article#9881 to start and https://github.com/broadinstitute/gatk/blob/master/README.md for details.

Base recalibration & realignment without known sites

gilgigilgi Member
edited October 2012 in Ask the GATK team

Dear GATK team,

Thanks a lot for the new GATK version and GATK forum!

I am trying to use GATK for yeast strains. I do not have files of known sites of SNPs/indels.
I understand that the BaseRecalibrator must get such a file.
Do you suggest to skip calibration and realignment, or is there another way to go here?

Post edited by Geraldine_VdAuwera on

Best Answer

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie
    edited September 2012 Accepted Answer

    Hi Gilgi,

    Glad you like the new digs!

    Skipping Base recalibration and realignment are never a good option. Theere is an article on Base Quality Score Recalibration in the Methods & Workflows section of the Guide that you need to read. Towards the end, in the "Troubleshooting" section, there is a discussion of what to do if you don't have a set of known variants for your organism.

    Good luck!

    Post edited by Geraldine_VdAuwera on

Answers

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie
    edited September 2012 Accepted Answer

    Hi Gilgi,

    Glad you like the new digs!

    Skipping Base recalibration and realignment are never a good option. Theere is an article on Base Quality Score Recalibration in the Methods & Workflows section of the Guide that you need to read. Towards the end, in the "Troubleshooting" section, there is a discussion of what to do if you don't have a set of known variants for your organism.

    Good luck!

    Post edited by Geraldine_VdAuwera on
  • gilgigilgi Member

    Thanks a lot!!!

  • The link is broken right now.

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie

    Link fixed, sorry about that.

  • mhtmht Member

    Hi, do you have any solution for the IndelRealigner as well? I have a non-model organism as reference and there are no known indel sites/document which I can use.

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie

    If you don't have any known indels you can run without.

  • mhtmht Member

    Hi, I ran it without any known indels, but some regions didn't seem to realign as expected. As a result, GATK is calling both a SNP and single base insertion at the same location (the true variant is the insert). Any idea why this is happening?

  • ebanksebanks Broad InstituteMember, Broadie, Dev

    What's your command line?

  • mhtmht Member

    Hi Eric,

    to identify intervals:
    java -Xmx1g -jar GenomeAnalysisTKLite.jar -T RealignerTargetCreator -R reference.fa -I sampleA_clean_paired.sorted.dedup.bam -o sampleA_clean_paired.sorted.dedup.intervals

    to realign:
    java -Xmx2g -jar GenomeAnalysisTKLite.jar -T IndelRealigner -R reference.fa -I sampleA_clean_paired.sorted.dedup.bam -targetIntervals sampleA_clean_paired.sorted.dedup.intervals -o sampleA_clean_paired.sorted.dedup.realigned.bam

    Thanks for your help

  • ebanksebanks Broad InstituteMember, Broadie, Dev

    Okay, good. Those look reasonable.
    Now why do you think some regions didn't realign properly? Could you please post an example IGV screenshot with an explanation? Thanks.

  • mhtmht Member

    Hi Eric,

    I've attached a screenshot of a region as an example. (see: http://imageshack.us/photo/my-images/803/realigned.png)
    Also, I noticed some messages during running GATK that goes like this:

    INFO 15:56:32,707 IndelRealigner - Not attempting realignment in interval contig_31331:6142-6275 because there are too many reads.

    I'm wondering if this is the reason why... and if so, any suggestions on how to deal with this?

    Thanks so much for your help.

  • ebanksebanks Broad InstituteMember, Broadie, Dev

    That is certainly a possible reason. Please take a look at the documentation for command-line arguments that allow you to increase the maximum number of reads allowed in a given realignment region.

  • mhtmht Member

    Hi Eric,

    I found the --maxReadsForRealignment option and I set it at a higher number (~ 100000) and am getting less of those messages. There are, however, still a few of these cases. There is a note in the documentation which says "For expert users only!". Is there other harms of increasing this number, other than the toll it takes on memory? How do I decide what to set this number at?

    Much thanks.

  • ebanksebanks Broad InstituteMember, Broadie, Dev

    Just memory. The value to use really depends on your average coverage and whether it's important to realign every single region that might need it.

  • Hi, I would like to use the BaseRecalibrator but I do not have a database of known SNPs for the organism I am studying. I have a VCF files with potential real SNPs and I want to use the best ones as known sites. How many SNPS at least do we need? Are 500 enough?

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie

    There is no hard number, it depends how divergent you expect your organism samples to be. Do you have any way to estimate this?

  • nanacnanac Member
    edited November 2012

    We will study 40 genotypes. So far, the genetic differentiation between 2 genotypes has been assessed from microsatellite markers (Fsc), and varies from 50% to 70%. We currently have GATK results on 7 genotypes from a previous study reporting at least 14500 polymorphic sites present on each genotype with a quality between 1329 and 4867. Maybe we can use these ones or a part of them?

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie

    Yes, you can use a high-confidence subset of those.

  • joanajoana Member

    Hi
    I also work on nonmodel organisms and don't have a set of known polymorphisms. I have 1-4 individuals each from several quite divergent species. Would you use high-confidence SNPs found with Unified Genotyper run over all individuals of all species together or would you define the set of high-confidence SNPs for each species separately, even though I only have few individuals each?
    Thanks a lot.

  • Geraldine_VdAuweraGeraldine_VdAuwera Cambridge, MAMember, Administrator, Broadie

    Hi Joana,

    If the species are quite divergent I would probably generate separate high-confidence sets for each. It may be interesting to then look at how the separate sets intersect, how much overlap there is.

Sign In or Register to comment.