get male-only var IDs from ped file and GTs not log p-value (#254)

epiercehoffman · web-flow · commit ba508ac3ffa2 · 2021-11-05T17:02:38.000-04:00
diff --git a/wdl/GATKSVPipelinePhase1.wdl b/wdl/GATKSVPipelinePhase1.wdl
@@ -223,6 +223,7 @@ workflow GATKSVPipelinePhase1 {
     RuntimeAttr? runtime_attr_exclude_outliers
     RuntimeAttr? runtime_attr_cat_outliers
     RuntimeAttr? runtime_attr_filter_samples
+    RuntimeAttr? runtime_attr_get_male_only
 
     ############################################################
     ## Module metrics parameters for GatherBatchEvidence, ClusterBatch, GenerateBatchMetrics, FilterBatch metrics
diff --git a/wdl/GenerateBatchMetrics.wdl b/wdl/GenerateBatchMetrics.wdl
@@ -65,6 +65,7 @@ workflow GenerateBatchMetrics {
     RuntimeAttr? runtime_attr_merge_allo
     RuntimeAttr? runtime_attr_merge_baf
     RuntimeAttr? runtime_attr_merge_stats
+    RuntimeAttr? runtime_attr_get_male_only
   }
 
   Array[String] algorithms = ["depth", "melt", "delly", "wham", "manta"]
@@ -101,6 +102,15 @@ workflow GenerateBatchMetrics {
       String algorithm = algorithms[i]
       File vcf = select_first([vcfs[i]])
 
+      call GetMaleOnlyVariantIDs {
+        input:
+          vcf = vcf,
+          female_samples = GetSampleLists.female_samples,
+          male_samples = GetSampleLists.male_samples,
+          sv_pipeline_docker = sv_pipeline_docker,
+          runtime_attr_override = runtime_attr_get_male_only
+      }
+
       if (algorithm != "melt") {
         call rdt.RDTest as RDTest {
           input:
@@ -118,6 +128,7 @@ workflow GenerateBatchMetrics {
             samples = GetSampleLists.samples_file,
             male_samples = GetSampleLists.male_samples,
             female_samples = GetSampleLists.female_samples,
+            male_only_variant_ids = GetMaleOnlyVariantIDs.male_only_variant_ids,
             sv_pipeline_docker = sv_pipeline_docker,
             sv_pipeline_rdtest_docker = sv_pipeline_rdtest_docker,
             linux_docker = linux_docker,
@@ -162,6 +173,7 @@ workflow GenerateBatchMetrics {
             samples = GetSampleLists.samples_file,
             male_samples = GetSampleLists.male_samples,
             female_samples = GetSampleLists.female_samples,
+            male_only_variant_ids = GetMaleOnlyVariantIDs.male_only_variant_ids,
             run_common = true,
             common_cnv_size_cutoff = common_cnv_size_cutoff,
             sv_base_mini_docker = sv_base_mini_docker,
@@ -190,6 +202,7 @@ workflow GenerateBatchMetrics {
             samples = GetSampleLists.samples_file,
             male_samples = GetSampleLists.male_samples,
             female_samples = GetSampleLists.female_samples,
+            male_only_variant_ids = GetMaleOnlyVariantIDs.male_only_variant_ids,
             common_cnv_size_cutoff = common_cnv_size_cutoff,
             sv_base_mini_docker = sv_base_mini_docker,
             linux_docker = linux_docker,
@@ -331,6 +344,46 @@ task GetSampleLists {
   }
 }
 
+task GetMaleOnlyVariantIDs {
+  input {
+    File vcf
+    File female_samples
+    File male_samples
+    String sv_pipeline_docker
+    RuntimeAttr? runtime_attr_override
+  }
+
+  RuntimeAttr default_attr = object {
+    cpu_cores: 1, 
+    mem_gb: 3.75,
+    disk_gb: 10,
+    boot_disk_gb: 10,
+    preemptible_tries: 3,
+    max_retries: 1
+  }
+  RuntimeAttr runtime_attr = select_first([runtime_attr_override, default_attr])
+
+  output {
+    File male_only_variant_ids = "male_only_variant_ids.txt"
+  }
+  command <<<
+    bcftools view -S ~{male_samples} ~{vcf} | bcftools view --min-ac 1 | bcftools query -f '%ID\n' > variant_ids_in_males.txt
+    bcftools view -S ~{female_samples} ~{vcf} | bcftools view --min-ac 1 | bcftools query -f '%ID\n' > variant_ids_in_females.txt
+    awk 'NR==FNR{a[$0];next} !($0 in a)' variant_ids_in_females.txt variant_ids_in_males.txt > male_only_variant_ids.txt
+    
+  >>>
+  runtime {
+    cpu: select_first([runtime_attr.cpu_cores, default_attr.cpu_cores])
+    memory: select_first([runtime_attr.mem_gb, default_attr.mem_gb]) + " GiB"
+    disks: "local-disk " + select_first([runtime_attr.disk_gb, default_attr.disk_gb]) + " HDD"
+    bootDiskSizeGb: select_first([runtime_attr.boot_disk_gb, default_attr.boot_disk_gb])
+    docker: sv_pipeline_docker
+    preemptible: select_first([runtime_attr.preemptible_tries, default_attr.preemptible_tries])
+    maxRetries: select_first([runtime_attr.max_retries, default_attr.max_retries])
+  }
+}
+
+
 task AggregateTests {
   input {
     File vcf
diff --git a/wdl/PETest.wdl b/wdl/PETest.wdl
@@ -17,6 +17,7 @@ workflow PETest {
     File ped_file
     File male_samples
     File female_samples
+    File male_only_variant_ids
     File samples
     Int common_cnv_size_cutoff
 
@@ -47,6 +48,7 @@ workflow PETest {
         samples = samples,
         male_samples = male_samples,
         female_samples = female_samples,
+        male_only_variant_ids = male_only_variant_ids,
         allosome = false,
         ref_dict = ref_dict,
         common_cnv_size_cutoff = common_cnv_size_cutoff,
@@ -75,6 +77,7 @@ workflow PETest {
         samples = samples,
         male_samples = male_samples,
         female_samples = female_samples,
+        male_only_variant_ids = male_only_variant_ids,
         allosome = true,
         ref_dict = ref_dict,
         common_cnv_size_cutoff = common_cnv_size_cutoff,
diff --git a/wdl/PETestChromosome.wdl b/wdl/PETestChromosome.wdl
@@ -16,6 +16,7 @@ workflow PETestChromosome {
     Int? suffix_len
     File male_samples
     File female_samples
+    File male_only_variant_ids
     File samples
     Boolean allosome
     Int common_cnv_size_cutoff
@@ -77,10 +78,11 @@ workflow PETestChromosome {
         input:
           male_test = PETestMale.stats,
           female_test = PETestFemale.stats,
+          male_only_ids_list = male_only_variant_ids,
           chrom = chrom,
           sv_pipeline_docker = sv_pipeline_docker,
           runtime_attr_override = runtime_attr_merge_allo,
-          male_only_expr = "females.log_pval == 0"
+          male_only_expr = "females.name.isin(male_only_ids)"
       }
     }
     if (!allosome) {
diff --git a/wdl/RDTest.wdl b/wdl/RDTest.wdl
@@ -17,6 +17,7 @@ workflow RDTest {
     File ped_file
     File male_samples
     File female_samples
+    File male_only_variant_ids
     File samples
     File ref_dict
 
@@ -47,6 +48,7 @@ workflow RDTest {
         samples = samples,
         male_samples = male_samples,
         female_samples = female_samples,
+        male_only_variant_ids = male_only_variant_ids,
         allosome = false,
         ref_dict = ref_dict,
         sv_pipeline_docker = sv_pipeline_docker,
@@ -74,6 +76,7 @@ workflow RDTest {
         samples = samples,
         male_samples = male_samples,
         female_samples = female_samples,
+        male_only_variant_ids = male_only_variant_ids,
         allosome = true,
         ref_dict = ref_dict,
         sv_pipeline_docker = sv_pipeline_docker,
diff --git a/wdl/RDTestChromosome.wdl b/wdl/RDTestChromosome.wdl
@@ -16,6 +16,7 @@ workflow RDTestChromosome {
     Int? suffix_len
     File male_samples
     File female_samples
+    File male_only_variant_ids
     File samples
     Boolean allosome
     File ref_dict
@@ -79,10 +80,11 @@ workflow RDTestChromosome {
         input:
           male_test = RDTestMale.stats,
           female_test = RDTestFemale.stats,
+          male_only_ids_list = male_only_variant_ids,
           chrom = chrom,
           sv_pipeline_docker = sv_pipeline_docker,
           runtime_attr_override = runtime_attr_merge_allo,
-          male_only_expr = "females.P.astype(str) == 'No_samples_for_analysis'"
+          male_only_expr = "females.CNVID.isin(male_only_ids)"
       }
     }
 
diff --git a/wdl/SRTest.wdl b/wdl/SRTest.wdl
@@ -16,6 +16,7 @@ workflow SRTest {
     String algorithm
     File male_samples
     File female_samples
+    File male_only_variant_ids
     File samples
     Boolean run_common
     Int? common_cnv_size_cutoff  # Required if run_common is true
@@ -49,6 +50,7 @@ workflow SRTest {
         samples = samples,
         male_samples = male_samples,
         female_samples = female_samples,
+        male_only_variant_ids = male_only_variant_ids,
         allosome = false,
         run_common = run_common,
         common_cnv_size_cutoff = common_cnv_size_cutoff,
@@ -78,6 +80,7 @@ workflow SRTest {
         samples = samples,
         male_samples = male_samples,
         female_samples = female_samples,
+        male_only_variant_ids = male_only_variant_ids,
         allosome = true,
         run_common = run_common,
         common_cnv_size_cutoff = common_cnv_size_cutoff,
diff --git a/wdl/SRTestChromosome.wdl b/wdl/SRTestChromosome.wdl
@@ -15,6 +15,7 @@ workflow SRTestChromosome {
     Int? suffix_len
     File male_samples
     File female_samples
+    File male_only_variant_ids
     File samples
     File ref_dict
     Boolean allosome
@@ -78,10 +79,11 @@ workflow SRTestChromosome {
         input:
           male_test = SRTestMale.stats,
           female_test = SRTestFemale.stats,
+          male_only_ids_list = male_only_variant_ids,
           chrom = chrom,
           runtime_attr_override = runtime_attr_merge_allo,
           sv_pipeline_docker = sv_pipeline_docker,
-          male_only_expr = "females.log_pval == 0"
+          male_only_expr = "females.name.isin(male_only_ids)"
       }
     } 
 
diff --git a/wdl/TasksGenerateBatchMetrics.wdl b/wdl/TasksGenerateBatchMetrics.wdl
@@ -117,6 +117,7 @@ task MergeAllosomes {
   input {
     File male_test
     File female_test
+    File male_only_ids_list
     String chrom
     String male_only_expr
     String sv_pipeline_docker
@@ -143,6 +144,10 @@ task MergeAllosomes {
     import pandas as pd
     males = pd.read_table("~{male_test}")
     females = pd.read_table("~{female_test}")
+    male_only_ids = set()
+    with open("~{male_only_ids_list}", 'r') as male_only_file:
+      for line in male_only_file:
+        male_only_ids.add(line.strip())
     if "~{chrom}" == 'Y' or "~{chrom}" == 'chrY':
       males.to_csv("~{basename(male_test)}.merged.csv", sep='\t', index=False, na_rep='NA')
     else:
