Skip to content
Snippets Groups Projects
Commit 0bca3f1e authored by Peter van 't Hof's avatar Peter van 't Hof
Browse files

Merge branch 'issue186' into 'develop'

Issue186

See #186. 

See merge request !209
parents 23bbd96e d55c2684
No related branches found
No related tags found
No related merge requests found
...@@ -17,7 +17,7 @@ package nl.lumc.sasc.biopet.tools ...@@ -17,7 +17,7 @@ package nl.lumc.sasc.biopet.tools
import java.io.File import java.io.File
import htsjdk.variant.variantcontext.VariantContext import htsjdk.variant.variantcontext.{ GenotypeType, VariantContext }
import htsjdk.variant.variantcontext.writer.{ AsyncVariantContextWriter, VariantContextWriterBuilder } import htsjdk.variant.variantcontext.writer.{ AsyncVariantContextWriter, VariantContextWriterBuilder }
import htsjdk.variant.vcf.VCFFileReader import htsjdk.variant.vcf.VCFFileReader
import nl.lumc.sasc.biopet.core.config.Configurable import nl.lumc.sasc.biopet.core.config.Configurable
...@@ -72,6 +72,7 @@ object VcfFilter extends ToolCommand { ...@@ -72,6 +72,7 @@ object VcfFilter extends ToolCommand {
minSamplesPass: Int = 1, minSamplesPass: Int = 1,
mustHaveVariant: List[String] = Nil, mustHaveVariant: List[String] = Nil,
calledIn: List[String] = Nil, calledIn: List[String] = Nil,
mustHaveGenotype: List[(String, GenotypeType)] = Nil,
deNovoInSample: String = null, deNovoInSample: String = null,
resToDom: List[Trio] = Nil, resToDom: List[Trio] = Nil,
trioCompound: List[Trio] = Nil, trioCompound: List[Trio] = Nil,
...@@ -126,6 +127,12 @@ object VcfFilter extends ToolCommand { ...@@ -126,6 +127,12 @@ object VcfFilter extends ToolCommand {
opt[String]("calledIn") unbounded () valueName "<sample>" action { (x, c) => opt[String]("calledIn") unbounded () valueName "<sample>" action { (x, c) =>
c.copy(calledIn = x :: c.calledIn) c.copy(calledIn = x :: c.calledIn)
} text "Must be called in this sample" } text "Must be called in this sample"
opt[String]("mustHaveGenotype") unbounded () valueName "<sample:genotype>" action { (x, c) =>
c.copy(mustHaveGenotype = (x.split(":")(0), GenotypeType.valueOf(x.split(":")(1))) :: c.mustHaveGenotype)
} validate { x =>
if (x.split(":").length == 2 && GenotypeType.values().contains(x.split(":")(1))) success
else failure("--mustHaveGenotype should be in this format: sample:genotype")
} text "Must have genotoype <genotype> for this sample. Genotype can be " + GenotypeType.values().mkString(", ")
opt[String]("diffGenotype") unbounded () valueName "<sample:sample>" action { (x, c) => opt[String]("diffGenotype") unbounded () valueName "<sample:sample>" action { (x, c) =>
c.copy(diffGenotype = (x.split(":")(0), x.split(":")(1)) :: c.diffGenotype) c.copy(diffGenotype = (x.split(":")(0), x.split(":")(1)) :: c.diffGenotype)
} validate { x => if (x.split(":").length == 2) success else failure("--notSameGenotype should be in this format: sample:sample") } validate { x => if (x.split(":").length == 2) success else failure("--notSameGenotype should be in this format: sample:sample")
...@@ -184,6 +191,7 @@ object VcfFilter extends ToolCommand { ...@@ -184,6 +191,7 @@ object VcfFilter extends ToolCommand {
minAlternateDepth(record, cmdArgs.minAlternateDepth, cmdArgs.minSamplesPass) && minAlternateDepth(record, cmdArgs.minAlternateDepth, cmdArgs.minSamplesPass) &&
(cmdArgs.mustHaveVariant.isEmpty || mustHaveVariant(record, cmdArgs.mustHaveVariant)) && (cmdArgs.mustHaveVariant.isEmpty || mustHaveVariant(record, cmdArgs.mustHaveVariant)) &&
calledIn(record, cmdArgs.calledIn) && calledIn(record, cmdArgs.calledIn) &&
hasGenotype(record, cmdArgs.mustHaveGenotype) &&
(cmdArgs.diffGenotype.isEmpty || cmdArgs.diffGenotype.forall(x => notSameGenotype(record, x._1, x._2))) && (cmdArgs.diffGenotype.isEmpty || cmdArgs.diffGenotype.forall(x => notSameGenotype(record, x._1, x._2))) &&
( (
cmdArgs.filterHetVarToHomVar.isEmpty || cmdArgs.filterHetVarToHomVar.isEmpty ||
...@@ -220,6 +228,16 @@ object VcfFilter extends ToolCommand { ...@@ -220,6 +228,16 @@ object VcfFilter extends ToolCommand {
else true else true
} }
/**
* Checks if given genotypes for given samples are there
* @param record VCF record
* @param samplesGenotypes samples and their associated genotypes to be checked (of format sample:genotype)
* @return false when filter fails
*/
def hasGenotype(record: VariantContext, samplesGenotypes: List[(String, GenotypeType)]): Boolean = {
samplesGenotypes.forall(x => record.getGenotype(x._1).getType == x._2)
}
/** /**
* Checks if record has atleast minQualScore * Checks if record has atleast minQualScore
* @param record VCF record * @param record VCF record
......
...@@ -15,8 +15,11 @@ ...@@ -15,8 +15,11 @@
*/ */
package nl.lumc.sasc.biopet.tools package nl.lumc.sasc.biopet.tools
import java.io.File
import java.nio.file.Paths import java.nio.file.Paths
import htsjdk.variant.variantcontext.GenotypeType
import htsjdk.variant.vcf.VCFFileReader
import org.scalatest.Matchers import org.scalatest.Matchers
import org.scalatest.mock.MockitoSugar import org.scalatest.mock.MockitoSugar
import org.scalatest.testng.TestNGSuite import org.scalatest.testng.TestNGSuite
...@@ -37,6 +40,7 @@ class VcfFilterTest extends TestNGSuite with MockitoSugar with Matchers { ...@@ -37,6 +40,7 @@ class VcfFilterTest extends TestNGSuite with MockitoSugar with Matchers {
} }
val vepped_path = resourcePath("/VEP_oneline.vcf") val vepped_path = resourcePath("/VEP_oneline.vcf")
val vepped = new File(vepped_path)
val rand = new Random() val rand = new Random()
@Test def testOutputTypeVcf() = { @Test def testOutputTypeVcf() = {
...@@ -57,4 +61,24 @@ class VcfFilterTest extends TestNGSuite with MockitoSugar with Matchers { ...@@ -57,4 +61,24 @@ class VcfFilterTest extends TestNGSuite with MockitoSugar with Matchers {
main(arguments) main(arguments)
} }
@Test def testHasGenotype() = {
val reader = new VCFFileReader(vepped, false)
val record = reader.iterator().next()
hasGenotype(record, List(("Child_7006504", GenotypeType.HET))) shouldBe true
hasGenotype(record, List(("Child_7006504", GenotypeType.HOM_VAR))) shouldBe false
hasGenotype(record, List(("Child_7006504", GenotypeType.HOM_REF))) shouldBe false
hasGenotype(record, List(("Child_7006504", GenotypeType.NO_CALL))) shouldBe false
hasGenotype(record, List(("Child_7006504", GenotypeType.MIXED))) shouldBe false
hasGenotype(record, List(("Mother_7006508", GenotypeType.HET))) shouldBe false
hasGenotype(record, List(("Mother_7006508", GenotypeType.HOM_VAR))) shouldBe false
hasGenotype(record, List(("Mother_7006508", GenotypeType.HOM_REF))) shouldBe true
hasGenotype(record, List(("Mother_7006508", GenotypeType.NO_CALL))) shouldBe false
hasGenotype(record, List(("Mother_7006508", GenotypeType.MIXED))) shouldBe false
hasGenotype(record, List(("Mother_7006508", GenotypeType.HOM_REF), ("Child_7006504", GenotypeType.HET))) shouldBe true
hasGenotype(record, List(("Mother_7006508", GenotypeType.HET), ("Child_7006504", GenotypeType.HOM_REF))) shouldBe false
}
} }
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment