diff --git a/codes/benchmark/bench_fcts.py b/codes/benchmark/bench_fcts.py index 14acd8a..265f7a3 100644 --- a/codes/benchmark/bench_fcts.py +++ b/codes/benchmark/bench_fcts.py @@ -1,3 +1,4 @@ +import os from contextlib import redirect_stdout from typing import Any @@ -378,7 +379,6 @@ def evaluate_iterative_predictions( # container for the piecewise predictions; seed t=0 with ground truth so errors # are computed only on actual predictions for t>=1 while keeping shape intact iterative_preds = np.zeros_like(targets) - iterative_preds[:, 0, :] = targets[:, 0, :] # number of chunks n_chunks = (n_timesteps + iter_interval - 1) // iter_interval @@ -433,14 +433,17 @@ def evaluate_iterative_predictions( ) # We predict steps 1..(chunk_len-1) relative to the provided init state (index 0). # Map these to global indices [start+1 .. end] inclusively. + if i == 0: + iterative_preds[:, start : end + 1, :] = preds_chunk[:, : model.n_timesteps, :].detach().cpu().numpy() iterative_preds[:, start + 1 : end + 1, :] = ( preds_chunk[:, 1 : model.n_timesteps, :].detach().cpu().numpy() ) iterative_preds_log = model.denormalize(iterative_preds, leave_log=True) + full_preds_log = model.denormalize(full_preds, leave_log=True) targets_log = model.denormalize(targets, leave_log=True) iterative_preds = model.denormalize(iterative_preds) - full_preds = model.denormalize(full_preds.detach().cpu().numpy()) + full_preds_real = model.denormalize(full_preds.detach().cpu().numpy()) targets = model.denormalize(targets) # compute error metrics @@ -466,7 +469,7 @@ def evaluate_iterative_predictions( surr_name, conf, iterative_preds, - full_preds, + full_preds_real, targets, timesteps, iter_interval=iter_interval, @@ -1583,7 +1586,7 @@ def compare_UQ(all_metrics: dict, config: dict) -> None: ensemble_errors, ensemble_std, config, - flag_fractions=(0.01, 0.05, 0.10, 0.20, 0.30, 0.40, 0.50), + flag_fractions=(0, 0.025, 0.05, 0.10, 0.20, 0.30, 0.40, 0.50), save=True, show_title=True, ) diff --git a/codes/benchmark/bench_plots.py b/codes/benchmark/bench_plots.py index fb5980c..f2804c4 100644 --- a/codes/benchmark/bench_plots.py +++ b/codes/benchmark/bench_plots.py @@ -1661,6 +1661,7 @@ def plot_errors_over_time( elif mode == "iterative": # Single backslash inside raw string to render the LaTeX Delta properly plt.ylabel(r"Log-MAE ($\Delta dex$)") + plt.ylim(bottom=0, top=min(np.max(list(mean_errors.values())) * 1.1, 5)) fname = "iterative_delta_dex_time.png" title = "Comparison of Δdex Errors Over Time for Iterative Predictions" # Add subtle dashed vertical lines at every n-th timestep if provided and valid @@ -2035,8 +2036,6 @@ def inference_time_bar_plot( # Calculate the upper y-limit to provide space for text max_bar = max(means[i] + stds[i] for i in range(len(means))) # min_bar = min(means[i] - stds[i] for i in range(len(means))) - # Temp! - # ax.set_ylim(min_bar * 0.3, max_bar * 2) # Set limits with some padding ax.set_ylim(0, max_bar * 1.2) # Set limits with some padding # Add inference time as text to the bars using the format_time function @@ -2053,8 +2052,6 @@ def inference_time_bar_plot( ax.set_xlabel("Surrogate Model") ax.set_ylabel("Mean Inference Time per Run") - # Temp! - # ax.set_yscale("log") if show_title: ax.set_title("Surrogate Mean Inference Time Comparison") @@ -2507,11 +2504,16 @@ def plot_catastrophic_detection_curves( xs, ys = [], [] for f in flag_fractions: - unc_thr = np.percentile(u, 100.0 * (1.0 - float(f))) - flagged = u >= unc_thr - recall = (flagged & is_cat).sum() / n_cat if n_cat > 0 else 0.0 - xs.append(100.0 * flagged.mean()) - ys.append(100.0 * recall) + if f <= 0.0: + xs.append(0.0) + ys.append(0.0) + recall = 0.0 + else: + unc_thr = np.percentile(u, 100.0 * (1.0 - float(f))) + flagged = u >= unc_thr + recall = (flagged & is_cat).sum() / n_cat if n_cat > 0 else 0.0 + xs.append(100.0 * flagged.mean()) + ys.append(100.0 * recall) ax.plot( xs, @@ -2988,8 +2990,6 @@ def rel_errors_and_uq( ax1.set_xlabel("Time") ax1.set_xlim(timesteps[0], timesteps[-1]) - # Temp! - # ax1.set_ylim(3e-4, 1) ax1.set_ylabel("Relative Error") ax1.set_yscale("log") ax1.set_title("Comparison of Relative Errors Over Time") @@ -3020,8 +3020,6 @@ def rel_errors_and_uq( ax2.set_xlabel("Time") ax2.set_xlim(timesteps[0], timesteps[-1]) - # Temp! - # ax2.set_ylim(0, 0.04) ax2.set_ylabel("Uncertainty/Absolute Error") if show_title: ax2.set_title("Comparison of Predictive Uncertainty Over Time") diff --git a/codes/surrogates/AbstractSurrogate/abstract_surrogate.py b/codes/surrogates/AbstractSurrogate/abstract_surrogate.py index 8c62d7b..4b2a735 100644 --- a/codes/surrogates/AbstractSurrogate/abstract_surrogate.py +++ b/codes/surrogates/AbstractSurrogate/abstract_surrogate.py @@ -449,8 +449,10 @@ def denormalize( Returns: Tensor | np.ndarray: The denormalized data. """ + data_type = None if self.normalisation is not None: if not leave_norm: + data_type = data.dtype if self.normalisation["mode"] == "disabled": ... elif self.normalisation["mode"] == "minmax": @@ -475,6 +477,13 @@ def denormalize( if self.normalisation["log10_transform"] and not leave_log: data = 10**data + # Conserve dtype + if data_type is not None: + if isinstance(data, Tensor): + return data.to(dtype=data_type) + if isinstance(data, np.ndarray): + return data.astype(data_type) + return data def denormalize_old(self, data: Tensor) -> Tensor: diff --git a/codes/tune/evaluate_tuning.py b/codes/tune/evaluate_tuning.py index a421835..c1ae48d 100644 --- a/codes/tune/evaluate_tuning.py +++ b/codes/tune/evaluate_tuning.py @@ -11,6 +11,8 @@ import optuna import psycopg2 import torch +from matplotlib.lines import Line2D +from matplotlib.ticker import FormatStrFormatter from optuna.trial import TrialState from psycopg2 import sql @@ -22,9 +24,35 @@ from codes.tune import load_yaml_config from codes.utils import nice_print - -def pareto_front(points: np.ndarray) -> np.ndarray: - # lower-is-better for both objectives +os.environ.setdefault("PGCONNECT_TIMEOUT", "3") + +TRIAL_COLOR = "#c7c7c7" +PARETO_COLOR = "#1f77b4" +BEST_F1_COLOR = "#d62728" +CHOSEN_COLOR = "#2ca02c" +EDGE_COLOR = "#000000" +PARETO_X_RANGE = (0, 6) +SURROGATE_ORDER = [ + "multionet", + "fullyconnected", + "latentpoly", + "latentneuralode", +] + + +def format_surrogate_title(name: str) -> str: + mapping = { + "multionet": "MultiONet", + "fullyconnected": "FullyConnected", + "latentneuralode": "LatentNeuralODE", + "latentpoly": "LatentPoly", + } + return mapping.get(name.lower(), name) + + +def pareto_front_mask(points: np.ndarray) -> np.ndarray: + if points.size == 0: + return np.array([], dtype=bool) is_efficient = np.ones(points.shape[0], dtype=bool) for i, p in enumerate(points): if not is_efficient[i]: @@ -34,7 +62,12 @@ def pareto_front(points: np.ndarray) -> np.ndarray: dominated = better & (np.arange(points.shape[0]) != i) if np.any(dominated): is_efficient[i] = False - return points[is_efficient] + return is_efficient + + +def pareto_front(points: np.ndarray) -> np.ndarray: + mask = pareto_front_mask(points) + return points[mask] def hypervolume_2d(pareto_points: np.ndarray, reference: np.ndarray) -> float: @@ -99,13 +132,413 @@ def compute_hypervolume_over_time( return hypervolumes, reference +def _pareto_legend_handles() -> list[Line2D]: + return [ + Line2D( + [0], + [0], + marker="o", + color="none", + markerfacecolor=TRIAL_COLOR, + markeredgecolor=TRIAL_COLOR, + label="Trial Outcome", + markersize=6, + ), + Line2D( + [0], + [0], + marker="o", + color="none", + markerfacecolor=PARETO_COLOR, + markeredgecolor=EDGE_COLOR, + label="Pareto front", + markersize=6, + ), + Line2D( + [0], + [0], + marker="v", + color="none", + markerfacecolor=BEST_F1_COLOR, + markeredgecolor=EDGE_COLOR, + label="Lowest Δdex99", + markersize=6, + ), + Line2D( + [0], + [0], + marker="^", + color="none", + markerfacecolor=CHOSEN_COLOR, + markeredgecolor=EDGE_COLOR, + label="Chosen Trial", + markersize=6, + ), + ] + + +def compute_pareto_plot_data( + study: optuna.Study, + suffix: str, + ignore_last_n: int = 10, + chosen_trial_number: int | None = None, +): + if len(study.directions) != 2: + return None + completed = [t for t in study.trials if t.state == TrialState.COMPLETE] + if not completed: + print(f"Skipping Pareto plot for {suffix}: no completed trials.") + return None + completed.sort(key=lambda t: t.datetime_complete or t.datetime_start) + cutoff = max(0, len(completed) - max(0, ignore_last_n)) + eligible_trials = completed[:cutoff] + filtered = [] + for t in eligible_trials: + if t.values is None: + continue + vals = np.array(t.values, dtype=float) + if vals.shape[0] != 2 or np.any(~np.isfinite(vals)): + continue + filtered.append((vals, t.number)) + if not filtered: + print( + f"Skipping Pareto plot for {suffix}: no eligible trials after ignoring last {ignore_last_n}." + ) + return None + points = np.vstack([vals for vals, _ in filtered]) + trial_numbers = np.array([num for _, num in filtered]) + mask = pareto_front_mask(points) + best_idx = int(np.argmin(points[:, 0])) + best_point = points[best_idx] + chosen_point = None + if chosen_trial_number is not None: + matches = np.where(trial_numbers == chosen_trial_number)[0] + if matches.size: + chosen_point = points[matches[0]] + else: + print( + f"Chosen trial {chosen_trial_number} not in eligible set for {suffix}; " + "skipping highlight." + ) + + y = points[:, 1] + y_min = 0 + y_max = best_point[1] * 1.3 + + return { + "suffix": suffix, + "points": points, + "mask": mask, + "best_point": best_point, + "chosen_point": chosen_point, + "xlim": PARETO_X_RANGE, + "ylim": (y_min, y_max), + "ignore_last_n": ignore_last_n, + } + + +def summarize_pareto_tradeoff(data: dict): + suffix = data["suffix"] + best_point = data["best_point"] + chosen_point = data["chosen_point"] + pts = data["points"][data["mask"]] + x_span = float(pts[:, 0].max() - pts[:, 0].min()) + y_span = float(pts[:, 1].max() - pts[:, 1].min()) + min_x = float(pts[:, 0].min()) + min_y = float(pts[:, 1].min()) + x_rel = x_span / max(abs(min_x), 1e-12) + y_rel = y_span / max(abs(min_y), 1e-12) + print( + f"{suffix}: Pareto span Δdex={x_span:.4f} " + f"(x{ x_rel:.2f} relative), inference time span={y_span:.4f}s " + f"(x{ y_rel:.2f} relative)" + ) + best_err, best_time = float(best_point[0]), float(best_point[1]) + print( + f"{suffix}: lowest-error trial Δdex={best_err:.4f}, inference time={best_time:.4f}s" + ) + if chosen_point is None: + print( + f"{suffix}: chosen trial not available after filtering; skipping summary." + ) + return + chosen_err, chosen_time = float(chosen_point[0]), float(chosen_point[1]) + print( + f"{suffix}: chosen trial Δdex={chosen_err:.4f}, inference time={chosen_time:.4f}s" + ) + err_reduction = 1.0 - best_err / chosen_err if chosen_err > 0 else float("nan") + time_ratio = best_time / chosen_time if chosen_time > 0 else float("nan") + if math.isnan(err_reduction) or math.isnan(time_ratio): + print(f"{suffix}: insufficient data to compute tradeoff summary.") + return + err_phrase = ( + f"{abs(err_reduction) * 100:.1f}% lower error" + if err_reduction >= 0 + else f"{abs(err_reduction) * 100:.1f}% higher error" + ) + ratio_phrase = ( + f"but {time_ratio:.2f}x higher inference time" + if time_ratio >= 1 + else f"and {time_ratio:.2f}x lower inference time" + ) + print(f"{suffix}: {err_phrase}, {ratio_phrase}.") + + +def _render_pareto_scatter( + ax, + data, + show_xlabel: bool, + show_ylabel: bool, + title: str | None, + hide_xticklabels: bool = False, +): + pts = data["points"] + mask = data["mask"] + best_point = data["best_point"] + chosen_point = data["chosen_point"] + + # Remove best point and chosen point from pareto points to avoid double-plotting + if mask.sum() > 0: + pareto_pts = pts[mask] + pareto_pts = pareto_pts[ + ~np.all(pareto_pts == best_point, axis=1) + ] # remove best point + if chosen_point is not None: + pareto_pts = pareto_pts[ + ~np.all(pareto_pts == chosen_point, axis=1) + ] # remove chosen point + # Recompute mask + new_mask = np.array( + [any(np.all(p == pp) for pp in pareto_pts) for p in pts], dtype=bool + ) + mask = new_mask + else: + mask = np.array([False] * pts.shape[0], dtype=bool) + + ax.scatter(pts[:, 0], pts[:, 1], color=TRIAL_COLOR, alpha=0.7, label=None) + ax.scatter( + pts[mask, 0], + pts[mask, 1], + marker="o", + color=PARETO_COLOR, + edgecolor=EDGE_COLOR, + linewidth=0.5, + label=None, + ) + ax.scatter( + best_point[0], + best_point[1], + marker="v", + color=BEST_F1_COLOR, + edgecolor=EDGE_COLOR, + linewidth=0.6, + label=None, + zorder=3, + ) + if chosen_point is not None: + ax.scatter( + chosen_point[0], + chosen_point[1], + marker="^", + color=CHOSEN_COLOR, + edgecolor=EDGE_COLOR, + linewidth=0.6, + label=None, + zorder=4, + ) + + ax.set_xlim(*data["xlim"]) + ax.set_ylim(*data["ylim"]) + if show_xlabel: + ax.set_xlabel(r"LAE$_{99}$ [dex]") + if show_ylabel: + ax.set_ylabel("Inference Time [s]") + if title: + ax.set_title(title) + ax.tick_params(labelbottom=not hide_xticklabels) + ax.grid(True, linestyle="--", alpha=0.3) + + ax.yaxis.set_major_formatter(FormatStrFormatter("%.3f")) + + +def save_individual_pareto_plot(data, out_dir: str): + fig, ax = plt.subplots(figsize=(5, 3)) + _render_pareto_scatter( + ax, data, show_xlabel=True, show_ylabel=True, title=data["suffix"] + ) + ax.legend(handles=_pareto_legend_handles(), loc="best", ncol=2) + fig.tight_layout() + os.makedirs(out_dir, exist_ok=True) + fig.savefig(os.path.join(out_dir, f"pareto_front_{data['suffix']}.png"), dpi=300) + plt.close(fig) + print( + f"Saved Pareto front plot for {data['suffix']}; " + f"ignored last {data['ignore_last_n']} completed trials." + ) + + +def save_pareto_front_grid(datasets: list[dict], out_dir: str): + if not datasets: + return + n_cols = 2 + n_rows = max(1, math.ceil(len(datasets) / n_cols)) + fig, axes = plt.subplots( + n_rows, n_cols, figsize=(5 * n_cols, 3 * n_rows), squeeze=False + ) + + for idx, data in enumerate(datasets): + row, col = divmod(idx, n_cols) + ax = axes[row][col] + show_xlabel = row == n_rows - 1 + show_ylabel = col == 0 + proper_title = format_surrogate_title(data["suffix"]) + _render_pareto_scatter( + ax, + data, + show_xlabel=show_xlabel, + show_ylabel=show_ylabel, + title=proper_title, + hide_xticklabels=not show_xlabel, + ) + for idx in range(len(datasets), n_rows * n_cols): + row, col = divmod(idx, n_cols) + axes[row][col].axis("off") + + fig.tight_layout(rect=(0.01, 0.035, 0.99, 0.99)) + fig.legend( + handles=_pareto_legend_handles(), + loc="lower center", + ncol=2, + bbox_to_anchor=(0.53, 0.0), + fontsize=10, + ) + os.makedirs(out_dir, exist_ok=True) + fig.savefig(os.path.join(out_dir, "pareto_front_grid.png"), dpi=300) + plt.close(fig) + print( + f"Saved combined Pareto front grid with {len(datasets)} subplot(s) " + f"to {os.path.join(out_dir, 'pareto_front_grid.png')}." + ) + + +def save_relative_hv_grid(datasets: list[dict], out_dir: str): + if not datasets: + return + n_cols = 2 + n_rows = max(1, math.ceil(len(datasets) / n_cols)) + fig, axes = plt.subplots( + n_rows, n_cols, figsize=(5 * n_cols, 3 * n_rows), squeeze=False + ) + + for idx, data in enumerate(datasets): + row, col = divmod(idx, n_cols) + ax = axes[row][col] + show_xlabel = row == n_rows - 1 + show_ylabel = col == 0 + ax.plot(data["x"], data["y"], color="#1f77b4") + ax.set_xlim(0, max(data["x"]) if len(data["x"]) else 1) + ax.set_ylim(0.6, 1.03) + ax.set_title(format_surrogate_title(data["suffix"])) + if show_xlabel: + ax.set_xlabel("Completed Trials") + if show_ylabel: + ax.set_ylabel("Normalized Hypervolume") + else: + ax.tick_params(labelleft=False) + ax.grid(True, linestyle="--", alpha=0.3) + + for idx in range(len(datasets), n_rows * n_cols): + row, col = divmod(idx, n_cols) + axes[row][col].axis("off") + + fig.tight_layout(rect=(0.01, 0.01, 0.99, 0.99)) + os.makedirs(out_dir, exist_ok=True) + fig.savefig(os.path.join(out_dir, "hypervolume_normalized_grid.png"), dpi=300) + plt.close(fig) + print( + f"Saved combined normalized-hypervolume grid with {len(datasets)} subplot(s) " + f"to {os.path.join(out_dir, 'hypervolume_normalized_grid.png')}." + ) + + +def save_relative_hv_combined(datasets: list[dict], out_dir: str): + if len(datasets) < 2: + return + fig, ax = plt.subplots(figsize=(6, 4)) + color_palette = plt.cm.viridis(np.linspace(0, 0.95, len(SURROGATE_ORDER))) + order_index = {name: idx for idx, name in enumerate(SURROGATE_ORDER)} + + def color_for(name: str): + key = name.lower() + idx = order_index.get(key) + if idx is not None: + return color_palette[idx] + return plt.cm.tab10(0) + + ordered = sorted( + datasets, + key=lambda d: order_index.get(d["suffix"].lower(), len(SURROGATE_ORDER)), + ) + + for data in ordered: + color = color_for(data["suffix"]) + ax.plot( + data["x"], + data["y"], + label=format_surrogate_title(data["suffix"]), + color=color, + ) + if len(data["x"]): + last_x = data["x"][-1] + last_y = data["y"][-1] + ax.vlines( + last_x, + max(0, last_y - 0.01), + min(1.05, last_y + 0.01), + colors=color, + linewidth=2, + ) + ax.set_xlim(0, max(max(d["x"]) for d in datasets if len(d["x"])) if datasets else 1) + ax.set_ylim(0.6, 1.02) + ax.set_xlabel("Completed Trials") + ax.set_ylabel("Normalized Hypervolume") + # ax.set_title("Normalized Hypervolume") + ax.grid(True, linestyle="--", alpha=0.3) + ax.legend(ncol=2, loc="lower right") + fig.tight_layout() + os.makedirs(out_dir, exist_ok=True) + fig.savefig(os.path.join(out_dir, "hypervolume_normalized_combined.png"), dpi=300) + plt.close(fig) + print( + "Saved combined normalized-hypervolume line plot to " + f"{os.path.join(out_dir, 'hypervolume_normalized_combined.png')}." + ) + + def load_loss_history(model_path: str) -> tuple[np.ndarray, np.ndarray, int]: """ - Load loss histories from a saved model file (.pth). - Returns (train_loss, test_loss, n_epochs). + Load loss histories from a saved .pth. + If the checkpoint contains pickled Optuna storage (which may try to connect + to a remote DB during unpickle), bound the timeout and fall back to skipping. """ - model_dict = torch.load(model_path, map_location="cpu", weights_only=False) - attributes = model_dict.get("attributes", {}) + # Prefer the Tensor-only safe path if your PyTorch supports it (PyTorch >= 2.0) + try: + obj = torch.load(model_path, map_location="cpu", weights_only=True) + # weights_only returns just state_dict; no custom 'attributes' available + # We can’t recover losses from state_dict → skip gracefully. + return None, None, 0 + except TypeError: + # weights_only not supported → careful unpickle with timeout + catch + pass + + try: + obj = torch.load(model_path, map_location="cpu") + except Exception as e: + # Any error (including timeout from unpickling Optuna storage) → skip losses + print(f"[warn] Could not safely load {model_path}: {e}. Skipping loss curves.") + return None, None, 0 + + attributes = obj.get("attributes", {}) if isinstance(obj, dict) else {} train_loss = ( np.array(attributes.get("train_loss")) if attributes.get("train_loss") is not None @@ -219,6 +652,7 @@ def evaluate_tuning( top_n: int = 10, storage_name: str = "optuna_db", ignore_last_n: int = 10, + chosen_indices: list[int] | None = None, ) -> None: """ For all surrogate studies named '_', @@ -243,6 +677,8 @@ def evaluate_tuning( connect_timeout=5, ).close() except Exception: + if str(pg.get("host", "localhost")) not in ("localhost", "127.0.0.1", "::1"): + raise # don’t try to start a local server if host isn't local pg_data = pg.get("data_dir", os.path.expanduser("~/postgres/data")) pg_ctl = os.path.join( pg.get("database_folder", os.path.expanduser("~/postgres")), "bin", "pg_ctl" @@ -276,10 +712,10 @@ def evaluate_tuning( from optuna.study import get_all_study_summaries summaries = get_all_study_summaries(storage=storage_url) - study_names = [ + study_names_all = [ s.study_name for s in summaries if s.study_name.startswith(f"{study_prefix}_") ] - if not study_names: + if not study_names_all: print(f"No studies found with prefix '{study_prefix}_' in {storage_url}") return @@ -287,8 +723,41 @@ def evaluate_tuning( save_dir = os.path.join("tuned", study_prefix) os.makedirs(save_dir, exist_ok=True) + # Derive deterministic ordering from config surrogate list if available + surrogate_entries = config.get("surrogates", []) + surrogate_names = [ + str(entry.get("name")) + for entry in surrogate_entries + if isinstance(entry, dict) and entry.get("name") + ] + study_names_ordered = [] + used: set[str] = set() + study_lookup = {name.lower(): name for name in study_names_all} + if surrogate_names: + for surrogate in surrogate_names: + candidate = f"{study_prefix}_{surrogate}".lower() + match = study_lookup.get(candidate) + if match and match not in used: + study_names_ordered.append(match) + used.add(match) + for name in sorted(study_names_all): + if name not in used: + study_names_ordered.append(name) + used.add(name) + study_names = study_names_ordered + + if chosen_indices is not None and len(chosen_indices) != len(study_names): + print( + "Provided chosen_indices length does not match number of studies; " + "skipping chosen-trial highlighting." + ) + chosen_indices = None + + pareto_datasets: list[dict] = [] + hypervolume_datasets: list[dict] = [] + # Loop over each surrogate study - for full_name in study_names: + for idx, full_name in enumerate(study_names): suffix = full_name[len(study_prefix) + 1 :] print(f"--- Evaluating study {full_name} -> surrogate '{suffix}' ---") try: @@ -324,11 +793,11 @@ def evaluate_tuning( plt.plot( np.arange(1, len(rel_hvs) + 1), rel_hvs, - label="Relative Hypervolume", + label="Normalized Hypervolume", ) plt.xlabel("Completed Trials") - plt.ylabel("Fraction of Final HV") - plt.title(f"{suffix} Relative Hypervolume") + plt.ylabel("Normalized Hypervolume") + plt.title(f"{suffix} Normalized Hypervolume") plt.grid(True) plt.tight_layout() plt.savefig( @@ -336,12 +805,33 @@ def evaluate_tuning( dpi=300, ) plt.close() + max_x = max(1, len(rel_hvs) - max(0, ignore_last_n)) + hypervolume_datasets.append( + { + "suffix": suffix, + "x": np.arange(1, len(rel_hvs) + 1)[:max_x], + "y": rel_hvs[:max_x], + } + ) print( f"Saved hypervolume plots for {suffix} (final HV={final_hv:.3e}); " f"ignored last {ignore_last_n} trials in Pareto/HV." ) else: print("No hypervolume computed (no complete trials).") + chosen_trial = None + if chosen_indices is not None: + chosen_trial = chosen_indices[idx] + plot_data = compute_pareto_plot_data( + study, + suffix, + ignore_last_n=ignore_last_n, + chosen_trial_number=chosen_trial, + ) + if plot_data: + save_individual_pareto_plot(plot_data, save_dir) + pareto_datasets.append(plot_data) + summarize_pareto_tradeoff(plot_data) else: print("Skipping hypervolume: study is not two-objective.") @@ -369,40 +859,44 @@ def evaluate_tuning( continue models_dir = os.path.join(models_root, match) - # Collect top trials' loss histories - records = [] - epochs = None - for fname in os.listdir(models_dir): - if not fname.endswith(".pth"): - continue - m = re.search(r"_(\d+)\.pth$", fname) - if not m: - continue - tnum = int(m.group(1)) - if tnum not in best: - continue - path = os.path.join(models_dir, fname) - train, test, ep = load_loss_history(path) - if test is None: - continue - records.append((tnum, test)) - epochs = ep - if not records: - print(f"No top-{top_n} trials for surrogate '{suffix}'") - continue - records.sort(key=lambda x: best.index(x[0])) - tnums, losses = zip(*records) - labels = [f"Trial {n}" for n in tnums] - out_dir = save_dir - plot_losses( - list(losses), - epochs, - labels, - title=f"{suffix} Top-{top_n}", - save=True, - out_dir=out_dir, - mode=suffix, - ) + # # Collect top trials' loss histories + # records = [] + # epochs = None + # for fname in os.listdir(models_dir): + # if not fname.endswith(".pth"): + # continue + # m = re.search(r"_(\d+)\.pth$", fname) + # if not m: + # continue + # tnum = int(m.group(1)) + # if tnum not in best: + # continue + # path = os.path.join(models_dir, fname) + # train, test, ep = load_loss_history(path) + # if test is None: + # continue + # records.append((tnum, test)) + # epochs = ep + # if not records: + # print(f"No top-{top_n} trials for surrogate '{suffix}'") + # continue + # records.sort(key=lambda x: best.index(x[0])) + # tnums, losses = zip(*records) + # labels = [f"Trial {n}" for n in tnums] + # out_dir = save_dir + # plot_losses( + # list(losses), + # epochs, + # labels, + # title=f"{suffix} Top-{top_n}", + # save=True, + # out_dir=out_dir, + # mode=suffix, + # ) + + save_pareto_front_grid(pareto_datasets, save_dir) + save_relative_hv_grid(hypervolume_datasets, save_dir) + save_relative_hv_combined(hypervolume_datasets, save_dir) def parse_args(): @@ -412,13 +906,13 @@ def parse_args(): p.add_argument( "--study_name", type=str, - default="primordial_parametric_final", + default="primordial_final", help="Main study prefix (e.g. lvparams5)", ) p.add_argument( "--storage_name", type=str, - default="primordial_parametric_final", + default="primordial_final", help="Main study prefix (e.g. lvparams5)", ) p.add_argument( @@ -435,16 +929,37 @@ def parse_args(): "Number of most-recent completed trials to exclude from Pareto/hypervolume" ), ) + p.add_argument( + "--chosen_indices", + type=str, + # default="171,114,135,237", # cloud_final + # default="27,61,13,299", # cloud_parametric_final + # default="18,1,16,234", # primordial_parametric_final + default="196,107,31,243", # primordial_final + help=( + "Comma-separated list of Optuna trial numbers chosen per study (order " + "matches sorted study names)." + ), + ) return p.parse_args() def main(): args = parse_args() + chosen_list = None + if args.chosen_indices: + chosen_list = [] + for item in args.chosen_indices.split(","): + item = item.strip() + if not item: + continue + chosen_list.append(int(item)) evaluate_tuning( args.study_name, args.top_n, args.storage_name, ignore_last_n=args.ignore_last_n, + chosen_indices=chosen_list, ) diff --git a/config.yaml b/config.yaml index 1642023..4d10e24 100644 --- a/config.yaml +++ b/config.yaml @@ -11,32 +11,32 @@ dataset: tolerance: 1e-25 normalise_per_species: True log_timesteps: True -devices: ["cuda:3", "cuda:2", "cuda:3", "cuda:5", "cuda:6", "cuda:7", "cuda:8", "cuda:9"] +devices: ["cuda:4", "cuda:2", "cuda:3", "cuda:5", "cuda:6", "cuda:7", "cuda:8", "cuda:9"] seed: 42 verbose: False checkpoint: True # Models to train interpolation: - enabled: False + enabled: True intervals: [2, 3, 4, 5, 6, 7, 8, 9, 10] extrapolation: - enabled: False + enabled: True cutoffs: [50, 60, 70, 80, 90] sparse: - enabled: False + enabled: True factors: [2, 4, 8, 16, 32] batch_scaling: - enabled: False + enabled: True sizes: [1/16, 1/8, 1/4, 1/2] uncertainty: - enabled: False + enabled: True ensemble_size: 5 # Number of models for deep ensemble # Evaluations during benchmark -iterative: False -losses: False -gradients: False +iterative: True +losses: True +gradients: True timing: True -compute: False +compute: True compare: True # Whether to compare the surrogates diff --git a/config_full.yaml b/config_full.yaml index e20ee71..b578f98 100644 --- a/config_full.yaml +++ b/config_full.yaml @@ -13,7 +13,7 @@ dataset: tolerance: 1e-25 subset_factor: 1 log_timesteps: True -devices: ["cuda:0", "cuda:1", "cuda:2", "cuda:3"] # ["cuda:0", "cuda:1", "cuda:2", "cuda:3", "cuda:4", "cuda:5", "cuda:6", "cuda:7", "cuda:8"] +devices: ["cuda:1", "cuda:1", "cuda:2", "cuda:3"] # ["cuda:0", "cuda:1", "cuda:2", "cuda:3", "cuda:4", "cuda:5", "cuda:6", "cuda:7", "cuda:8"] seed: 42 verbose: False relative_error_threshold: 1e-10