Fix core size application with freebayes. Fix test according.

metagenopolis · Sep 2, 2024 · f924ba5 · f924ba5
1 parent d820b2a
commit f924ba5
Show file tree

Hide file tree

Showing 10 changed files with 30,039 additions and 20 deletions.
diff --git a/meteor/strain.py b/meteor/strain.py
@@ -36,7 +36,7 @@ class Strain(Session):
     MIN_MIN_SNP_DEPTH: ClassVar[int] = 1
     MAX_MIN_SNP_DEPTH: ClassVar[int] = 10000
     DEFAULT_MIN_SNP_DEPTH: ClassVar[int] = 3
-    DEFAULT_MIN_FREQUENCY: ClassVar[float] = 0.8
+    DEFAULT_MIN_FREQUENCY: ClassVar[float] = 0.01
     DEFAULT_PLOIDY: ClassVar[int] = 1
     MIN_MIN_MSP_COVERAGE: ClassVar[int] = 1
     MAX_MIN_MSP_COVERAGE: ClassVar[int] = 100
@@ -79,7 +79,7 @@ def filter_coverage(
                 "startpos",
                 "endpos",
             ],
-            header=1,
+            # header=1,
         )
         cov_df = pd.read_csv(
             StringIO(
@@ -92,7 +92,7 @@ def filter_coverage(
                 )  # type: ignore[attr-defined]
             ),
             sep="\t",
-            header=1,
+            header=0,
             names=[
                 "gene_id",
                 "startpos",
@@ -146,8 +146,8 @@ def get_msp_variant(
         msp_content = pd.read_csv(
             msp_file,
             sep="\t",
-            names=["msp_name", "gene_id", "gene_name", "gene_category"],
-            header=1,
+            # names=["msp_name", "gene_id", "gene_name", "gene_category"],
+            header=0,
         )
         msp_content = msp_content.loc[msp_content["gene_category"] == "core"]
         msp_content = msp_content.groupby("msp_name").head(self.core_size).reset_index()
@@ -258,7 +258,6 @@ def execute(self) -> None:
             self.json_data["directory"] = stage3_dir
             self.json_data["census"] = census_json
             self.json_data["reference"] = ref_json
-            self.json_data["reference"] = ref_json
             self.json_data["Stage3FileName"] = (
                 stage3_dir / f"{sample_info['sample_name']}_census_stage_3.json"
             )
@@ -274,6 +273,7 @@ def execute(self) -> None:
                 self.min_snp_depth,
                 self.min_frequency,
                 self.ploidy,
+                self.core_size,
             )
             if self.json_data["Stage3FileName"].exists():
                 logging.info(

diff --git a/meteor/tests/test_strain.py b/meteor/tests/test_strain.py
@@ -37,7 +37,7 @@ def strain_builder(datadir: Path, tmp_path: Path) -> Strain:
     }
     # meteor.mapping_dir = tmp_path
     # return Strain(meteor, 100, 3, 3, 0.5, 1, 3, 0.8, True, False, json_data=json_data)
-    return Strain(meteor, 100, 3, 3, 0.5, 1, 1, 0.8, 10, True, json_data=json_data)
+    return Strain(meteor, 100, 3, 3, 0.5, 1, 1, 0.8, 100, True, json_data=json_data)
 
 
 def test_filter_coverage(
@@ -49,7 +49,11 @@ def test_filter_coverage(
         strain_builder.meteor.ref_dir / "database" / "mock.bed",
         strain_builder.meteor.ref_dir / "fasta" / "mock.fasta.gz",
     )
+    filtered_cov.to_csv("/pasteur/zeus/projets/p01/BioIT/amine/filtered.tsv", sep="\t")
     expected_output = pd.read_table(datadir / "expected_output" / "filtered_cov.tsv")
+    expected_output.to_csv(
+        "/pasteur/zeus/projets/p01/BioIT/amine/expected.tsv", sep="\t"
+    )
     assert filtered_cov.reset_index(drop=True).equals(expected_output)
 
 
@@ -87,15 +91,15 @@ def test_get_msp_variant(strain_builder, datadir: Path, tmp_path: Path):
     assert PA.exists()
     with BS.open("rb") as out:
         # assert md5(out.read()).hexdigest() == "958199bf4afa9adb56c9e0100a0cc23c"
-        assert md5(out.read()).hexdigest() == "b6b87174554ca33cc4a8923be4e780f7"
+        assert md5(out.read()).hexdigest() == "c0f943768ca8629b8be4a995bb3af341"
     with PA.open("rb") as out:
         # assert md5(out.read()).hexdigest() == "cec667aad0449853b3bd7db689dd7ed3"
-        assert md5(out.read()).hexdigest() == "b4eadbd604c765a2579bf742ae27386a"
+        assert md5(out.read()).hexdigest() == "cec667aad0449853b3bd7db689dd7ed3"
 
 
 def test_execute(strain_builder, tmp_path: Path) -> None:
     strain_builder.execute()
     BS = tmp_path / "strain" / "test" / "BS.fasta.xz"
     assert BS.exists()
-    # with BS.open("rb") as out:
-    #     assert md5(out.read()).hexdigest() == "83d3a5f89b89df438dc1751009801fc6"
+    with BS.open("rb") as out:
+        assert md5(out.read()).hexdigest() == "1da5a4bf1efd32b58e12e42ae971fe3f"