机器学习中空间和时间自相关的分析：从理论基础到实践应用

import numpy as np
 import pandas as pd
 import matplotlib.pyplot as plt
 from sklearn.model_selection import train_test_split, GridSearchCV, TimeSeriesSplit
 from sklearn.metrics import mean_squared_error, r2_score
 from sklearn.ensemble import RandomForestRegressor
 from statsmodels.tsa.stattools import acf, pacf
 import geopandas as gpd
 from sklearn.preprocessing import StandardScaler
 
 # 设置随机种子以确保可重现性
 np.random.seed(42)

# 生成空间数据点（纬度、经度）
 n_points = 500
 latitudes = np.random.uniform(low=49, high=60, size=n_points)
 longitudes = np.random.uniform(low=-125, high=-100, size=n_points)
 
 # 生成时间序列数据（30天模拟）
 time_series_length = 30
 latitudes = np.repeat(latitudes, time_series_length)
 longitudes = np.repeat(longitudes, time_series_length)
 days = np.tile(np.arange(1, time_series_length + 1), n_points)
 
 # 模拟环境特征：温度、风速、湿度
 temperature = 20 + 10 * np.sin(0.1 * days) + np.random.normal(scale=2, size=len(days))
 wind_speed = 10 + 3 * np.cos(0.1 * days) + np.random.normal(scale=1, size=len(days))
 humidity = 50 + 10 * np.sin(0.05 * days) + np.random.normal(scale=5, size=len(days))
 
 # 生成具有时空自相关的火灾风险
 fire_risk = 0.3 * temperature + 0.4 * wind_speed - 0.2 * humidity + np.random.normal(scale=1, size=len(days))
 
 # 整合数据到DataFrame
 data = pd.DataFrame({
     'latitude': latitudes,
     'longitude': longitudes,
     'day': days,
     'temperature': temperature,
     'wind_speed': wind_speed,
     'humidity': humidity,
     'fire_risk': fire_risk
 })

# 添加时空滞后特征
 data['temp_lag1'] = data.groupby(['latitude', 'longitude'])['temperature'].shift(1).fillna(data['temperature'].mean())
 data['wind_lag1'] = data.groupby(['latitude', 'longitude'])['wind_speed'].shift(1).fillna(data['wind_speed'].mean())
 data['humidity_lag1'] = data.groupby(['latitude', 'longitude'])['humidity'].shift(1).fillna(data['humidity'].mean())
 
 # 移除包含NaN的行
 data.dropna(inplace=True)
 
 # 定义特征集和目标变量
 X = data[['temperature', 'wind_speed', 'humidity', 'temp_lag1', 'wind_lag1', 'humidity_lag1']]
 y = data['fire_risk']

# 划分训练集和测试集
 X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, shuffle=False)
 
 # 定义超参数网格
 param_grid = {
     'n_estimators': [50, 100, 200],
     'max_depth': [5, 10, 20],
     'min_samples_split': [2, 5, 10]
 }
 
 # 初始化随机森林回归器和时间序列交叉验证
 rf = RandomForestRegressor(random_state=42)
 tscv = TimeSeriesSplit(n_splits=5)
 
 # 执行网格搜索
 grid_search = GridSearchCV(rf, param_grid, cv=tscv, scoring='neg_mean_squared_error')
 grid_search.fit(X_train, y_train)
 
 # 获取最佳模型
 best_rf = grid_search.best_estimator_
 
 # 在测试集上进行预测
 y_pred = best_rf.predict(X_test)

# 计算RMSE和R²
 rmse = np.sqrt(mean_squared_error(y_test, y_pred))
 r2 = r2_score(y_test, y_pred)
 
 print(f'Root Mean Squared Error (RMSE): {rmse:.4f}')
 print(f'R-squared: {r2:.4f}')

plt.figure(figsize=(10, 6))
 plt.plot(y_test.values, label='Actual Fire Risk')
 plt.plot(y_pred, label='Predicted Fire Risk', linestyle='--')
 plt.legend()
 plt.title('Actual vs Predicted Fire Risk')
 plt.xlabel('Time')
 plt.ylabel('Fire Risk')
 plt.show()

# 计算残差
 residuals = y_test - y_pred
 lag_acf = acf(residuals, nlags=20)
 lag_pacf = pacf(residuals, nlags=20)
 
 # 绘制ACF和PACF
 plt.figure(figsize=(12, 6))
 plt.subplot(121)
 plt.stem(lag_acf)
 plt.title('ACF of Residuals')
 
 plt.subplot(122)
 plt.stem(lag_pacf)
 plt.title('PACF of Residuals')
 plt.show()

# 创建GeoDataFrame
 gdf = gpd.GeoDataFrame(data, geometry=gpd.points_from_xy(data.longitude, data.latitude))
 
 # 绘制火灾风险的空间分布
 gdf.plot(column='fire_risk', cmap='OrRd', legend=True)
 plt.title('Spatial Distribution of Fire Risk')
 plt.show()