### Mann-Kendall 趋势分析与突变检验的 Python 示例代码 以下是基于 `pymannkendall` 库实现 Mann-Kendall 趋势分析和突变检测的完整示例代码： #### 安装依赖库 首先需要安装必要的库，可以通过以下命令完成： ```bash pip install pymannkendall numpy pandas matplotlib scipy statsmodels ``` #### 示例代码 下面是一个完整的 Python 实现，展示如何使用 `pymannkendall` 进行趋势分析并绘制结果。 ```python import numpy as np import pandas as pd import pymannkendall as mk import matplotlib.pyplot as plt # 创建模拟时间序列数据 np.random.seed(42) dates = pd.date_range('2000-01-01', '2020-12-31', freq='Y') data = np.cumsum(np.random.randn(len(dates))) # 随机累积噪声作为示例数据 df = pd.DataFrame({'date': dates, 'value': data}) df.set_index('date', inplace=True) # 执行 Mann-Kendall 趋势测试 result = mk.original_test(df['value'], alpha=0.05) print("Mann-Kendall Test Result:") print(f"Trend: {result.trend}") print(f"Slope: {result.slope}") print(f"P-value: {result.p}") print(f"Z-statistic: {result.z}") # 绘制时间序列图及趋势线 plt.figure(figsize=(10, 6)) plt.plot(df.index, df['value'], label="Original Data", color="blue") if result.trend != "no trend": slope = result.slope intercept = df['value'][0] - slope * (df.index[0].year - 1970) fitted_values = [slope * (yr.year - 1970) + intercept for yr in df.index] plt.plot(df.index, fitted_values, label=f"Fitted Trend Line ({result.trend})", color="red", linestyle="--") plt.title("Time Series with Mann-Kendall Trend Analysis") plt.xlabel("Year") plt.ylabel("Value") plt.legend() plt.grid(True) plt.show() ``` 此代码实现了以下几个功能： 1. **生成随机时间序列数据**：创建一个具有年频率的时间索引，并填充随机累积噪声作为示例数据。 2. **执行 Mann-Kendall 测试**：调用 `mk.original_test()` 函数来评估是否存在显著趋势及其方向[^2]。 3. **可视化结果**：如果存在明显趋势，则拟合一条直线表示该趋势，并将其叠加到原始数据图表上。 --- #### 关于突变点检测 对于突变点检测（Change Point Detection），可以结合其他方法如 Pettitt 检验或二分法。这里给出一种简单的方法——利用滑动窗口技术寻找潜在的变化位置。 ```python from scipy.stats import rankdata def pettitt_test(data): n = len(data) k = range(n) r = rankdata(data, method='average') U = [(2*(i+1)-n)*np.sum(r[:i+1])/(n*(n+1)) for i in k] tau = np.argmax(U)+1 K = max(U)*(6/(n*(n+1)))**(0.5) p_value = 2*np.exp(-6*K**2/n) return {"change_point": tau, "statistic": K, "p_value": p_value} pettitt_result = pettitt_test(df['value']) print("\nPettitt Change Point Test Result:") print(f"Potential change point at year: {str(pettitt_result['change_point']+df.index[0].year)}") print(f"Statistic value: {pettitt_result['statistic']}, P-value: {pettitt_result['p_value']}") # 可视化突变点 cp_year = str(pettitt_result['change_point']+df.index[0].year) plt.axvline(pd.to_datetime(cp_year), color='green', linestyle='--', label=f'Pettitt Change Point ({cp_year})') plt.legend() plt.show() ``` 以上代码展示了如何应用 Pettitt 方法识别可能存在的突变点，并在图形中标记出来。 --- #### 结果解释 - 如果 Mann-Kendall 的 p 值小于指定阈值（通常为 0.05），则认为有显著趋势。 - 斜率正值表明上升趋势；负值代表下降趋势。 - 对于 Pettitt 检测，较低的 p 值意味着更有可能存在突变点。 ---