Search in sources :

Example 1 with ScatterChart

use of com.github.lindenb.jvarkit.chart.ScatterChart in project jvarkit by lindenb.

the class CaseControlJfx method doWork.

@Override
public int doWork(final List<String> args) {
    final VariantPartition partition;
    Pedigree pedigree = null;
    VCFIterator in = null;
    try {
        switch(this.partitionType) {
            case variantType:
                partition = new VariantTypePartition();
                break;
            case chromosome:
                partition = new ChromosomePartition();
                break;
            case autosomes:
                partition = new SexualContigPartition();
                break;
            case qual:
                partition = new QualPartition();
                break;
            case vqslod:
                partition = new VQSLODPartition();
                break;
            case typeFilter:
                partition = new TypeAndFilterPartiton();
                break;
            case distance:
                partition = new DisanceToDiagonalPartiton();
                break;
            case n_alts:
                partition = new NAltsPartition();
                break;
            default:
                throw new IllegalStateException(this.partitionType.name());
        }
        in = openVCFIterator(oneFileOrNull(args));
        if (this.pedigreeFile != null) {
            pedigree = Pedigree.newParser().parse(this.pedigreeFile);
        } else {
            pedigree = Pedigree.newParser().parse(in.getHeader());
        }
        final AFExtractor controlAFExtractor;
        if (!StringUtil.isBlank(this.controlFields)) {
            final List<AFExtractor> extractors = new AFExtractorFactory().parseFieldExtractors(this.controlFields);
            if (extractors.size() != 1) {
                LOG.error("extractor list should have size==1 . got " + extractors);
                return -1;
            }
            controlAFExtractor = extractors.get(0);
            if (!controlAFExtractor.validateHeader(in.getHeader())) {
                LOG.error("Invalid : " + controlAFExtractor);
                return -1;
            }
        } else {
            controlAFExtractor = null;
        }
        int count = 0;
        final ProgressFactory.Watcher<VariantContext> progress = ProgressFactory.newInstance().dictionary(in.getHeader()).logger(LOG).build();
        while (in.hasNext() && (this.limit_to_N_variants < 0 || count < this.limit_to_N_variants)) {
            final VariantContext ctx = progress.apply(in.next());
            if (this.ignore_ctx_filtered && ctx.isFiltered())
                continue;
            ++count;
            final List<Allele> alternates = ctx.getAlternateAlleles();
            for (int alt_idx = 0; alt_idx < alternates.size(); ++alt_idx) {
                final Allele alt = alternates.get(alt_idx);
                final Double[] mafs = { null, null };
                for (int i = 0; i < 2; ++i) {
                    if (i == 1 && controlAFExtractor != null) {
                        final List<Double> dvals = controlAFExtractor.parse(ctx);
                        if (alt_idx < dvals.size() && dvals.get(alt_idx) != null) {
                            final double d = dvals.get(alt_idx);
                            if (!Double.isNaN(d) && d >= 0 && d <= 1.0)
                                mafs[1] = d;
                        }
                    } else {
                        final MafCalculator mafCalculator = new MafCalculator(alt, ctx.getContig());
                        mafCalculator.setNoCallIsHomRef(no_call_is_homref);
                        for (Pedigree.Person person : (i == 0 ? pedigree.getAffected() : pedigree.getUnaffected())) {
                            if (selectSamples.equals(SelectSamples.males) && !person.isMale())
                                continue;
                            if (selectSamples.equals(SelectSamples.females) && !person.isFemale())
                                continue;
                            final Genotype genotype = ctx.getGenotype(person.getId());
                            if (genotype == null)
                                continue;
                            if (ignore_gt_filtered && genotype.isFiltered())
                                continue;
                            mafCalculator.add(genotype, person.isMale());
                        }
                        if (!mafCalculator.isEmpty()) {
                            mafs[i] = mafCalculator.getMaf();
                        }
                    }
                }
                if (mafs[0] == null || mafs[1] == null)
                    continue;
                final XYChart.Data<Number, Number> data = new XYChart.Data<Number, Number>(mafs[0], mafs[1]);
                partition.add(ctx, pedigree, data);
            }
        }
        progress.close();
        in.close();
        in = null;
        final NumberAxis xAxis = new NumberAxis(0.0, 1.0, 0.1);
        xAxis.setLabel("Cases");
        final NumberAxis yAxis = new NumberAxis(0.0, 1.0, 0.1);
        yAxis.setLabel("Controls" + (StringUtil.isBlank(this.controlFields) ? "" : "[" + this.controlFields + "]"));
        final ScatterChart<Number, Number> chart = new ScatterChart<>(xAxis, yAxis);
        for (final XYChart.Series<Number, Number> series : partition.getSeries()) {
            chart.getData().add(series);
        }
        String title = "Case/Control";
        if (!args.isEmpty()) {
            title = args.get(0);
            int slash = title.lastIndexOf("/");
            if (slash != -1)
                title = title.substring(slash + 1);
            if (title.endsWith(".vcf.gz"))
                title = title.substring(0, title.length() - 7);
            if (title.endsWith(".vcf"))
                title = title.substring(0, title.length() - 4);
        }
        if (userTitle != null)
            title = userTitle;
        chart.setTitle(title);
        chart.setAnimated(false);
        // chart.setLegendSide(this.legendSide);
        final RExporter rExporter = new RExporter();
        final PrintWriter pw = super.openFileOrStdoutAsPrintWriter(this.outputFile);
        rExporter.exportToR(pw, chart);
        pw.flush();
        pw.close();
    } catch (final Exception err) {
        LOG.error(err);
        return -1;
    } finally {
        CloserUtil.close(in);
    }
    return 0;
}
Also used : NumberAxis(com.github.lindenb.jvarkit.chart.NumberAxis) ProgressFactory(com.github.lindenb.jvarkit.util.log.ProgressFactory) ScatterChart(com.github.lindenb.jvarkit.chart.ScatterChart) VariantContext(htsjdk.variant.variantcontext.VariantContext) AFExtractorFactory(com.github.lindenb.jvarkit.util.vcf.AFExtractorFactory) VCFIterator(htsjdk.variant.vcf.VCFIterator) PrintWriter(java.io.PrintWriter) Genotype(htsjdk.variant.variantcontext.Genotype) RExporter(com.github.lindenb.jvarkit.chart.RExporter) Allele(htsjdk.variant.variantcontext.Allele) Pedigree(com.github.lindenb.jvarkit.util.Pedigree) AFExtractor(com.github.lindenb.jvarkit.util.vcf.AFExtractorFactory.AFExtractor) XYChart(com.github.lindenb.jvarkit.chart.XYChart)

Aggregations

NumberAxis (com.github.lindenb.jvarkit.chart.NumberAxis)1 RExporter (com.github.lindenb.jvarkit.chart.RExporter)1 ScatterChart (com.github.lindenb.jvarkit.chart.ScatterChart)1 XYChart (com.github.lindenb.jvarkit.chart.XYChart)1 Pedigree (com.github.lindenb.jvarkit.util.Pedigree)1 ProgressFactory (com.github.lindenb.jvarkit.util.log.ProgressFactory)1 AFExtractorFactory (com.github.lindenb.jvarkit.util.vcf.AFExtractorFactory)1 AFExtractor (com.github.lindenb.jvarkit.util.vcf.AFExtractorFactory.AFExtractor)1 Allele (htsjdk.variant.variantcontext.Allele)1 Genotype (htsjdk.variant.variantcontext.Genotype)1 VariantContext (htsjdk.variant.variantcontext.VariantContext)1 VCFIterator (htsjdk.variant.vcf.VCFIterator)1 PrintWriter (java.io.PrintWriter)1