evaluate_model.py 7.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138
  1. #!/usr/bin/env python3
  2. """Evaluate a fixed PKS assistant model against the sealed test directory."""
  3. from __future__ import annotations
  4. import argparse
  5. import csv
  6. import importlib
  7. import json
  8. from collections import Counter, defaultdict
  9. from pathlib import Path
  10. import joblib
  11. import numpy as np
  12. from sklearn.metrics import average_precision_score, confusion_matrix, precision_recall_fscore_support
  13. def feature_module(model_type: str):
  14. """Pressure and vibration use independent feature modules."""
  15. return importlib.import_module(
  16. "generalized_features" if model_type == "pressure" else "generalized_features_vibration"
  17. )
  18. def load_samples(directory: Path) -> list[dict]:
  19. return [json.loads(path.read_text(encoding="utf-8")) for path in sorted(directory.glob("*.json"))]
  20. def matrix(feature_dicts: list[dict], names: list[str]) -> np.ndarray:
  21. return np.asarray([[item.get(name) for name in names] for item in feature_dicts], dtype=float)
  22. def metrics(y: np.ndarray, predicted: np.ndarray, probabilities: np.ndarray) -> dict:
  23. labels = [0, 1]
  24. per_precision, per_recall, per_f1, per_support = precision_recall_fscore_support(
  25. y, predicted, labels=labels, zero_division=0
  26. )
  27. tn, fp, fn, tp = confusion_matrix(y, predicted, labels=labels).ravel()
  28. samples = int(len(y))
  29. normal_support, anomaly_support = int(per_support[0]), int(per_support[1])
  30. specificity = float(tn / (tn + fp)) if (tn + fp) else 0.0
  31. accuracy = float((tp + tn) / samples) if samples else 0.0
  32. return {
  33. "samples": samples, "normal_samples": normal_support, "anomaly_samples": anomaly_support,
  34. "tp": int(tp), "fp": int(fp), "fn": int(fn), "tn": int(tn),
  35. "precision": float(per_precision[1]), "recall": float(per_recall[1]), "f1": float(per_f1[1]),
  36. "pr_auc": float(average_precision_score(y, probabilities)) if len(set(y)) == 2 else None,
  37. "specificity": specificity, "accuracy": accuracy,
  38. "balanced_accuracy": float((per_recall[0] + per_recall[1]) / 2),
  39. "macro_f1": float((per_f1[0] + per_f1[1]) / 2),
  40. "weighted_f1": float((per_f1[0] * normal_support + per_f1[1] * anomaly_support) / samples) if samples else 0.0,
  41. "class_metrics": {
  42. "normal": {"precision": float(per_precision[0]), "recall": float(per_recall[0]), "f1": float(per_f1[0]), "support": normal_support},
  43. "anomaly": {"precision": float(per_precision[1]), "recall": float(per_recall[1]), "f1": float(per_f1[1]), "support": anomaly_support},
  44. },
  45. "confusion_matrix_counts": {"tn": int(tn), "fp": int(fp), "fn": int(fn), "tp": int(tp)},
  46. "confusion_matrix_normalized": {
  47. "normal": {"predicted_normal": float(tn / normal_support) if normal_support else 0.0, "predicted_anomaly": float(fp / normal_support) if normal_support else 0.0},
  48. "anomaly": {"predicted_normal": float(fn / anomaly_support) if anomaly_support else 0.0, "predicted_anomaly": float(tp / anomaly_support) if anomaly_support else 0.0},
  49. },
  50. }
  51. def event_metrics(samples: list[dict], predicted: np.ndarray) -> dict:
  52. groups = defaultdict(list)
  53. for sample, value in zip(samples, predicted):
  54. groups[sample["event_group_id"]].append((sample, int(value)))
  55. positives = [values for values in groups.values() if values[0][0]["label"] == 1]
  56. negatives = [values for values in groups.values() if values[0][0]["label"] == 0]
  57. tp = sum(any(value for _, value in values) for values in positives)
  58. fn = len(positives) - tp
  59. fp = sum(any(value for _, value in values) for values in negatives)
  60. return {
  61. "anomaly_events": len(positives), "normal_events": len(negatives),
  62. "detected_anomaly_events": tp, "missed_anomaly_events": fn, "false_positive_normal_events": fp,
  63. "event_recall": tp / len(positives) if positives else None,
  64. }
  65. def main() -> None:
  66. parser = argparse.ArgumentParser(description="仅对封存PKS测试集做最终模型评估")
  67. parser.add_argument("model_type", choices=("pressure", "vibration"))
  68. parser.add_argument("--datasets-dir", default=str(Path(__file__).resolve().parent / "datasets"))
  69. parser.add_argument("--models-dir", default=str(Path(__file__).resolve().parent / "models"))
  70. parser.add_argument("--output-dir", default=str(Path(__file__).resolve().parents[2] / "cache" / "pks_scanner" / "model_assistant_generalized"))
  71. args = parser.parse_args()
  72. gf = feature_module(args.model_type)
  73. bundle = joblib.load(Path(args.models_dir) / f"{args.model_type}_model.joblib")
  74. metadata = json.loads((Path(args.models_dir) / f"{args.model_type}_metadata.json").read_text(encoding="utf-8"))
  75. samples = load_samples(Path(args.datasets_dir) / args.model_type / "test")
  76. if not samples or not any(sample["label"] == 1 for sample in samples):
  77. parser.error("封存测试集必须包含至少一个异常样本")
  78. x = bundle["imputer"].transform(matrix(gf.feature_matrix(samples, bundle["unit_baselines"]), bundle["feature_names"]))
  79. y = np.asarray([sample["label"] for sample in samples])
  80. probabilities = bundle["model"].predict_proba(x)[:, 1]
  81. predicted = (probabilities >= metadata["threshold"]).astype(int)
  82. distances, neighbors = bundle["knn"].named_steps["knn"].kneighbors(bundle["knn"].named_steps["scale"].transform(x), n_neighbors=1)
  83. rows = []
  84. for index, sample in enumerate(samples):
  85. rows.append({
  86. "sample_id": sample["sample_id"], "unit": sample["unit"], "batch": sample["batch"],
  87. "event_group_id": sample["event_group_id"], "start_hour": sample["start_hour"], "end_hour": sample["end_hour"],
  88. "true_label": sample["label"], "anomaly_probability": probabilities[index], "predicted_label": int(predicted[index]),
  89. "nearest_training_sample_index": int(neighbors[index][0]), "nearest_training_distance": float(distances[index][0]),
  90. "nearest_training_similarity": float(1 / (1 + distances[index][0])),
  91. })
  92. by_unit = {}
  93. for unit in sorted({sample["unit"] for sample in samples}):
  94. indices = [index for index, sample in enumerate(samples) if sample["unit"] == unit]
  95. by_unit[unit] = metrics(y[indices], predicted[indices], probabilities[indices])
  96. report = {
  97. "model_type": args.model_type, "evaluation_set": "sealed_test", "threshold": metadata["threshold"],
  98. "window_metrics": metrics(y, predicted, probabilities), "event_metrics": event_metrics(samples, predicted),
  99. "by_unit": by_unit,
  100. "warning": "Labels originate from current scanner rules unless manually updated; metrics measure held-out rule-label agreement, not confirmed physical-fault diagnosis.",
  101. }
  102. output = Path(args.output_dir)
  103. output.mkdir(parents=True, exist_ok=True)
  104. with (output / f"{args.model_type}_test_predictions.csv").open("w", encoding="utf-8-sig", newline="") as handle:
  105. writer = csv.DictWriter(handle, fieldnames=list(rows[0]))
  106. writer.writeheader()
  107. writer.writerows(rows)
  108. (output / f"{args.model_type}_test_evaluation.json").write_text(json.dumps(report, ensure_ascii=False, indent=2), encoding="utf-8")
  109. summary = report["window_metrics"]
  110. print(f"测试集:{summary['samples']} 窗口(正常 {summary['normal_samples']}、异常 {summary['anomaly_samples']}),阈值 {report['threshold']:.6f}")
  111. print(f"{'类别':<6}{'查准率':>10}{'召回率':>10}{'F1':>10}{'样本数':>10}")
  112. for name, label in (("正常", "normal"), ("异常", "anomaly")):
  113. item = summary["class_metrics"][label]
  114. print(f"{name:<6}{item['precision']:>10.4f}{item['recall']:>10.4f}{item['f1']:>10.4f}{item['support']:>10}")
  115. print(
  116. f"准确率 {summary['accuracy']:.4f} 平衡准确率 {summary['balanced_accuracy']:.4f} "
  117. f"macro-F1 {summary['macro_f1']:.4f} 加权F1 {summary['weighted_f1']:.4f} "
  118. f"特异度 {summary['specificity']:.4f} PR-AUC {summary['pr_auc']:.4f}"
  119. )
  120. print(json.dumps(report, ensure_ascii=False, indent=2))
  121. if __name__ == "__main__":
  122. main()