From 90d6c2286d76c41a1a0ec781a40363168b8352c1 Mon Sep 17 00:00:00 2001
From: Leonardo Schwarz <leonardo.schwarz@fgcz.ethz.ch>
Date: Fri, 7 Jun 2024 12:40:17 +0200
Subject: [PATCH] more clearly label the plots

---
 .../example_compare/run_compare.py            |  4 ++-
 .../workflow/exp/plot_compare_peak_density.py | 26 +++++++++++++++++++
 .../workflow/experimental.smk                 | 13 +++++++++-
 .../workflow/qc/plot_peak_density.py          |  3 ++-
 .../workflow/rules/rules_proc.smk             |  3 ++-
 5 files changed, 45 insertions(+), 4 deletions(-)
 create mode 100644 src/depiction_targeted_preproc/workflow/exp/plot_compare_peak_density.py

diff --git a/src/depiction_targeted_preproc/example_compare/run_compare.py b/src/depiction_targeted_preproc/example_compare/run_compare.py
index 1de19b9..32f18b1 100644
--- a/src/depiction_targeted_preproc/example_compare/run_compare.py
+++ b/src/depiction_targeted_preproc/example_compare/run_compare.py
@@ -16,7 +16,9 @@ def prepare_tasks(input_imzml_path: Path, work_dir: Path) -> list[Path]:
     folders = set_up_work_dir(work_dir, input_imzml_path, input_mass_list)
     requested_files = get_all_output_files(folders)
 
-    exp_files = [work_dir / input_imzml_path.stem / "exp_compare_cluster_stats.pdf"]
+    combined_dir = work_dir / input_imzml_path.stem
+    exp_files = [combined_dir / "exp_compare_cluster_stats.pdf",
+                 combined_dir / "exp_plot_compare_peak_density.pdf"]
     return requested_files + exp_files
 
 
diff --git a/src/depiction_targeted_preproc/workflow/exp/plot_compare_peak_density.py b/src/depiction_targeted_preproc/workflow/exp/plot_compare_peak_density.py
new file mode 100644
index 0000000..aaba086
--- /dev/null
+++ b/src/depiction_targeted_preproc/workflow/exp/plot_compare_peak_density.py
@@ -0,0 +1,26 @@
+from pathlib import Path
+from typing import Annotated
+
+import polars as pl
+import typer
+import vegafusion
+from typer import Option, Argument
+
+from depiction_targeted_preproc.workflow.qc.plot_peak_density import plot_density_combined
+
+
+def exp_plot_compare_peak_density(
+    tables_marker_distances_calib: Annotated[list[Path], Argument()],
+    output_pdf: Annotated[Path, Option()],
+) -> None:
+    vegafusion.enable()
+
+    table = pl.concat(
+        [pl.read_parquet(path).with_columns(variant=pl.lit(path.parents[1].name)) for path in
+         tables_marker_distances_calib])
+
+    plot_density_combined(df_peak_dist=table, out_pdf=output_pdf)
+
+
+if __name__ == "__main__":
+    typer.run(exp_plot_compare_peak_density)
diff --git a/src/depiction_targeted_preproc/workflow/experimental.smk b/src/depiction_targeted_preproc/workflow/experimental.smk
index 41e3c68..ed90050 100644
--- a/src/depiction_targeted_preproc/workflow/experimental.smk
+++ b/src/depiction_targeted_preproc/workflow/experimental.smk
@@ -29,4 +29,15 @@ rule exp_mass_list_preparation:
         " --input-csv-path {input.csv}"
         " --out-calibration-csv-path {output.calibration_csv}"
         " --out-standards-csv-path {output.standards_csv}"
-        " --out-visualization-csv-path {output.visualization_csv}"
\ No newline at end of file
+        " --out-visualization-csv-path {output.visualization_csv}"
+
+
+rule exp_plot_compare_peak_density:
+    input:
+        tables_marker_distance=expand("{{sample}}/{exp_variant}/qc/table_marker_distances_calib.parquet", exp_variant=exp_variants)
+    output:
+        pdf="{sample}/exp_plot_compare_peak_density.pdf"
+    shell:
+        "python -m depiction_targeted_preproc.workflow.exp.plot_compare_peak_density"
+        " {input.tables_marker_distance}"
+        " --output-pdf {output.pdf}"
\ No newline at end of file
diff --git a/src/depiction_targeted_preproc/workflow/qc/plot_peak_density.py b/src/depiction_targeted_preproc/workflow/qc/plot_peak_density.py
index 5816ece..de5f6b9 100644
--- a/src/depiction_targeted_preproc/workflow/qc/plot_peak_density.py
+++ b/src/depiction_targeted_preproc/workflow/qc/plot_peak_density.py
@@ -16,6 +16,7 @@ def subsample_dataframe(df: pl.DataFrame) -> pl.DataFrame:
 
 
 def plot_density_combined(df_peak_dist: pl.DataFrame, out_pdf: Path) -> None:
+    n_tot = len(df_peak_dist)
     df_peak_dist = subsample_dataframe(df_peak_dist)
     chart = (
         (
@@ -30,7 +31,7 @@ def plot_density_combined(df_peak_dist: pl.DataFrame, out_pdf: Path) -> None:
         .encode(y=alt.Y("density:Q"))
         .properties(title="Linear scale")
     )
-    chart = chart.properties(title="Density of target-surrounding peak distances")
+    chart = chart.properties(title=f"Density of target-surrounding peak distances (n_tot = {n_tot} sampled to n = {len(df_peak_dist)})")
     chart.save(out_pdf)
 
 
diff --git a/src/depiction_targeted_preproc/workflow/rules/rules_proc.smk b/src/depiction_targeted_preproc/workflow/rules/rules_proc.smk
index 00168bc..ce820b2 100644
--- a/src/depiction_targeted_preproc/workflow/rules/rules_proc.smk
+++ b/src/depiction_targeted_preproc/workflow/rules/rules_proc.smk
@@ -27,7 +27,8 @@ rule proc_calibrate_remove_global_shift:
     input:
         imzml=multiext("{sample}/corrected.peaks",".imzML",".ibd"),
         config="{sample}/pipeline_params.yml",
-        mass_list="{sample}/mass_list.standards.csv",
+        #mass_list="{sample}/mass_list.standards.csv",
+        mass_list="{sample}/mass_list.calibration.csv",
     output:
         imzml=temp(multiext("{sample}/calibrated.tmp",".imzML",".ibd")),
     shell: