Convert original ngsrelate rule to one that is triggered by asking fo…

…r ngsrelate with ibsrelate method only. This also pars down the output to only the usable columns with the frequency free method.
zjnolen · Sep 30, 2024 · 7b731a9 · 7b731a9
1 parent 775d130
commit 7b731a9
Show file tree

Hide file tree

Showing 4 changed files with 70 additions and 21 deletions.
diff --git a/.test/config/config.yaml b/.test/config/config.yaml
@@ -55,10 +55,11 @@ analyses:
   ld_decay: true
   pca_pcangsd: true
   admix_ngsadmix: true
-  relatedness: 
-    ngsrelate: true
+  relatedness:
     ibsrelate_ibs: true
     ibsrelate_sfs: true
+    ngsrelate_ibsrelate-only: true
+    ngsrelate_freqbased: false
   1dsfs: true
   1dsfs_boot: true
   2dsfs: true
@@ -85,9 +86,10 @@ subsample_analyses:
   pca_pcangsd: false
   admix_ngsadmix: true
   relatedness:
-    ngsrelate: true
     ibsrelate_ibs: true
     ibsrelate_sfs: true
+    ngsrelate_ibsrelate-only: true
+    ngsrelate_freqbased: false
   1dsfs: true
   1dsfs_boot: true
   2dsfs: true

diff --git a/config/config.yaml b/config/config.yaml
@@ -55,10 +55,11 @@ analyses:
   ld_decay: false
   pca_pcangsd: false
   admix_ngsadmix: false
-  relatedness: 
-    ngsrelate: false
+  relatedness:
     ibsrelate_ibs: false
     ibsrelate_sfs: false
+    ngsrelate_ibsrelate-only: false
+    ngsrelate_freqbased: false
   1dsfs: false
   1dsfs_boot: false
   2dsfs: false
@@ -85,9 +86,10 @@ subsample_analyses:
   pca_pcangsd: false
   admix_ngsadmix: false
   relatedness:
-    ngsrelate: false
     ibsrelate_ibs: false
     ibsrelate_sfs: false
+    ngsrelate_ibsrelate-only: false
+    ngsrelate_freqbased: false
   1dsfs: false
   1dsfs_boot: false
   2dsfs: false

diff --git a/workflow/Snakefile b/workflow/Snakefile
@@ -179,9 +179,9 @@ if config["analyses"]["2dsfs_boot"]:
         )
     )
 
-if config["analyses"]["relatedness"]["ngsrelate"]:
+if config["analyses"]["relatedness"]["ngsrelate_ibsrelate-only"]:
     all_outputs.append(
-        "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_all_{sites}-filts_relate.html"
+        "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_all_{sites}-filts_ibsrelate-nofreq.html"
     )
 
 if config["analyses"]["relatedness"]["ibsrelate_sfs"]:
@@ -346,10 +346,10 @@ if config["subsample_analyses"]["2dsfs_boot"]:
         )
     )
 
-if config["subsample_analyses"]["relatedness"]["ngsrelate"]:
+if config["subsample_analyses"]["relatedness"]["ngsrelate_ibsrelate-only"]:
     all_outputs.extend(
         expand(
-            "results/datasets/{{dataset}}/analyses/kinship/ngsrelate/{{dataset}}.{{ref}}_all{dp}_{{sites}}-filts_relate.html",
+            "results/datasets/{{dataset}}/analyses/kinship/ngsrelate/{{dataset}}.{{ref}}_all{dp}_{{sites}}-filts_ibsrelate-nofreq.html",
             dp=subdp,
         )
     )

diff --git a/workflow/rules/5.0_relatedness.smk b/workflow/rules/5.0_relatedness.smk
@@ -157,20 +157,25 @@ rule kinship_table_html:
         "../scripts/tsv2html.Rmd"
 
 
-rule ngsrelate:
+rule ngsrelate_ibsrelate_only:
     """
-    Estimates inbreeding and relatedness measures using NGSrelate.
+    Estimates inbreeding and relatedness measures using NGSrelate. Only the
+    coefficients for the IBS relate method will be estimated and kept, which
+    does not require allele frequencies.
     """
     input:
         unpack(get_ngsrelate_input),
         inds="results/datasets/{dataset}/poplists/{dataset}_{population}{dp}.indiv.list",
     output:
-        relate="results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_relate.tsv",
-        samples="results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_samples.list",
+        ngsrelate=temp(
+            "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_ngsrelate-nofreq.tsv"
+        ),
+        ibsrelate="results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_ibsrelate-nofreq.tsv",
+        samples="results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_samples.ibsrelate-nofreq.list",
     wildcard_constraints:
         population="all",
     log:
-        "logs/{dataset}/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts.log",
+        "logs/{dataset}/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts.ibsrelate-nofreq.log",
     container:
         ngsrelate_container
     threads: lambda wildcards, attempt: attempt * 4
@@ -184,26 +189,66 @@ rule ngsrelate:
         echo "nsites nind"
         echo $nsites {params.nind}
         cut -f1 {input.inds} | tail -n +2 > {output.samples}
-        ngsRelate -G {input.beagle} -n {params.nind} -L $nsites -O {output.relate} \
-            -z {output.samples}) &> {log}
+        ngsRelate -G {input.beagle} -n {params.nind} -L $nsites \
+            -O {output.relate} -z {output.samples}
+        cut -f3-5,30-35 {output.ngsrelate} > {output.ibsrelate}) &> {log}
         """
 
 
+# rule ngsrelate_freqbased:
+#     """
+#     Estimates inbreeding and relatedness measures using NGSrelate. This will use
+#     allele frequencies, enabling all coefficients, including IBS relate ones.
+#     All will be kept.
+#     """
+#     input:
+#         unpack(get_ngsrelate_input),
+#         inds="results/datasets/{dataset}/poplists/{dataset}_{population}{dp}.indiv.list",
+#     output:
+#         relate="results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_relate.tsv",
+#         samples="results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts_samples.list",
+#     wildcard_constraints:
+#         population="all",
+#     log:
+#         "logs/{dataset}/kinship/ngsrelate/{dataset}.{ref}_{population}{dp}_{sites}-filts.log",
+#     container:
+#         ngsrelate_container
+#     threads: lambda wildcards, attempt: attempt * 4
+#     params:
+#         nind=get_nind,
+#     resources:
+#         runtime=lambda wildcards, attempt: attempt * 360,
+#     shell:
+#         r"""
+#         (nsites=$(zcat {input.beagle} | tail -n +2 | wc -l)
+#         echo "nsites nind"
+#         echo $nsites {params.nind}
+#         cut -f1 {input.inds} | tail -n +2 > {output.samples}
+#         ngsRelate -G {input.beagle} -n {params.nind} -L $nsites -O {output.relate} \
+#             -z {output.samples}) &> {log}
+#         """
+
+
 rule ngsrelate_summary:
     """
     Converts NGSrelate table to html.
     """
     input:
-        "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_all{dp}_{sites}-filts_relate.tsv",
+        "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_all{dp}_{sites}-filts_{method}.tsv",
     output:
         report(
-            "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_all{dp}_{sites}-filts_relate.html",
+            "results/datasets/{dataset}/analyses/kinship/ngsrelate/{dataset}.{ref}_all{dp}_{sites}-filts_{method}.html",
             category="02 Relatedness",
             subcategory="NgsRelate",
-            labels=lambda w: {"Filter": "{sites}", **dp_report(w), "Type": "Table"},
+            labels=lambda w: {
+                "Filter": "{sites}",
+                **dp_report(w),
+                "Method": "{method}",
+                "Type": "Table",
+            },
         ),
     log:
-        "logs/{dataset}/kinship/ngsrelate/{dataset}.{ref}_all{dp}_{sites}-filts_tsv2html.log",
+        "logs/{dataset}/kinship/ngsrelate/{dataset}.{ref}_all{dp}_{sites}-filts_{method}.tsv2html.log",
     container:
         r_container
     script: