From 61bc334df1eff953d93866c1b107bf5fbf38a2f1 Mon Sep 17 00:00:00 2001
From: Chris Hartl <chartl@broadinstitute.org>
Date: Mon, 14 Jan 2013 09:21:30 -0500
Subject: [PATCH] Ensure output table formatting does not contain NaNs. For (0
 eval ref calls)/(0 comp ref calls), set the proportion to 0.00.

Added integration tests (checked against manual tabulation)
---
 .../variantutils/GenotypeConcordance.java     | 25 +++++---
 .../GenotypeConcordanceIntegrationTest.java   | 63 +++++++++++++++++++
 2 files changed, 79 insertions(+), 9 deletions(-)
 create mode 100644 protected/java/test/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordanceIntegrationTest.java

diff --git a/protected/java/src/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordance.java b/protected/java/src/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordance.java
index de7b14ddb..ab137d4d5 100644
--- a/protected/java/src/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordance.java
+++ b/protected/java/src/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordance.java
@@ -85,6 +85,13 @@ public class GenotypeConcordance extends RodWalker<Pair<VariantContext,VariantCo
         return metrics;
     }
 
+    private static double repairNaN(double d) {
+     if ( Double.isNaN(d) ) {
+      return 0.0;
+     }
+     return d;
+    }
+
     public void onTraversalDone(ConcordanceMetrics metrics) {
         GATKReport report = new GATKReport();
         GATKReportTable concordanceCounts = new GATKReportTable("GenotypeConcordance_Counts","Per-sample concordance tables: comparison counts",2+GenotypeType.values().length*GenotypeType.values().length);
@@ -126,13 +133,13 @@ public class GenotypeConcordance extends RodWalker<Pair<VariantContext,VariantCo
                     int count = table.get(evalType, compType);
                     concordanceCounts.set(entry.getKey(),colKey,count);
                     if ( evalType == GenotypeType.HET || evalType == GenotypeType.HOM_REF || evalType == GenotypeType.HOM_VAR)
-                        concordanceEvalProportions.set(entry.getKey(),colKey,( (double) count)/table.getnEvalGenotypes(evalType));
+                        concordanceEvalProportions.set(entry.getKey(),colKey,repairNaN(( (double) count)/table.getnEvalGenotypes(evalType)));
                     if ( compType == GenotypeType.HET || compType == GenotypeType.HOM_VAR || compType == GenotypeType.HOM_REF )
-                        concordanceCompProportions.set(entry.getKey(),colKey,( (double) count)/table.getnCompGenotypes(compType));
+                        concordanceCompProportions.set(entry.getKey(),colKey,repairNaN(( (double) count)/table.getnCompGenotypes(compType)));
                 }
             }
-            concordanceEvalProportions.set(entry.getKey(),"Mismatching_Alleles", ( (double) table.getnMismatchingAlt() )/table.getnCalledEvalGenotypes());
-            concordanceCompProportions.set(entry.getKey(),"Mismatching_Alleles", ( (double) table.getnMismatchingAlt() )/table.getnCalledCompGenotypes());
+            concordanceEvalProportions.set(entry.getKey(),"Mismatching_Alleles", repairNaN(( (double) table.getnMismatchingAlt() )/table.getnCalledEvalGenotypes()));
+            concordanceCompProportions.set(entry.getKey(),"Mismatching_Alleles", repairNaN(( (double) table.getnMismatchingAlt() )/table.getnCalledCompGenotypes()));
             concordanceCounts.set(entry.getKey(),"Mismatching_Alleles",table.getnMismatchingAlt());
         }
 
@@ -147,13 +154,13 @@ public class GenotypeConcordance extends RodWalker<Pair<VariantContext,VariantCo
                 int count = table.get(evalType,compType);
                 concordanceCounts.set(rowKey,colKey,count);
                 if ( evalType == GenotypeType.HET || evalType == GenotypeType.HOM_REF || evalType == GenotypeType.HOM_VAR)
-                    concordanceEvalProportions.set(rowKey,colKey,( (double) count)/table.getnEvalGenotypes(evalType));
+                    concordanceEvalProportions.set(rowKey,colKey,repairNaN(( (double) count)/table.getnEvalGenotypes(evalType)));
                 if ( compType == GenotypeType.HET || compType == GenotypeType.HOM_VAR || compType == GenotypeType.HOM_REF )
