2025-03-24 09:57:14 +08:00
{
"cells": [
{
"metadata": {},
"cell_type": "markdown",
"source": [
"# 预测建模\n",
"北京市空气质量指数预测( 推荐难度系数10) \n",
"\n",
"这个数据集是北京市2022年11月1日至2023年10月31日期间空气质量相关数据。\n",
"根据这个数据集,回答以下问题"
],
"id": "b610f839dca4877"
},
{
"metadata": {
"collapsed": true,
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T07:01:00.767221Z",
"start_time": "2025-03-24T07:00:54.883547Z"
2025-03-24 09:57:14 +08:00
}
},
2025-03-24 15:19:11 +08:00
"cell_type": "code",
2025-03-24 09:57:14 +08:00
"source": [
2025-03-24 15:19:11 +08:00
"import os\n",
"import sys\n",
"\n",
2025-03-24 09:57:14 +08:00
"#导入基础包\n",
"import pandas as pd\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
2025-03-24 15:19:11 +08:00
"from calculate import *\n",
"from heatmap import *\n",
"from statsmodels.graphics.tsaplots import plot_acf\n",
"import matplotlib.font_manager as fm\n",
"font_path = '/System/Library/Fonts/STHeiti Light.ttc' # 替换为实际可用的字体文件路径\n",
"font_prop = fm.FontProperties(fname=font_path)\n",
"plt.rcParams['font.family'] = font_prop.get_name()\n",
2025-03-24 09:57:14 +08:00
"\n",
2025-03-24 15:19:11 +08:00
"# 导入主成分分析相关包\n",
2025-03-24 09:57:14 +08:00
"from factor_analyzer import Rotator\n",
"from factor_analyzer.factor_analyzer import calculate_bartlett_sphericity, calculate_kmo\n",
2025-03-24 15:19:11 +08:00
"from sklearn.decomposition import PCA\n",
"from sklearn.preprocessing import StandardScaler\n",
2025-03-24 09:57:14 +08:00
"\n",
2025-03-24 15:19:11 +08:00
"# 导入SARIMA相关包\n",
"from statsmodels.tsa.statespace.sarimax import SARIMAX\n",
"from sklearn.metrics import mean_squared_error, mean_absolute_error\n",
"\n",
"# 导入XGBOOST相关包\n",
2025-03-24 09:57:14 +08:00
"from xgboost import XGBRegressor\n",
"from scipy.stats import randint, uniform\n",
"from sklearn.model_selection import RandomizedSearchCV\n",
"from matplotlib.dates import DateFormatter, HourLocator\n",
"\n",
2025-03-24 15:19:11 +08:00
"# 导入单独写的函数\n",
2025-03-24 09:57:14 +08:00
"from calculate import *\n",
"from heatmap import *\n",
"from sort_matrix import *"
],
2025-03-24 15:19:11 +08:00
"id": "initial_id",
2025-03-24 09:57:14 +08:00
"outputs": [],
2025-03-24 15:19:11 +08:00
"execution_count": 2
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:16.367026Z",
"start_time": "2025-03-24T03:36:15.877757Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
2025-03-24 15:19:11 +08:00
"# 设置字体\n",
"if sys.platform == 'darwin': # macOS\n",
" font_path = '/System/Library/Fonts/STHeiti Light.ttc'\n",
"elif sys.platform == 'win32': # Windows\n",
" plt.rcParams['font.sans-serif'] = ['SimHei'] # Windows系统自带黑体\n",
" plt.rcParams['axes.unicode_minus'] = False # 解决负号显示问题\n",
"else: # Linux/其他系统\n",
" font_path = '/usr/share/fonts/truetype/wqy/wqy-zenhei.ttc' # 文泉驿字体\n",
"\n",
"# 仅非Windows系统需要加载字体文件\n",
"if sys.platform != 'win32':\n",
" try:\n",
" font_prop = fm.FontProperties(fname=font_path)\n",
" plt.rcParams['font.family'] = font_prop.get_name()\n",
" except:\n",
" print(f\"警告:{font_path} 字体加载失败,请检查路径有效性\")\n",
"# 读取数据\n",
2025-03-24 09:57:14 +08:00
"data=pd.read_excel('北京市空气质量指数与气象数据.xlsx')\n",
2025-03-24 15:19:11 +08:00
"data.head()\n",
"\n",
"try:\n",
" os.mkdir('./images')\n",
"except FileExistsError:\n",
" pass"
2025-03-24 09:57:14 +08:00
],
"id": "92ea7ba1218799cd",
2025-03-24 15:19:11 +08:00
"outputs": [],
"execution_count": 16
2025-03-24 09:57:14 +08:00
},
{
"metadata": {},
"cell_type": "markdown",
"source": [
"## 题目1\n",
"研究单日内空气质量指数与各项指标的变化趋势,这种趋势是否具有周期性?"
],
"id": "bca65e544d8bef55"
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:16.615153Z",
"start_time": "2025-03-24T03:36:16.534329Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
2025-03-24 15:19:11 +08:00
"# 数据预处理:将数据按小时分组,计算每个小时各指标的平均值\n",
"# 转换Excel日期序列值为实际日期并分组\n",
"data['datetime'] = pd.to_datetime(data['date']) + pd.to_timedelta(data['hour'], unit='h')\n",
"valid_hours = sorted(data['hour'].unique())\n",
"hourly_data = data.groupby('hour').mean().loc[valid_hours]\n",
"plt.figure(figsize=(12, 8))\n",
"indicators = ['AQI', 'PM2.5', 'PM10', 'CO', 'NO2', 'O3','SO2']\n",
"colors = ['#2d87bb', '#ff7f0e', '#2ca02c', '#d62728', '#9467bd', '#8c564b', '#e377c2', '#7f7f7f', '#bcbd22', '#17becf', '#1f77b4', '#ffbb78', '#98df8a', '#d62728', '#9467bd', '#8c564b', '#e377c2', '#7f7f7f', '#bcbd22', '#17becf', '#1f77b4', '#ffbb78', '#98df8a', '#d62728', '#9467bd', '#8c564b', '#e377c2', '#7f7f7f', '#bcbd22', '#17becf', '#1f77b4', '#ffbb78', '#98df8a', '#d62728',]\n",
2025-03-24 09:57:14 +08:00
"\n",
2025-03-24 15:19:11 +08:00
"normalized = (hourly_data[indicators] - hourly_data[indicators].mean(axis=0)) / hourly_data[indicators].std(axis=0)\n"
2025-03-24 09:57:14 +08:00
],
2025-03-24 15:19:11 +08:00
"id": "118b1b48e798a7ba",
"outputs": [
{
"data": {
"text/plain": [
"<Figure size 1200x800 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"execution_count": 17
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:19.304829Z",
"start_time": "2025-03-24T03:36:16.655941Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
2025-03-24 15:19:11 +08:00
"source": [
"# 绘制各指标小时均值变化趋势(标准化后)折线图\n",
"for i, indicator in enumerate(indicators):\n",
" plt.plot(normalized.index, normalized[indicator], \n",
" marker='o',label=indicator, color=colors[i], linewidth=2)\n",
"\n",
"plt.title('各指标小时均值变化趋势(标准化后)', fontsize=14)\n",
"plt.xlabel('小时', fontsize=12)\n",
"plt.ylabel('标准化值', fontsize=12)\n",
"plt.xticks(range(0, 24))\n",
"plt.grid(alpha=0.3)\n",
"plt.legend(bbox_to_anchor=(1.05, 1), loc='upper left')\n",
"plt.tight_layout()\n",
"\n",
"# 新增保存代码\n",
"plt.savefig('images/hourly_trends_combined.png', dpi=300, bbox_inches='tight') # 保存组合大图\n",
"plt.show()\n",
"\n",
"# 新增保存子图代码\n",
"for i, indicator in enumerate(indicators):\n",
" plt.figure(figsize=(8, 5))\n",
" plt.plot(normalized.index, normalized[indicator], \n",
" marker='o', color=colors[i], linewidth=2)\n",
" plt.title(f'{indicator}小时均值变化趋势(标准化后)')\n",
" plt.xlabel('小时')\n",
" plt.ylabel('标准化值')\n",
" plt.xticks(range(0, 24))\n",
" plt.grid(alpha=0.3)\n",
" plt.tight_layout()\n",
" plt.savefig(f'images/hourly_{indicator}.png', dpi=300) # 保存单个指标子图\n",
" plt.close()"
],
"id": "57dedbd9b7bbe12d",
"outputs": [
{
"data": {
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
],
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAnIAAAHWCAYAAADzS2TwAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOydd3xN5//A3/fe3Oy9SYgIoqi9ZxRFaatoVWu22v761aIturToUKOt6tCp6FSlirbUnrGpLUgiyN47dz6/P25yiezkZuB5v15ecp7znOd8zl3ncz5TIYQQSCQSiUQikUhuO5S1LYBEIpFIJBKJpHJIRU4ikUgkEonkNkUqchKJRCKRSCS3KVKRk0gkEolEIrlNkYqcRCKRSCQSyW2KVOQkEolEIpFIblOkIieRSCQSiURymyIVOYlEIpFIJJLbFKnISSQSiUQikdymSEXuLkM28pBIap7Tp09z4cKF2hbjjuHkyZOEh4fXthgSSZ1AKnJ3Ef/++y8tWrTgm2++qdTxUVFRpKWlYTQaAdBqtSxfvpxNmzZVeK3du3czceJEfv7550rJUh62bdvGP//8w+nTp6vtHH5+frRu3ZrMzMxKHb9582YOHTpEbm5uuY+ZM2cOK1euLPacgwYNIiQkhEuXLlVKnry8PDZv3sxff/1VqeOBIjfYjIwM8vLyKrVWZGQk48aN4/vvvy+yb+HChcyfP5/o6OhyrbVz504WL15MREREkX379+9n+/btxMfHV0rO0oiLi6N3794kJycXu7+6vwsajYbJkyezefNmcnNz0Wq1/PPPP+zcubPca1jyPS0Per2eQ4cOERkZWez+c+fOMWDAALKzs6tNBonktkFI7hr69+8vPD09RVhYWKWO79mzp1Cr1eLEiRNCCCGMRqOwsrIS/fv3LzI3PT1dvP322+Ljjz8udq2NGzcKQGzevNk8duXKFTFx4kQxZcoUMWPGDPHqq6+KV199VUyfPl28+OKLYuLEiRWSt1u3bgIQO3fuLNf8pKQkkZGRITQajdDpdEKr1YrMzEyRkJBQ4jEBAQGiS5cupa47YcIEERwcXOzr3q5dOwGIWbNmlUtGjUYjnJychLe3tzhz5oxITU0VOp3OvD8wMFA0b9680DFGo1FkZ2eLrKysMtdPS0sTnp6eolWrVuWS51aysrKEUqkUHTp0EAcPHhRCCLFq1SoBiJEjR1Z4vcjISAGIr7/+utC4wWAQfn5+onXr1iIzM7Nca3300UcCKHb+o48+KgBx9erVCstYFqNGjRIDBgwocX91fxf2798vANGyZUuRkZEh0tLShLOzs+jatWu5jrfUe2o0Gs3nT05OFtHR0eLkyZNix44d4rfffhPz5s0T48ePF926dRP29vYCECNGjCh2LYPBIJo0aSJmzpxZ7vNLJHcqVrWkP0pqmPXr17Nt2zYAgoODS5xnZ2fH5cuXqV+/fpF9Xl5e1K9fn7Zt2wKgUCiws7PDw8MDAKPRyIEDB/j111/54YcfyMzMxN7ensGDB9O8efNCazk7OwNgY2NjHktNTWX58uUolUqsra1RqVQIITAYDGi1WoQQxVpmSsLR0RGA9u3bl2v+vHnz+Pjjj4uMjxgxgjVr1hR7jL29fZnr5uTkEBYWRuPGjQuNp6WlcerUKRo3bszcuXPLJeOGDRvIzMzknXfeYePGjbz++uuA6b2wsrJCp9OhVCqxtbUFTO+JXq9HCMHs2bOZM2cOAMuXL+fKlSs4OTlhZ2eHjY0NVlamn4Pg4GD279/P0qVLzden1WrJy8sjJyeH9PR0Ro0aZf4c3Mz+/fsxGo0YjUY6dOhQ6DUaPHhwua4RIDk5GTc3N/N1WFtbI4QgJSUFDw8PVq9eTXR0NJ6enjz++OMAZGVlERISYr5GgOvXr5OVlYWdnZ3Z6hkfH09kZCT33HMPe/bsQa1Wk5qaip2dHdHR0Vy7do3MzEx69+6NnZ1duWUujsOHD/Pbb7+xZcuWEudY6ruwf/9+Zs2axfr1681rAubv/YcffoiTkxNgei/WrVtHWFhYqb8HBeta4j1NSEjA19e3xP1qtRo3Nzfc3Nxo3bo1rq6u6PV6rl+/jr+/f6G5SqWSl156iZdeeonnn3+eRo0alVsOieROQypydwHJycm8+OKLAIwbN67YH981a9awdu1apk+fXqwSB4VvNAUolUr+++8/HnroIQ4cOEBSUhL169dn+PDhhISE0LZtW/z8/IocV5wC1KpVK9LS0nBxcSmyTwhBenp6mdd6M2q12ixjeXBwcABMSk4B06dPN48XR4HyUxolzdmwYQMGg4EBAwZw/Phx87heryc3NxdbW1u6detW6Jhly5bh6OjIM888Q3h4OH5+ftjZ2Zlv9g899BDDhw/n6aefRqFQoNfr0ev1aDQamjVrZl7n77//Zu3atWb5rKys0Ov1GAwG7O3tcXBwYMqUKQghsLGxwWg0mhUIgA4dOhSryK1fvx4wuT0Lrruk9yE8PJzAwMBi35+WLVsSHx9vXmPSpElMnDgRlUpFZmYmb7/9NsOHD6dTp04ArF27lsuXL7Nq1apC63z44YcsWbIEGxsbs+xNmjQBIDExkaeeeorr168jhEChUNC/f3/zdebk5BT7vlWERYsW0bhxY/r371/iHEt9F3777Td27drF2LFjWbdunfl1Xbt2Lb6+vuj1erPLvG/fvjzyyCOEhYVx9uxZdDod7du3p2nTpkXWtdR76unpab62L774Ajs7O5ycnHBxccHZ2bnU71lxjBkzhhkzZjB//ny++uqrCh0rkdxR1IodUFJj6HQ60b9/f1GvXj3x4osvCi8vryKuwr179wobGxvRqVMnodVqC+1LSUkR3333nVizZo3o06eP8Pb2FmvXrhXfffedyMvLEy4uLsLFxUWMHDlSLFmyRJw5c6bQ8VeuXCm0/cEHH4hvvvlGvPPOOwIQr732mvj666/FV199ZfFrHzJkSImutOKYPXu2uPUrERAQIMaPHy9mzpwpJk+eLF5++WUxc+ZMs6vL29tb1KtXT0yePFk8+eSTIikpqci648aNE0AhF6gQJlc1UOK/hx9+uND8o0ePCkAEBASYx7Kzs0Vubq4QQoj4+HgBiA8//NC832AwFOtSjYuLE/Hx8SIhIUFoNBohhBBPP/20CAoKMs/p16+f2W2s1+uFECY3W3R0dLGvaXZ2tnB3dxedO3cuNL5p0yYBiOXLlxea26hRI9G3b19x8eLFImt98cUX4vvvvxfLli0TgJg8ebJYsWKF+Oijj8SkSZNEw4YNzdcVGRkp7O3ti7hfhRAiMzPTfH1ffvmlAIRerxfJycnCYDCI7OxsIYQQzZs3F2PHjhVCmD4H9erVK7JWRYmNjRVKpVK88sorxe639HfBYDCIYcOGCUBMnTpVCHHjM1Oefze/PwVY8j0VQgg7OzvRo0ePSl1fcYwYMULY29uXK2xAIrlTkYrcHUxeXp54+OGHhUKhEJs2bRJ5eXkiODhY3HfffeYb2MGDB4W7u7to3rx5sbFg58+fL/GHPzU1Vbi4uIhRo0aVKEOrVq3MNxUhRIlr3axAVIZt27aJZcuWiQsXLpjHblXkcnNzxffffy+uXbtW7BqlKXItW7YUCoVC2NjYCCcnJ7MCq1QqhUKhEAqFQgCF1n7//ffF7NmzRfv27QUg5s6dK2bOnCnOnz8vdu3aJQDxwgsviNOnTxf6999//wlAPP7444VkGThwYBFF7s033xSAUCqVwsbGRgDCxsZG2Nvbm7dHjx5d4us2fPhw0bZtW3H58uUSFbnU1FTRtWtXMXXqVLNCVxyffPKJAMSTTz5ZaLy4m/7kyZMFIBo2bGiOuyqOxMTEQseeOHFCBAYGiiVLloioqCgRExMjhgwZItq0aSOuXbsmLl68KFJSUoqsk5mZaX5/33jjDfHAAw+ItLQ08fzzz4tXXnlFKBQK8e233wohTJ+Dm1/jyvLjjz8KQPz666/F7q+O70JKSopo1aqV2LJlixBCiPvvv18A4uOPPxbXrl0T165dEwcOHDCPRUVFicuXL4tz584V+xBi6ffUw8OjwoqcTqcTaWl
},
"metadata": {},
"output_type": "display_data"
}
],
"execution_count": 18
},
{
"metadata": {
"ExecuteTime": {
"end_time": "2025-03-24T03:36:25.435504Z",
"start_time": "2025-03-24T03:36:19.402487Z"
}
},
"cell_type": "code",
"source": [
"# ACF检验周期性\n",
"# 创建完整时间序列(每小时一个样本,缺失值用线性插值)\n",
"full_idx = pd.date_range(start=data['datetime'].min(),\n",
" end=data['datetime'].max(),\n",
" freq='h')\n",
"full_series = data.set_index('datetime').reindex(full_idx)\n",
"interpolated = full_series[indicators].interpolate(method='time')\n",
"\n",
"# 绘制ACF图, 检验3天周期( 24*3) \n",
"plt.figure(figsize=(60, 20)) # 调整整体画布尺寸\n",
"for i, indicator in enumerate(indicators):\n",
" ax = plt.subplot(2, 4, i+1) # 创建2行4列的子图布局\n",
" plot_acf(interpolated[indicator].dropna(),\n",
" lags=72,\n",
" alpha=0.05,\n",
" title=f'{indicator}',\n",
" color=colors[i],\n",
" ax=ax)\n",
" plt.xticks(np.arange(0, 73, 12))\n",
"plt.tight_layout()\n",
"plt.savefig('./images/all_acf_subplots.png', dpi=200, bbox_inches='tight')\n",
"plt.show()\n",
"\n",
"for i, indicator in enumerate(indicators):\n",
" plt.figure(figsize=(12, 6))\n",
" plot_acf(interpolated[indicator].dropna(),\n",
" lags=72,\n",
" alpha=0.05,\n",
" title=f'{indicator} ACF',\n",
" color=colors[i])\n",
" plt.xticks(np.arange(0, 73, 12))\n",
" plt.savefig(f'./images/acf_{indicator}.png', dpi=200, bbox_inches='tight')\n",
" plt.close()\n"
],
"id": "5f8e89a8d1561e4f",
"outputs": [
{
"data": {
"text/plain": [
"<Figure size 6000x2000 with 7 Axes>"
],
"image/png": "iVBORw0KGgoAAAANSUhEUgAAF2YAAAfFCAYAAADkPWjkAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdfZhWBZ038O+5gTEmxAE0lxmZIcFNcQMxc43K1XJDzOXVjBZdfHq2rNRerK18FpJJeirTdd1l07b1ZVfbXlZewjbb2qfdnlLzqZZkV2wjiqEYKhcaiACBmfP8MSwty6AgM3PPDJ/PdXGdi/M79833kPbHb47fU5RlWQYAAAAAAAAAAAAAAAAAAAAAAAAAAAAAYACrVDsAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAEBPU8wOAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAx4itkBAAAAAAAAAAAAAAAAAAAAAAAAAAAAgAFPMTsAAAAAAAAAAAAAAAAAAAAAAAAAAAAAMOApZgcAAAAAAAAAAAAAAAAAAAAAAAAAAAAABjzF7AAAAAAAAAAAAAAAAAAAAAAAAAAAAADAgKeYHQAAAAAAAAAAAAAAAAAAAAAAAAAAAAAY8BSzAwAAAAAAAAAAAAAAAAAAAAAAAAAAAAADnmJ2AAAAAAAAAAAAAAAAAAAAAAAAAAAAAGDAU8wOAAADzFe/+tW84Q1veNbrvva1r2Xu3LkZM2ZMjjvuuIwePToXXXRR7r333uzdu/eg68eOHZsFCxb0RGQAAAAAAOh2V111VYqiSFEUGTRoUGprazNu3LgsWLAg7e3t+6/7oz/6o/zgBz94xu8qyzKvf/3rDzrf0tKSN7zhDRk5cmSGDh2al7/85fnKV77yrNna29tTW1u7P99//vrIRz5y5DcKAAAAAAA9oDf27Fu3bs2YMWNSFEVOOeWULj+7adOmXHHFFRk1alSe//zn51WvelW+/e1vH93NAQAAAABAH7Bz587ceeed+Z3f+Z2cfPLJGTp0aMaOHZvf//3fzz/8wz90+Zknn3wyb3/72/OiF70otbW1Ofnkk3P++efn4x//eHbs2NHLdwAAAP2XYnYAABhgPvShD+Xv/u7vsm7dui7n7e3teetb35oLLrggP/rRj3L99dfn05/+dD70oQ/l5JNPzpve9Ka87GUvy49+9KNeTg4AAAAAAN1rwoQJefDBB/P5z38+n/jEJ/Ka17wmH/nIR/LOd75z/zVPPvlkFi9e/Izfs2LFiixbtuyAF5v+5Cc/yW//9m/nO9/5Tj70oQ/lE5/4RIYNG5aLL744Dz744DN+35NPPpmdO3fmz//8z/Pggw/u//W6173uqO4XAAAAAAC6U0/u2ZNk27ZtueOOO/KGN7yhy8/t3LkzF154Yb72ta/lpptuyic/+cmUZZkLLrgg3/3ud4/29gAAAAAAoGoef/zxTJ48Oe94xzvS1NSUj3zkI7n//vvzrne9Kz//+c9z8cUXZ86cOdm5c+f+z9x2220566yz8s///M954xvfmPvuuy+33HJLzjjjjLznPe/J5MmT8/jjj1fxrgAAoP8oyrIsqx0CAADoHt/85jfzspe9LEny5je/OZ/4xCcOuuYd73hH/vzP/zy33357rrvuuoPm3/3udzN9+vQMGzYs3/zmNzN8+PAkydixY3PFFVc860PzAAAAAADQF1x11VX5wQ9+kG984xsHnH/ve9+bP/mTP0lbW1uGDRuWF73oRVm3bl2efPLJnHbaaQd9T1mWOfvss/Pd7343a9euzfjx45Mk11xzTT73uc/lySefzIknnrj/2qlTp6a1tTX/9m//dshsf/3Xf503velN+eUvf5njjjuuG+8aAAAAAAC6R0/v2f+rRYsW5a/+6q/yk5/85IDzf/Znf5Z3v/vdWb16dc4444wkyZ49e/Kyl70sw4cPz1e/+tVuvGMAAAAAAOgdP/7xj/PSl740dXV1+fznP58XvehFB13z0EMPZe7cubnpppvy9re/PXfddVf+8A//MAsXLkxzc3OKojjg+vXr12f69On5j//4j3zrW99KQ0NDb90OAAD0S5VqBwAAALrPhz70oYwePTrvete78td//ddpbW09YP7444/nz//8z3PDDTd0WcqeJGeddVa++MUv5gc/+EFuueWW3ogNAAAAAAC9ZsSIEWlvb8+uXbvS3t6eH/3oRymKIjfddFOX169YsSKrV69Oknz/+9/ff/7rX/96ZsyYsb+UPUmKosj06dPz7//+78+Y4Tvf+U7OPPPM1NTUpL29vRvuCgAAAAAAekd37dkPx+c///lMnTp1fyl7kgwZMiR/9Ed/lH/+53/OT3/60+d+IwAAAAAAUCXvf//7s3PnznzlK1/pspQ9SaZNm5a1a9fm7W9/e7Zt25Z3vetdmT9/fj74wQ8eVMqeJGPHjs0//MM/ZNeuXfnjP/7jnr4FAADo9xSzAwDAALF69ep84QtfyPXXX5/3vve9SZLbbrvtgGs+9alP5bjjjssNN9zwjN/1W7/1W7n88stz11139VheAAAAAADoaR0dHdm+fXu2b9+en/70p1m6dGk+9rGPZcqUKTnxxBOzfv367NmzJ+95z3vyt3/7twcVwpRlmQ9+8IO5/PLLM3LkyKxdu3b/7Mtf/vJBe/gkefLJJzN27NhnzPWd73wn3/ve9zJ06NAcd9xx+d3f/d1873vf65Z7BgAAAACA7tKTe/bD8cQTT+Tcc8896PyUKVNSlmUef/zxo7o/AAAAAADobb/85S/zuc99Lm95y1syZsyYZ7z2BS94QZLkc5/7XH75y19m4cKFz3j96NGj85a3vCWf/vSn86tf/arbMgMAwECkmB0AAAaID33oQxkxYkSuvvrq/MZv/Ebmz5+fO++8M7/4xS/2X7NmzZqceeaZGTZs2LN+33nnnZfW1tZs3bq1J2MDAAAAAECPefTRR3P88cfn+OOPz+jRo/O6170u55xzTj772c8myf4CmPe9732ZMGFCPvjBDx7w+RUrVmT16tX5wAc+kHHjxh1QGPMbv/EbOf744w+4/jvf+U7uvvvuXHvttYfMVJZlTj755Nxyyy1ZuXJl7r777mzatCnnn39+Nm3a1F23DgAAAAAAR60n9+yHY8uWLTnppJMOOv8bv/EbSZKnnnrqudwWAAAAAABUzbp167J379684hWvOOzPrFmzJieeeGLGjRv3rNe+7GUvy+7du7Nu3bqjiQkAAAOeYnYAABgAvv/97+eBBx7I3Llzs3Xr1vzkJz/J3Llzs3379ixZsuSAa4uiOOjzmzZtyr/927/l6aefPmjW0dHRY7kBAAAAAKAnvfjFL843vvGNPPLII3n88cezefPmfPnLX84pp5ySpLMw5sQTT0xdXV1uvPHGfPrTn86///u/J+ksUP/gBz+Y17/+9TnjjDMyfvz4fP/73z/kn/WDH/wgl156aV75ylfmuuuuO+R1RVFk2bJlueaaa/Ka17wmf/AHf5Cvfe1rKcsyt99+e/f+BQAAAAAAwFHozT37oXT1/HtZlgccAQAAAACgv2hvb0+SDBo06LA/cyTdL5VKZ71kV/t1AADg1xSzAwDAAPDhD384HR0dueOOOzJmzJiMGTMmr3rVq5Ikt99+e3bs2JEkOeOMM/LEE09k+/btB3x++/bt+Z3f+Z286U1v2n/uW9/6VhoaGjJixIjeuxEAAAAAAOhGw4cPz8tf/vK87GUvy8SJEw/aea9duzbjxo1LksyePTsvfvGL88EPfjBJsmLFiqxevTof+MAHkiTjxo3L2rVru/xzWlpa8upXvzqjR4/OAw88sP9h9sM1atSovOIVr8h3vvOdI71FAAAAAADoMb21Zz+UkSNH5j/+4z8OOv+zn/0sSXLSSScd8T0BAAAAAEA1jR8/PoMHD863vvWtw/7MGWeckf/4j//ID3/4w2e99tvf/naOO+64/ft7AACga4rZAQCgn2tpacmnPvWpzJs3Lw8++OABv+66665s3rw5n/zkJ5MkV1xxRXbt2pUPf/jDB3zHaaedlqVLl+Yzn/lMbr311qxbty6f/exn88Y
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x600 with 0 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"execution_count": 19
2025-03-24 09:57:14 +08:00
},
{
"metadata": {},
"cell_type": "markdown",
"source": [
2025-03-24 15:19:11 +08:00
"\n",
"\n",
2025-03-24 09:57:14 +08:00
"## 题目2\n",
"简述各项指标间的相互关系。"
],
"id": "59e20f3463e819a6"
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:26.513462Z",
"start_time": "2025-03-24T03:36:25.565039Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"#计算相关系数矩阵\n",
"correlation_matrix = data.iloc[:, 1:].corr()\n",
"#绘制热力图\n",
2025-03-24 15:19:11 +08:00
"plot_heatmap(correlation_matrix,20,16,title=\"Correlation Matrix Heatmap\",save_path=\"./images/correlation_heatmap.png\")"
2025-03-24 09:57:14 +08:00
],
"id": "c917d14115569bcd",
"outputs": [
{
"data": {
"text/plain": [
"<Figure size 2000x1600 with 2 Axes>"
],
2025-03-24 15:19:11 +08:00
"image/png": "iVBORw0KGgoAAAANSUhEUgAABxUAAAY1CAYAAADkUIAZAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd1gTSR8H8C9I7yACitIFG4piFxXUs4uenr1X7GevZy/YsfeO7ezC6enZe0GxF7qg0hUVQgmBvH9wRkMCogcEeL+f58nzyOzs5jfj7mR2Z3dWSSwWi0FERERERERERERERERElANlRQdAREREREREREREREREREUbBxWJiIiIiIiIiIiIiIiIKFccVCQiIiIiIiIiIiIiIiKiXHFQkYiIiIiIiIiIiIiIiIhyxUFFIiIiIiIiIiIiIiIiIsoVBxWJiIiIiIiIiIiIiIiIKFccVCQiIiIiIiIiIiIiIiKiXHFQkYiIiIiIiIiIiIiIiIhyxUFFIiIiIiIiIiIiIiIiIsoVBxWJiIiI6P9GSkoKNm/ejKZNm8LU1BSampqwsrJCr169cO7cOUWHJ2XAgAFwcXH56fVdXV0xYMAAyd87d+6Enp4ezp49mw/R5d3jx48xb9687+bz8vLCrVu3/tN3ZS+zPP+1Hnbv3g0lJSWIRCK5y5WUlLB79+6f2vaPsrKywty5cwvlu4iIiIiIiIg4qEhERERE/xceP36MmjVr4vfff4elpSWWLFmCffv2Yfz48YiNjUXr1q3RpUsXpKSkKDrUApGSkoKUlBQIhcJC/d6goCB4eXnh48ePOeZJSUnBsmXLEBQUVODxKKoeiIiIiIiIiIo7DioSERERUYn35s0btGrVCgDw5MkT7N27FwMHDkSXLl3w+++/48KFCzhz5gwuXLiAbdu2KTjagjFq1CikpqbC3d29UL83KCgInz59gpeXV455Nm/ejOjo6EIZVFRUPRAREREREREVdxxUJCIiIqISb9q0aUhJScH58+fh4OAgN0+bNm0QFBSEsWPHSqX7+PigadOm0NXVhY6ODpo0aQJfX1+pPEpKSpgyZQqqV68OfX19PH/+XG7aF4cPH0aDBg2gra0NXV1duLi44PDhw98tx6NHj9C5c2eULl0ampqaqFu3Lo4cOYLVq1fDysoKwNfpOa9evYo9e/ZIpuO8cuUKVFRU8Pr1a8n2BAIBZs+eDXt7e6irq6Ns2bIYNGgQIiIipL7X1dUVffr0wb1799CiRQtoa2ujbNmymDp16nef+AsKCoKKigpWr16NhIQEmeVfnlJUVlZGYGBgvpbZ1dUVnTt3RseOHaGlpYUdO3ZI1cP169ehrKyMtWvXSr7vxo0bMmn/1YMHD+Du7g5DQ0NoaGigatWq8PT0RFpamkzeK1euoFWrVtDX14e2tjZcXV1x+fJljBs3Dq6urgCAuXPnQklJCeHh4Zg3bx6UlJRw5coVXLlyBUpKSrh79y4GDRoEIyMjGBoaYtiwYRAIBPj48SOGDh2K0qVLw9TUFKNGjUJSUpLU9yclJWHGjBmws7ODuro6LCwsMGXKFLx9+1byPV9iKF++PCIjIzFw4ECYmJhAQ0MDjo6OWLVqVY7TwxIREREREVHxxUFFIiIiIirREhMTcfjwYQwfPhwVKlTINa+JiYnU33PnzkXHjh2hp6cHLy8veHl5QUdHB+7u7liwYIFU3hUrVqB+/frYu3cvqlatmmPatGnTMHLkSLRo0QK7d+/Gpk2b4OjoiH79+mHChAk5xvbXX3+hfv36CA0NxaxZs7Bjxw40aNAAffr0kXq6skWLFvD19UW1atXQvHlz+Pr6okWLFjLb+/jxIxo1aoSVK1fi119/xZ49ezBu3DicP38eNWvWxOPHj6Xy+/n5wdXVFVWrVsXu3bsxePBgrFmzBiNGjMi1ToOCgtC/f3+oq6tj1apVMss3b96MpKQk9OnTR+ZJxfwo88mTJ/Hhwwds3boVPXr0kNp+48aNMXbsWMyYMQOvX79GWloahgwZAhcXF4wZMybXcgkEAiQlJcl8svPx8UGTJk1gZmaGdevWYe/evejSpQvWrVuHZs2aITU1VZJ3y5YtaN68OVJTU+Hp6YlNmzbBxsYGLVu2lHoHZK9eveDr6wsTExP07NkTvr6+cHR0lCzv3Lkz3r59iw0bNmDmzJk4cOAABg4ciObNm+Pdu3dYv349Jk+eDG9vb6lyJiQkoEGDBli7di1+++037NmzB2PHjsXx48clT/p+KzExEXXr1kVwcDDmzJmD7du3o3Hjxpg6dSratm3LKWaJiIiIiIhKGjERERERUQn28OFDMQCxj4/PD6139epVMQDxggULZJYtXbpUDEB8/fp1sVgsFgMQV69eXZyZmSnJIy/t0qVLYi0tLfGjR4/EiYmJUp8zZ86IlZSUxH5+fmKxWCzu37+/uFGjRmKxWCz+8OGD2NDQUNy+fXuxUCiUiuXOnTtiNTU1saWlpVR606ZNxf3795f8ffnyZTEAcVhYmFgsFosHDhwo1tXVFb948UJqvU+fPonr1asndnBwEItEIsm2AIj37dsnldfT01OsqqoqTkxMzLEeTU1NxWvWrBEvX75crKurK37//r1kWXJystjMzEw8ffp08YYNG8Ta2tqSZflR5qZNm4pVVVXFUVFROdaDQCAQ29nZiX/55Rfx9OnTxVpaWuKgoKAcy7Nr1y4xgFw/u3btkpTBwMBAvHPnTpn/77CwMHGZMmXEy5cvF4vFYnFgYKBYVVVVPHz4cKl9RiwWi48dOyYGIG7atKlUuqWlpXjOnDkyZWvbtq3UNpYsWSI3fcWKFWItLS1xRkaGWCzO2uf09PTEr169kvoegUAgrl+/vhiA+PLly2KxWCyeM2eOGIDYw8NDJt5r166JVVRUxPPmzcuxHomIiIiIiKj44ZOKRERERFSiZWRkAABKlSr1Q+vt2bMHNjY2mDFjhsyyKVOmSJ7Y+6Jdu3ZQUlKSypc9zdvbG8nJyXBycoKurq7Up23bthCLxTh37pzM9/n6+iIhIQEbN26Eqqqq1LJ69eph0KBBP1Q2oVCIQ4cOYeLEiahcubLUMj09Paxbtw4BAQG4ffu2JL1q1aro3bu3VN5WrVohPT0dISEhcr8nMTERMTExsLW1xciRI6GlpYWVK1dKln95SnHixImws7ODQCBAZGRkvpa5du3aMDMzy3G5lpYWdu7ciQsXLsDT0xOLFy+GnZ3dd7d75coVXL9+XebzLV9fX3z8+BGDBg2S+f+2trZGXFyc5P/74MGDUFNTg5eXl8x+1LlzZ7Ru3TpP5QUADw8PqW18eYpxyJAhUunVqlVDcnIyoqOjIRQK8eeff2LSpEkyUwRraWlh9erVMt+jrq4uN97GjRtj4MCB2Lt3b55jJiIiIiIioqJPRdEBEBEREREVJDs7O6ioqMDPzw9t27bN83pv3ryBo6MjlJXl34dXtWpVvHnzRvJ32bJlZfJkT3v79i0aNmyIpUuX5vi95ubmMmlv376Fvr5+jtO3Vq1aFX///XeO28wuPj4eKSkpqFGjhtzl1apVAwCp8lWqVEkmn76+PgDg8+fPcrfzZTpTW1tbaGlpYcqUKZg7dy4mTJgALS0tLFu2DGPHjkXp0qVha2srWadcuXL5VmZ5/y/Z1alTB+bm5nj37l2eB+8aNWoEFZXcT6fevn0LTU1N/PPPPznm0dHRkeS1s7ODhoaG3HxVq1bF/fv38xTbl3dNfqGpqQkAsLGxkUpXV1cHkPVey/j4eKSmpkpNo5r9++V9z5dtZ+fk5IQ9e/bkKV4iIiIiIiIqHvikIhERERGVaPr6+mjfvj02bNiA+Pj4XPMKBALJvytUqIDAwMAc8z579gzly5eX/C1vMCh7mrm5OaKiotCoUSO4uLhIferWrYvIyEhYWlrKbKdChQr49OmT1CDft168eJFrubIzNjaGhoYGAgIC5C5/+vQpAHy3fF8GXL88DZpdUFAQlJSUYG1tDQAYMWIEtLW1sXLlSmzevBkCgQATJ04EAFhaWkJFRUUyEJlfZc5pkO5bs2bNwqdPn2Bvb4/BgwcjMzM
2025-03-24 09:57:14 +08:00
},
"metadata": {},
"output_type": "display_data"
}
],
2025-03-24 15:19:11 +08:00
"execution_count": 20
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:27.534852Z",
"start_time": "2025-03-24T03:36:26.568836Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"#主成分分析( PCA)\n",
2025-03-24 15:19:11 +08:00
"data=pd.read_excel('北京市空气质量指数与气象数据.xlsx')\n",
2025-03-24 09:57:14 +08:00
"PCA_data=data.iloc[:,2:]#去除日期列\n",
"\n",
"# 计算KMO值\n",
"kmo_all, kmo_model = calculate_kmo(PCA_data)\n",
"print(f\"KMO值: {kmo_model.round(3)}\")\n",
"# 进行巴赫利特检验\n",
"chi_square_value, p_value = calculate_bartlett_sphericity(PCA_data)\n",
"print(f\"巴赫利特检验卡方值: {chi_square_value.round(3)}, p值: {p_value}\")\n",
"#判断\n",
"if kmo_model>0.7 and p_value<0.05:\n",
" print(\"数据适合进行主成分分析\",'\\n')\n",
"else:\n",
" print(\"数据不适合进行主成分分析\",'\\n')\n",
"\n",
"# 数据标准化\n",
"scaled_data = (PCA_data - PCA_data.mean()) / PCA_data.std()\n",
"scaled_data = scaled_data.dropna()#去除空值\n",
"\n",
"# 计算协方差矩阵\n",
"cov_matrix = np.cov(scaled_data, rowvar=False)\n",
"\n",
"# 计算特征值和特征向量\n",
"eigenvalues, eigenvectors = np.linalg.eigh(cov_matrix)\n",
"sorted_indices = np.argsort(eigenvalues)[::-1]\n",
"sorted_eigenvalues = eigenvalues[sorted_indices]\n",
"sorted_eigenvectors = eigenvectors[:, sorted_indices]\n",
"\n",
"# 绘制累计方差解释比例图\n",
"explained_variance_ratio = sorted_eigenvalues / np.sum(sorted_eigenvalues)\n",
"cumulative_explained_variance = np.cumsum(explained_variance_ratio)\n",
"print(\"累计方差解释比例:\", [f\"{cum * 100:.2f}%\" for cum in cumulative_explained_variance])\n",
"\n",
"plt.plot(range(1, len(cumulative_explained_variance) + 1), cumulative_explained_variance, marker='o')\n",
"plt.xlabel('主成分数量')\n",
"plt.ylabel('累计方差解释比例')\n",
"plt.title('PCA 累计方差解释比例')\n",
2025-03-24 15:19:11 +08:00
"plt.savefig('./images/PCA_cumulative_explained_variance.png', dpi=200, bbox_inches='tight')\n",
2025-03-24 09:57:14 +08:00
"plt.show()\n",
"\n",
"# 选择特征值大于1的作为主成分\n",
"mask = sorted_eigenvalues > 1\n",
"selected_eigenvectors = sorted_eigenvectors[:, mask]\n",
"\n",
"# 计算因子载荷矩阵\n",
"loadings = selected_eigenvectors * np.sqrt(sorted_eigenvalues[mask])\n",
"\n",
"# 使用Varimax旋转载荷矩阵\n",
"rotator = Rotator(method='varimax')\n",
"rotated_loadings = rotator.fit_transform(loadings)\n",
"\n",
"# 输出旋转后的成分矩阵\n",
"rotated_components_df = pd.DataFrame(rotated_loadings,\n",
" index=PCA_data.columns,\n",
" columns=[f'Factor{i+1}' for i in range(rotated_loadings.shape[1])])\n",
"rotated_components_df = rotated_components_df.round(3)\n",
"\n",
"# 输出排序后的载荷矩阵\n",
"rotated_components_df=sort_matrix_by_diag(rotated_components_df)\n",
"print(\"旋转后的载荷矩阵(排序后):\\n\", rotated_components_df)\n",
2025-03-24 15:19:11 +08:00
"plot_heatmap(rotated_components_df, 4, 8,save_path=\"./images/components_heatmap.png\")"
2025-03-24 09:57:14 +08:00
],
"id": "509d783a82bbdcb2",
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"KMO值: 0.762\n",
"巴赫利特检验卡方值: 90424.712, p值: 0.0\n",
"数据适合进行主成分分析 \n",
"\n",
"累计方差解释比例: ['31.41%', '54.60%', '66.53%', '73.02%', '78.89%', '84.04%', '88.27%', '91.46%', '93.59%', '95.70%', '97.14%', '98.29%', '98.91%', '99.26%', '99.55%', '99.79%', '99.96%', '100.00%', '100.00%']\n"
]
},
{
"data": {
"text/plain": [
"<Figure size 640x480 with 1 Axes>"
],
2025-03-24 15:19:11 +08:00
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAjQAAAHGCAYAAABjORGMAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAABvDklEQVR4nO3deVxUVf8H8M9szMAAA6goCiLuOyruGqVZmmUlrmVmz6/FzNQ000xzeyxzfVpMs8y0NC2tNM2tKNMKxH3BpRQUZFUEhmX2Ob8/kMlxAEGBYeDzfr3mJXPuufd+71xivp17FokQQoCIiIjIhUmdHQARERHRvWJCQ0RERC6PCQ0RERG5PCY0RERE5PKY0BAREZHLY0JDRERELo8JDREREbk8JjRERETk8pjQEFGFMBgMuHTpUrHblyxZgpUrV0Kn05X6mGazGVFRUWWOxWq1YsuWLTh48KBd+ccff4xPPvkEmZmZZT5mVfDPP//g1KlTyM3NLVX9tLS0Co6IyHmY0BBVsueeew4SiQQSiQQymQweHh5o0qQJZs2aBYvF4lBfp9Phk08+wf3334+6devC3d0djRo1wtNPP429e/fe8XwPPPAAzp07V6YY9Xo9srKyoNPpoNfri33pdDpkZWXBYDDY7S+EQO/evfHYY4/BZDIVeY73338fmzdvhru7e6njevzxx9GzZ0/Ex8eX6Xq0Wi2mTZuGsWPH2mLVarWYNWsW/vzzT2g0mmL3FUJg7Nixd5VIFdLpdDh58iQ+//xzTJs2DRkZGYiJicHff/+N+Ph4XL58GZcvX8alS5dw+vRp3DqBe3x8PJ555hn8/PPPDsddtWoVunbtWuxnfKsrV64gKCgIU6ZMwZ9//omzZ8/i/PnzOHfuHI4cOYKMjAxb3SlTpuCxxx676+slcga5swMgqolat26NRYsWAQAyMzPx119/4b333kN2djY++ugjW72TJ09ixIgRiI+Px4gRI/Dcc8/B29sbV69exY4dOzBgwABERERgw4YNRSYGf/zxB37//XcsWrQI69atK3V8P/74I0aMGFHq+jt27LD7ApRIJJg6dSpGjhyJd955B3PnznXYx9PTEx4eHqU+BwD83//9H3bv3o0vvvgC8+fPL7FuRkYG1q1bB29vb3h6eiIiIgIGgwE//vgjJBIJtm3bBgDo06cP1q9fj5ycHDRt2hQDBw60O87mzZvx1VdfYfz48TAajdDpdHBzc8OIESPwzz//QCaT2epaLBa0bt0a3333HQwGAx566CGcP38e165dg4eHBx544AG0bNkSu3fvxujRoyGVSu3uW15eHh566CHs27fPViaXy7Fx40aMGTPG4Rrd3NxQt25d+Pr63vGzO3DgAEwmEwYOHIjevXs7bN+9ezcGDBgAoKAlrLStPkRVBRMaIifw9fW1SwBGjx4NLy8vLF++HAsXLoSnpycSExPRv39/+Pj44NSpU2jRooXdMSZNmoTdu3dj5MiR+OyzzzBx4kSH8yxYsAAAsHHjRsyfPx8NGzYsVXyFicatLQXbtm3D4MGDcfvybxKJBHK545+SESNG4P3338eJEydgtVohldo3CCsUimLP/9tvv+H8+fNQqVRQKBSQSqWQSqXQ6XT4v//7P7Ro0QKbN2+GEAJWqxUmkwkGgwF16tRBREQEgIIWmI0bN0KtVtuOI5PJMG/ePMTFxeGhhx5C79698c0339gSlQEDBtglNNevX8frr78OnU6H0NBQAEBISAhOnjyJSZMmITc3F08++SQ++ugjtG3bFrNnz4ZarQYAKJVKzJw5EyqVCr/99hv27NmDn376CbNnz0bz5s1x9OhR9O7dG+fOnUNQUBD++OMP3H///Vi4cGGRn5NEIinysyqu/HY7duxAaGgo+vXrh8TERLz44ovw8vLC2rVrkZeXh9zcXAwePBjTp093uFdELkEQUaUaM2aM6NWrl0P5u+++KwCIa9euCSGEePrpp4W3t7dISEgo8XhpaWlFlh85ckQAEHPmzBFqtVpMmDCh1DHu2LFDABAymcz2kkqlDmUymUwAELt37y7yOFlZWcWeo02bNqJ///5Fbhs/frwAINzc3ISXl5fw9fUVbm5uQiqVCj8/P+Hj4yM0Go3w9PQUarVauLm5CYlEIh588EG74+Tn54sbN24Ig8FgK1uyZIkIDg62q2c0GkVGRoawWq22shs3bohevXqJ4cOHC7PZLFJSUkT9+vXFH3/8Yavz22+/CYVCIW7cuCGEEKJ169Zi1qxZtu3PP/+8eOutt8TChQtFt27dRFxcnPD09BTp6elCCCFeeOEF8fjjj4tLly6JoKAgMW3aNIfPIiUlRQAQP//8s8O26dOn211LamqqGDFihNBqtXb1srOzhVqtFmPGjLGVBQcHi6VLl9reJyUlCQDi8OHDYtKkSeL+++93OB9RVcYWGiInsFqttib93Nxc/Pnnn1iyZAl69uyJ2rVrIycnB99++y2mTJmCoKCgEo/l7+9fZPk777yDOnXqYNq0acjKysKnn36Kt99+G3Xq1Cl1nLc+dtixYweGDx/u8Cji9kddO3fuREZGBqRSKfLy8jBkyBCkpqZi7dq18Pb2hlKphFwuR0ZGBsxmMxYtWgSTyYScnBx06NABTz31FP73v//ho48+smt9eOSRR5CSkoITJ04UG6/RaLR7/91332H06NG294UtSWazGSqVClarFWaz2dbqlJOTA09PTwDA4cOHYTAYsHPnTmg0GpjNZhiNRjzwwAOoXbs2oqOj8eqrr+Kll17CF198gffffx+JiYl45513bOd74okn8NZbb2HUqFEACh7r/Oc//7Hdg48++gjdunVD69at0bNnT7z77rvFXtuhQ4eg0+lgNpuRl5eHRx991G57dHQ0hg8fjsTERNSpU8fu0eXKlSuRl5dne5+YmIgrV67gxx9/xOHDh+Hr64s5c+YUe24il+DsjIqophkzZowAYPeSSCTioYceEomJiUIIIY4fPy4AiB9//PGuzhEbGyskEolYsGCBEEKIK1euCLlcLmbOnFmq/b/99ltbXIUtMbf+XPi+MP7t27fb9n3iiSfsru3atWu2Fh+FQiG8vLyEj4+PkMlkQi6XCy8vL6FUKgUAMXbs2GJjuu+++8TAgQOL3KbX64ssT0tLEzExMeL8+fPiypUrIiUlRcydO1cEBQWJtLQ0kZKSIq5cuSLOnTsnDh06JMxms93+GRkZovDP5BNPPCG++OIL8fPPP4vg4GAxZswYMXDgQKHT6cSJEyfEhAkTxJYtW2z7Wq1WER0dLY4dOybeeust0alTJ3HlyhVx4sQJodfrRWRkpIiIiBAymUz07dtXyGQycf/994vZs2eLzz//XOTm5goh/m2hUavVwtPTU7i5uQkA4vTp02L69OkiMDBQzJ07V8jlctGsWTO7eyGEEJcuXRKenp5CpVLZWmhWrVol6tWrJ/773/+Kbt26iU6dOtnOwxYaclVsoSFygnbt2mHVqlWQSqVQq9UICgqy69hZONrp1g6nZfHuu+/Czc0NTzzxBK5evQqpVIpHHnkEH3/8MaZNmwZvb+8S93/yySdhMBigUCiwc+dOzJs3D82bN8fq1avh7u6On376CcuWLUPfvn0xdepUKJVK276ff/451q9fjy+//BITJ06Et7c3BgwYAIPBADc3N1u9tm3bIjAwEHv27LFds9lsLjamtLQ0/PHHHw79daxWK15++WWsXLnSrvyzzz5DfHy8rf+MXC6HTCbDsWPHkJubi6+//trWB8dsNsNgMOCHH37AwIEDcd999wGALd6srCyYTCbk5+fbWjo+++wz/Pnnn7BYLNi9ezfOnz+PDz74wHZ+g8GA7t272+6hxWJBcHAwAGDq1KlYunQp+vbti3Xr1kGtVmPatGn44YcfsGnTJtSpUwf/93//Z3c927ZtQ79+/QAUtDAVHvfq1atYtGgR5s2bh6lTp9p9xkBBf6T69evb+gAVxj5u3DjMmjULp0+fdmjZInJFTGiInMDb2xu9evUqdnvTpk0hl8tx+PBhh1E3dxIXF4fNmzfDYrGgXbt2Dts/+eQTTJs
2025-03-24 09:57:14 +08:00
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"旋转后的载荷矩阵(排序后):\n",
" Factor1 Factor3 Factor2 Factor4 Factor5\n",
2025-03-24 15:19:11 +08:00
"Tn -0.963 0.035 0.071 -0.079 0.035\n",
"T -0.958 -0.138 0.033 0.074 -0.028\n",
"Tx -0.954 -0.014 0.045 -0.052 -0.063\n",
"P 0.924 -0.029 0.071 -0.032 -0.000\n",
"Po 0.921 -0.029 0.073 -0.033 -0.000\n",
"Td -0.898 0.366 0.043 -0.013 0.076\n",
"O3 -0.637 -0.529 -0.030 0.239 -0.084\n",
"U -0.322 0.824 -0.008 -0.156 0.229\n",
"Ff -0.045 -0.772 -0.126 0.024 0.172\n",
"NO2 0.300 0.728 -0.290 0.110 -0.202\n",
"CO -0.101 0.695 -0.449 0.298 -0.007\n",
"VV 0.153 -0.667 0.531 -0.093 -0.175\n",
"AQI -0.017 0.038 -0.967 0.025 -0.029\n",
"PM10 0.037 -0.060 -0.933 -0.092 0.003\n",
"PM2.5 0.049 0.359 -0.879 0.149 -0.007\n",
"Pa 0.006 0.055 -0.147 -0.747 -0.130\n",
"SO2 -0.035 0.099 -0.208 0.694 -0.065\n",
"RRR -0.139 0.094 0.103 -0.077 0.819\n",
"tR 0.163 -0.120 -0.087 0.131 0.512\n"
2025-03-24 09:57:14 +08:00
]
},
{
"data": {
"text/plain": [
"<Figure size 400x800 with 2 Axes>"
],
2025-03-24 15:19:11 +08:00
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYwAAAMWCAYAAADiW7FGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd1hT1//A8XfYIwwZAirbrSiKe+OorXu2rlbrQmutdVTFvXHUXRVH3ba27lH7ax11K+JGZA9RNm5WQoDfH2g0JrR8lVSk5/U893nIuZ9788khyck55w5Jfn5+PoIgCILwD3TedwKCIAjCh0E0GIIgCEKRiAZDEARBKBLRYAiCIAhFIhoMQRAEoUhEgyEIgiAUiWgwBEEQhCIRDYYgCIJQJKLBEARBEIpENBjC/2TQoEE0a9ZM47rTp08jkUiIjY3Vag6xsbFIJBJOnz6t1ecRBEGVaDAEQRCEIhENhiAIglAkosEQtOrXX3+lcePGmJqaYmZmRrNmzfj111/V4vLz89m4cSNeXl4YGxtjbW1N3759iY6OxtPTk1mzZgHQqlUrXF1dAfD29kYikQAwa9YsKlSowO3bt2nbti0mJiY4OTmxbNkyAIKDg/n444+RSqW4u7uzatUqtRyio6MZOHAgdnZ2GBoa4uHhgb+/PwcOHFA+z8scBgwYQEBAAO3bt8fMzAypVErLli05fPhwcVehIJQYeu87AeHDk5eXR3p6ulp5VlaWyuPJkyezadMmRo4cybhx45DJZFy4cIEvvviCy5cvK7/M8/Pz+eKLL9i1axcDBgxgzJgxZGRksGfPHho0aICurq5yn35+foSEhDBkyBAWLFiAh4eHct2TJ09o2bIlffv2xcfHh7NnzzJ+/HhSU1NZv349AwYMYMiQIZw/f54xY8ZgYWHBwIEDAbh69Spt27alTJkyjB8/HicnJ4KCgpg4cSKOjo5qr/Xy5cvs3buX3r17s3r1auRyOYcPH6Zr167MnDlT2cAJQqmSLwj/g4EDB+YDf7vExMTknzp1Kt/ExCT/5s2b+c+fP1dZjh07li+RSPIDAwPz8/Pz87ds2ZIP5P/0009qzzdmzJh8IH/mzJnKspiYmHwg/6+//lKWzZw5Mx/IX7hwocr2H3/8cT6Qv2jRIpXyTp065bdv3z4/Pz8/PycnJ79ixYr5Xl5e+c+fP1eJi4qKyi9Tpkz+6x+Vli1b5gP5u3btUsvXz88vH8g/c+ZM0SpUED4gYkhK+J95eHhw7tw5teX1YZ4dO3aQmZmJp6cnZmZmKkuHDh3Iz8/njz/+AGD79u20bt2avn37qj3X4sWLMTc3L1JeEokEHx8ftVwBhg0bplJes2ZNoqOjAbh48SKRkZEsX74cqVSqEufm5sakSZPUnqt58+b069dPrXzSpElUrlyZ7du3FylnQfiQiCEp4X9mbm6u8dBahUKh/PvBgwc0adKERYsWFbqf8uXLK2M7dOigMcbAwIBKlSoVKS9LS0ssLS1VyoyNjbG0tKRMmTIq5YaGhsohtAcPHgCoDG+9rkaNGmpl1atX1xgrkUioXbs29+/fL1LOgvAhEQ2GoBXly5fnzJkzNG3aVGXCGEAul3Pw4EGaNGkCgKOjI3fu3NG4n5ycHCIiIor0nCYmJmplEolEYzkUzJ28fH6AO3fuaGwI7969W6Syl8LDw6lXr16RchaED4kYkhK0YsCAAcTExLB69WqV8vz8fMaNG8fQoUNJSEgACk4GPHnyJD///LPafnx9fXn27Nk75fJmg/WmJk2aUKlSJb799lu1yfzY2FgWL16sts25c+fYuXOnWvmqVau4deuWcjJdEEoT0cMQtKJNmzZ8++23jBkzhsDAQNq3b09WVha7du3i/Pnz7NixgwoVKgAFjcuff/5J//79+f3332nbti1ZWVns2bOHsLAwHBwcVPatr68PwJkzZ7C1tdU4ZPS/0NXV5aeffqJt27bUrFmTr776CmdnZ4KCgli9ejVVq1bl4cOHKttUr16dkSNH8ttvv/HRRx+Rl5fHwYMHOXr0KDNnzqR58+bvlJMglESiwRC0Zvny5Xh5ebF69Wr27duHkZERXl5enDhxglatWinjJBIJ27dvp1mzZqxbt45ff/0VCwsLunbtyq5du2jfvr3KfsuVK0f79u2ZP38+T548Yfny5e+ca7169bh27RozZ85kyZIlpKenU7NmTTZs2IChoSHdu3dXia9Tpw579+5l2rRpjB07FrlcTp06ddi/f79arCCUFpL8lwO5giAUSatWrahQoYLGISlBKM3EHIYgCIJQJKLBEARBEIpENBiC8D86ffq0GI4StOblNdL69Onzt3HHjh2jRo0aGBkZ0aBBA65cuaKyftGiRZQvXx5TU1N69OhBcnLyO+cmGgxBEIQSQqFQsG7dOkaNGvW3cXfv3qVHjx7079+fy5cv06pVK9q3b098fDwAmzZtYuHChSxbtoxTp06RlZVFly5deNcpazHpLQiCUMLMmjWL0NBQdu/erXG9j48PycnJHDx4UFnm7e1N/fr1Wbx4MVWqVGHkyJF8++23ADx//pzy5cuze/fuQq+qUBSihyEIgvCBOXfuHD179lQpGzBgAMePHyclJYXw8HCV9WZmZnTt2pXjx4+/0/OKBkMQBOEDExcXp7wvzEsVK1YkOjqauLg49PT01C7L/3L9uxAn7gmCIGiRTCZDJpOplBkaGmJoaPjW+8zKylK7oGaZMmXIyMggKytL7SKcr69/F6LBEASh1PlNv8r7TkEpcGpfZs+erVL2rjfZMjY2Jjs7W6XsyZMnmJiYaFz3+vp38Z9sMFr1uvS+U3grp/c2Ji4i5H2n8VacKlUjMfTm+07jrThU9SQo8t0PSXwfPCrafdDvmdLA19eXcePGqZS9S+8CwNnZmejoaLy8vJRlkZGRuLm54ezsTHp6Oqmpqdja2qqtfxdiDkMQBEGLDA0NMTc3V1netcFo0aIFR44cUSnbvXs3bdu2xdbWlqpVq6qsz8zM5MiRI7Rt2/adnvc/2cMQBKF0k+j//SXtPzRZWVl4eXnRo0cP5s2bx+jRo/Hy8qJ27dq0bduWPXv2EBAQwObNmwGYMGECkyZNwtLSEkdHR+bOnUvFihXf6ZBaEA2GIAjCB0EulyOXy4GCy+vv27ePCRMmMHXqVGrVqsWff/6pvGXAkCFDSE1NZdSoUTx79oyPPvqII0eOoKPzboNK/8kT98Qcxr9PzGG8H//VOYxjJlWLMZN30yEz9H2nUGxED0MQhFJHR690DUmVFGLSWxAEQSgS0cMQBKHUkeiL38LaIGpVEARBKBLRYAiCIAhFIoak3iNjIx28m9hw7FTK+05FEEoVMemtHSWuwXBxceHevXsqZS1btuT06dPvJyGgYR1LRnzhTHk7I6LuZbLyxxhCI9M1xhoa6DC0nyPeTWzQ0YEzlx/hv/0eMnmeWuzgPo7Ur23Jn2dTUSiK9+jmx4+fsHKtP9eu38DU1JRuXTrRp1cPjbEBgVfZuGUbCYlJuLm48PWIYVStUlljbFR0DOcvXWZg/77Fmu/rHj15wvJ1mwi8cRupqQk9On1Cv55dNcZevnoD/607iU9Kxt3FiTHDB1OtckUAxkydza07d9W26d6xPWOGD9Za/gDXAy+xffM6khMTcHZ1Y8iIb6lUpbrG2NxcBbu2buCvE7+Tq1DQqGlLvvT5BmPjguv+PHn8iK0bf+DmtQAMDI1o074jvfoMRFdXt9jz/pDfN4L2lbghqZMnTxISEkK3bt0YNmwYISEhbN++/b3l41zBmLnfVeHEuTS+mhLEzeCnLJlWDRsrA43x3w51pU4NC/xWRzJ3RQTVK0kZO9xVLa6SqyndP7Zn6froYm8sAOb4LSIjI4NF8+cwymcoe/cf5NDRY2px9+LuM3vBIlq3asGqpYuo7VED3xmzSUt7qBYrz8lh4dLl5OTkFHu+r5u5aDnpGZksnTuN0cO+5JeDRzjw2x9qcbFxD5i+cCltWzZj7eJ5eNasznezFpD68BEAU74dxbY1y5TLlLFfo6+vT7+e3bSa//24WJbMn07zVu1YsHQdNTz
2025-03-24 09:57:14 +08:00
},
"metadata": {},
"output_type": "display_data"
}
],
2025-03-24 15:19:11 +08:00
"execution_count": 21
2025-03-24 09:57:14 +08:00
},
{
"metadata": {},
"cell_type": "markdown",
"source": [
"## 题目3\n",
"令2022年11月1日至2023年9月30日的空气质量数据为训练集, 剩余数据为测试集。基于训练集, 尝试使用两种不同的方法构建空气质量指数预测模型, 并在测试集上测试。比较所选模型的预测效果。"
],
"id": "3f89fa62a897a3e3"
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:27.837635Z",
"start_time": "2025-03-24T03:36:27.550331Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
2025-03-24 15:19:11 +08:00
"outputs": [],
"execution_count": 22,
2025-03-24 09:57:14 +08:00
"source": [
"#重新读取数据\n",
"data=pd.read_excel('北京市空气质量指数与气象数据.xlsx')\n",
"# 合并 date 和 hour 为新的 data_hour 列\n",
"data['data_hour'] = pd.to_datetime(data['date']) + pd.to_timedelta(data['hour'], unit='h')\n",
"# 设置 data_hour 为索引列\n",
"data = data[['data_hour', 'AQI']].set_index('data_hour') # 仅保留时间和AQI"
],
2025-03-24 15:19:11 +08:00
"id": "d1bdac1e4e1562f2"
2025-03-24 09:57:14 +08:00
},
{
"metadata": {},
"cell_type": "markdown",
2025-03-24 15:19:11 +08:00
"source": "### (1)SARIMA模型\n",
"id": "1fc53937767d55fd"
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T07:17:24.611373Z",
"start_time": "2025-03-24T07:17:22.170632Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"\"\"\"\n",
"该模型在假设不知道测试集其他指标的情况下, 仅使用AQI历史数据预测未来AQI\n",
"\"\"\"\n",
"\n",
2025-03-24 15:19:11 +08:00
"data=pd.read_excel('北京市空气质量指数与气象数据.xlsx')\n",
"data['date'] = pd.to_datetime(data['date'])\n",
"data = data.set_index('date')\n",
"# 确保数据按时间排序\n",
"ts = data['AQI'].resample('D').mean() # 按天重采样\n",
"ts = ts.bfill() # 处理缺失值\n",
"\n",
"# SARIMA模型训练\n",
"model = SARIMAX(ts,\n",
" order=(1, 1, 1),\n",
" seasonal_order=(1, 1, 1, 7))\n",
"results = model.fit(disp=False)\n",
2025-03-24 09:57:14 +08:00
"\n",
2025-03-24 15:19:11 +08:00
"# 生成预测结果(包含样本内拟合)\n",
"forecast = results.get_prediction(start=pd.to_datetime('2023-01-01'), dynamic=False)\n",
"forecast_mean = forecast.predicted_mean\n",
"confidence_int = forecast.conf_int()\n",
"\n",
"# 可视化对比\n",
"fig, ax = plt.subplots(figsize=(15,6))\n",
"ts.plot(ax=ax, label='实际观测值', alpha=0.7)\n",
"forecast_mean.plot(ax=ax, style='--', label='模型预测值', color='red')\n",
"ax.fill_between(confidence_int.index,\n",
" confidence_int.iloc[:, 0],\n",
" confidence_int.iloc[:, 1], color='pink', alpha=0.3, label='置信区间')\n",
"\n",
"ax.set_title('空气质量实际观测 vs 模型预测', fontsize=14)\n",
"ax.set_xlabel('日期', fontsize=12)\n",
"ax.set_ylabel('AQI数值', fontsize=12)\n",
"plt.legend()\n",
"plt.show()\n",
2025-03-24 09:57:14 +08:00
"\n",
2025-03-24 15:19:11 +08:00
"# 计算拟合度指标\n",
"from sklearn.metrics import r2_score, mean_squared_error\n",
"\n",
"# 获取实际值和预测值的交集\n",
"y_actual = ts[forecast_mean.index]\n",
"y_pred = forecast_mean\n",
"\n",
"# 计算评估指标\n",
"r2 = r2_score(y_actual, y_pred)\n",
"rmse = np.sqrt(mean_squared_error(y_actual, y_pred))\n",
"mae = mean_absolute_error(y_actual, y_pred)\n",
"\n",
"print(f'''模型拟合度评估:\n",
"R² = {r2:.3f}\n",
"RMSE = {rmse:.3f}\n",
"MAE = {mae:.3f}''')\n",
"\n"
2025-03-24 09:57:14 +08:00
],
"id": "24996a0c06820cdc",
"outputs": [
2025-03-24 15:19:11 +08:00
{
"name": "stdout",
"output_type": "stream",
"text": [
"<class 'pandas.core.indexes.datetimes.DatetimeIndex'>\n"
]
},
2025-03-24 09:57:14 +08:00
{
"data": {
"text/plain": [
2025-03-24 15:19:11 +08:00
"<Figure size 1500x600 with 1 Axes>"
],
"image/png": "iVBORw0KGgoAAAANSUhEUgAABN0AAAI3CAYAAACvXS1VAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd5wTdfoH8M9M6nbKUqVL8+wCinh2EMTyw0MRBCmCgmDhQLCcCuqd2AXLIVgQBRXRE6WogKAiolKlSS/S+/a0mfn+/piSursJrGzJ5/166bLJJJlMdrPJJ8/zfSQhhAARERERERERERGVGbm8d4CIiIiIiIiIiKiqYehGRERERERERERUxhi6ERERERERERERlTGGbkRERERERERERGWMoRsREREREREREVEZY+hGRERERERERERUxhi6ERERERERERERlTGGbkRERERERERERGWMoRsRERGdkmXLluH//u//8Oabb57U5ZcvX44vv/wSmqaV8Z4VT9M0rFy5EoFA4LTdZqT9+/fjmWeewejRo8ttH8rasWPHEr7M3r178cgjj2D+/Pl/wR4l5vjx4/jll1+wadOmqPN+/vlnLFu2DABw+PBhLFq0CPv37y/T2//oo4/w0UcfQVXVMr3eSCfzOHk8HowePRozZsz4C/aIiIioapKEEKK8d4KIiIgqr40bN+Lss89G+/btrVAiET179sSMGTPQpUsXzJs3D5IkndR+HDlyBKqqQpZleDweFBUV4dixYzh69CgOHjyIPXv24M8//8TmzZuxceNGFBYWYubMmbj11ltjXt+AAQNw7Ngx2O12yHLszymFEAgEAigqKsLIkSNx/fXXx72/+fn5qFOnDrxeL/bs2YMzzjijxO3Xrl2Lt99+G5mZmbDb7dZtCyEwbtw4TJw4EVlZWXA6ndYxFELA5/NB0zT07ds37n3Ly8vDxx9/DLfbDZvNFnZ/VVVFUVERhgwZEnZc1q5diw4dOqB///549tlnMXv2bNjtdqSlpVnbmI/Jvffea522a9cuNG3aFAsWLEDHjh3j3se/wsqVK9G2bVu8/PLLGDFiRNh5N954IwoLC7F48WL88ssvuPTSS/HTTz/hsssui7qeTZs24amnnkJqaiocDkfYeUII+P1+tG3bFgMHDoTH40FWVhZkWcZFF10EVVXx+++/W9sfPHgQ3377rfW4S5JkPQYdOnTAjBkzIEkSMjIyYLPZIEmS9Th5vV6cOHECI0aMQFZWFgA97L3gggtwzTXXYPz48fjtt99QUFCAzMxM6za9Xi/279+PwYMHw+VyAdCDaofDgUmTJmHQoEFldsyJiIiqMnt57wARERFVfDk5Ofjqq6+QlpYGt9sNu91uhVGSJCE7Oxtr1qzB999/b73hDwQC8Pv9KCoqwjnnnINzzz036no3bdqE//3vf3C5XLjpppswffr0YvdB0zT4fD4cOXIE999/PzIyMsLO/89//oMJEyaEnWaz2ZCamgpFURAIBPD3v/8dzZo1w4UXXgin0xmzoil037KysvDYY4+VeGxyc3PRpUuXYoOIp59+Gg6HAw6HA7IshwVVderUwa5duzBy5Ei0b98eAKCqqnXcrrzySiuIys3NxRtvvBF1/f3794eqqhg2bFix+9i4ceOo0G3nzp14++234Xa7Icuy9bjVrl0b//jHPzBkyBDY7faw0A0AFEUBAAwdOjTs9GbNmmHQoEF47bXXcOjQIRw4cABLly6N2pcGDRqge/fueOmll9C1a1f87W9/AwDY7eX/sjQ9PR0AUL9+fQB6Fd7evXuRmZkJVVUhSRK2bduGo0ePAtAr47Zs2QKv1wubzYazzz4bAOD3+3HgwAHIsozFixfj/PPPR4MGDbBz505s3LgRV199NZo2bYp9+/ahefPmeOmllzBy5EhkZmZax9e0a9cu9O/fP2pfL7roIvz666/473//i927d8Nms0EIYYVjmqZZFXN33323FbrVqlULw4cPx9ixY7Fjxw60bNmy2N+7wYMHY9SoUejQoQNuueUWuFyuCvE4ERERVRb8q0lERESlOnjwIPr161fqdldffXXM01988cWo0E0IgQcffBCKosDpdEZVFkVSVdUKJO66666o0M0MA3bs2IEaNWogJSUFTqcTgB4Qffnll1i8eHGp98HkdDrxyy+/lFrVU1or4AsvvIDCwsISt5kxY0bMtj2bzWaFbqmpqQCA119/Hffddx8++ugj9O7dGz179oTNZsO8efOQmZlpVaeFVjtFBjkAsG/fPowbNy4sVFNVFbfccgt69eoFAJg4cWLU/R8yZAjef//9qOtLT0/H+PHjccUVV1hVXL1798Zbb71lbfPss8/i888/h9PpxIsvvoiLLroIF1xwQYnH5nSKDJQWLlyIAQMGhJ3WokUL698333yz9e/evXtj2rRpAIDzzjsP33//Pfbu3YuGDRti+PDh6N+/P8aPH49//vOfWLRoEQBY4V2dOnUAAC6XK+qxuuCCC7B3717rcf31119x00034b333oPdbseqVauQnp4Op9OJ22+/HV9//TWOHTsGh8MBv9+PvLw8VK9e3bo+h8OBxx57DJdddhkyMjIwceJE/P3vf8fXX39tbTN9+nQ88MADcLlcePvtt1GzZk0AsH6fiIiIKD4M3YiIiChu48ePx8CBAwEAPp8Pb7zxBu6//37rzXibNm3Qrl07K2j5888/reqfSP/5z38wf/58q8onHkII5ObmhrXCmczwqFq1alZLoxlgFBYWIiUlBYqiWGGUoijw+XzIz89Hdna2VeVkCgQCuPHGG60gpThHjx5FrVq1Styme/fu+Oyzz2Kep2laVPvqwoUL0alTp7DTzDY/0//+9z9kZ2fj2muvhSRJCbW2AkD79u0RCASsoEkIAbfbjSZNmkS1+Mbax+L84x//sPbXbrcjNzcXCxYswOWXX261vppVV/EaMmQI3n33XezcuRMNGjSIOn/QoEH44IMPsHv3btSrVw9ff/01HnvsMWzZsgVnnXUW/v3vf6NLly7FXv+2bdswa9Ys+Hw+AMB3332HQ4cOoWnTpti7dy/S09Nxzz33oLCwEHPmzMH69etx7rnnYvny5WjWrJlVXZaoyJAvVmu12+3GsWPHcM899+Ctt97CiBEjMHLkSJx//vkAgBo1aljb7ty5E1dffbW1L06nE9nZ2TFv+8orrwSgP05mSDtz5ky0adMGLpfL2pdq1aolfL+IiIhIx9CNiIiI4uZyuZCeng4hBO677z5MnToVK1euxJdffglJkiBJEux2uxVgmdVZkT788EM8+eST6NKlC8466yy8++67UaGSSdM0KIqCvLw8DB48uNgQwAwJQkOISMUFI4sXL8ZVV10VdlphYSE+/fRTzJkzB5IkWaGTuWaWpmkQQsBcHtcMbEry559/Ytq0aRgxYgTcbjcAoG3btti+fTt+/PFHK0gpjcfjwTfffIM777wTPp8PhYWFSE9Pj2oFNfn9fhQUFIQdm8jA5+jRo/D7/WGVXKb09HSceeaZYWuNRVq9ejXq169vVW2Z9uzZgwEDBmDJkiXWaZIkFbuvsfzzn//E5MmT8cYbb+C5556L2u/p06ejZ8+eqFevHjZu3Ihu3bqhWbNm6NevH7755ht07doVP//8s9XCG6mgoAB//PEHvF4vAL0K0GazoWbNmhg3blzYkJDQYKxdu3YAgPnz54eFpHPmzIHT6UROTg4AvVX5+++/x7Zt2wAACxYsQGFhYdTPXHHM9f7atWuHxo0b48knn7TOu/fee5GTkwO3241Vq1bh3HPPxcCBAxEIBJCbm4t77703LHDcsmUL7HY7mjVrFnYbhYWFGDhwID788MOw09lOSkREdPL4V5SIiIhKlZ2djTFjxqBt27bweDzo168fZs6ciXvvvRevvfZascMPGjRogK1bt4ZV27z11lsYOnQozjzzTEyfPh0jRozABx98UGzrmlmZpqoq+vTpg5SUlGK3AxAVGgB6e2thYSHGjBljBWaKosDr9aKoqAhNmzYFoAdna9euRUpKCqZMmWKtYRcrEBRCQFEU+P1+eDweKIq
2025-03-24 09:57:14 +08:00
},
"metadata": {},
2025-03-24 15:19:11 +08:00
"output_type": "display_data"
},
{
"name": "stdout",
"output_type": "stream",
"text": [
"模型拟合度评估:\n",
"R² = 0.354\n",
"RMSE = 39.860\n",
"MAE = 25.894\n"
]
2025-03-24 09:57:14 +08:00
}
],
2025-03-24 15:19:11 +08:00
"execution_count": 12
2025-03-24 09:57:14 +08:00
},
{
"metadata": {},
"cell_type": "markdown",
"source": "### (2)XGBOOST模型",
"id": "ebe88094b6c13e0c"
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:27.945366Z",
"start_time": "2025-03-24T03:36:27.921468Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"\"\"\"\n",
"该模型在假设不考虑测试集其他指标的情况下, 仅使用AQI数据对未来AQI进行<单步预测>, 即每次预测都是根据之前时间点的真实AQI值进行的。\n",
"整体运行时间约为20s, 请耐心等待。\n",
"\"\"\"\n",
"#特征工程\n",
"data_processed = data.copy()\n",
"\n",
"#时间分解特征\n",
"# 基础特征\n",
"data_processed['hour'] = data_processed.index.hour\n",
"data_processed['day_of_week'] = data_processed.index.dayofweek\n",
"data_processed['month'] = data_processed.index.month\n",
"\n",
"# 周期性编码\n",
"data_processed['hour_sin'] = np.sin(2 * np.pi * data_processed['hour'] / 24)\n",
"data_processed['hour_cos'] = np.cos(2 * np.pi * data_processed['hour'] / 24)\n",
"data_processed['week_sin'] = np.sin(2 * np.pi * data_processed['day_of_week'] / 7)\n",
"data_processed['week_cos'] = np.cos(2 * np.pi * data_processed['day_of_week'] / 7)\n",
"\n",
"#滞后特征\n",
"# 生成3小时粒度的滞后特征( 最多7天) \n",
"lags = [i for i in range(1, 7 * 8 + 1)] # 7天*每天8个时间点( 3小时间隔) \n",
"for lag in lags:\n",
" data_processed[f'AQI_lag_{lag}'] = data_processed['AQI'].shift(lag)\n",
"\n",
"# 划分数据集\n",
"train_data = data_processed.loc['2022-11-01':'2023-09-30']\n",
"test_data = data_processed.loc['2023-10-01':]\n",
"\n",
"# 特征选择\n",
"features = [col for col in train_data.columns if col != 'AQI']\n",
"X_train, y_train = train_data[features], train_data['AQI']\n",
"X_test, y_test = test_data[features], test_data['AQI']"
],
"id": "66f104e110aba36",
"outputs": [],
2025-03-24 15:19:11 +08:00
"execution_count": 24
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:57.925705Z",
"start_time": "2025-03-24T03:36:27.984878Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"#随机搜索法参数调优(这里耗时较长,请耐心等待)\n",
"param_dist = {\n",
" 'n_estimators': [100, 200, 300],\n",
" 'max_depth': randint(5, 10),\n",
" 'learning_rate': uniform(0.01, 0.2),\n",
" 'subsample': uniform(0.7, 0.3),\n",
" 'colsample_bytree': uniform(0.7, 0.3),\n",
" 'gamma': uniform(0, 0.3)\n",
"}\n",
"\n",
"search = RandomizedSearchCV(\n",
" XGBRegressor(n_jobs=-1, random_state=42),\n",
" param_distributions=param_dist,\n",
" n_iter=10,\n",
" cv=3,\n",
" scoring='neg_mean_absolute_error',\n",
" verbose=1\n",
")\n",
"search.fit(X_train, y_train)"
],
"id": "199aa487e826c1ac",
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"Fitting 3 folds for each of 10 candidates, totalling 30 fits\n"
]
},
{
"data": {
"text/plain": [
"RandomizedSearchCV(cv=3,\n",
" estimator=XGBRegressor(base_score=None, booster=None,\n",
" callbacks=None,\n",
" colsample_bylevel=None,\n",
" colsample_bynode=None,\n",
" colsample_bytree=None, device=None,\n",
" early_stopping_rounds=None,\n",
" enable_categorical=False,\n",
" eval_metric=None, feature_types=None,\n",
2025-03-24 15:19:11 +08:00
" feature_weights=None, gamma=None,\n",
" grow_policy=None,\n",
2025-03-24 09:57:14 +08:00
" importance_type=None,\n",
2025-03-24 15:19:11 +08:00
" interaction_constraint...\n",
" 'gamma': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x142527890>,\n",
" 'learning_rate': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x144aa2de0>,\n",
" 'max_depth': <scipy.stats._distn_infrastructure.rv_discrete_frozen object at 0x144a6e690>,\n",
2025-03-24 09:57:14 +08:00
" 'n_estimators': [100, 200, 300],\n",
2025-03-24 15:19:11 +08:00
" 'subsample': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x144aa1d90>},\n",
2025-03-24 09:57:14 +08:00
" scoring='neg_mean_absolute_error', verbose=1)"
],
"text/html": [
2025-03-24 15:19:11 +08:00
"<style>#sk-container-id-2 {\n",
2025-03-24 09:57:14 +08:00
" /* Definition of color scheme common for light and dark mode */\n",
" --sklearn-color-text: black;\n",
" --sklearn-color-line: gray;\n",
" /* Definition of color scheme for unfitted estimators */\n",
" --sklearn-color-unfitted-level-0: #fff5e6;\n",
" --sklearn-color-unfitted-level-1: #f6e4d2;\n",
" --sklearn-color-unfitted-level-2: #ffe0b3;\n",
" --sklearn-color-unfitted-level-3: chocolate;\n",
" /* Definition of color scheme for fitted estimators */\n",
" --sklearn-color-fitted-level-0: #f0f8ff;\n",
" --sklearn-color-fitted-level-1: #d4ebff;\n",
" --sklearn-color-fitted-level-2: #b3dbfd;\n",
" --sklearn-color-fitted-level-3: cornflowerblue;\n",
"\n",
" /* Specific color for light theme */\n",
" --sklearn-color-text-on-default-background: var(--sg-text-color, var(--theme-code-foreground, var(--jp-content-font-color1, black)));\n",
" --sklearn-color-background: var(--sg-background-color, var(--theme-background, var(--jp-layout-color0, white)));\n",
" --sklearn-color-border-box: var(--sg-text-color, var(--theme-code-foreground, var(--jp-content-font-color1, black)));\n",
" --sklearn-color-icon: #696969;\n",
"\n",
" @media (prefers-color-scheme: dark) {\n",
" /* Redefinition of color scheme for dark theme */\n",
" --sklearn-color-text-on-default-background: var(--sg-text-color, var(--theme-code-foreground, var(--jp-content-font-color1, white)));\n",
" --sklearn-color-background: var(--sg-background-color, var(--theme-background, var(--jp-layout-color0, #111)));\n",
" --sklearn-color-border-box: var(--sg-text-color, var(--theme-code-foreground, var(--jp-content-font-color1, white)));\n",
" --sklearn-color-icon: #878787;\n",
" }\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 {\n",
2025-03-24 09:57:14 +08:00
" color: var(--sklearn-color-text);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 pre {\n",
2025-03-24 09:57:14 +08:00
" padding: 0;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 input.sk-hidden--visually {\n",
2025-03-24 09:57:14 +08:00
" border: 0;\n",
" clip: rect(1px 1px 1px 1px);\n",
" clip: rect(1px, 1px, 1px, 1px);\n",
" height: 1px;\n",
" margin: -1px;\n",
" overflow: hidden;\n",
" padding: 0;\n",
" position: absolute;\n",
" width: 1px;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-dashed-wrapped {\n",
2025-03-24 09:57:14 +08:00
" border: 1px dashed var(--sklearn-color-line);\n",
" margin: 0 0.4em 0.5em 0.4em;\n",
" box-sizing: border-box;\n",
" padding-bottom: 0.4em;\n",
" background-color: var(--sklearn-color-background);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-container {\n",
2025-03-24 09:57:14 +08:00
" /* jupyter's `normalize.less` sets `[hidden] { display: none; }`\n",
" but bootstrap.min.css set `[hidden] { display: none !important; }`\n",
" so we also need the `!important` here to be able to override the\n",
" default hidden behavior on the sphinx rendered scikit-learn.org.\n",
" See: https://github.com/scikit-learn/scikit-learn/issues/21755 */\n",
" display: inline-block !important;\n",
" position: relative;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-text-repr-fallback {\n",
2025-03-24 09:57:14 +08:00
" display: none;\n",
"}\n",
"\n",
"div.sk-parallel-item,\n",
"div.sk-serial,\n",
"div.sk-item {\n",
" /* draw centered vertical line to link estimators */\n",
" background-image: linear-gradient(var(--sklearn-color-text-on-default-background), var(--sklearn-color-text-on-default-background));\n",
" background-size: 2px 100%;\n",
" background-repeat: no-repeat;\n",
" background-position: center center;\n",
"}\n",
"\n",
"/* Parallel-specific style estimator block */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-parallel-item::after {\n",
2025-03-24 09:57:14 +08:00
" content: \"\";\n",
" width: 100%;\n",
" border-bottom: 2px solid var(--sklearn-color-text-on-default-background);\n",
" flex-grow: 1;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-parallel {\n",
2025-03-24 09:57:14 +08:00
" display: flex;\n",
" align-items: stretch;\n",
" justify-content: center;\n",
" background-color: var(--sklearn-color-background);\n",
" position: relative;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-parallel-item {\n",
2025-03-24 09:57:14 +08:00
" display: flex;\n",
" flex-direction: column;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-parallel-item:first-child::after {\n",
2025-03-24 09:57:14 +08:00
" align-self: flex-end;\n",
" width: 50%;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-parallel-item:last-child::after {\n",
2025-03-24 09:57:14 +08:00
" align-self: flex-start;\n",
" width: 50%;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-parallel-item:only-child::after {\n",
2025-03-24 09:57:14 +08:00
" width: 0;\n",
"}\n",
"\n",
"/* Serial-specific style estimator block */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-serial {\n",
2025-03-24 09:57:14 +08:00
" display: flex;\n",
" flex-direction: column;\n",
" align-items: center;\n",
" background-color: var(--sklearn-color-background);\n",
" padding-right: 1em;\n",
" padding-left: 1em;\n",
"}\n",
"\n",
"\n",
"/* Toggleable style: style used for estimator/Pipeline/ColumnTransformer box that is\n",
"clickable and can be expanded/collapsed.\n",
"- Pipeline and ColumnTransformer use this feature and define the default style\n",
"- Estimators will overwrite some part of the style using the `sk-estimator` class\n",
"*/\n",
"\n",
"/* Pipeline and ColumnTransformer style (default) */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-toggleable {\n",
2025-03-24 09:57:14 +08:00
" /* Default theme specific background. It is overwritten whether we have a\n",
" specific estimator or a Pipeline/ColumnTransformer */\n",
" background-color: var(--sklearn-color-background);\n",
"}\n",
"\n",
"/* Toggleable label */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 label.sk-toggleable__label {\n",
2025-03-24 09:57:14 +08:00
" cursor: pointer;\n",
" display: block;\n",
" width: 100%;\n",
" margin-bottom: 0;\n",
" padding: 0.5em;\n",
" box-sizing: border-box;\n",
" text-align: center;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 label.sk-toggleable__label-arrow:before {\n",
2025-03-24 09:57:14 +08:00
" /* Arrow on the left of the label */\n",
" content: \"▸\";\n",
" float: left;\n",
" margin-right: 0.25em;\n",
" color: var(--sklearn-color-icon);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 label.sk-toggleable__label-arrow:hover:before {\n",
2025-03-24 09:57:14 +08:00
" color: var(--sklearn-color-text);\n",
"}\n",
"\n",
"/* Toggleable content - dropdown */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-toggleable__content {\n",
2025-03-24 09:57:14 +08:00
" max-height: 0;\n",
" max-width: 0;\n",
" overflow: hidden;\n",
" text-align: left;\n",
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-0);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-toggleable__content.fitted {\n",
2025-03-24 09:57:14 +08:00
" /* fitted */\n",
" background-color: var(--sklearn-color-fitted-level-0);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-toggleable__content pre {\n",
2025-03-24 09:57:14 +08:00
" margin: 0.2em;\n",
" border-radius: 0.25em;\n",
" color: var(--sklearn-color-text);\n",
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-0);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-toggleable__content.fitted pre {\n",
2025-03-24 09:57:14 +08:00
" /* unfitted */\n",
" background-color: var(--sklearn-color-fitted-level-0);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 input.sk-toggleable__control:checked~div.sk-toggleable__content {\n",
2025-03-24 09:57:14 +08:00
" /* Expand drop-down */\n",
" max-height: 200px;\n",
" max-width: 100%;\n",
" overflow: auto;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 input.sk-toggleable__control:checked~label.sk-toggleable__label-arrow:before {\n",
2025-03-24 09:57:14 +08:00
" content: \"▾\";\n",
"}\n",
"\n",
"/* Pipeline/ColumnTransformer-specific style */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label input.sk-toggleable__control:checked~label.sk-toggleable__label {\n",
2025-03-24 09:57:14 +08:00
" color: var(--sklearn-color-text);\n",
" background-color: var(--sklearn-color-unfitted-level-2);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label.fitted input.sk-toggleable__control:checked~label.sk-toggleable__label {\n",
2025-03-24 09:57:14 +08:00
" background-color: var(--sklearn-color-fitted-level-2);\n",
"}\n",
"\n",
"/* Estimator-specific style */\n",
"\n",
"/* Colorize estimator box */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-estimator input.sk-toggleable__control:checked~label.sk-toggleable__label {\n",
2025-03-24 09:57:14 +08:00
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-2);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-estimator.fitted input.sk-toggleable__control:checked~label.sk-toggleable__label {\n",
2025-03-24 09:57:14 +08:00
" /* fitted */\n",
" background-color: var(--sklearn-color-fitted-level-2);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label label.sk-toggleable__label,\n",
"#sk-container-id-2 div.sk-label label {\n",
2025-03-24 09:57:14 +08:00
" /* The background is the default theme color */\n",
" color: var(--sklearn-color-text-on-default-background);\n",
"}\n",
"\n",
"/* On hover, darken the color of the background */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label:hover label.sk-toggleable__label {\n",
2025-03-24 09:57:14 +08:00
" color: var(--sklearn-color-text);\n",
" background-color: var(--sklearn-color-unfitted-level-2);\n",
"}\n",
"\n",
"/* Label box, darken color on hover, fitted */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label.fitted:hover label.sk-toggleable__label.fitted {\n",
2025-03-24 09:57:14 +08:00
" color: var(--sklearn-color-text);\n",
" background-color: var(--sklearn-color-fitted-level-2);\n",
"}\n",
"\n",
"/* Estimator label */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label label {\n",
2025-03-24 09:57:14 +08:00
" font-family: monospace;\n",
" font-weight: bold;\n",
" display: inline-block;\n",
" line-height: 1.2em;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-label-container {\n",
2025-03-24 09:57:14 +08:00
" text-align: center;\n",
"}\n",
"\n",
"/* Estimator-specific */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-estimator {\n",
2025-03-24 09:57:14 +08:00
" font-family: monospace;\n",
" border: 1px dotted var(--sklearn-color-border-box);\n",
" border-radius: 0.25em;\n",
" box-sizing: border-box;\n",
" margin-bottom: 0.5em;\n",
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-0);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-estimator.fitted {\n",
2025-03-24 09:57:14 +08:00
" /* fitted */\n",
" background-color: var(--sklearn-color-fitted-level-0);\n",
"}\n",
"\n",
"/* on hover */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-estimator:hover {\n",
2025-03-24 09:57:14 +08:00
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-2);\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 div.sk-estimator.fitted:hover {\n",
2025-03-24 09:57:14 +08:00
" /* fitted */\n",
" background-color: var(--sklearn-color-fitted-level-2);\n",
"}\n",
"\n",
"/* Specification for estimator info (e.g. \"i\" and \"?\") */\n",
"\n",
"/* Common style for \"i\" and \"?\" */\n",
"\n",
".sk-estimator-doc-link,\n",
"a:link.sk-estimator-doc-link,\n",
"a:visited.sk-estimator-doc-link {\n",
" float: right;\n",
" font-size: smaller;\n",
" line-height: 1em;\n",
" font-family: monospace;\n",
" background-color: var(--sklearn-color-background);\n",
" border-radius: 1em;\n",
" height: 1em;\n",
" width: 1em;\n",
" text-decoration: none !important;\n",
" margin-left: 1ex;\n",
" /* unfitted */\n",
" border: var(--sklearn-color-unfitted-level-1) 1pt solid;\n",
" color: var(--sklearn-color-unfitted-level-1);\n",
"}\n",
"\n",
".sk-estimator-doc-link.fitted,\n",
"a:link.sk-estimator-doc-link.fitted,\n",
"a:visited.sk-estimator-doc-link.fitted {\n",
" /* fitted */\n",
" border: var(--sklearn-color-fitted-level-1) 1pt solid;\n",
" color: var(--sklearn-color-fitted-level-1);\n",
"}\n",
"\n",
"/* On hover */\n",
"div.sk-estimator:hover .sk-estimator-doc-link:hover,\n",
".sk-estimator-doc-link:hover,\n",
"div.sk-label-container:hover .sk-estimator-doc-link:hover,\n",
".sk-estimator-doc-link:hover {\n",
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-3);\n",
" color: var(--sklearn-color-background);\n",
" text-decoration: none;\n",
"}\n",
"\n",
"div.sk-estimator.fitted:hover .sk-estimator-doc-link.fitted:hover,\n",
".sk-estimator-doc-link.fitted:hover,\n",
"div.sk-label-container:hover .sk-estimator-doc-link.fitted:hover,\n",
".sk-estimator-doc-link.fitted:hover {\n",
" /* fitted */\n",
" background-color: var(--sklearn-color-fitted-level-3);\n",
" color: var(--sklearn-color-background);\n",
" text-decoration: none;\n",
"}\n",
"\n",
"/* Span, style for the box shown on hovering the info icon */\n",
".sk-estimator-doc-link span {\n",
" display: none;\n",
" z-index: 9999;\n",
" position: relative;\n",
" font-weight: normal;\n",
" right: .2ex;\n",
" padding: .5ex;\n",
" margin: .5ex;\n",
" width: min-content;\n",
" min-width: 20ex;\n",
" max-width: 50ex;\n",
" color: var(--sklearn-color-text);\n",
" box-shadow: 2pt 2pt 4pt #999;\n",
" /* unfitted */\n",
" background: var(--sklearn-color-unfitted-level-0);\n",
" border: .5pt solid var(--sklearn-color-unfitted-level-3);\n",
"}\n",
"\n",
".sk-estimator-doc-link.fitted span {\n",
" /* fitted */\n",
" background: var(--sklearn-color-fitted-level-0);\n",
" border: var(--sklearn-color-fitted-level-3);\n",
"}\n",
"\n",
".sk-estimator-doc-link:hover span {\n",
" display: block;\n",
"}\n",
"\n",
"/* \"?\"-specific style due to the `<a>` HTML tag */\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 a.estimator_doc_link {\n",
2025-03-24 09:57:14 +08:00
" float: right;\n",
" font-size: 1rem;\n",
" line-height: 1em;\n",
" font-family: monospace;\n",
" background-color: var(--sklearn-color-background);\n",
" border-radius: 1rem;\n",
" height: 1rem;\n",
" width: 1rem;\n",
" text-decoration: none;\n",
" /* unfitted */\n",
" color: var(--sklearn-color-unfitted-level-1);\n",
" border: var(--sklearn-color-unfitted-level-1) 1pt solid;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 a.estimator_doc_link.fitted {\n",
2025-03-24 09:57:14 +08:00
" /* fitted */\n",
" border: var(--sklearn-color-fitted-level-1) 1pt solid;\n",
" color: var(--sklearn-color-fitted-level-1);\n",
"}\n",
"\n",
"/* On hover */\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 a.estimator_doc_link:hover {\n",
2025-03-24 09:57:14 +08:00
" /* unfitted */\n",
" background-color: var(--sklearn-color-unfitted-level-3);\n",
" color: var(--sklearn-color-background);\n",
" text-decoration: none;\n",
"}\n",
"\n",
2025-03-24 15:19:11 +08:00
"#sk-container-id-2 a.estimator_doc_link.fitted:hover {\n",
2025-03-24 09:57:14 +08:00
" /* fitted */\n",
" background-color: var(--sklearn-color-fitted-level-3);\n",
"}\n",
2025-03-24 15:19:11 +08:00
"</style><div id=\"sk-container-id-2\" class=\"sk-top-container\"><div class=\"sk-text-repr-fallback\"><pre>RandomizedSearchCV(cv=3,\n",
2025-03-24 09:57:14 +08:00
" estimator=XGBRegressor(base_score=None, booster=None,\n",
" callbacks=None,\n",
" colsample_bylevel=None,\n",
" colsample_bynode=None,\n",
" colsample_bytree=None, device=None,\n",
" early_stopping_rounds=None,\n",
" enable_categorical=False,\n",
" eval_metric=None, feature_types=None,\n",
2025-03-24 15:19:11 +08:00
" feature_weights=None, gamma=None,\n",
" grow_policy=None,\n",
2025-03-24 09:57:14 +08:00
" importance_type=None,\n",
2025-03-24 15:19:11 +08:00
" interaction_constraint...\n",
" 'gamma': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x142527890>,\n",
" 'learning_rate': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x144aa2de0>,\n",
" 'max_depth': <scipy.stats._distn_infrastructure.rv_discrete_frozen object at 0x144a6e690>,\n",
2025-03-24 09:57:14 +08:00
" 'n_estimators': [100, 200, 300],\n",
2025-03-24 15:19:11 +08:00
" 'subsample': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x144aa1d90>},\n",
" scoring='neg_mean_absolute_error', verbose=1)</pre><b>In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook. <br />On GitHub, the HTML representation is unable to render, please try loading this page with nbviewer.org.</b></div><div class=\"sk-container\" hidden><div class=\"sk-item sk-dashed-wrapped\"><div class=\"sk-label-container\"><div class=\"sk-label fitted sk-toggleable\"><input class=\"sk-toggleable__control sk-hidden--visually\" id=\"sk-estimator-id-4\" type=\"checkbox\" ><label for=\"sk-estimator-id-4\" class=\"sk-toggleable__label fitted sk-toggleable__label-arrow fitted\"> RandomizedSearchCV<a class=\"sk-estimator-doc-link fitted\" rel=\"noreferrer\" target=\"_blank\" href=\"https://scikit-learn.org/1.4/modules/generated/sklearn.model_selection.RandomizedSearchCV.html\">?<span>Documentation for RandomizedSearchCV</span></a><span class=\"sk-estimator-doc-link fitted\">i<span>Fitted</span></span></label><div class=\"sk-toggleable__content fitted\"><pre>RandomizedSearchCV(cv=3,\n",
2025-03-24 09:57:14 +08:00
" estimator=XGBRegressor(base_score=None, booster=None,\n",
" callbacks=None,\n",
" colsample_bylevel=None,\n",
" colsample_bynode=None,\n",
" colsample_bytree=None, device=None,\n",
" early_stopping_rounds=None,\n",
" enable_categorical=False,\n",
" eval_metric=None, feature_types=None,\n",
2025-03-24 15:19:11 +08:00
" feature_weights=None, gamma=None,\n",
" grow_policy=None,\n",
2025-03-24 09:57:14 +08:00
" importance_type=None,\n",
2025-03-24 15:19:11 +08:00
" interaction_constraint...\n",
" 'gamma': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x142527890>,\n",
" 'learning_rate': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x144aa2de0>,\n",
" 'max_depth': <scipy.stats._distn_infrastructure.rv_discrete_frozen object at 0x144a6e690>,\n",
2025-03-24 09:57:14 +08:00
" 'n_estimators': [100, 200, 300],\n",
2025-03-24 15:19:11 +08:00
" 'subsample': <scipy.stats._distn_infrastructure.rv_continuous_frozen object at 0x144aa1d90>},\n",
" scoring='neg_mean_absolute_error', verbose=1)</pre></div> </div></div><div class=\"sk-parallel\"><div class=\"sk-parallel-item\"><div class=\"sk-item\"><div class=\"sk-label-container\"><div class=\"sk-label fitted sk-toggleable\"><input class=\"sk-toggleable__control sk-hidden--visually\" id=\"sk-estimator-id-5\" type=\"checkbox\" ><label for=\"sk-estimator-id-5\" class=\"sk-toggleable__label fitted sk-toggleable__label-arrow fitted\">estimator: XGBRegressor</label><div class=\"sk-toggleable__content fitted\"><pre>XGBRegressor(base_score=None, booster=None, callbacks=None,\n",
2025-03-24 09:57:14 +08:00
" colsample_bylevel=None, colsample_bynode=None,\n",
2025-03-24 15:19:11 +08:00
" colsample_bytree=None, device=None, early_stopping_rounds=None,\n",
" enable_categorical=False, eval_metric=None, feature_types=None,\n",
" feature_weights=None, gamma=None, grow_policy=None,\n",
" importance_type=None, interaction_constraints=None,\n",
" learning_rate=None, max_bin=None, max_cat_threshold=None,\n",
" max_cat_to_onehot=None, max_delta_step=None, max_depth=None,\n",
" max_leaves=None, min_child_weight=None, missing=nan,\n",
" monotone_constraints=None, multi_strategy=None, n_estimators=None,\n",
" n_jobs=-1, num_parallel_tree=None, ...)</pre></div> </div></div><div class=\"sk-serial\"><div class=\"sk-item\"><div class=\"sk-estimator fitted sk-toggleable\"><input class=\"sk-toggleable__control sk-hidden--visually\" id=\"sk-estimator-id-6\" type=\"checkbox\" ><label for=\"sk-estimator-id-6\" class=\"sk-toggleable__label fitted sk-toggleable__label-arrow fitted\"> XGBRegressor<a class=\"sk-estimator-doc-link fitted\" rel=\"noreferrer\" target=\"_blank\" href=\"https://xgboost.readthedocs.io/en/release_3.0.0/python/python_api.html#xgboost.XGBRegressor\">?<span>Documentation for XGBRegressor</span></a></label><div class=\"sk-toggleable__content fitted\"><pre>XGBRegressor(base_score=None, booster=None, callbacks=None,\n",
2025-03-24 09:57:14 +08:00
" colsample_bylevel=None, colsample_bynode=None,\n",
2025-03-24 15:19:11 +08:00
" colsample_bytree=None, device=None, early_stopping_rounds=None,\n",
" enable_categorical=False, eval_metric=None, feature_types=None,\n",
" feature_weights=None, gamma=None, grow_policy=None,\n",
" importance_type=None, interaction_constraints=None,\n",
" learning_rate=None, max_bin=None, max_cat_threshold=None,\n",
" max_cat_to_onehot=None, max_delta_step=None, max_depth=None,\n",
" max_leaves=None, min_child_weight=None, missing=nan,\n",
" monotone_constraints=None, multi_strategy=None, n_estimators=None,\n",
" n_jobs=-1, num_parallel_tree=None, ...)</pre></div> </div></div></div></div></div></div></div></div></div>"
2025-03-24 09:57:14 +08:00
]
},
2025-03-24 15:19:11 +08:00
"execution_count": 25,
2025-03-24 09:57:14 +08:00
"metadata": {},
"output_type": "execute_result"
}
],
2025-03-24 15:19:11 +08:00
"execution_count": 25
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:36:58.084316Z",
"start_time": "2025-03-24T03:36:58.067915Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"#模型预测\n",
"best_model = search.best_estimator_\n",
"y_pred = best_model.predict(X_test)\n",
"#评估指标\n",
"metrics=cal_metrics(y_pred, y_test)\n",
"#输出结果\n",
"print(\"最佳参数组合:\", search.best_params_)\n",
"print(\"评估指标:\")\n",
"for k, v in metrics.items():\n",
" print(f\"{k}: {v:.2f}\")"
],
"id": "fe076794bae89ccb",
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
2025-03-24 15:19:11 +08:00
"最佳参数组合: {'colsample_bytree': 0.901544465464015, 'gamma': 0.179916517642253, 'learning_rate': 0.0342954642051878, 'max_depth': 5, 'n_estimators': 300, 'subsample': 0.8476270714159938}\n",
2025-03-24 09:57:14 +08:00
"评估指标:\n",
2025-03-24 15:19:11 +08:00
"RMSE: 11.54\n",
"R-squared: 0.93\n",
"MAE: 7.42\n"
2025-03-24 09:57:14 +08:00
]
}
],
2025-03-24 15:19:11 +08:00
"execution_count": 26
2025-03-24 09:57:14 +08:00
},
{
"metadata": {
"ExecuteTime": {
2025-03-24 15:19:11 +08:00
"end_time": "2025-03-24T03:37:00.507030Z",
"start_time": "2025-03-24T03:36:58.146270Z"
2025-03-24 09:57:14 +08:00
}
},
"cell_type": "code",
"source": [
"#预测结果可视化\n",
"def plot_results(y_true, y_pred, timestamps):\n",
" plt.figure(figsize=(18, 8))\n",
" ax = plt.gca()\n",
"\n",
" # 绘制预测曲线\n",
" ax.plot(timestamps, y_true, label='真实值',\n",
" marker='o', markersize=4, linewidth=1, alpha=0.8)\n",
" ax.plot(timestamps, y_pred, label='预测值',\n",
" linestyle='--', marker='x', markersize=5, alpha=0.9)\n",
"\n",
" # 设置时间轴格式\n",
" ax.xaxis.set_major_locator(HourLocator(interval=12))\n",
" ax.xaxis.set_minor_locator(HourLocator(interval=3))\n",
" ax.xaxis.set_major_formatter(DateFormatter(\"%m-%d %H:%M\"))\n",
"\n",
" # 增强可视化元素\n",
" plt.title(f'AQI预测效果对比( MAE={metrics[\"MAE\"]:.2f}, R-squared={metrics[\"R-squared\"]:.2f}) ',\n",
" fontsize=14, pad=20)\n",
" plt.xlabel('时间', fontsize=12)\n",
" plt.ylabel('AQI', fontsize=12)\n",
" plt.grid(True, which='both', linestyle='--', alpha=0.5)\n",
" plt.legend()\n",
"\n",
" # 自动调整标签\n",
" plt.xticks(rotation=45, ha='right')\n",
" plt.tight_layout()\n",
2025-03-24 15:19:11 +08:00
" plt.savefig('./images/xg_by_step.png', dpi=200, bbox_inches='tight')\n",
2025-03-24 09:57:14 +08:00
" plt.show()\n",
"\n",
"plot_results(y_test, y_pred, test_data.index)\n",
"\n",
"#特征重要性可视化\n",
"def plot_importance(model, features, top_n=20):\n",
" importance = pd.Series(model.feature_importances_, index=features)\n",
" top_features = importance.sort_values(ascending=False)[:top_n]\n",
"\n",
" plt.figure(figsize=(12, 8))\n",
" ax = top_features.sort_values().plot.barh()\n",
"\n",
" # 添加数据标签\n",
" for i in ax.patches:\n",
" ax.text(i.get_width() + 0.02, i.get_y() + 0.2,\n",
" f'{i.get_width():.2f}',\n",
" fontsize=10, color='dimgrey')\n",
"\n",
" plt.title('Top {} 重要特征'.format(top_n), fontsize=14)\n",
" plt.xlabel('特征重要性', fontsize=12)\n",
" plt.tight_layout()\n",
2025-03-24 15:19:11 +08:00
" plt.savefig('./images/xg_feature_importance.png', dpi=200, bbox_inches='tight')\n",
2025-03-24 09:57:14 +08:00
" plt.show()\n",
"\n",
"\n",
"plot_importance(best_model, features)"
],
"id": "2551eec52baeb4cb",
"outputs": [
{
"data": {
"text/plain": [
"<Figure size 1800x800 with 1 Axes>"
],
2025-03-24 15:19:11 +08:00
"image/png": "iVBORw0KGgoAAAANSUhEUgAABv0AAAMVCAYAAABUfzjNAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAEAAElEQVR4nOzdd3gU1f4G8Hd2N5tNYTchCYRO6AgoSC+Wa8VegWsXC1IsV8F6VbBjBwS7VxEQbOBV8IfgVWyggIJIkRo6AdJ2k022zvf3R9whm91NNpjNTuD9PE8e5Tuz57xzznK8l8PMKCIiICIiIiIiIiIiIiIiIqIGyxDvAERERERERERERERERET093DTj4iIiIiIiIiIiIiIiKiB46YfERERERERERERERERUQPHTT8iIiIiIiIiIiIiIiKiBo6bfkREREREREREREREREQNHDf9iIiIiIiIiIiIiIiIiBo4bvoRERERERERERERERERNXDc9CMiIiIiIiIiIiIiIiJq4LjpR0RERERERERERERERNTAcdOPiIiIiIhIZ/bt24fly5fHOwYR1aMtW7Zg/fr18Y5BRERERA0YN/2IiIiIiCr5/vvvceGFF8Y7Rp0pLCys9WeWLVuGf//739i6dWvY4/v27cPUqVOxaNGisMffffdd/Pbbb/D5fNX2U1ZWhsLCQqiqGlUup9OJgoKCqM6tybPPPov7778fv/76a8ix8vJyeL3eiJ998skncf311+Obb76pkyxViQiGDRuGb7/9Nibth7Nz50707Nmz3vqrL7/88gtGjx4d7xhUh2bPng1FUbBs2bJ4R6lzBQUFOOWUU3DgwIF4RyEiIiKiBoqbfkRERERElbz00kv48ssvsXnz5qg/83//93+4+uqr0bZtWyQlJSEjIwMnn3wyHn74YWzfvj3sZ5YtWwZFUTB79myt9v3332PevHn473//i4ULF+KLL77AggULMG/ePGzbtg0TJ07EtGnT8MYbb+Dtt9/G22+/jbfeeguvvvoqnn/+ebzzzjtBfbz11lvIycnB9OnTUV5ejunTp2P+/PlYuHCh1v7s2bPx2WefBX1u586dePrpp5GTkxM2++7du/Gvf/0LTzzxRMix8vJyjBo1ChdffHGNG47r1q1DRkYGjEYjLBYLrFYrMjIykJmZiczMTGRkZCAtLQ3JyckwGo1ITU3FpEmTqm0zWj/99BOee+459OjRI+TYd999B7PZjHvvvTfsZ7dv345Zs2ahY8eOdZKlqjfeeAN//vkn7rrrrrDHhw4dGvVG6aRJk/D999/XeN6UKVPw+++/Y+nSpdWe17ZtWyiKEvFnypQpUeWqbbvJycnYv39/rdosKyvDddddF3bzef369Rg2bBgyMzNhsVjQvXt3vPjii/D7/UeV/8Ybb9SyGgwGJCYmIjMzE4MHD8Ybb7wR9XxR/Zk7dy769++PlJQUNG7cGJdccgn++OOPqD+/ePFinH766UhJSUGjRo1wzjnnYMWKFSHnuVwuPPvss+jWrRuSk5PRrFkzjB49GocPHw45d+DAgTj55JPxr3/96+9cGhEREREdx0zxDkBEREREpBdbt27FF198ARHB1KlT8eqrr1Z7/uHDh3H11Vfj66+/RlJSEk477TRceOGFcLvd+P333/H000/j2WefxcSJE/Hggw/CaDRW294HH3yAN954I6T+0ksvYe/evXj88ccBACaTCT6fD4qiaP8uIjjzzDNx8803a587++yz0b9/f9xxxx3w+/0R/yD52muvRZMmTfDFF19g7NixSE5O1voJZ+DAgWjevDmKi4tx8OBBNG3aVDv2/fffw+fz4aGHHkKTJk2qvV6r1QoAePrpp3HZZZdp9SeffBILFiwIugvP7/eje/fu2mcq+/nnn6EoCsxmMxRFAVBxt1xZWRmSk5PRq1evkM8Ers1sNoccS0tLAwB069YtbO7APGZlZVV7fUejrKwMjzzyCMaMGYPU1NSQ4+Xl5ViyZAk+/PBDXHXVVdW2ZbfbMXXqVLRq1QqnnnpqxPOKi4u1DeOXX34ZZ599drXttm7dGsOGDQt7rHfv3tV+NpJRo0aF3SRev349vvrqK9x+++1o3rx5rdqcMGECtm7dGnLtX3/9NS6++GIAwPnnn4/09HR89913mDBhAn744QfMnz8fBkPt/36soii45557AABerxcHDx7E119/jdGjR+P333+vcT2h+vPMM8/goYceQosWLXDttdeitLQU//3vf7FkyRIsWbIEp5xySrWfnz59Ou644w40adIEV155JdxuNxYtWoTTTjsNCxYswAUXXACgYt268sorsWjRIvTr1w8jR47E7t278dZbb+Grr77Czz//HLR+AsD48eNxwQUX4K677sKgQYNiNgZEREREdIwSIiIiIiISEZExY8YIAGnWrJkkJydLQUFBxHOLioqka9euAkBuvfVWyc/PDzln3bp1MmDAAAEgt9xyS9Cxb7/9VgDIrFmztFp+fr7k5eVJQUGBFBcXyw033CD9+/cXv98vbrdbiouLRVVVKSsrEwBy2223aZ91Op1SVFQUksHv98urr74qXq9XTCaTvP7661JSUqL9DBkyREaNGiULFiwQAHLgwAHt3yubMmWKzJgxQ2bOnCkffPCBvP322zJr1ix555135MUXX5Q333xTRERuuukmsdlsUlpaqn22rKxM9u3bF5Jt06ZNAkDefffdoPq4ceMkJSUlqOb1egWA/Pvf/w5pp3fv3gIg7M8dd9wRcr6IyOWXXx5yjQGrVq0Kmyvg5ptvFgBSXl4e9vjfMX36dFEURXJzc8Me/+OPPwSAdOnSRfx+f7VtTZw4UQDI/fffX+15zz77rPa9VxRFNm7cGPHcNm3ayJlnnlnjddSVIUOGiNVqDfv7qzqLFy/WvgM333yzVne5XNK8eXNp3LixbNq0Sat7PB4ZOXJkyO/JaN1www1iNBpD6mVlZdK3b18BIFu3bq11uxRq1qxZAkC+/fbbo/r85s2bxWQySd++fcVutwfVMzMzJScnR1wuV8TP7969W8xms5x88slSWFio1XNzc6VZs2aSnZ2trQ0ffvihAJC77rorqI1PP/1UAMjdd98d0r7f75fWrVvLueeee1TXR0RERETHNz7ek4iIiIgIFe9SmjlzJvr06YNp06ahrKws7F13AePGjcOmTZvwxBNP4M0330RGRkbIOT169MB3332Hc845B2+//TZmzZpVbYaMjAzMnj0bY8eOxa+//oq5c+firbfegsFggNlshs1mg6IoyM3NBQCce+652meTk5O1O9QqMxgMGDNmDEwmExITE5GYmIhNmzZh/vz5UBQFRqMRiqKE/WxlL7/8MsaNG4cbbrgBV199NW655RZcd911uPnmmzF+/HjMnDkTBQUF+PDDD1FSUoKmTZtqj+VMTk5Gu3btIrZ94MAB/Pnnn9pPcXExRCSoVt3jVl999VX8+OOPWLVqlfYTeAxo165d8cQTT2Dq1KmYPn269rNjxw4AFXfsvPDCC1i1alW111+Tbdu24bHHHvtbbQAV13LyySejbdu2YY8H3rP4559/4oMPPojYTuAuv8qfCcfr9WLatGlo3rw5PvnkE4jIUT+is7LqHtUZ+Il0jQFffvklfvzxR4wfPz7s769IioqKcNNNN4W9K3TZsmXYv38/HnroIXTp0kWrJyQk4JVXXoHJZMIXX3wRdV81SUpKwtVXXw0AYd8fWdX+/fsxfPhwNG7cGM2bN8ejjz4Kt9sddI7f78fbb7+NgQMHolGjRrBarTjjjDOwaNEi7N27F4qiBD0GV1EU3HLLLWH7i/R+PFVV8e6772Lw4MFaH6eeeirmz5+P/Pz8oD4mTZoERVHw+eefo1u3bmjUqFHQY4uBikdpnnbaabBarUhOTkbPnj0xefJklJeXh81lt9vx6KOPomvXrrBYLMjKysLw4cOxYcOGkHPfe++9qL5vgbzvvvsufD4fZsyYEfQd6dSpE55++mnk5ubi888/D5sLABYtWgSPx4PJkycjPT1dq7dt2xb//ve/kZeXh+XLl2u1uXPnYvLkyUFtBO7
2025-03-24 09:57:14 +08:00
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/plain": [
"<Figure size 1200x800 with 1 Axes>"
],
2025-03-24 15:19:11 +08:00
"image/png": "iVBORw0KGgoAAAANSUhEUgAABKUAAAMWCAYAAAAgRDUeAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguNCwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8fJSN1AAAACXBIWXMAAA9hAAAPYQGoP6dpAAD5uElEQVR4nOzde1iVVf7//9cW5LAVUJE+IqKQkA6iRkMmVgpoiaJTap7SKDXL1LI002pKRWdkPo4/5NCM4SExQ820clBHTeiAOkWOYZqIlodMTDwBAoLg/v3R1/1xt1GRYIP5fFzXfV261rrv9b6VP/R1rbVug8lkMgkAAAAAAACwoQZ1XQAAAAAAAABuP4RSAAAAAAAAsDlCKQAAAAAAANgcoRQAAAAAAABsjlAKAAAAAAAANkcoBQAAAAAAAJsjlAIAAAAAAIDNEUoBAAAAAADA5gilAAAAfmeOHz+uoqKiSvsKCwtVXl5ereeWlpYqPz+/0r5Lly4pPT1du3fvrrR/9uzZev7553X27Nkqz3fx4sVq1QkAAG4NhFIAAKBe8vHxkcFguOZVH5w5c0aTJk1S69at5eDgoNatW+uFF17QmTNnrMbu2rVLkZGRatq0qRo1aqRu3brpww8/rNI8eXl5ev/997Vx40Zt2rRJGzZs0Lp167Rs2TItWLBAr732mkaPHq3Q0FB5eXnJ29tbr776aqXPmjp1qho2bCh7e3s5OTmpcePG172MRqMcHBzUoEEDOTk5acyYMZU+187OTj179lRcXFyl/f/973/18ccfq1mzZlV6588++0zt27dXWlpalcYDAIBbj31dFwAAAFCZZ555xryqZvXq1Tp+/LimTJlSx1X9n59++kn333+/jh49qvDwcEVGRurbb79VQkKCNm7cqB07duiOO+6Q9EvA8vDDD6thw4Z65JFH1LhxY23cuFEDBw7UvHnz9PLLL193Ljc3N40bN07nzp0ztxmNRjVp0kTnzp1To0aN1LVrV7Vo0UJ33nmn7OzslJ+fr6NHj6pNmzYWz3rllVc0fvx4NW7cWA4ODmrYsKH27Nmjhx9+WP/973/l7e2tLVu2aMSIETp8+LCMRqMuXbqksrIylZSUyN7e8p+Pe/bsUbNmzdSyZUsZjUa1atXK3Pftt9/KYDDIxcVFJpNJJpNJx48f18WLF1VQUKAmTZrozjvvtHrfoqIiPf300yooKJC7u7vOnDmjGTNmqGnTpjIajbK3tzcHk5cvX9alS5d04cIFtWzZUs8///zN/UUCAIA6YzCZTKa6LgIAAOB6QkND9dlnn6k+/bNl0KBBWrdunVavXq0hQ4aY25csWaKnn35aY8aM0eLFi1VeXq727dvrzJkz2rFjh/7whz9Iki5cuKCIiAh9+eWXysrKUkBAQKXz5Ofn6/Lly2rQ4JcF7vb29nJ2djb/PjQ0VI0bN1ZqaqrVvWVlZSouLpYkNWnSxNxuMplUUlIiJycnNWjQQF9//bXuvfdeHT58WD4+PkpNTVX//v1VWFioxo0b6/Lly7p48aKcnZ2tVqkNGDBAGzdu1NGjR9W+fXtNnDhRc+bMkSRFRERo8+bN1/wzfP31181jr65t+PDhWrt2rf7973+rZ8+eOnLkiHx9fSVJjo6Osre3V1FRkbme0tJSVVRUqHfv3vr3v/99zfkAAED9wvY9AABwy/vhhx/09NNPy9vbW46OjmrVqpXGjBmjH374wWpsaGio/Pz8VF5ergULFqhTp05ydnZW8+bN1b9/f3322Wc3nO/ChQtav369+vXrZxFISdKYMWMUEhKif/3rX5Kkbdu26fvvv9drr71mDqQkqXHjxlqyZIkqKir01ltvXXOuGTNmqFmzZmrSpImaNm0qd3d3ubq6ys3NTU2aNFFGRoY2b96sJk2aqEmTJnJzc5Ozs7Ps7Ozk6Oiopk2b6umnn7aqv1GjRrKzs5O9vb26desmSWrXrp2cnJw0cOBASZK7u7vs7e1lZ2enRo0a6cKFC1b1nTp1SiEhIWrRooVV31NPPaVFixZp9erV6tGjh+644w6lpqbqgw8+0JIlSxQZGWkx/vLly3r22We1evVqLVy4UD169JAkeXt768yZM6qoqNDFixe1d+9eSdI//vEPFRUVqby8XEVFRVq1atU1/xwBAED9w/Y9AABwS9u8ebMee+wxlZaWKiIiQj4+Pvrhhx+0fPlyrVmzRuvWrVOvXr0s7iktLVWfPn20bds2hYWFqUePHjpz5ow2bdqk1NRUzZkzR6+//vo15zx9+rQ6dOigAQMGVNrfokULHThwQJLMIdfgwYOtxrVr10733nuvtmzZcs25Jk6cqMcff1yNGjWSs7OzGjZsqIYNG6pBgwYyGAx69NFH1ahRI61YsUImk0mXL19WeXm5Ll26pJKSEhUWFsrR0dHimS4uLmrQoIHeeOMNjRw5Unv37tWAAQO0ZcsWeXl5KT09Xc8884wyMzNlNBq1bNky/eUvf7F6jiQdPXpU58+fl4+PjwoKCpSQkKC3335bkydPtjjXav369Tpy5IhVEHW13bt366OPPtL8+fM1ZswYTZo0SYWFhfr73/9+w7OojEajjEbjdccAAID6hVAKAADcsr7//ns99thjcnd314YNG9ShQwdz3549e9SvXz8NHDhQWVlZ5u1fksznGmVkZJhXCUnSuXPnNGTIEP35z3+Wv7+/1SqoK3x8fPTNN99U2ldUVKSMjAz98Y9/lCQdPHhQjRs3lo+PT6Xj7777bi1evFiXLl1Sw4YNrfr9/PzMv/7qq6/k5OQkJycn82Hlly9f1uXLl1VWVibpl9VGpaWlOn/+vO666y65ublVOq+9vb2aN28uPz8/nT9/XpLUpk0b+fj4KDs7W5J05513qnHjxmrevLkkmbcMXnH69Gn99NNPevjhh+Xv769FixbJz89PnTp1UllZmV599VU1bdpUzs7OOnjwoAoLC7Vw4ULzOVDnzp2Tl5eXxo4dK0n64x//qJycHDVp0kS7d+/WW2+9pR49esjJyanSdwAAALc2QikAAHDLmjt3roqKivTZZ59ZBFKS1KlTJ61du1b33Xef/vrXv2rRokUW/UuXLrUIpCSpadOm+vDDDxUQEKDXXnvtmqHU9UyaNEl5eXmaPn26pF/OhLreKp/mzZvr8uXLys/PN4c/1/LAAw/o0qVLlfZ5e3tbtX377beVhlImk0llZWX6+eeflZ2draNHj0r6JeS7ePGijh8/LknKycmR0WjUqVOnzPddbevWrZKkv/71r/rjH/+oFStWqE+fPpozZ475XCrplzDL0dFRdnZ2evHFF1VaWioHBweVlZUpMjLSHEpJv5x9VVZWptGjR8vFxUXvvvuujEajpk6dqqysLK1cudLqfSoqKnTs2DH9z//8D6ulAAC4hXCmFAAAuGX9+9//1v3336977rmn0v57771X3bp1szr82tPT0xyY/Frjxo311FNP6fvvv9f3339/U/XMnj1bS5Ys0eTJkxUeHm5u//Xh4Fe7fPmyJOvApzLOzs56/fXXVVJSYr4efPBB9enTx6ItKSlJkq4Z0BQUFEiS5syZo8DAQA0fPlyOjo6KjIzU3XffrRdffFGOjo7q2rWrOnTooLlz50qS+dD0K4KCgrR8+XJ17tzZao7w8HCdPHlSxcXFqqioUHFxsQoLC81h3c8//6yysrJKQ6bnnntO33zzjV555RW1bNlSknT27Fnt27dP7u7u5nHjxo1To0aN1LBhQ915552VniEGAADqL1ZKAQCAW9apU6fUs2fP64658847lZmZadFW2aqiq13Zavfzzz+rbdu2VaolPj5eb775pgYMGKB58+aZ293c3HTu3Llr3nf69GkZDIZrbrO7mr29vf7yl7/oL3/5i1Wfs7OzVduvt9tdXVNZWZnFdsHc3FzNmzdPTzzxhIKCgizGX9kS+OttdO3bt1f79u0t2rKyshQVFaVnn31Wfn5+atCggUwmk1UwV1ZWpkuXLlmdU/Xyyy9r6dKlkn5ZuXaFo6Oj1fbG7t2766677jIfdF6VP0MAAFB/EEoBAIBbloeHh3766afrjvnxxx/
2025-03-24 09:57:14 +08:00
},
"metadata": {},
"output_type": "display_data"
}
],
2025-03-24 15:19:11 +08:00
"execution_count": 27
2025-03-24 09:57:14 +08:00
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3",
"language": "python",
"name": "python3"
},
"language_info": {
"codemirror_mode": {
"name": "ipython",
"version": 2
},
"file_extension": ".py",
"mimetype": "text/x-python",
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython2",
"version": "2.7.6"
}
},
"nbformat": 4,
"nbformat_minor": 5
}