optimization hyperparameter run on specific timestep

rules/plots.smk

@@ -1,15 +1,19 @@
 import os
 
 plot_config = config["plots"]
+
 input_replicates = ""
 for replicate in HMM["replicates"]:
     input_replicates += replicate + ","
+
 input_timesteps = ""
 for timestep in HMM["timesteps"]:
     input_timesteps += timestep + ","
+
 input_references = ""
 for reference in HMM["references"]:
     input_references += reference + ","
+
 for replicate in HMM["replicates"]:
     os.makedirs(f"{out_fld}/genomewide_recombination/{replicate}/", exist_ok=True)
     os.makedirs(f"{out_fld}/coverage_arrays/{replicate}/", exist_ok=True)
@@ -92,6 +96,9 @@ rule unique_plot:
         """
 
 HMM_config = config["HMM_parameters"]
+optimization_timesteps = ""
+for timestep in config["optimization_recombination_parameter"]["timesteps"]:
+    optimization_timesteps += timestep + ","
 
 rule optimize_recombination_parameter:
     input:
@@ -104,7 +111,7 @@ rule optimize_recombination_parameter:
         '../conda_envs/sci_py.yml'
     params:
         replicates = input_replicates,
-        timesteps = input_timesteps,
+        timesteps = optimization_timesteps,
         cores = HMM_config["cores"],
         initial_probability = HMM_config["initial_probability"]["A"]+","+HMM_config["initial_probability"]["B"],
         transition_probability = config["optimization_recombination_parameter"]["values"],

run_config.yml

@@ -44,8 +44,10 @@ HMM_parameters:
 
 optimization_recombination_parameter:
   flag: true
+  timesteps: #timestep data to use for optimization
+    - 'D10'
   values: '0.00001,0.00002,0.00003,0.00004,0.00005,0.00006,0.00008,0.0001' #comma separated values
-  subsample: 10000
+  subsample: 1000
 
 plots:
   coverage_threshold: 0

scripts/optimize_recombination_parameter.py

@@ -110,7 +110,8 @@ def get_evidence_arrays(evidences_file):
                 evidence_arrays = [evidence_arrays[i] for i in idx]
                 mapping_starts = [mapping_starts[i] for i in idx]
                 mapping_ends = [mapping_ends[i] for i in idx]
-                tot_reads+=subsample
+
+            tot_reads+=len(evidence_arrays)
 
             for prob in transition_probabilities: