rodts28/Inverse_3DCP_Evolutionary_Algorithm
0
1# -*- coding: utf-8 -*-
2"""
3Created on Tue Dec 9 18:45:06 2025
4
5@author: rod-t
6"""
7
8# analysis_quality.py
9"""
10Quality & sensitivity analysis for inverse-designed mixes.
11
12Works with any inverse-design engine (evolutionary search, RL, or
13nearest-neighbor) as long as a *_REPORT.csv and *_ELITES.csv exist.
14
15Includes:
16- analyze_inverse_and_sensitivity: per-output error + kNN manifold distance +
17 central-difference local sensitivity (around the inverse mix).
18- local_sensitivity_with_constraints: constraint-aware perturbation of mix
19 variables and loss evaluation (range-normalized RMSE + engineering + manifold).
20"""
21
22import os
23import math
24import numpy as np
25import pandas as pd
26
27from core_inverse import CFG, load_and_fit, engineering_constraints_penalty
28
29
30# ---------------------------------------------------------
31# Helpers
32# ---------------------------------------------------------
33
34def _extract_targets_preds_mix_from_report(report_df: pd.DataFrame):
35 """
36 Return dicts: targets, preds, mix (canonical names), and ordered lists.
37 Expects columns like Target::Pump-Speed, Pred::Pump-Speed, Mix::Cement, ...
38 """
39 targ_cols = [c for c in report_df.columns if c.startswith("Target::")]
40 pred_cols = [c for c in report_df.columns if c.startswith("Pred::")]
41 mix_cols = [c for c in report_df.columns if c.startswith("Mix::")]
42
43 # Strip prefixes for canonical keys
44 targets = {
45 c.replace("Target::", ""): float(report_df[c].iloc[0])
46 for c in targ_cols
47 }
48 preds = {
49 c.replace("Pred::", ""): float(report_df[c].iloc[0])
50 for c in pred_cols
51 }
52 mix = {
53 c.replace("Mix::", ""): float(report_df[c].iloc[0])
54 for c in mix_cols
55 }
56
57 # Orders for nice printing
58 out_order = [c.replace("Target::", "") for c in targ_cols] # keep report order
59 mix_order = [c.replace("Mix::", "") for c in mix_cols]
60
61 return targets, preds, mix, out_order, mix_order
62
63
64def _per_output_metrics(targets, preds, y_rng_map):
65 """
66 Build per-output error table:
67 - AbsErr
68 - PctErr(%)
69 - RangeNormErr (abs error normalized by dataset output range)
70 """
71 rows = []
72 for o in targets.keys():
73 t = targets[o]
74 p = preds.get(o, np.nan)
75
76 abs_err = abs(p - t)
77 pct_err = abs_err / (abs(t) + 1e-9) * 100.0
78
79 rng = y_rng_map.get(o, np.nan)
80 if rng is None or (isinstance(rng, float) and rng <= 0):
81 rn_err = abs_err
82 else:
83 rn_err = abs_err / rng
84
85 rows.append(
86 {
87 "Output": o,
88 "Target": t,
89 "Pred": p,
90 "AbsErr": abs_err,
91 "PctErr(%)": pct_err,
92 "RangeNormErr": rn_err,
93 }
94 )
95
96 return pd.DataFrame(rows)
97
98
99def _elite_variability_tables(elites_df: pd.DataFrame):
100 """
101 Return two Series: variability of mix cols and pred cols among top elites.
102 - var_mix: variance by mix variable (top 10 elites)
103 - var_pred: variance by predicted outputs (top 10 elites)
104 """
105 if elites_df is None or elites_df.empty:
106 return pd.Series(dtype=float), pd.Series(dtype=float)
107
108 mix_cols = [
109 c
110 for c in elites_df.columns
111 if not c.startswith("Pred::") and not c.startswith("Score")
112 ]
113 pred_cols = [c for c in elites_df.columns if c.startswith("Pred::")]
114
115 top = elites_df.head(min(10, len(elites_df))) # top 10
116
117 var_mix = (
118 top[mix_cols].var(numeric_only=True).sort_values(ascending=False)
119 if mix_cols
120 else pd.Series(dtype=float)
121 )
122 var_pred = (
123 top[pred_cols].var(numeric_only=True).sort_values(ascending=False)
124 if pred_cols
125 else pd.Series(dtype=float)
126 )
127
128 # Strip prefixes for readability
129 if not var_pred.empty:
130 var_pred.index = [c.replace("Pred::", "") for c in var_pred.index]
131
132 return var_mix, var_pred
133
134
135def _engineering_checks(mix: dict):
136 """
137 Quick engineering plausibility checks on the inverse mix.
138 Uses global CFG for:
139 - enforce_sum_100, sum_target, sum_tolerance
140 - wc_min, wc_max
141 - cap_sp_percent, sp_max
142 """
143 notes = []
144
145 # Sum ~ 100 (if enabled in CFG)
146 if CFG.get("enforce_sum_100", False):
147 s = sum(mix.values())
148 tol = CFG.get("sum_tolerance", 1.0)
149 tgt = CFG.get("sum_target", 100.0)
150 if abs(s - tgt) > tol:
151 notes.append(f"Sum check: {s:.2f} (target {tgt}±{tol})")
152
153 # w/c ratio band
154 if "Water" in mix and "Cement" in mix and mix["Cement"] > 0:
155 wc = mix["Water"] / max(mix["Cement"], 1e-9)
156 if not (CFG.get("wc_min", 0.0) <= wc <= CFG.get("wc_max", 999.0)):
157 notes.append(
158 f"w/c={wc:.3f} outside "
159 f"[{CFG.get('wc_min', 0.0)},{CFG.get('wc_max', 999.0)}]"
160 )
161
162 # SP cap (if enabled)
163 if CFG.get("cap_sp_percent", False) and ("SP" in mix):
164 sp_max = CFG.get("sp_max", np.inf)
165 if mix["SP"] > sp_max:
166 notes.append(f"SP={mix['SP']:.2f} exceeds cap {sp_max}")
167
168 return notes
169
170
171# ---------------------------------------------------------
172# Global quality + local sensitivity (central difference)
173# ---------------------------------------------------------
174
175def analyze_inverse_and_sensitivity(
176 report_csv_path: str,
177 elites_csv_path: str,
178 perturb_frac: float = 0.05,
179):
180 """
181 Full analysis around an inverse-designed solution (EA, RL, etc.):
182 - per-output errors (abs, %, range-normalized)
183 - elites variability (mix & preds)
184 - manifold proximity score (scaled kNN distance)
185 - local sensitivities (central difference, ±perturb_frac in each mix var)
186 """
187
188 # Load artifacts from inverse run
189 report_df = pd.read_csv(report_csv_path)
190 elites_df = pd.read_csv(elites_csv_path)
191
192 # Unified loader (matches core_inverse.load_and_fit)
193 (
194 df,
195 model,
196 mix_map,
197 out_map,
198 X,
199 Y,
200 lb,
201 ub,
202 y_min,
203 y_max,
204 y_rng,
205 nn,
206 nn_scale,
207 ) = load_and_fit()
208
209 # Maps for ranges by canonical output name
210 out_actual_to_canon = {v: k for k, v in out_map.items()}
211 y_rng_map = {
212 out_actual_to_canon[col]: float(y_rng[col]) for col in y_rng.index
213 }
214
215 # Parse report row
216 targets, preds, mix, out_order, mix_order = _extract_targets_preds_mix_from_report(
217 report_df
218 )
219
220 # ---------- Per-output error table ----------
221 per_out = _per_output_metrics(targets, preds, y_rng_map)
222 print("\n=== PER-OUTPUT ERRORS ===")
223 print(per_out.to_string(index=False))
224
225 # ---------- Aggregate diagnostics ----------
226 mean_pct = per_out["PctErr(%)"].mean()
227 mean_rng = per_out["RangeNormErr"].mean()
228 print(f"\nMean % error: {mean_pct:.2f}%")
229 print(f"Mean range-normalized abs error: {mean_rng:.4f}")
230
231 # ---------- Elite variability ----------
232 var_mix, var_pred = _elite_variability_tables(elites_df)
233 print("\n=== ELITE VARIABILITY — MIX (top 10) ===")
234 print(var_mix)
235 print("\n=== ELITE VARIABILITY — PREDICTIONS (top 10) ===")
236 print(var_pred)
237
238 # ---------- Engineering plausibility ----------
239 eng_notes = _engineering_checks(mix)
240 if eng_notes:
241 print("\n[ENGINEERING CHECKS] Issues:")
242 for n in eng_notes:
243 print(" -", n)
244 else:
245 print("\n[ENGINEERING CHECKS] OK")
246
247 # ---------- Manifold proximity ----------
248 # Build DF with actual training column names and order for model/NN
249 mix_cols_order = [mix_map[c] for c in CFG["canon_mix_cols"] if c in mix_map]
250 x_vec = [mix.get(c, 0.0) for c in CFG["canon_mix_cols"] if c in mix_map]
251 x_df = pd.DataFrame([x_vec], columns=mix_cols_order)
252
253 d, _ = nn.kneighbors(x_df.values)
254 d0 = float(np.mean(d[0]))
255 d_scaled = d0 / nn_scale
256 print(
257 f"\nManifold proximity (kNN mean distance / scale): {d_scaled:.3f} "
258 f"(closer to 1 or below is better)"
259 )
260
261 # ---------- Local sensitivity (central difference) ----------
262 base = x_df.values.copy()
263 y_base = model.predict(x_df).flatten()
264 out_cols_order = list(out_map.values())
265
266 sens_rows = []
267 for j, canon in enumerate([c for c in CFG["canon_mix_cols"] if c in mix_map]):
268 x0 = base.copy()
269 step = max(abs(x0[0, j]) * perturb_frac, 1e-6)
270
271 x_up = x0.copy()
272 x_up[0, j] += step
273
274 x_dn = x0.copy()
275 x_dn[0, j] -= step
276 x_dn[0, j] = max(x_dn[0, j], 0.0)
277
278 y_up = model.predict(pd.DataFrame(x_up, columns=mix_cols_order)).flatten()
279 y_dn = model.predict(pd.DataFrame(x_dn, columns=mix_cols_order)).flatten()
280
281 # Central difference sensitivity dy/dx_j
282 dy = (y_up - y_dn) / (2.0 * step)
283
284 # Map to canonical output names
285 sens_rows.append(
286 {
287 out_actual_to_canon[col]: dy[i]
288 for i, col in enumerate(out_cols_order)
289 }
290 | {"_var": canon}
291 )
292
293 sens_df = (
294 pd.DataFrame(sens_rows)
295 .set_index("_var")
296 .reindex(CFG["canon_mix_cols"], fill_value=np.nan)
297 )
298
299 # Rank variables by average absolute influence
300 sens_rank = sens_df.abs().mean(axis=1).sort_values(ascending=False)
301 print(
302 "\n=== LOCAL SENSITIVITY (central diff, ±{:.0f}%) ===".format(
303 perturb_frac * 100
304 )
305 )
306 print(sens_rank)
307
308 return {
309 "per_output_table": per_out,
310 "mean_pct_error": float(mean_pct),
311 "mean_range_norm_abs_error": float(mean_rng),
312 "elite_var_mix": var_mix,
313 "elite_var_pred": var_pred,
314 "eng_issues": eng_notes,
315 "manifold_scaled_distance": float(d_scaled),
316 "sensitivity_matrix": sens_df,
317 "sensitivity_rank": sens_rank,
318 }
319
320
321# ---------------------------------------------------------
322# Constraint-aware local sensitivity (simplified)
323# ---------------------------------------------------------
324
325def local_sensitivity_with_constraints(
326 report_csv_path: str,
327 elites_csv_path: str,
328 delta_frac_of_sum: float = 0.05,
329 outputs_to_match=None,
330 save_excel_path: str | None = None,
331):
332 """
333 Perturb each ingredient by ± delta (fraction of total mix),
334 enforce bounds & simple constraints, and recompute:
335 - outputs
336 - a scalar loss (range-normalized RMSE + engineering penalties + manifold).
337 Agnostic to how the inverse mix was obtained (EA, RL, etc.).
338 """
339 report_df = pd.read_csv(report_csv_path)
340
341 # Load data/model
342 (
343 df,
344 model,
345 mix_map,
346 out_map,
347 X,
348 Y,
349 lb,
350 ub,
351 y_min,
352 y_max,
353 y_rng,
354 nn,
355 nn_scale,
356 ) = load_and_fit()
357
358 out_actual_to_canon = {v: k for k, v in out_map.items()}
359 out_cols_order = list(out_map.values())
360 mix_canon_order = [c for c in CFG["canon_mix_cols"] if c in mix_map]
361 mix_cols_order = [mix_map[c] for c in mix_canon_order]
362
363 # Base targets, mix
364 targ_cols = [c for c in report_df.columns if c.startswith("Target::")]
365 mix_cols = [c for c in report_df.columns if c.startswith("Mix::")]
366
367 base_targets = {
368 c.replace("Target::", ""): float(report_df[c].iloc[0]) for c in targ_cols
369 }
370 base_mix = {c.replace("Mix::", ""): float(report_df[c].iloc[0]) for c in mix_cols}
371
372 if outputs_to_match is None:
373 outputs_to_match = list(base_targets.keys())
374
375 # Build base X and Y
376 x_base = np.array([base_mix.get(c, 0.0) for c in mix_canon_order], dtype=float)
377 x_df = pd.DataFrame([x_base], columns=mix_cols_order)
378 y_all = model.predict(x_df.values)[0]
379
380 y_rng_vec = y_rng[out_cols_order].values
381 y_rng_vec = np.where(y_rng_vec > 0, y_rng_vec, 1.0)
382
383 # Range-normalized RMSE loss
384 def loss_fn(x_vec: np.ndarray) -> tuple[float, np.ndarray]:
385 x_vec = np.maximum(x_vec, 0.0)
386 x_df_ = pd.DataFrame([x_vec], columns=mix_cols_order)
387 y_ = model.predict(x_df_.values)[0]
388
389 errs = []
390 for i, col in enumerate(out_cols_order):
391 canon = out_actual_to_canon[col]
392 if canon in outputs_to_match and canon in base_targets:
393 t = base_targets[canon]
394 e = (y_[i] - t) / y_rng_vec[i]
395 errs.append(e)
396 rmse = math.sqrt(np.mean(np.square(errs))) if errs else 0.0
397
398 mix_dict = {c: float(x_vec[j]) for j, c in enumerate(mix_canon_order)}
399 eng_pen = engineering_constraints_penalty(mix_dict)
400
401 # Manifold penalty
402 d, _ = nn.kneighbors(x_df_.values)
403 d0 = float(np.mean(d[0]))
404 man_pen = CFG["manifold_lambda"] * (d0 / (nn_scale if nn_scale > 0 else 1.0))
405
406 total = rmse + eng_pen + man_pen
407 return float(total), y_
408
409 base_loss, base_y = loss_fn(x_base)
410 base_sum = float(np.sum(x_base))
411 base_y_map = {
412 out_actual_to_canon[col]: float(base_y[i])
413 for i, col in enumerate(out_cols_order)
414 }
415
416 delta_abs = delta_frac_of_sum * base_sum
417 rows = []
418
419 for j, var in enumerate(mix_canon_order):
420 for direction in ["+", "-"]:
421 x_vec = x_base.copy()
422 x_vec[j] = max(
423 x_vec[j] + (delta_abs if direction == "+" else -delta_abs), 0.0
424 )
425
426 loss, y_new = loss_fn(x_vec)
427 y_map = {
428 out_actual_to_canon[col]: float(y_new[i])
429 for i, col in enumerate(out_cols_order)
430 }
431
432 row = {
433 "Var": var,
434 "Dir": direction,
435 "Delta_abs_of_sum": float(delta_abs),
436 "Loss_new": loss,
437 "Loss_delta": float(loss - base_loss),
438 "Sum_after": float(np.sum(x_vec)),
439 }
440
441 for o in base_y_map.keys():
442 row[f"{o}::Pred_new"] = y_map[o]
443 row[f"{o}::dPred"] = y_map[o] - base_y_map[o]
444
445 for m in mix_canon_order:
446 idx = mix_canon_order.index(m)
447 row[f"Mix::{m}"] = float(x_vec[idx])
448
449 rows.append(row)
450
451 sens_table = pd.DataFrame(rows)
452
453 saved = None
454 if save_excel_path:
455 with pd.ExcelWriter(save_excel_path, engine="xlsxwriter") as xw:
456 pd.DataFrame(
457 [{"Base_Loss": base_loss, "Base_Sum": base_sum}]
458 ).to_excel(xw, sheet_name="Base", index=False)
459 pd.DataFrame(
460 [
461 {
462 "Output": k,
463 "Base_Pred": v,
464 "Target": base_targets.get(k, np.nan),
465 }
466 for k, v in base_y_map.items()
467 ]
468 ).to_excel(xw, sheet_name="Base_Preds", index=False)
469 sens_table.to_excel(xw, sheet_name="LocalSensitivity", index=False)
470 saved = save_excel_path
471
472 return {
473 "table": sens_table,
474 "base": {
475 "mix": base_mix,
476 "preds": base_y_map,
477 "targets": base_targets,
478 "loss": float(base_loss),
479 "sum": float(base_sum),
480 },
481 "paths": {"excel": saved},
482 }
483