Revise member selection (#54)

- clarify naming convenction, for `member_num` use either `member_id`,
`member_ids`, or `member_count`
- revise the member selection algorithm, instead of checking random it
uses this deterministic approach

At each step, every ensemble member is evaluated to determine how many
of the remaining members satisfy the given tolerances when compared to
it. The ensemble member that results in the highest number of other
members passing the tolerances is selected. This process is repeated
until either all members are successfully validated using the current
tolerances, or a maximum number of attempts is reached. In the latter
case, the tolerance factor is increased, and the process continues.

---------

Co-authored-by: Annika Lauber <annika.lauber@c2sm.ethz.ch>

Author: huppd

README.md

@@ -191,7 +191,7 @@ python ../externals/probtest/probtest.py tolerance
 These commands will generate a number of files:
 
 - `stats_ref.csv`: contains the post-processed output from the unperturbed reference run
-- `stats_{member_num}.csv`: contain the post-processed output from the perturbed reference runs (only needed temporarily to generate the tolerance file)
+- `stats_{member_id}.csv`: contain the post-processed output from the perturbed reference runs (only needed temporarily to generate the tolerance file)
 - `exp_name_tolerance.csv`: contains tolerance ranges computed from the stats-files
 
 These can then be used to compare against the output of a test binary (usually a

engine/cdo_table.py

@@ -14,7 +14,7 @@
 import xarray as xr
 
 from util import model_output_parser
-from util.click_util import CommaSeperatedInts, cli_help
+from util.click_util import cli_help
 from util.constants import cdo_bins
 from util.dataframe_ops import df_from_file_ids
 from util.file_system import file_names_from_pattern
@@ -93,10 +93,10 @@ def rel_diff_stats(
     help=cli_help["file_id"],
 )
 @click.option(
-    "--member-num",
-    type=CommaSeperatedInts(),
-    default="10",
-    help=cli_help["member_num"],
+    "--member-id",
+    type=int,
+    default=1,
+    help=cli_help["member_id"],
 )
 @click.option(
     "--member-type",
@@ -120,20 +120,17 @@ def rel_diff_stats(
 def cdo_table(
     model_output_dir,
     file_id,
-    member_num,
+    member_id: int,
     member_type,
     perturbed_model_output_dir,
     cdo_table_file,
     file_specification,
 ):  # pylint: disable=too-many-positional-arguments
-    # TODO: A single perturbed run provides enough data to make proper statistics.
-    #       refactor cdo_table interface to reflect that
-    if len(member_num) == 1:
-        member_num = list(range(1, member_num[0] + 1))
+
     if member_type:
-        member_id = member_type + "_" + str(member_num[0])
+        member_id_str = member_type + "_" + str(member_id)
     else:
-        member_id = str(member_num[0])
+        member_id_str = str(member_id)
 
     file_specification = file_specification[0]  # can't store dicts as defaults in click
     assert isinstance(file_specification, dict), "must be dict"
@@ -154,7 +151,7 @@ def cdo_table(
                 continue
             ref_files.sort()
             perturb_files, err = file_names_from_pattern(
-                perturbed_model_output_dir.format(member_id=member_id), file_pattern
+                perturbed_model_output_dir.format(member_id=member_id_str), file_pattern
             )
             if err > 0:
                 logger.info(
@@ -168,7 +165,7 @@ def cdo_table(
                     continue
                 ref_data = xr.open_dataset(f"{model_output_dir}/{rf}")
                 perturb_data = xr.open_dataset(
-                    f"{perturbed_model_output_dir.format(member_id=member_id)}/{pf}"
+                    f"{perturbed_model_output_dir.format(member_id=member_id_str)}/{pf}"
                 )
                 diff_data = ref_data.copy()
                 varnames = [

engine/check.py

@@ -12,7 +12,7 @@
 import click
 
 from util.click_util import cli_help
-from util.dataframe_ops import compute_div_dataframe, test_stats_file_with_tolerances
+from util.dataframe_ops import check_stats_file_with_tolerances, compute_div_dataframe
 from util.log_handler import logger
 
 
@@ -36,7 +36,7 @@
 )
 def check(input_file_ref, input_file_cur, tolerance_file_name, factor):
 
-    out, err, tol = test_stats_file_with_tolerances(
+    out, err, tol = check_stats_file_with_tolerances(
         tolerance_file_name, input_file_ref, input_file_cur, factor
     )
 

engine/init.py

@@ -43,10 +43,10 @@
     help=cli_help["config"],
 )
 @click.option(
-    "--member-num",
+    "--member-ids",
     type=CommaSeperatedInts(),
-    default="10",
-    help=cli_help["member_num"],
+    default="1,2,3,4,5,6,7,8,9,10",
+    help=cli_help["member_ids"],
 )
 @click.option(
     "--member-type",
@@ -89,14 +89,15 @@ def init(
     reference,
     config,
     template_name,
-    member_num,
+    member_ids,
     member_type,
     factor,
     perturb_amplitude,
     timing_current,
     timing_reference,
     append_time,
 ):  # pylint: disable=too-many-positional-arguments
+
     template_partition = str(template_name).rpartition("/")
     env = Environment(
         loader=FileSystemLoader(template_partition[0]), undefined=StrictUndefined
@@ -115,8 +116,8 @@ def init(
         logger.warning(warn_template, "file_id", "")
     if not reference:
         logger.warning(warn_template, "reference", "")
-    if not member_num:
-        logger.warning(warn_template, "member_num", member_num)
+    if not member_ids:
+        logger.warning(warn_template, "member_ids", member_ids)
     if not member_type:
         logger.warning(warn_template, "member_type", member_type)
     if not factor:
@@ -135,7 +136,7 @@ def init(
     render_dict["experiment_name"] = experiment_name
     render_dict["codebase_install"] = Path(codebase_install).resolve()
     render_dict["reference"] = Path(reference).resolve()
-    render_dict["member_num"] = member_num
+    render_dict["member_ids"] = member_ids
     render_dict["member_type"] = member_type
     render_dict["factor"] = factor
     render_dict["perturb_amplitude"] = perturb_amplitude
@@ -148,7 +149,7 @@ def init(
     # append file_id via json
     json_dict = json.loads(rendered)
     json_dict["default"]["file_id"] = file_id
-    json_dict["default"]["member_num"] = member_num
+    json_dict["default"]["member_ids"] = member_ids
     json_dict["default"]["factor"] = factor
     rendered = json.dumps(json_dict, indent=2)
     # print file

engine/perturb.py

@@ -17,7 +17,7 @@
 from util.click_util import CommaSeperatedInts, CommaSeperatedStrings, cli_help
 from util.log_handler import logger
 from util.netcdf_io import nc4_get_copy
-from util.utils import get_seed_from_member_number, process_member_num
+from util.utils import get_seed_from_member_id
 
 
 def create_perturb_files(in_path, in_files, out_path, copy_all_files=False):
@@ -64,10 +64,10 @@ def perturb_array(array, s, a):
     help=cli_help["files"],
 )
 @click.option(
-    "--member-num",
+    "--member-ids",
     type=CommaSeperatedInts(),
-    default="10",
-    help=cli_help["member_num"],
+    default="1,2,3,4,5,6,7,8,9,10",
+    help=cli_help["member_ids"],
 )
 @click.option(
     "--member-type",
@@ -95,33 +95,36 @@ def perturb(
     model_input_dir,
     perturbed_model_input_dir,
     files,
-    member_num,
+    member_ids,
     member_type,
     variable_names,
     perturb_amplitude,
     copy_all_files,
 ):  # pylint: disable=unused-argument, too-many-positional-arguments
 
-    processed_member_num = process_member_num(member_num)
-
-    for m_num, m_id in processed_member_num:
+    for m_id in member_ids:
 
         if member_type:
-            m_id = member_type + "_" + m_id
+            m_id_str = member_type + "_" + str(m_id)
+        else:
+            m_id_str = str(m_id)
+
         perturbed_model_input_dir_member_id = perturbed_model_input_dir.format(
-            member_id=m_id
+            member_id=m_id_str
         )
+
         data = create_perturb_files(
             model_input_dir,
             files,
             perturbed_model_input_dir_member_id,
             copy_all_files,
         )
+
         for d in data:
             for vn in variable_names:
                 d.variables[vn][:] = perturb_array(
                     d.variables[vn][:],
-                    get_seed_from_member_number(m_num),
+                    get_seed_from_member_id(m_id),
                     perturb_amplitude,
                 )
             d.close()

engine/run_ensemble.py

@@ -17,7 +17,7 @@
 
 from util.click_util import CommaSeperatedInts, CommaSeperatedStrings, cli_help
 from util.log_handler import logger
-from util.utils import get_seed_from_member_number, process_member_num
+from util.utils import get_seed_from_member_id
 
 
 def is_float(string):
@@ -114,7 +114,7 @@ def append_job(job, job_list, parallel):
         try:
             time.sleep(5)
             p.wait()
-            test_job_returncode(p)
+            check_job_returncode(p)
         finally:
             p.kill()
     else:
@@ -127,15 +127,15 @@ def finalize_jobs(job_list, dry, parallel):
         for job in job_list:
             job.communicate()
             try:
-                test_job_returncode(job)
+                check_job_returncode(job)
             except subprocess.CalledProcessError as e:
                 logger.error(e)
                 last_exception = e
         if last_exception:
             raise last_exception
 
 
-def test_job_returncode(job):
+def check_job_returncode(job):
     """Test job return code."""
     if job.returncode != 0:
         raise subprocess.CalledProcessError(returncode=job.returncode, cmd=job.args)
@@ -171,10 +171,10 @@ def test_job_returncode(job):
     help=cli_help["submit_command"],
 )
 @click.option(
-    "--member-num",
-    default="10",
+    "--member-ids",
+    default="1,2,3,4,5,6,7,8,9,10",
     type=CommaSeperatedInts(),
-    help=cli_help["member_num"],
+    help=cli_help["member_ids"],
 )
 @click.option(
     "--member-type",
@@ -215,7 +215,7 @@ def run_ensemble(
     experiment_name,
     perturbed_experiment_name,
     submit_command,
-    member_num,
+    member_ids,
     member_type,
     parallel,
     dry,
@@ -234,31 +234,33 @@ def run_ensemble(
         append_job(job, job_list, parallel)
 
     # run the ensemble
-    processed_member_num = process_member_num(member_num)
+    for m_id in member_ids:
 
-    for m_num, m_id in processed_member_num:
-
-        Path(perturbed_run_dir.format(member_id=m_id)).mkdir(
+        Path(perturbed_run_dir.format(member_id=str(m_id))).mkdir(
             exist_ok=True, parents=True
         )
-        os.chdir(perturbed_run_dir.format(member_id=m_id))
+        os.chdir(perturbed_run_dir.format(member_id=str(m_id)))
+
         if member_type:
-            m_id = member_type + "_" + m_id
+            m_id_str = member_type + "_" + str(m_id)
+        else:
+            m_id_str = str(m_id)
+
         runscript = f"{run_dir}/{run_script_name}"
 
-        perturbed_run_dir_path = perturbed_run_dir.format(member_id=m_id)
-        perturbed_run_script_path = perturbed_run_script_name.format(member_id=m_id)
+        perturbed_run_dir_path = perturbed_run_dir.format(member_id=m_id_str)
+        perturbed_run_script_path = perturbed_run_script_name.format(member_id=m_id_str)
         perturbed_runscript = f"{perturbed_run_dir_path}/{perturbed_run_script_path}"
 
         prepare_perturbed_run_script(
             runscript,
             perturbed_runscript,
             experiment_name,
-            perturbed_experiment_name.format(member_id=m_id),
+            perturbed_experiment_name.format(member_id=m_id_str),
             lhs,
             rhs_new,
             rhs_old,
-            get_seed_from_member_number(m_num),
+            get_seed_from_member_id(m_id),
         )
 
         if not dry: