The current GATK version is 3.6-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!

Powered by Vanilla. Made with Bootstrap.

Using the GATK API to annotate my VCF (2)

lindenblindenb FrancePosts: 37Member

So I started to play with the GATK API to annotate my VCF. :-)

My first simple aim is to translate mt VCFBigWig to the GATK . Here is my code so far (see below)

I compiled it with:

javac -cp /home/lindenb/package/GenomeAnalysisTK-2.4-3-g2a7af43/GenomeAnalysisTK.jar:/home/lindenb/package/bigwig-read-only/dist/BigWig.jar:/home/lindenb/package/bigwig-read-only/lib/log4j-1.2.15.jar -sourcepath src/main/java \
        ./src/main/java/com/github/lindenb/gatk/tools/vcfbigwig/VCFBigWig.java

(it is compiled but there is a bunch of warnings like "/home/lindenb/package/GenomeAnalysisTK-2.4-3-g2a7af43/GenomeAnalysisTK.jar(org/broadinstitute/sting/gatk/GenomeAnalysisEngine.class): warning: Cannot find annotation method 'value()' in type 'Ensures': class file for com.google.java.contract.Ensures not found"

And I tried to run my code using various invocations of java, like

java  -cp /home/lindenb/package/GenomeAnalysisTK-2.4-3-g2a7af43/GenomeAnalysisTK.jar:/home/lindenb/package/bigwig-read-only/dist/BigWig.jar:/home/lindenb/package/bigwig-read-only/lib/log4j-1.2.15.jar:. org.broadinstitute.sting.gatk.CommandLineGATK \
            -T VariantAnnotator \
            -A com.github.lindenb.gatk.tools.vcfbigwig.VCFBigWig \
            -R /test.fa test.vcf.gz

but it raises an exeception

java.lang.ExceptionInInitializerError
    at org.broadinstitute.sting.gatk.walkers.annotator.VariantAnnotatorEngine.initializeAnnotations(VariantAnnotatorEngine.java:116)
(...)
Caused by: java.lang.NullPointerException
    at org.broadinstitute.sting.utils.classloader.JVMUtils.isAnonymous(JVMUtils.java:91)
    at org.broadinstitute.sting.utils.classloader.PluginManager.<init>(PluginManager.java:158)
    at org.broadinstitute.sting.utils.classloader.PluginManager.<init>(PluginManager.java:107)
    at org.broadinstitute.sting.gatk.walkers.annotator.interfaces.AnnotationInterfaceManager.<clinit>(AnnotationInterfaceManager.java:34)
    ... 9 more

Questions:

  • What the proper to declare my class VCFBigWig , how can I make it available for the VariantAnnotatorEngine ?

  • I used the @Argument annotation, is it the right way to catch an argument from the command line?

  • Can I use the same Annotator for one GATK invocation but with different arguments ?
  • there is a "initialize" method; Is there a "dispose" method to release the resources ?
  • what is the interface "RodRequiringAnnotation" ?
  • is "UserException" the right way to send an error ?
  • Should an InfoFieldAnnotation or a Walker be thread free (like a HttpServlet) ? I mean can I safely use class members ?
  • I saw in VariantAnnotatorEngine that an annotator can either change the ID or the INFO. How should I design my class in order to alter both fields ? For example my tool VCFTabix with use another VCF file (e.g: the VCF for 1KGenomes) to update the ID and the INFO fields.

Thank you for your help,

Pierre

.

package com.github.lindenb.gatk.tools.vcfbigwig;

import java.io.IOException;
import java.util.Arrays;
import java.util.Collections;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
import java.util.Set;

import org.broad.igv.bbfile.BBFileReader;
import org.broad.igv.bbfile.BigWigIterator;
import org.broad.igv.bbfile.WigItem;
import org.broadinstitute.sting.commandline.Argument;
import org.broadinstitute.sting.gatk.GenomeAnalysisEngine;
import org.broadinstitute.sting.gatk.contexts.AlignmentContext;
import org.broadinstitute.sting.gatk.contexts.ReferenceContext;
import org.broadinstitute.sting.gatk.refdata.RefMetaDataTracker;
import org.broadinstitute.sting.gatk.walkers.annotator.interfaces.AnnotatorCompatible;
import org.broadinstitute.sting.gatk.walkers.annotator.interfaces.InfoFieldAnnotation;
import org.broadinstitute.sting.gatk.walkers.annotator.interfaces.RodRequiringAnnotation;
import org.broadinstitute.sting.utils.exceptions.UserException;
import org.broadinstitute.sting.utils.genotyper.PerReadAlleleLikelihoodMap;
import org.broadinstitute.variant.variantcontext.VariantContext;
import org.broadinstitute.variant.vcf.VCFHeaderLine;
import org.broadinstitute.variant.vcf.VCFHeaderLineType;
import org.broadinstitute.variant.vcf.VCFInfoHeaderLine;


public class VCFBigWig
extends InfoFieldAnnotation 
implements RodRequiringAnnotation
{
/** PATH to the bigwig file */
@Argument(fullName="bigWigFile",shortName="bw",doc="BigWig to process",required=true)
private String bigWigFile=null;
/** the ID in the INFO field */
@Argument(fullName = "bigwigid", shortName = "bwid", doc="bigwig ID ", required=true)
private String bigWigId=null;
/** BBFileReader to the bigwig file */
private BBFileReader bbFileReader=null;

public VCFBigWig()
    {
    }

@Override
public void initialize(AnnotatorCompatible walker,
            GenomeAnalysisEngine toolkit, Set<VCFHeaderLine> headerLines) {
    super.initialize(walker, toolkit, headerLines);
    /* open the BIGFile, check it is a bigwig */
    try {
        this.bbFileReader=new BBFileReader(bigWigFile);
        }
    catch (IOException e) {
        throw new UserException("Cannot read "+bigWigFile,e);
        }
    if(!this.bbFileReader.isBigWigFile())
        {
        throw new UserException("The following file is not a BigWig File : "+bigWigFile);
        }
    }


@Override
public Map<String, Object> annotate(
        RefMetaDataTracker tracker,
        AnnotatorCompatible walker,
        ReferenceContext ref,
        Map<String, AlignmentContext> stratifiedContexts,
        VariantContext vc,
        Map<String, PerReadAlleleLikelihoodMap> stratifiedPerReadAlleleLikelihoodMap)
    {
    return annotate(tracker, walker,ref,stratifiedContexts,vc);
    }



@Override
public Map<String, Object> annotate(
        RefMetaDataTracker tracker,
        AnnotatorCompatible walker,
        ReferenceContext ref,
        Map<String, AlignmentContext> stratifiedContexts,
        VariantContext vc)
    {
    final boolean contained=true;
    double sum=0;
    int count=0;
    BigWigIterator iter=this.bbFileReader.getBigWigIterator(
            vc.getChr(), vc.getStart()-1,
            vc.getChr(), vc.getEnd(),
            contained
            );

    while(iter.hasNext())
        {
        WigItem item=iter.next();
        float v=item.getWigValue();
        sum+=v;
        count++;
        }

    if(count==0)
        {
        return Collections.emptyMap();
        }
    Map<String,Object> m=new HashMap<String, Object>(1);
    m.put(this.bigWigId, (float)(sum/count));
    return m ;
    }
@Override
public List<VCFInfoHeaderLine> getDescriptions()
    {
    return Arrays.asList(new VCFInfoHeaderLine[]{
        new VCFInfoHeaderLine(bigWigId,1, VCFHeaderLineType.Float,"Mean values for bigwig file "+this.bigWigFile)   
        });
    }
@Override
public List<String> getKeyNames() {
    return Arrays.asList(new String[]{
            bigWigId    
            });
    }

}

Best Answer

Answers

Sign In or Register to comment.