-                    concordanceCompProportions.set(rowKey,colKey,( (double) count)/table.getnCompGenotypes(compType));
+                    concordanceCompProportions.set(rowKey,colKey,repairNaN(( (double) count)/table.getnCompGenotypes(compType)));
             }
         }
-        concordanceEvalProportions.set(rowKey,"Mismatching_Alleles", ( (double) table.getnMismatchingAlt() )/table.getnCalledEvalGenotypes());
-        concordanceCompProportions.set(rowKey,"Mismatching_Alleles", ( (double) table.getnMismatchingAlt() )/table.getnCalledCompGenotypes());
+        concordanceEvalProportions.set(rowKey,"Mismatching_Alleles", repairNaN(( (double) table.getnMismatchingAlt() )/table.getnCalledEvalGenotypes()));
+        concordanceCompProportions.set(rowKey,"Mismatching_Alleles", repairNaN(( (double) table.getnMismatchingAlt() )/table.getnCalledCompGenotypes()));
         concordanceCounts.set(rowKey,"Mismatching_Alleles",table.getnMismatchingAlt());
 
         for ( Map.Entry<String,Double> nrsEntry : metrics.getPerSampleNRS().entrySet() ) {
@@ -210,4 +217,4 @@ public class GenotypeConcordance extends RodWalker<Pair<VariantContext,VariantCo
         builder.genotypes(newGeno);
         return builder.make();
     }
-}
\ No newline at end of file
+}
diff --git a/protected/java/test/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordanceIntegrationTest.java b/protected/java/test/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordanceIntegrationTest.java
new file mode 100644
index 000000000..d4d8d6f8c
--- /dev/null
+++ b/protected/java/test/org/broadinstitute/sting/gatk/walkers/variantutils/GenotypeConcordanceIntegrationTest.java
@@ -0,0 +1,63 @@
+/*
+ * Copyright (c) 2010.
+ *
+ * Permission is hereby granted, free of charge, to any person
+ * obtaining a copy of this software and associated documentation
+ * files (the "Software"), to deal in the Software without
+ * restriction, including without limitation the rights to use,
+ * copy, modify, merge, publish, distribute, sublicense, and/or sell
+ * copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following
+ * conditions:
+ *
+ * The above copyright notice and this permission notice shall be
+ * included in all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
+ * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
+ * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
+ * NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
+ * HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
+ * WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR
+ * THE USE OR OTHER DEALINGS IN THE SOFTWARE.
+ */
+
+package org.broadinstitute.sting.gatk.walkers.variantutils;
+
+import org.broadinstitute.sting.WalkerTest;
+import org.broadinstitute.sting.utils.exceptions.UserException;
+import org.testng.annotations.Test;
+
+import java.util.Arrays;
+
+public class GenotypeConcordanceIntegrationTest extends WalkerTest {
+
+    protected static final String emptyMd5 = "d41d8cd98f00b204e9800998ecf8427e";
+
+    public static String baseTestString(String eval, String comp) {
+        return "-T GenotypeConcordance -R " + b37KGReference + " --eval " + validationDataLocation + eval + " --comp " + validationDataLocation + comp + " -o %s";
+    }
+
+    @Test
+    public void testIndelConcordance() {
+        WalkerTestSpec spec = new WalkerTestSpec(
+                baseTestString("NA12878.Jan2013.haplotypeCaller.subset.indels.vcf", "NA12878.Jan2013.bestPractices.subset.indels.vcf"),
+                0,
+                Arrays.asList("0f29a0c6dc44066228c8cb204fd53ec0")
+        );
+
+        executeTest("test indel concordance", spec);
+    }
+    
+    @Test
+    public void testNonoverlapingSamples() {
+        WalkerTestSpec spec = new WalkerTestSpec(
+                baseTestString("GenotypeConcordanceNonOverlapTest_Eval.vcf", "GenotypeConcordanceNonOverlapTest_Comp.vcf"),
+                0,
+                Arrays.asList("fc725022d47b4b5f8a6ef87f0f1ffe89")
+        );
+
+        executeTest("test non-overlapping samples", spec);
+    }
+}