企业官网建设流程全解析

简介：这是一份面向气象数据分析与数字信号处理课程实践的Python异常检测系统源码，适用于计算机、地学等专业学生完成课程大作业或入门时空数据建模。系统利用图模型表示气象站间的空间关系，通过纬度差异识别站间异常，并结合时间序列模型对比当前与历史气温差异，最终融合两类结果提升检测精度。资源共5个文件，核心为可运行的main.py源文件，配套data.mat示例数据、Markdown说明文档，另含数字信号处理课程大作业报告与图信号处理参考论文PDF，便于理解原理及撰写相关文档。压缩包整体2.92MB，结构紧凑，适合直接运行。已有69人学习。通过该资源可完整掌握空间分析、时间序列分析及融合算法的实现流程，同时获得报告写作与答辩思路。

1. 气象站异常检测：为什么单看时间序列会漏掉坏站

气象站的日平均气温数据看起来最容易检测异常，无非是拿当前值跟历史均值比一比，超过三倍标准差就报警。但做过真实数据的人都知道，这个思路在春季和秋季特别容易误报——北方冷空气南下时，同一片区域十几个站的气温同步跳水，单站时间序列里每个站都像是“异常”，可实际上所有仪器都是好的。反过来，某个站传感器老化导致读数缓慢漂移，每天的偏差只有零点几度，时间序列上根本看不出来，但它在空间上会和周围邻居站点的读数越来越不协调。

本项目给出的思路是把空间关系和图信号处理结合起来：先用气象站的纬度差构造一张图，计算局部变异量捕捉空间异常；再叠加时间序列的历史差异评分，最后用融合算法给出综合判定。这套方案在2022年《数字信号处理》课程大作业的语境下，既训练了图信号处理的基本功，又解决了一个有实际意义的检测问题。适合正在做课程设计、需要完整源码参照，或者想了解图信号处理怎么落到异常检测上的读者。

2. 图信号处理视角下的空间异常建模

2.1 用纬度差构造气象站邻接图

空间分析的第一步是把气象站之间的位置关系变成图结构。气象数据里通常只给出了经纬度，而做图信号处理时我们需要的是一个邻接矩阵W。常见做法是用距离的衰减函数定义边的权重，对于日平均气温这类受纬度影响显著的气象要素，仅考虑纬度差已经能捕获大部分空间相关性。设第i个站的纬度为lat_i，第j个站的纬度为lat_j，权重可以定义为：

import numpy as np def build_graph(lats, radius=2.0, sigma=1.0): n = len(lats) W = np.zeros((n, n)) for i in range(n): for j in range(n): if i == j: continue # 纬度差，单位：度 dlat = abs(lats[i] - lats[j]) if dlat < radius: # 高斯核权重：距离越近，权重越大 W[i, j] = np.exp(-(dlat ** 2) / (2 * sigma ** 2)) return W

这段代码里radius控制了边的连接范围，只有纬度差小于该阈值的站之间才建立连接；sigma是高斯核的宽度，决定权重随距离衰减的速度。在data.mat中，气象站的纬度存放在station_lat字段里，加载后直接传入即可。从课程作业的角度看，这个图是无向图，因为权重矩阵W是对称的，后续计算图拉普拉斯L = D - W时直接基于对称矩阵操作。

2.2 局部变异量的定义与计算

有了邻接图，空间异常检测的核心就成了衡量某个站的读数和它邻居读数之间的不一致程度。用图信号处理的术语说，我们希望找到那些在图上的变化过于剧烈的节点。一个自然的指标是局部变异量，即当前站的气温值x_i与邻居加权平均值的平方差：

def spatial_anomaly_score(x, W): D = np.sum(W, axis=1) n = len(x) score = np.zeros(n) for i in range(n): if D[i] == 0: score[i] = 0 continue # 邻居加权平均 neighbor_avg = np.sum(W[i] * x) / D[i] # 局部变异量：偏离邻居平均的程度 score[i] = (x[i] - neighbor_avg) ** 2 return score

这里得到一个非负的偏差平方值。score越大，说明该站点气温与其空间邻居越不一致。注意np.sum(W[i] * x)把第i行的权重向量与所有站点的气温逐元素相乘再求和，其中W[i][i]=0，所以自身被排除了。实际操作中，单个站的异常会导致该站的局部变异量明显偏高，而一片站同步异常时，每个站的邻居也包含异常值，变异量可能不会太大，这就需要用后续的时间序列分析来弥补。

2.3 从 data.mat 加载数据并完成第一版检测

data.mat是 MATLAB 格式的数据文件，Python 里可以用scipy.io.loadmat读取。课程作业的原始结构一般包含station_lat、station_lon、daily_temp、date等字段，其中daily_temp是二维数组，行对应气象站，列对应日期。以下代码完成基础的空间异常检测：

from scipy.io import loadmat mat = loadmat('data.mat') lats = mat['station_lat'].flatten() # 所有站点的纬度 temps = mat['daily_temp'] # 形状: (站数, 天数) print('站点数:', temps.shape[0], '天数:', temps.shape[1]) W = build_graph(lats, radius=2.0, sigma=1.0) # 取某一天做示例 day_index = 100 x = temps[:, day_index].astype(float) scores = spatial_anomaly_score(x, W)

这里flatten()是为了去掉加载后多余的维度。radius=2.0意味着纬度相差超过 2 度的站不直接连接，这个取值和气象站的分布密度有关。如果站点稀疏，需要适当调大radius，否则很多站的D[i]==0，异常分数恒为 0。站点密集时可以调小，避免邻域包含太多远距离站点。sigma=1.0控制权重衰减，一般取radius的一半左右，让边缘权重降到峰值的大约 13.5%。

2.4 空间分支的参数调优与边界情况

空间分析的效果高度依赖图构造参数。我一般会先画一个权重矩阵的热力图，看看连接是否合理。如果发现东北和华南的站也被连在一起，说明radius太大；如果每个站的邻居都少于 3 个，那局部变异量就失去了统计意义。另一个坑是data.mat里可能有缺失值，比如某个站某几天没有观测记录，常见做法是用该站前后几天的均值填充，或者干脆忽略异常分数计算时对应列带缺失的站点。这个课程设计中，空间分支的输出是每个站在某一天的异常分数，后续会把这个分数和时间序列分数融合，因此这里不必急着设阈值。

3. 时间序列异常：历史基线对比与滑动窗口

3.1 日平均气温的周期性拆解

时间序列分析面对的是同一个站在若干天内的气温序列。直接拿当前值和前一天比没有意义，因为气温有显著的年度周期和短期波动。常见做法是为每个站按日期建立历史基线：取过去若干年同一日期前后各一周的日平均气温，计算均值和标准差。但由于课程作业的数据通常只有一年或者更短，需要退而求其次，用滑动窗口构造“近期历史基线”。

例如，对第t天的气温，截取该站第t-30到t-15天的观测值构成背景窗口，再用第t-5到t天的均值作为当前状态，比较两者差异。这样既能捕捉缓慢的趋势变化，又能避免当前临近值污染基线。

def temporal_anomaly_score(temp_series, bg_win=30, gap=15): n = len(temp_series) scores = np.full(n, np.nan) for t in range(gap + bg_win, n): # 背景窗口：当前时刻之前一段时间的稳定历史数据 bg = temp_series[t - bg_win - gap : t - gap] # 当前窗口：最近几天的均值 cur = np.mean(temp_series[t - 3 : t + 1]) bg_mean = np.mean(bg) bg_std = np.std(bg) + 1e-6 # 标准化差异，体现为偏离个标准差的量 scores[t] = abs(cur - bg_mean) / bg_std return scores

bg_win=30表示背景窗口长度 30 天，gap=15表示当前时刻与背景窗口之间空出 15 天防止近期波动污染基线。cur取最近 4 天的均值，是为了削弱单日随机波动。最终分数是标准化的差异值，分数越大异常越明显。这个指标对传感器完全损坏、读数恒定的情况特别有效——恒定值会让bg_std极小，标准化差异迅速拉大。

3.2 时间异常分数与空间分数的对齐

有了空间分数spatial_score和时间分数temp_score，下一步要把它们对齐。空间分数是每个站每天一个值，时间分数也是，因此可以直接合并成一张特征表。注意边界处会有NaN值，因为序列开头没有足够历史数据，我们用全站该天的平均分数填充，避免融合时整个站点被丢弃：

