diff --git a/Snake.config.json b/Snake.config.json index 6c41edd84776579cc1e7f9ae6f9c3df489cfb7d5..73afc042334f6e9dd892dd2806daccc55dbded9c 100644 --- a/Snake.config.json +++ b/Snake.config.json @@ -19,10 +19,9 @@ "R_reference" : "BSgenome.Hsapiens.UCSC.hg38", "bp_density" : { - "few" : 0.1, - "medium" : 0.25, - "more" : 0.33, - "many" : 0.5 + "few" : 0.2, + "medium" : 0.4, + "many" : 0.6 }, "paired_end" : true, diff --git a/Snakefile b/Snakefile index ef51fa31165512330bc013f90168031d476abd09..805bf3c33f0159ff72b118bc6a69ffdc44da3740 100644 --- a/Snakefile +++ b/Snakefile @@ -18,7 +18,7 @@ import os.path # * calculate a segmentation into potential SVs using Mosaicatcher -METHODS = ["simpleCalls", "biAllelic"] +METHODS = ["simpleCalls_llr1", "simpleCalls_llr4", "biAllelic_llr1", "biAllelic_llr4"] rule all: input: @@ -27,7 +27,7 @@ rule all: sample = SAMPLE, chrom = config["chromosomes"], window = [50000, 100000], - bpdens = ["few","medium","more","many"], + bpdens = ["few","medium","many"], method = METHODS) @@ -349,9 +349,9 @@ rule mosaiClassifier_make_call: input: probs = "sv_probabilities/{sample}/{windows}.{bpdens}/probabilities.Rdata" output: - "sv_calls/{sample}/{windows}.{bpdens}/simpleCalls.txt" + "sv_calls/{sample}/{windows}.{bpdens}/simpleCalls_llr{llr}.txt" log: - "log/{sample}/mosaiClassifier_make_call.{windows}.{bpdens}.txt" + "log/{sample}/mosaiClassifier_make_call.{windows}.{bpdens}.{llr}.txt" script: "utils/mosaiClassifier_call.snakemake.R" @@ -379,9 +379,9 @@ rule mosaiClassifier_make_call_biallelic: input: probs = "sv_probabilities/{sample}/{windows}.{bpdens}/probabilities.Rdata" output: - "sv_calls/{sample}/{windows}.{bpdens}/biAllelic.txt" + "sv_calls/{sample}/{windows}.{bpdens}/biAllelic_llr{llr}.txt" log: - "log/{sample}/mosaiClassifier_make_call_biallelic.{windows}.{bpdens}.txt" + "log/{sample}/mosaiClassifier_make_call_biallelic.{windows}.{bpdens}.{llr}.txt" script: "utils/mosaiClassifier_call_biallelic.snakemake.R" diff --git a/utils/mosaiClassifier.snakemake.R b/utils/mosaiClassifier.snakemake.R index 8471cbb609146f44fc3bebb6072a676c7a617730..60a716cecd666561ecec03b25094c28153fc9c48 100644 --- a/utils/mosaiClassifier.snakemake.R +++ b/utils/mosaiClassifier.snakemake.R @@ -12,7 +12,8 @@ segs = fread(snakemake@input[["bp"]]) # is there a normalization file given? -if ("norm" %in% snakemake@input) { +if ("norm" %in% names(snakemake@input)) { + message("[MosaiClassifier] Read normalization from ", snakemake@input[["norm"]]) normalization = fread(snakemake@input[["norm"]]) } else { normalization = NULL diff --git a/utils/mosaiClassifier/mosaiClassifier.R b/utils/mosaiClassifier/mosaiClassifier.R index d03ef9f4e3f8ea62e23f3b41490cf8875c7a50ac..547cb83277aaaf02e537d92d3d0af910dfdfab91 100644 --- a/utils/mosaiClassifier/mosaiClassifier.R +++ b/utils/mosaiClassifier/mosaiClassifier.R @@ -127,6 +127,7 @@ mosaiClassifierPrepare <- function(counts, info, strand, segs, normVector = NULL # Add normalization factors to the expected counts ("scalar") if (!is.null(normVector)) { + message("[MosaiClassifier] Normalize coverage expectation") addNormalizationScalar(probs, counts, normVector) } else { probs[, scalar := 1.0] diff --git a/utils/mosaiClassifier_call.snakemake.R b/utils/mosaiClassifier_call.snakemake.R index 19c3a35682b4b67491d2449da964417fde1c8aea..1cbd445dfeafa0ac3d2933761b1fe377b66814a7 100644 --- a/utils/mosaiClassifier_call.snakemake.R +++ b/utils/mosaiClassifier_call.snakemake.R @@ -3,9 +3,9 @@ library(data.table) source("utils/mosaiClassifier/makeSVcalls.R") probs = readRDS(snakemake@input[["probs"]]) -llr = as.numeric(snakemake@params[["llr"]]) +llr = as.numeric(snakemake@wildcards[["llr"]]) probs <- mosaiClassifierPostProcessing(probs) -tab <- makeSVCallSimple(probs) +tab <- makeSVCallSimple(probs, llr_thr = llr) write.table(tab, file = snakemake@output[[1]], sep = "\t", quote=F, row.names = F, col.names = T) diff --git a/utils/mosaiClassifier_call_biallelic.snakemake.R b/utils/mosaiClassifier_call_biallelic.snakemake.R index da228389974286ecc758b7144a70c8083562fb20..fa42cd821eb1953c5ca5bdf56b333da6e67087dd 100644 --- a/utils/mosaiClassifier_call_biallelic.snakemake.R +++ b/utils/mosaiClassifier_call_biallelic.snakemake.R @@ -3,10 +3,10 @@ library(data.table) source("utils/mosaiClassifier/makeSVcalls.R") probs = readRDS(snakemake@input[["probs"]]) -llr = as.numeric(snakemake@params[["llr"]]) +llr = as.numeric(snakemake@wildcards[["llr"]]) probs <- mosaiClassifierPostProcessing(probs) probs <- forceBiallelic(probs) -tab <- makeSVCallSimple(probs) +tab <- makeSVCallSimple(probs, llr_thr = llr) write.table(tab, file = snakemake@output[[1]], sep = "\t", quote=F, row.names = F, col.names = T)