jiang
5966d039de
Reorganize algorithm packages by business responsibility, move orchestration into services, and keep database access behind pooled repositories. Harden analysis API validation, remove unsafe legacy simulation endpoints, and add regression and architecture boundary coverage. BREAKING CHANGE: legacy algorithm module paths and obsolete simulation endpoints are removed.
260 lines
9.8 KiB
Python
260 lines
9.8 KiB
Python
from __future__ import annotations
|
|
|
|
from typing import Any
|
|
|
|
import numpy as np
|
|
import pandas as pd
|
|
from scipy.fft import fft, ifft
|
|
from sklearn.ensemble import IsolationForest
|
|
|
|
PressureDataInput = (
|
|
pd.DataFrame
|
|
| dict[str, list[Any]]
|
|
| list[dict[str, Any]]
|
|
| list[list[Any]]
|
|
| np.ndarray
|
|
)
|
|
IGNORED_OBSERVATION_COLUMNS = {"time", "timestamp", "datetime", "date"}
|
|
|
|
|
|
class PressureAnomalyDetector:
|
|
"""FFT + IsolationForest based burst detection for daily aligned pressure data."""
|
|
|
|
def __init__(
|
|
self,
|
|
*,
|
|
mu: int = 100,
|
|
points_per_day: int = 1440,
|
|
iforest_params: dict[str, Any] | None = None,
|
|
) -> None:
|
|
if points_per_day <= 0:
|
|
raise ValueError("points_per_day 必须大于 0。")
|
|
if mu <= 0:
|
|
raise ValueError("mu 必须大于 0。")
|
|
|
|
self.mu = int(mu)
|
|
self.points_per_day = int(points_per_day)
|
|
self.iforest_params = {
|
|
"n_estimators": 50,
|
|
"random_state": 42,
|
|
"contamination": "auto",
|
|
}
|
|
if iforest_params:
|
|
self.iforest_params.update(iforest_params)
|
|
|
|
self.data: np.ndarray | None = None
|
|
self.sensor_names: list[str] = []
|
|
self.high_freq_features: np.ndarray | None = None
|
|
|
|
def load_data(
|
|
self,
|
|
data_source: PressureDataInput,
|
|
*,
|
|
sensor_nodes: list[str] | None = None,
|
|
) -> pd.DataFrame:
|
|
"""
|
|
标准化输入观测数据为 DataFrame。
|
|
|
|
支持的 `data_source` 格式：
|
|
- `pd.DataFrame`
|
|
每一列代表一个传感器，每一行代表一个时间点。
|
|
- `dict[str, list[Any]]`
|
|
键为传感器 ID，值为该传感器按时间顺序排列的压力序列。
|
|
例如：`{"J1": [101.2, 101.0], "J2": [99.8, 99.7]}`。
|
|
- `list[dict[str, Any]]`
|
|
每个字典代表一个时间点，键为传感器 ID，值为该时刻压力。
|
|
例如：`[{"J1": 101.2, "J2": 99.8}, {"J1": 101.0, "J2": 99.7}]`。
|
|
- `list[list[Any]]`
|
|
二维列表，格式为 `(时间点数, 传感器数)`。
|
|
例如：`[[101.2, 99.8], [101.0, 99.7]]`。
|
|
- `np.ndarray`
|
|
二维数组，形状必须为 `(时间点数, 传感器数)`。
|
|
|
|
参数：
|
|
- `sensor_nodes`:
|
|
可选的传感器列筛选列表。传入后，数据中必须包含这些列名。
|
|
|
|
返回：
|
|
- 标准化后的 `pd.DataFrame`，列为传感器，行为时间点。
|
|
"""
|
|
if isinstance(data_source, np.ndarray):
|
|
observation_df = pd.DataFrame(data_source)
|
|
elif isinstance(data_source, pd.DataFrame):
|
|
observation_df = data_source.copy()
|
|
else:
|
|
observation_df = pd.DataFrame(data_source)
|
|
|
|
return self._normalize_observation_frame(
|
|
observation_df=observation_df, sensor_nodes=sensor_nodes
|
|
)
|
|
|
|
def process(
|
|
self,
|
|
observed_pressure_data: PressureDataInput,
|
|
*,
|
|
sensor_nodes: list[str] | None = None,
|
|
) -> np.ndarray:
|
|
"""
|
|
对输入压力序列按天切片，并提取每天末时刻的高频特征。
|
|
|
|
`observed_pressure_data` 的格式与 `load_data()` 一致，统一要求：
|
|
- 数据必须表示为“行=时间点、列=传感器”。
|
|
- 总行数必须是 `points_per_day` 的整数倍。
|
|
- 至少需要 2 天数据，即总行数 `>= 2 * points_per_day`。
|
|
|
|
例如：
|
|
- 当 `points_per_day=1440` 时，15 天数据的形状通常为 `(21600, 传感器数)`。
|
|
- 若传入 `sensor_nodes=["J1", "J2"]`，则输入中必须存在 `J1/J2` 两列。
|
|
|
|
返回：
|
|
- `np.ndarray`，形状为 `(天数, 传感器数)`，
|
|
每个值表示对应传感器在当天末时刻提取出的高频分量。
|
|
"""
|
|
observation_df = self.load_data(
|
|
observed_pressure_data,
|
|
sensor_nodes=sensor_nodes,
|
|
)
|
|
matrix = observation_df.to_numpy(dtype=float)
|
|
total_points, sensor_count = matrix.shape
|
|
if sensor_count == 0:
|
|
raise ValueError("压力观测数据中未找到可用传感器列。")
|
|
if total_points < self.points_per_day * 2:
|
|
raise ValueError("至少需要 2 天的观测数据才能执行爆管侦测。")
|
|
if total_points % self.points_per_day != 0:
|
|
raise ValueError("观测数据长度必须能被每日采样点数整除，以便按天切分。")
|
|
|
|
day_count = total_points // self.points_per_day
|
|
high_freq_features = np.zeros((day_count, sensor_count), dtype=float)
|
|
|
|
for sensor_idx in range(sensor_count):
|
|
sensor_series = matrix[:, sensor_idx]
|
|
for day_idx in range(day_count):
|
|
start = day_idx * self.points_per_day
|
|
end = (day_idx + 1) * self.points_per_day
|
|
day_data = sensor_series[start:end]
|
|
mirrored_data = np.concatenate([day_data, day_data[::-1]])
|
|
transformed = fft(mirrored_data)
|
|
transformed[self.mu : len(mirrored_data) - self.mu + 1] = 0
|
|
low_freq = ifft(transformed).real
|
|
high_freq = day_data - low_freq[: self.points_per_day]
|
|
high_freq_features[day_idx, sensor_idx] = float(high_freq[-1])
|
|
|
|
self.data = matrix
|
|
self.sensor_names = [str(column) for column in observation_df.columns]
|
|
self.high_freq_features = high_freq_features
|
|
return high_freq_features
|
|
|
|
def detect(self) -> pd.DataFrame:
|
|
if self.high_freq_features is None:
|
|
raise ValueError("特征未提取。请先调用 process()。")
|
|
|
|
day_count = self.high_freq_features.shape[0]
|
|
if day_count < 2:
|
|
raise ValueError("孤立森林至少需要 2 天特征数据。")
|
|
|
|
clf = IsolationForest(
|
|
n_estimators=self.iforest_params.get("n_estimators", 50),
|
|
max_samples=day_count,
|
|
random_state=self.iforest_params.get("random_state", 42),
|
|
contamination=self.iforest_params.get("contamination", "auto"),
|
|
**{
|
|
key: value
|
|
for key, value in self.iforest_params.items()
|
|
if key not in {"n_estimators", "random_state", "contamination"}
|
|
},
|
|
)
|
|
clf.fit(self.high_freq_features)
|
|
|
|
scores = clf.decision_function(self.high_freq_features)
|
|
predictions = clf.predict(self.high_freq_features)
|
|
result_df = pd.DataFrame(
|
|
{
|
|
"Day": range(1, day_count + 1),
|
|
"Score": scores.astype(float),
|
|
"Prediction": predictions.astype(int),
|
|
}
|
|
)
|
|
result_df["IsBurst"] = result_df["Prediction"].eq(-1)
|
|
result_df.attrs["sensor_nodes"] = self.sensor_names.copy()
|
|
result_df.attrs["high_freq_features"] = self.high_freq_features.copy()
|
|
result_df.attrs["day_count"] = day_count
|
|
result_df.attrs["points_per_day"] = self.points_per_day
|
|
result_df.attrs["sample_count"] = (
|
|
int(self.data.shape[0]) if self.data is not None else 0
|
|
)
|
|
return result_df
|
|
|
|
def run_detection(
|
|
self,
|
|
observed_pressure_data: PressureDataInput,
|
|
*,
|
|
sensor_nodes: list[str] | None = None,
|
|
) -> pd.DataFrame:
|
|
"""
|
|
执行完整爆管侦测流程。
|
|
|
|
输入格式与 `process()` 相同：
|
|
- `DataFrame` / `dict[str, list[Any]]` / `list[dict[str, Any]]` / `list[list[Any]]` / `np.ndarray`
|
|
- 行表示时间点，列表示传感器
|
|
- 总行数必须能被 `points_per_day` 整除
|
|
|
|
返回结果包含列：
|
|
- `Day`: 第几天（从 1 开始）
|
|
- `Score`: IsolationForest 异常分数，越小越异常
|
|
- `Prediction`: `-1` 表示异常，`1` 表示正常
|
|
- `IsBurst`: 是否判定为异常日
|
|
"""
|
|
self.process(observed_pressure_data, sensor_nodes=sensor_nodes)
|
|
return self.detect()
|
|
|
|
@staticmethod
|
|
def _normalize_observation_frame(
|
|
*,
|
|
observation_df: pd.DataFrame,
|
|
sensor_nodes: list[str] | None,
|
|
) -> pd.DataFrame:
|
|
if observation_df.empty:
|
|
raise ValueError("压力观测数据为空。")
|
|
|
|
normalized_df = observation_df.copy()
|
|
normalized_df.columns = [str(column) for column in normalized_df.columns]
|
|
normalized_df = normalized_df.drop(
|
|
columns=[
|
|
column
|
|
for column in normalized_df.columns
|
|
if column.lower() in IGNORED_OBSERVATION_COLUMNS
|
|
or column.lower().startswith("unnamed:")
|
|
],
|
|
errors="ignore",
|
|
)
|
|
|
|
if sensor_nodes:
|
|
selected_columns = [str(node) for node in sensor_nodes]
|
|
missing_columns = [
|
|
column
|
|
for column in selected_columns
|
|
if column not in normalized_df.columns
|
|
]
|
|
if missing_columns:
|
|
preview = ", ".join(missing_columns[:10])
|
|
raise ValueError(f"观测数据缺少传感器列: {preview}")
|
|
normalized_df = normalized_df.loc[:, selected_columns]
|
|
else:
|
|
candidate_df = normalized_df.apply(pd.to_numeric, errors="coerce")
|
|
normalized_df = candidate_df.loc[:, candidate_df.notna().any(axis=0)]
|
|
|
|
if normalized_df.empty:
|
|
raise ValueError("未识别到可用的数值型压力观测列。")
|
|
|
|
normalized_df = normalized_df.apply(pd.to_numeric, errors="coerce")
|
|
invalid_columns = [
|
|
column
|
|
for column in normalized_df.columns
|
|
if normalized_df[column].isna().any()
|
|
]
|
|
if invalid_columns:
|
|
preview = ", ".join(invalid_columns[:10])
|
|
raise ValueError(f"压力观测数据包含非数值或缺失值: {preview}")
|
|
|
|
return normalized_df.reset_index(drop=True)
|