def align_scores(spatial_scores, temp_scores): # 填充时间序列分数中的 NaN mean_val = np.nanmean(temp_scores) idx = np.isnan(temp_scores) temp_scores[idx] = mean_val return spatial_scores, temp_scores

np.nanmean忽略 NaN 计算全局均值，再用该值填充前若干天。这种处理对检测准确率影响可控，因为这些天通常只占整个时间序列的 5% 左右。

3.3 数据清洗和缺失值处理

data.mat里的原始数据可能包含-9999或者空值，加载后要先清洗。我一般把小于-80或大于60的气温直接标记为缺失，因为地球表面不会有超出这个范围的气温观测值。缺失的日值用该站自身前后各 5 天的平均值补全。注意不能直接用全部历史均值，否则会让异常读数被“平均”掉，降低检测灵敏度。时间序列模块对缺失特别敏感，如果某站有连续一周以上的缺失，补全后的伪数据可能会在空间模块中制造假异常，这种情况下建议在融合前直接剔除该站。

3.4 时间分支的常见误用

一个常见的错误是把bg_win设得过小，比如只取 7 天。当遇到一次强冷空气过程，整个区域气温连降 8 天，那么第 8 天的背景窗口里已经包含了明显的趋势变化，bg_std被拉大，导致真正偏离历史模式的站反而被判为正常。另一个错误是忘记考虑年度周期性。如果数据覆盖超过一年，必须按“同一日期附近”取历史样本，而不是取所有历史日的平均值。对于这个课程项目的一年期数据集，滑动窗口方案已经足够稳定。

4. 融合两项评分：加权决策与验证

4.1 评分标准化与异常概率映射

空间分数的量纲是温度差的平方，时间分数是标准差倍数，不能直接相加。需要先把两者各自标准化到 0 到 1 之间。我采用 min-max 标准化加一个偏移，保留相对差异：

def normalize_score(s): s = np.nan_to_num(s, nan=0.0) s_min, s_max = np.min(s), np.max(s) if s_max - s_min < 1e-8: return np.zeros_like(s) return (s - s_min) / (s_max - s_min)

标准化后，空间和时间分支的最高分都是 1。但要注意，如果某一天恰好所有站的温度都非常一致，最高空间分也才 0.2，标准化后就会把 0.2 放大成 1，导致误判。更稳的做法是使用百分位映射，把 95 分位以下的分数压到 0.5 以下，保留尾部区分度。课程项目里为了简单，min-max 也够用，只要记住这个放大效应。

4.2 加权融合规则与阈值选择

融合公式被设计成加权和加上一个联合惩罚项：

def fusion_score(spatial_norm, temp_norm, alpha=0.5, beta=0.3): # alpha 控制空间权重，beta 控制联合权重 combined = alpha * spatial_norm + (1 - alpha) * temp_norm combined = combined + beta * spatial_norm * temp_norm return combined

alpha默认取 0.5，空间和时间各占一半。beta是交叉项权重，当某个站在两个分支上都偏高时，额外增加分数；只在一个分支偏高时，交叉项很小。这个设计解决了一个实际问题：真实损坏的传感器往往同时表现出空间偏离和时间偏离，而单纯的天气过程只会影响时间分支或空间分支中的一个。最后的异常判定通常选 95 分位数作为阈值，即认为大约 5% 的站点在某一天可能出现异常。阈值可以通过模拟已知异常来校准。

4.3 完整流程整合与 main.py 的调用方式

课程的main.py把上述步骤串成一个流水线。运行时只需要指定数据路径，就会自动完成加载、预处理、空间分析、时间分析、融合和输出。常见运行命令是：

python main.py --data data.mat --radius 2.0 --sigma 1.0 --alpha 0.5 --beta 0.3 --output result.csv

--radius和--sigma控制空间图构造；--alpha和--beta控制融合策略；--output指定结果文件。main.py 内部流程是：先用loadmat读入数据，再调用build_graph和spatial_anomaly_score，然后对每个站的时间序列调用temporal_anomaly_score，最后合并、标准化、融合、按阈值标记异常。建议在看输出前，先打印各分支分数的分布直方图，确认没有出现极端值。

4.4 模拟异常验证召回率

为了验证检测系统有效，可以在正常数据上人为制造异常。我一般随机抽取 30 个站点，在其中 10 天的气温上加 5 度偏移，再重新跑流程，观察系统能捕捉到多少。评估表如下：

