用PythonNumPy手把手实现一个马尔可夫链天气预测模型附完整代码天气预报总是让人又爱又恨——明明说今天晴天结果出门就被淋成落汤鸡。但你知道吗那些看似随机的阴晴变化其实可以用数学来捕捉规律。今天我们就用Python和NumPy从零开始构建一个能预测明天天气的马尔可夫链模型。不需要高深的数学知识跟着代码一步步来你会发现概率模型原来如此有趣又实用1. 马尔可夫链的核心概念想象你正在观察一个善变的天气系统。今天的天气会影响明天的天气但不会影响后天——这就是马尔可夫链的精髓无记忆性。用技术术语来说就是系统下一个状态的概率分布只取决于当前状态。三个关键要素构成了我们的天气预测模型状态集合晴天、雨天、阴天三种天气状态转移概率矩阵描述从今天天气到明天天气的转换概率初始状态向量起始日的天气概率分布让我们用具体数字来说明。假设今天的转移概率矩阵是这样的# 行代表今天天气列代表明天天气 # 顺序为[晴天, 雨天, 阴天] transition_matrix np.array([ [0.6, 0.2, 0.2], # 今天晴天 → 明天 [0.3, 0.4, 0.3], # 今天雨天 → 明天 [0.1, 0.5, 0.4] # 今天阴天 → 明天 ])这个矩阵告诉我们如果今天是晴天那么明天有60%概率继续晴天20%转雨20%转阴。其他行的解读方式相同。2. 环境准备与数据建模2.1 安装必要库确保你的Python环境已安装以下库pip install numpy matplotlib2.2 构建天气状态模型我们将天气系统抽象为离散的状态空间import numpy as np # 定义天气状态 STATES [晴天, 雨天, 阴天] NUM_STATES len(STATES) # 创建状态映射字典 state_to_idx {state: i for i, state in enumerate(STATES)} idx_to_state {i: state for i, state in enumerate(STATES)}转移矩阵的验证非常重要——每行概率之和必须为1def validate_transition_matrix(matrix): for row in matrix: if not np.isclose(row.sum(), 1.0): raise ValueError(每行概率之和必须等于1) return True3. 预测模型实现3.1 单步预测实现给定今天的天气预测明天的天气分布def predict_next_day(current_state, transition_matrix): current_idx state_to_idx[current_state] return transition_matrix[current_idx]试试看预测晴天后一天的天气weather_probs predict_next_day(晴天, transition_matrix) print(dict(zip(STATES, weather_probs))) # 输出{晴天: 0.6, 雨天: 0.2, 阴天: 0.2}3.2 多步预测与收敛想知道一周后的天气趋势连续应用转移矩阵即可def predict_n_days(initial_state, transition_matrix, n_days): current_vec np.zeros(NUM_STATES) current_vec[state_to_idx[initial_state]] 1.0 predictions [] for _ in range(n_days): current_vec np.dot(current_vec, transition_matrix) predictions.append(current_vec.copy()) return np.array(predictions)可视化7天预测结果import matplotlib.pyplot as plt def plot_predictions(predictions, days7): plt.figure(figsize(10, 6)) for i, state in enumerate(STATES): plt.plot(range(days), predictions[:days, i], labelstate) plt.xlabel(天数) plt.ylabel(概率) plt.title(未来天气趋势预测) plt.legend() plt.grid(True) plt.show() # 假设今天是雨天 predictions predict_n_days(雨天, transition_matrix, 50) plot_predictions(predictions, 20)你会观察到概率逐渐收敛到稳定值这就是马尔可夫链的稳态分布。4. 从理论到实践完整案例4.1 基于历史数据的参数估计真实的转移概率如何获得从历史数据统计# 示例历史数据30天的天气记录 historical_data [ 晴天, 晴天, 阴天, 雨天, 阴天, 阴天, 雨天, 雨天, 晴天, 晴天, 晴天, 阴天, 阴天, 雨天, 晴天, 晴天, 晴天, 阴天, 雨天, 雨天, 阴天, 晴天, 晴天, 晴天, 阴天, 雨天, 阴天, 阴天, 晴天, 晴天 ] def estimate_transition_matrix(data): matrix np.zeros((NUM_STATES, NUM_STATES)) for today, tomorrow in zip(data[:-1], data[1:]): today_idx state_to_idx[today] tomorrow_idx state_to_idx[tomorrow] matrix[today_idx][tomorrow_idx] 1 # 归一化 row_sums matrix.sum(axis1, keepdimsTrue) return matrix / row_sums estimated_matrix estimate_transition_matrix(historical_data) print(从历史数据估计的转移矩阵) print(estimated_matrix)4.2 完整预测流程整合所有步骤的完整示例class WeatherPredictor: def __init__(self, states): self.states states self.state_to_idx {s:i for i,s in enumerate(states)} self.idx_to_state {i:s for i,s in enumerate(states)} self.transition_matrix None def fit(self, historical_data): n len(self.states) matrix np.zeros((n, n)) for today, tomorrow in zip(historical_data[:-1], historical_data[1:]): i self.state_to_idx[today] j self.state_to_idx[tomorrow] matrix[i,j] 1 row_sums matrix.sum(axis1, keepdimsTrue) self.transition_matrix matrix / row_sums return self def predict(self, current_state, n_days1): if self.transition_matrix is None: raise ValueError(请先使用fit方法训练模型) current_vec np.zeros(len(self.states)) current_vec[self.state_to_idx[current_state]] 1 predictions [] for _ in range(n_days): current_vec current_vec self.transition_matrix predictions.append(current_vec.copy()) return predictions def simulate(self, initial_state, days30): current_state initial_state sequence [current_state] for _ in range(days-1): probs self.transition_matrix[self.state_to_idx[current_state]] next_state np.random.choice(self.states, pprobs) sequence.append(next_state) current_state next_state return sequence # 使用示例 predictor WeatherPredictor([晴天, 雨天, 阴天]) predictor.fit(historical_data) # 预测未来5天的概率分布 forecast predictor.predict(晴天, 5) for day, probs in enumerate(forecast, 1): print(f第{day}天预测:, dict(zip(STATES, probs))) # 模拟未来30天的天气序列 simulated_weather predictor.simulate(晴天, 30) print(\n模拟天气序列:, simulated_weather)5. 模型优化与扩展5.1 处理稀疏数据问题当某些状态转换在历史数据中未出现时可以添加平滑处理def estimate_with_smoothing(data, alpha0.1): 使用加一平滑(Laplace smoothing)估计转移矩阵 matrix np.ones((NUM_STATES, NUM_STATES)) * alpha for today, tomorrow in zip(data[:-1], data[1:]): i state_to_idx[today] j state_to_idx[tomorrow] matrix[i,j] 1 row_sums matrix.sum(axis1, keepdimsTrue) return matrix / row_sums5.2 季节性因素整合天气往往具有季节性模式我们可以为不同季节建立不同的转移矩阵class SeasonalWeatherPredictor: def __init__(self, states, seasons): self.states states self.seasons seasons self.models { season: WeatherPredictor(states) for season in seasons } def fit(self, historical_data, seasonal_labels): from collections import defaultdict season_data defaultdict(list) for season, weather in zip(seasonal_labels, historical_data): season_data[season].append(weather) for season, data in season_data.items(): self.models[season].fit(data) return self def predict(self, current_state, current_season, n_days1): return self.models[current_season].predict(current_state, n_days)5.3 可视化工具增强使用更丰富的可视化展示预测结果def plot_weather_simulation(sequence): plt.figure(figsize(12, 4)) colors {晴天: gold, 雨天: dodgerblue, 阴天: lightgray} for i, state in enumerate(sequence): plt.bar(i, 1, colorcolors[state], edgecolorwhite) plt.yticks([]) plt.xlabel(天数) plt.title(天气模拟序列) # 创建图例 handles [plt.Rectangle((0,0),1,1, colorcolors[s]) for s in colors] plt.legend(handles, colors.keys()) plt.show() # 使用示例 simulated predictor.simulate(晴天, 60) plot_weather_simulation(simulated)这个马尔可夫链天气预测模型虽然简单但已经包含了构建概率预测系统的核心要素。在实际项目中你可以通过增加更多天气状态如雪天、雾天、引入更复杂的时间依赖关系或者结合其他气象数据来进一步提升预测准确度。