Package org.broadinstitute.gatk.utils.codecs.beagle

Examples of org.broadinstitute.gatk.utils.codecs.beagle.BeagleFeature


        if (vc_input.isFiltered()) {
            vcfWriter.add(vc_input);
            return 1;
        }

        BeagleFeature beagleR2Feature = tracker.getFirstValue(beagleR2);
        BeagleFeature beagleProbsFeature = tracker.getFirstValue(beagleProbs);
        BeagleFeature beaglePhasedFeature = tracker.getFirstValue(beaglePhased);

        // ignore places where we don't have a variant
        if ( beagleR2Feature == null || beagleProbsFeature == null ||  beaglePhasedFeature == null)
        {
            vcfWriter.add(vc_input);
            return 1;
        }


        // get reference base for current position
        byte refByte = ref.getBase();

        // make new Genotypes based on Beagle results
        GenotypesContext genotypes = GenotypesContext.create(vc_input.getGenotypes().size());

        // for each genotype, create a new object with Beagle information on it

        int numGenotypesChangedByBeagle = 0;
        Integer alleleCountH = 0, chrCountH = 0;
        Double alleleFrequencyH = 0.0;
        int beagleVarCounts = 0;

        GenotypesContext hapmapGenotypes = null;

        if (vc_comp != null) {
            hapmapGenotypes = vc_comp.getGenotypes();
        }

        for ( final Genotype g : vc_input.getGenotypes() ) {
            boolean genotypeIsPhased = true;
            String sample = g.getSampleName();

            // If we have  a Hapmap (comp) ROD, compute Hapmap AC, AN and AF
            // use sample as key into genotypes structure
            if (vc_comp != null) {

                if (vc_input.getGenotypes().containsSample(sample) && hapmapGenotypes.containsSample(sample))  {

                    Genotype hapmapGenotype = hapmapGenotypes.get(sample);
                    if (hapmapGenotype.isCalled()){
                        chrCountH += 2;
                        if (hapmapGenotype.isHet()) {
                            alleleCountH += 1;
                        }    else if (hapmapGenotype.isHomVar()) {
                            alleleCountH += 2;
                        }
                    }
                }
            }

            ArrayList<String> beagleProbabilities = beagleProbsFeature.getProbLikelihoods().get(sample);
            ArrayList<String> beagleGenotypePairs = beaglePhasedFeature.getGenotypes().get(sample);

            // original alleles at this genotype
            Allele originalAlleleA = g.getAllele(0);

            Allele originalAlleleB = (g.getAlleles().size() == 2) ? g.getAllele(1) : g.getAllele(0); // hack to deal with no-call genotypes
View Full Code Here

TOP

Related Classes of org.broadinstitute.gatk.utils.codecs.beagle.BeagleFeature

Copyright © 2018 www.massapicom. All rights reserved.
All source code are property of their respective owners. Java is a trademark of Sun Microsystems, Inc and owned by ORACLE Inc. Contact coftware#gmail.com.