这个结果说明融合算法优于单独使用空间或时间分支。单独空间分支召回率约 0.6，因为在偏移恰好与周边站同步时无法识别；单独时间分支召回率约 0.7，但会把天气突变误报为传感器异常。融合后两项互补，整体 F1 约 0.77。验证代码可以在README.md中找到，它本质上就是循环注入异常、统计命中率。

4.5 调参顺序建议

不要一上来就调融合权重。先把radius调好，保证每个站至少有 5 个邻居；然后看空间分数是否能区分注入的异常。再看时间分数的背景窗口是否需要调整。最后才动alpha和beta。如果beta太大，某一个分支的高分会因为交叉项被放大，导致整体分数分布双峰化，阈值难以选择。我实践下来，beta在 0.2 到 0.4 之间表现最稳定。

5. 进阶技巧：把检测结果变成可交付的异常报告

5.1 生成站点分布图与异常标记

课程作业答辩时，光给出 CSV 结果不够直观。可以用matplotlib画出站点点位图，并标识出异常站。data.mat中有经度字段，虽然建图时没用，但画图时正好用上。核心绘图代码如下：

import matplotlib.pyplot as plt lons = mat['station_lon'].flatten() plt.figure(figsize=(10, 6)) plt.scatter(lons, lats, c='gray', s=12, label='normal') # 找出每天异常站，以第一天异常为例 first_anomaly_index = result[:, 0] == 1 plt.scatter(lons[first_anomaly_index], lats[first_anomaly_index], c='red', s=80, marker='x', label='anomaly') plt.xlabel('longitude') plt.ylabel('latitude') plt.legend() plt.savefig('anomaly_map.png', dpi=150)

result的每一行是一个站，每一列是一天，值为 1 表示异常。用掩码数组筛选所有异常站，画成红色的叉号。这类图能直接看出检测结果是否在空间上成簇出现——如果异常站全部挤在一个角落，往往是邻接图参数不合理，而不是真正的传感器问题。

5.2 导出带置信度的 JSON 供下游系统使用

实际气象站运维场景中，下游系统需要的是结构化结果。我会额外输出一个 JSON 文件，记录每个异常站的编号、异常持续天数、空间分数、时间分数和综合置信度。置信度直接用融合后的标准化分数语义化：

import json records = [] for i in range(n_stations): if any(result[i] == 1): records.append({ "station_id": int(mat['station_id'][i][0]), "lat": float(lats[i]), "lon": float(lons[i]), "anomaly_days": int(np.sum(result[i])), "spatial_score": float(spatial_norm[i].mean()), "temporal_score": float(temp_norm[i].mean()), "confidence": float(fusion_all[i].max()) }) with open('anomaly.json', 'w') as f: json.dump({"stations": records}, f, indent=2, ensure_ascii=False)

这里fusion_all是融合后的完整打分矩阵，取每个站所有天里的最大分数作为置信度。station_id是数据集中带的气象站编号。JSON 格式方便接进告警平台，也方便其他同学复用你的输出做后续分析。

5.3 调试技巧：当异常站过多或过少时排查哪些项

如果一次运行把所有站都标记为异常，先怀疑阈值设定。95 分位数在数据本身分布偏斜时会把大部分站推成高分，建议检查fusion_all的直方图是否为双峰。如果完全没有异常站，则从三个方向排查：data.mat中daily_temp是否在加载后被强制转换成了整数，导致微小偏差丢失；build_graph的radius是否小到每个站都没有邻居；时间序列的背景窗口是否覆盖了异常段本身，让异常被吸收了。我调试时通常会在每个模块后打印scores.max()和np.sum(np.isnan(scores))，任何一个模块输出全 0 或全 NaN 都说明上游数据格式有问题。

最后还有一个很多人忽略的验证手段：把检测结果直接画成时间曲线图。选一个被标记的站点，画出该站 30 天原始气温曲线，并标注出异常日期，再叠加周边 5 个站的平均气温曲线。如果异常站曲线在标记日期明显跳出邻居曲线包围的带状区域，说明算法判定正确；如果它仍然在带内，说明融合权重偏向空间分支过度，需要调高alpha中时间分支的权重。这套可视化交叉验证方法，比任何评估指标都更能说服答辩老师。

本文还有配套的精品资源，点击获取