Add initialise command and test

Author: jeromekelleher

sc2ts/cli.py

@@ -181,6 +181,38 @@ def add_provenance(ts, output_file):
     logger.info(f"Wrote {output_file}")
 
 
+@click.command()
+@click.argument("ts", type=click.Path(dir_okay=False))
+@click.argument("match_db", type=click.Path(dir_okay=False))
+@click.option(
+    "--additional-problematic-sites",
+    default=None,
+    type=click.Path(exists=True, dir_okay=False),
+    help="File containing the list of additional problematic sites to exclude.",
+)
+@click.option("-v", "--verbose", count=True)
+@click.option("-l", "--log-file", default=None, type=click.Path(dir_okay=False))
+def initialise(ts, match_db, additional_problematic_sites, verbose, log_file):
+    """
+    Initialise a new base tree sequence to begin inference.
+    """
+    setup_logging(verbose, log_file)
+
+    additional_problematic = []
+    if additional_problematic_sites is not None:
+        additional_problematic = (
+            np.loadtxt(additional_problematic_sites).astype(int).tolist()
+        )
+        logger.info(
+            f"Excluding additional {len(additional_problematic)} problematic sites"
+        )
+
+    base_ts = inference.initial_ts(additional_problematic)
+    base_ts.dump(ts)
+    logger.info(f"New base ts at {ts}")
+    inference.MatchDb.initialise(match_db)
+
+
 @click.command()
 @click.argument("alignments", type=click.Path(exists=True, dir_okay=False))
 @click.argument("metadata", type=click.Path(exists=True, dir_okay=False))
@@ -202,7 +234,6 @@ def add_provenance(ts, output_file):
     default=10,
     type=int,
     help="Minimum size of groups of reconsidered samples",
-    show_default=True
 )
 @click.option(
     "--num-past-days",
@@ -477,6 +508,7 @@ def cli():
 cli.add_command(export_alignments)
 cli.add_command(export_metadata)
 
+cli.add_command(initialise)
 cli.add_command(daily_extend)
 cli.add_command(validate)
 cli.add_command(annotate_recombinants)

sc2ts/inference.py

@@ -165,6 +165,7 @@ def initialise(db_path):
             conn.execute(
                 "CREATE INDEX [ix_samples_match_date] on 'samples' " "([match_date]);"
             )
+        logger.info(f"Created new MatchDb at {db_path}")
         return MatchDb(db_path)
 
     def print_all(self):

tests/test_cli.py

@@ -0,0 +1,26 @@
+
+import click.testing as ct
+import pytest
+import tskit
+
+import sc2ts
+from sc2ts import __main__ as main
+from sc2ts import cli
+
+class TestInitialise:
+
+    def test_defaults(self, tmp_path):
+        ts_path = tmp_path / "trees.ts"
+        match_db_path = tmp_path / "match.db"
+        runner = ct.CliRunner(mix_stderr=False)
+        result = runner.invoke(
+            cli.cli,
+            f"initialise {ts_path} {match_db_path}",
+            catch_exceptions=False,
+        )
+        assert result.exit_code == 0
+        ts = tskit.load(ts_path)
+        other_ts = sc2ts.initial_ts()
+        other_ts.tables.assert_equals(ts.tables)
+        match_db = sc2ts.MatchDb(match_db_path)
+        assert len(match_db) == 0

tests/test_inference.py

@@ -486,6 +486,7 @@ def test_high_recomb_mutation(self):
 
 
 class TestRealData:
+    # def test_first_day(self, tmp_path, alignments_store, metadata_db):
     def test_first_day(self, tmp_path, alignments_store, metadata_db):
         ts = sc2ts.extend(
             alignment_store=alignments_store,