空气质量仿真软件：CALPUFF (California PUFF Model)_（10）.模式验证与评估

于 2025-05-05 07:22:11 发布
阅读量769 收藏 16
点赞数 20
分类专栏： 环境仿真 文章标签： 仿真模拟 环境仿真 模拟仿真
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/weixin_42749425/article/details/147707840
版权

模式验证与评估

在空气质量仿真软件的开发过程中，模式验证与评估是确保模型结果准确性和可靠性的关键步骤。CALPUFF 是一种广泛使用的三维大气扩散模型，用于模拟污染物在大气中的传输、扩散和沉降过程。为了确保模型的预测结果能够准确反映实际情况，模式验证与评估是必不可少的。本节将详细介绍模式验证与评估的基本概念、方法和具体操作步骤，包括如何准备验证数据、评估指标的计算方法以及如何使用 CALPUFF 的输出数据进行验证和评估。

在这里插入图片描述

1. 验证数据的准备

模式验证的第一步是准备验证数据。验证数据通常包括观测数据和模型预测数据。观测数据可以通过现场监测站获取，而模型预测数据则是通过运行 CALPUFF 模型得到的。以下是一些具体的步骤和注意事项：

1.1 获取观测数据

  1. 选择监测站：选择具有代表性的监测站，确保这些站点能够覆盖模型预测的区域。

  2. 数据质量控制：对观测数据进行质量控制，剔除异常值和缺失数据。

  3. 数据格式转换：将观测数据转换为 CALPUFF 可以读取的格式，如 ASCII 文本文件。

代码示例：数据质量控制

import pandas as pd



# 读取观测数据

obs_data = pd.read_csv('observation_data.csv')



# 检查数据缺失值

missing_values = obs_data.isnull().sum()

print(missing_values)



# 剔除缺失值

obs_data = obs_data.dropna()



# 检查异常值

def detect_outliers(df, column):

    Q1 = df[column].quantile(0.25)

    Q3 = df[column].quantile(0.75)

    IQR = Q3 - Q1

    lower_bound = Q1 - 1.5 * IQR

    upper_bound = Q3 + 1.5 * IQR

    outliers = df[(df[column] < lower_bound) | (df[column] > upper_bound)]

    return outliers



outliers = detect_outliers(obs_data, 'PM2.5')

print(outliers)



# 剔除异常值

obs_data = obs_data[(obs_data['PM2.5'] >= lower_bound) & (obs_data['PM2.5'] <= upper_bound)]



# 保存处理后的数据

obs_data.to_csv('cleaned_observation_data.csv', index=False)

1.2 准备模型预测数据

  1. 运行 CALPUFF 模型：根据研究区域和污染物类型，设置模型参数并运行模型。

  2. 提取预测数据：从模型输出文件中提取所需的预测数据。

  3. 数据格式转换：将预测数据转换为与观测数据相同的格式，以便进行对比分析。

代码示例：提取模型预测数据

import numpy as np



# 读取模型输出文件

model_output = np.loadtxt('model_output.txt')



# 提取特定站点的预测数据

site_id = 1  # 假设观测站ID为1

predicted_pm25 = model_output[:, site_id]



# 保存预测数据

np.savetxt('predicted_pm25.txt', predicted_pm25, fmt='%.2f')

2. 评估指标的计算方法

模式验证与评估通常使用一系列评估指标来衡量模型的预测效果。常见的评估指标包括均方根误差（RMSE）、平均绝对误差（MAE）、相关系数（R）和统计显著性检验等。

2.1 均方根误差（RMSE）

均方根误差是衡量模型预测值与观测值之间差异的常用指标。计算公式如下：

RMSE = 1 n ∑ i = 1 n ( y i − y ^ i ) 2 \text{RMSE} = \sqrt{\frac{1}{n} \sum_{i=1}^{n} (y_i - \hat{y}_i)^2} RMSE=n1​i=1∑n​(yi​−y^​i​)2 ​

其中， y i y_i yi​ 是观测值， y ^ i \hat{y}_i y^​i​ 是预测值， n n n 是数据点的数量。

代码示例：计算 RMSE

import numpy as np



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算 RMSE

rmse = np.sqrt(np.mean((obs_pm25 - pred_pm25) ** 2))

print(f'RMSE: {rmse:.2f}')

2.2 平均绝对误差（MAE）

平均绝对误差是衡量模型预测值与观测值之间差异的另一个常用指标。计算公式如下：

MAE = 1 n ∑ i = 1 n ∣ y i − y ^ i ∣ \text{MAE} = \frac{1}{n} \sum_{i=1}^{n} |y_i - \hat{y}_i| MAE=n1​i=1∑n​∣yi​−y^​i​∣

代码示例：计算 MAE

import numpy as np



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算 MAE

mae = np.mean(np.abs(obs_pm25 - pred_pm25))

print(f'MAE: {mae:.2f}')

2.3 相关系数（R）

相关系数用于衡量模型预测值与观测值之间的线性相关性。计算公式如下：

R = ∑ i = 1 n ( y i − y ˉ ) ( y ^ i − y ^ ˉ ) ∑ i = 1 n ( y i − y ˉ ) 2 ∑ i = 1 n ( y ^ i − y ^ ˉ ) 2 R = \frac{\sum_{i=1}^{n} (y_i - \bar{y})(\hat{y}_i - \bar{\hat{y}})}{\sqrt{\sum_{i=1}^{n} (y_i - \bar{y})^2 \sum_{i=1}^{n} (\hat{y}_i - \bar{\hat{y}})^2}} R=∑i=1n​(yi​−yˉ​)2∑i=1n​(y^​i​−y^​ˉ​)2 ​∑i=1n​(yi​−yˉ​)(y^​i​−y^​ˉ​)​

其中， y ˉ \bar{y} yˉ​ 和 y ^ ˉ \bar{\hat{y}} y^​ˉ​ 分别是观测值和预测值的均值。

代码示例：计算相关系数

import numpy as np

from scipy.stats import pearsonr



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算相关系数

r, _ = pearsonr(obs_pm25, pred_pm25)

print(f'Correlation Coefficient (R): {r:.2f}')

2.4 统计显著性检验

统计显著性检验用于评估模型预测值与观测值之间的差异是否具有统计学意义。常见的方法包括 t 检验和 F 检验等。

代码示例：进行 t 检验

import numpy as np

from scipy.stats import ttest_ind



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 进行 t 检验

t_stat, p_value = ttest_ind(obs_pm25, pred_pm25)

print(f't-statistic: {t_stat:.2f}')

print(f'p-value: {p_value:.2f}')

3. 模型结果的可视化

模型结果的可视化是验证与评估过程中不可或缺的一部分。通过可视化，可以直观地观察模型预测值与观测值之间的差异，识别潜在的问题和改进方向。

3.1 时间序列图

时间序列图可以展示观测值和预测值随时间的变化情况。以下是一个使用 Matplotlib 绘制时间序列图的示例：

代码示例：绘制时间序列图

import matplotlib.pyplot as plt

import pandas as pd



# 读取观测数据和预测数据

obs_data = pd.read_csv('cleaned_observation_data.csv', parse_dates=['date'], index_col='date')

pred_data = pd.read_csv('predicted_pm25.csv', parse_dates=['date'], index_col='date')



# 绘制时间序列图

plt.figure(figsize=(12, 6))

plt.plot(obs_data.index, obs_data['PM2.5'], label='Observation', color='blue')

plt.plot(pred_data.index, pred_data['PM2.5'], label='Prediction', color='red', linestyle='--')

plt.xlabel('Date')

plt.ylabel('PM2.5 Concentration (μg/m³)')

plt.title('PM2.5 Concentration Time Series')

plt.legend()

plt.grid(True)

plt.show()

3.2 散点图

散点图可以展示观测值与预测值之间的关系，帮助识别模型的预测误差。以下是一个使用 Matplotlib 绘制散点图的示例：

代码示例：绘制散点图

import matplotlib.pyplot as plt

import numpy as np



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 绘制散点图

plt.figure(figsize=(8, 8))

plt.scatter(obs_pm25, pred_pm25, color='blue', alpha=0.5)

plt.xlabel('Observation (μg/m³)')

plt.ylabel('Prediction (μg/m³)')

plt.title('PM2.5 Concentration Scatter Plot')

plt.grid(True)

plt.axline((0, 0), slope=1, color='red', linestyle='--', label='1:1 Line')

plt.legend()

plt.show()

3.3 残差图

残差图可以展示模型预测值与观测值之间的误差分布，帮助识别模型的系统性偏差。以下是一个使用 Matplotlib 绘制残差图的示例：

代码示例：绘制残差图

import matplotlib.pyplot as plt

import numpy as np



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算残差

residuals = obs_pm25 - pred_pm25



# 绘制残差图

plt.figure(figsize=(10, 6))

plt.scatter(pred_pm25, residuals, color='blue', alpha=0.5)

plt.axhline(y=0, color='red', linestyle='--', label='Zero Line')

plt.xlabel('Prediction (μg/m³)')

plt.ylabel('Residual (Observation - Prediction)')

plt.title('PM2.5 Concentration Residual Plot')

plt.grid(True)

plt.legend()

plt.show()

4. 模型性能的综合评估

综合评估模型性能时，通常需要考虑多个评估指标和可视化结果。以下是一些综合评估的方法和步骤：

4.1 多指标综合评估

  1. 计算多个评估指标：计算 RMSE、MAE、R 等多个评估指标。

  2. 对比分析：将计算得到的评估指标与行业标准或基准模型进行对比，评估模型的性能。

代码示例：计算多个评估指标

import numpy as np

from scipy.stats import pearsonr

from sklearn.metrics import mean_squared_error, mean_absolute_error



# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算 RMSE

rmse = np.sqrt(mean_squared_error(obs_pm25, pred_pm25))



# 计算 MAE

mae = mean_absolute_error(obs_pm25, pred_pm25)



# 计算相关系数

r, _ = pearsonr(obs_pm25, pred_pm25)



# 输出评估指标

print(f'RMSE: {rmse:.2f}')

print(f'MAE: {mae:.2f}')

print(f'Correlation Coefficient (R): {r:.2f}')

4.2 模型改进与优化

  1. 识别问题：通过评估指标和可视化结果，识别模型的主要问题。

  2. 调整参数：根据问题调整模型参数，如气象数据、源强、地形数据等。

  3. 重新验证：重新运行模型并进行验证，评估改进后的性能。

代码示例：调整气象数据参数

假设我们发现模型在特定气象条件下的预测误差较大，可以调整气象数据的输入参数。以下是一个调整气象数据参数的示例：


# 读取气象数据

met_data = pd.read_csv('met_data.csv')



# 调整气象参数，例如风速

met_data['WindSpeed'] = met_data['WindSpeed'] * 1.1  # 增加10%的风速



# 保存调整后的气象数据

met_data.to_csv('adjusted_met_data.csv', index=False)

5. 实例分析

为了更好地理解模式验证与评估的过程，我们通过一个具体的实例进行分析。假设我们有一个研究区域，包含多个监测站的 PM2.5 观测数据，并使用 CALPUFF 模型进行预测。

5.1 数据准备

  1. 获取观测数据：从监测站下载 PM2.5 观测数据，并进行质量控制。

  2. 运行模型：设置 CALPUFF 模型参数并运行模型，提取预测数据。

5.2 评估指标计算

  1. 计算 RMSE、MAE 和 R：使用上述代码示例计算评估指标。

  2. 进行 t 检验：使用上述代码示例进行 t 检验。

5.3 模型结果可视化

  1. 绘制时间序列图：使用上述代码示例绘制时间序列图。

  2. 绘制散点图：使用上述代码示例绘制散点图。

  3. 绘制残差图：使用上述代码示例绘制残差图。

5.4 模型改进与优化

  1. 识别问题：通过评估指标和可视化结果，识别模型的主要问题。

  2. 调整参数：根据问题调整模型参数，重新运行模型。

  3. 重新验证：重新运行模型并进行验证，评估改进后的性能。

代码示例：识别问题并调整参数

# 读取观测数据和预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算 RMSE

rmse = np.sqrt(mean_squared_error(obs_pm25, pred_pm25))



# 计算 MAE

mae = mean_absolute_error(obs_pm25, pred_pm25)



# 计算相关系数

r, _ = pearsonr(obs_pm25, pred_pm25)



# 输出评估指标

print(f'Initial RMSE: {rmse:.2f}')

print(f'Initial MAE: {mae:.2f}')

print(f'Initial Correlation Coefficient (R): {r:.2f}')



# 调整气象参数

met_data = pd.read_csv('met_data.csv')

met_data['WindSpeed'] = met_data['WindSpeed'] * 1.1  # 增加10%的风速

met_data.to_csv('adjusted_met_data.csv', index=False)



# 重新运行模型并提取预测数据

# 假设重新运行模型后得到的新预测数据文件为 'new_predicted_pm25.txt'

new_pred_pm25 = np.loadtxt('new_predicted_pm25.txt')



# 重新计算评估指标

new_rmse = np.sqrt(mean_squared_error(obs_pm25, new_pred_pm25))

new_mae = mean_absolute_error(obs_pm25, new_pred_pm25)

new_r, _ = pearsonr(obs_pm25, new_pred_pm25)



# 输出改进后的评估指标

print(f'Improved RMSE: {new_rmse:.2f}')

print(f'Improved MAE: {new_mae:.2f}')

print(f'Improved Correlation Coefficient (R): {new_r:.2f}')

5.5 结果讨论

  1. 评估指标变化：比较初始和改进后的评估指标，分析模型性能的变化。

  2. 可视化结果：重新绘制时间序列图、散点图和残差图，观察改进后的结果。

  3. 问题总结：总结模型的主要问题和改进措施，为后续研究提供参考。

代码示例：重新绘制时间序列图

import matplotlib.pyplot as plt

import pandas as pd



# 读取观测数据和改进后的预测数据

obs_data = pd.read_csv('cleaned_observation_data.csv', parse_dates=['date'], index_col='date')

new_pred_data = pd.read_csv('new_predicted_pm25.csv', parse_dates=['date'], index_col='date')



# 绘制时间序列图

plt.figure(figsize=(12, 6))

plt.plot(obs_data.index, obs_data['PM2.5'], label='Observation', color='blue')

plt.plot(new_pred_data.index, new_pred_data['PM2.5'], label='Improved Prediction', color='green', linestyle='--')

plt.xlabel('Date')

plt.ylabel('PM2.5 Concentration (μg/m³)')

plt.title('PM2.5 Concentration Time Series (Initial vs Improved)')

plt.legend()

plt.grid(True)

plt.show()

6. 验证报告的编写

编写验证报告是模式验证与评估的最后一步。验证报告应包括以下内容：

  1. 研究背景：简要介绍研究区域和污染物类型。

  2. 数据准备：描述观测数据和预测数据的准备过程。

  3. 评估方法：介绍使用的评估指标和统计方法。

  4. 结果分析：展示评估指标和可视化结果，分析模型性能。

  5. 改进措施：总结模型的主要问题和改进措施。

  6. 结论：提出模型的进一步改进方向和建议。

6.1 研究背景

本研究区域位于某一城市，包含多个监测站的 PM2.5 观测数据。我们使用 CALPUFF 模型对该区域的 PM2.5 浓度进行预测，以评估模型的准确性和可靠性。

6.2 数据准备

6.2.1 获取观测数据

  1. 选择监测站：选择具有代表性的监测站，确保这些站点能够覆盖模型预测的区域。

  2. 数据质量控制：对观测数据进行质量控制，剔除异常值和缺失数据。

  3. 数据格式转换：将观测数据转换为 CALPUFF 可以读取的格式，如 ASCII 文本文件。

6.2.2 准备模型预测数据

  1. 运行 CALPUFF 模型：根据研究区域和污染物类型，设置模型参数并运行模型。

  2. 提取预测数据：从模型输出文件中提取所需的预测数据。

  3. 数据格式转换：将预测数据转换为与观测数据相同的格式，以便进行对比分析。

6.3 评估方法

6.3.1 评估指标

  1. 均方根误差（RMSE）：衡量模型预测值与观测值之间差异的常用指标。

  2. 平均绝对误差（MAE）：衡量模型预测值与观测值之间差异的另一个常用指标。

  3. 相关系数（R）：衡量模型预测值与观测值之间的线性相关性。

  4. 统计显著性检验：评估模型预测值与观测值之间的差异是否具有统计学意义，常用方法包括 t 检验和 F 检验。

6.3.2 可视化方法

  1. 时间序列图：展示观测值和预测值随时间的变化情况。

  2. 散点图：展示观测值与预测值之间的关系，帮助识别模型的预测误差。

  3. 残差图：展示模型预测值与观测值之间的误差分布，帮助识别模型的系统性偏差。

6.4 结果分析

6.4.1 初始评估指标

初始评估指标如下：

  • RMSE：10.5

  • MAE：7.2

  • 相关系数（R）：0.75

  • t 检验 p 值：0.03

6.4.2 初始可视化结果

初始时间序列图、散点图和残差图如下：

  1. 时间序列图：

    
import matplotlib.pyplot as plt

import pandas as pd



# 读取观测数据和初始预测数据

obs_data = pd.read_csv('cleaned_observation_data.csv', parse_dates=['date'], index_col='date')

pred_data = pd.read_csv('predicted_pm25.csv', parse_dates=['date'], index_col='date')



# 绘制时间序列图

plt.figure(figsize=(12, 6))

plt.plot(obs_data.index, obs_data['PM2.5'], label='Observation', color='blue')

plt.plot(pred_data.index, pred_data['PM2.5'], label='Prediction', color='red', linestyle='--')

plt.xlabel('Date')

plt.ylabel('PM2.5 Concentration (μg/m³)')

plt.title('PM2.5 Concentration Time Series (Initial)')

plt.legend()

plt.grid(True)

plt.show()

  2. 散点图：

    
import matplotlib.pyplot as plt

import numpy as np



# 读取观测数据和初始预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 绘制散点图

plt.figure(figsize=(8, 8))

plt.scatter(obs_pm25, pred_pm25, color='blue', alpha=0.5)

plt.xlabel('Observation (μg/m³)')

plt.ylabel('Prediction (μg/m³)')

plt.title('PM2.5 Concentration Scatter Plot (Initial)')

plt.grid(True)

plt.axline((0, 0), slope=1, color='red', linestyle='--', label='1:1 Line')

plt.legend()

plt.show()

  3. 残差图：

    
import matplotlib.pyplot as plt

import numpy as np



# 读取观测数据和初始预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

pred_pm25 = np.loadtxt('predicted_pm25.txt')



# 计算残差

residuals = obs_pm25 - pred_pm25



# 绘制残差图

plt.figure(figsize=(10, 6))

plt.scatter(pred_pm25, residuals, color='blue', alpha=0.5)

plt.axhline(y=0, color='red', linestyle='--', label='Zero Line')

plt.xlabel('Prediction (μg/m³)')

plt.ylabel('Residual (Observation - Prediction)')

plt.title('PM2.5 Concentration Residual Plot (Initial)')

plt.grid(True)

plt.legend()

plt.show()

6.5 改进措施

6.5.1 识别问题

通过初始评估指标和可视化结果，我们发现模型在某些气象条件下的预测误差较大。具体表现为：

  • RMSE：10.5

  • MAE：7.2

  • 相关系数（R）：0.75

  • t 检验 p 值：0.03

6.5.2 调整气象参数

根据问题，我们调整了气象数据中的风速参数，增加10%的风速，重新运行模型并提取预测数据。


import pandas as pd



# 读取气象数据

met_data = pd.read_csv('met_data.csv')



# 调整气象参数，例如风速

met_data['WindSpeed'] = met_data['WindSpeed'] * 1.1  # 增加10%的风速



# 保存调整后的气象数据

met_data.to_csv('adjusted_met_data.csv', index=False)
6.5.3 重新验证

重新运行模型并提取预测数据后，我们重新计算评估指标并绘制可视化结果。


import numpy as np

from scipy.stats import pearsonr

from sklearn.metrics import mean_squared_error, mean_absolute_error



# 读取观测数据和改进后的预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

new_pred_pm25 = np.loadtxt('new_predicted_pm25.txt')



# 计算改进后的评估指标

new_rmse = np.sqrt(mean_squared_error(obs_pm25, new_pred_pm25))

new_mae = mean_absolute_error(obs_pm25, new_pred_pm25)

new_r, _ = pearsonr(obs_pm25, new_pred_pm25)



# 输出改进后的评估指标

print(f'Improved RMSE: {new_rmse:.2f}')

print(f'Improved MAE: {new_mae:.2f}')

print(f'Improved Correlation Coefficient (R): {new_r:.2f}')

改进后的评估指标如下：

  • RMSE：8.3

  • MAE：5.4

  • 相关系数（R）：0.85

  • t 检验 p 值：0.01

6.5.4 重新绘制可视化结果

  1. 时间序列图：

    
import matplotlib.pyplot as plt

import pandas as pd



# 读取观测数据和改进后的预测数据

obs_data = pd.read_csv('cleaned_observation_data.csv', parse_dates=['date'], index_col='date')

new_pred_data = pd.read_csv('new_predicted_pm25.csv', parse_dates=['date'], index_col='date')



# 绘制时间序列图

plt.figure(figsize=(12, 6))

plt.plot(obs_data.index, obs_data['PM2.5'], label='Observation', color='blue')

plt.plot(new_pred_data.index, new_pred_data['PM2.5'], label='Improved Prediction', color='green', linestyle='--')

plt.xlabel('Date')

plt.ylabel('PM2.5 Concentration (μg/m³)')

plt.title('PM2.5 Concentration Time Series (Initial vs Improved)')

plt.legend()

plt.grid(True)

plt.show()

  2. 散点图：

    
import matplotlib.pyplot as plt

import numpy as np



# 读取观测数据和改进后的预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

new_pred_pm25 = np.loadtxt('new_predicted_pm25.txt')



# 绘制散点图

plt.figure(figsize=(8, 8))

plt.scatter(obs_pm25, new_pred_pm25, color='green', alpha=0.5)

plt.xlabel('Observation (μg/m³)')

plt.ylabel('Prediction (μg/m³)')

plt.title('PM2.5 Concentration Scatter Plot (Improved)')

plt.grid(True)

plt.axline((0, 0), slope=1, color='red', linestyle='--', label='1:1 Line')

plt.legend()

plt.show()

  3. 残差图：

    
import matplotlib.pyplot as plt

import numpy as np



# 读取观测数据和改进后的预测数据

obs_pm25 = np.loadtxt('cleaned_observation_data.csv', delimiter=',', usecols=[1])

new_pred_pm25 = np.loadtxt('new_predicted_pm25.txt')



# 计算改进后的残差

new_residuals = obs_pm25 - new_pred_pm25



# 绘制改进后的残差图

plt.figure(figsize=(10, 6))

plt.scatter(new_pred_pm25, new_residuals, color='green', alpha=0.5)

plt.axhline(y=0, color='red', linestyle='--', label='Zero Line')

plt.xlabel('Prediction (μg/m³)')

plt.ylabel('Residual (Observation - Prediction)')

plt.title('PM2.5 Concentration Residual Plot (Improved)')

plt.grid(True)

plt.legend()

plt.show()

6.6 结论

通过模式验证与评估，我们发现初始模型在某些气象条件下的预测误差较大。通过调整气象参数（增加10%的风速），模型的预测性能得到了显著提升。改进后的评估指标如下：

  • RMSE：从 10.5 降低到 8.3

  • MAE：从 7.2 降低到 5.4

  • 相关系数（R）：从 0.75 提高到 0.85

  • t 检验 p 值：从 0.03 降低到 0.01

这些结果表明，模型在调整参数后，预测值与观测值之间的差异显著减小，线性相关性提高，统计显著性更强。然而，模型仍存在一些系统性偏差，未来的研究可以进一步调整其他参数，如源强和地形数据，以进一步提高模型的准确性。

代码示例：生成验证报告

import pandas as pd

from jinja2 import Environment, FileSystemLoader



# 读取评估结果

initial_rmse = 10.5

initial_mae = 7.2

initial_r = 0.75

initial_p_value = 0.03



new_rmse = 8.3

new_mae = 5.4

new_r = 0.85

new_p_value = 0.01



# 设置 Jinja2 模板环境

env = Environment(loader=FileSystemLoader('.'))

template = env.get_template('report_template.md')



# 渲染模板

report = template.render(

    initial_rmse=initial_rmse,

    initial_mae=initial_mae,

    initial_r=initial_r,

    initial_p_value=initial_p_value,

    new_rmse=new_rmse,

    new_mae=new_mae,

    new_r=new_r,

    new_p_value=new_p_value

)



# 保存验证报告

with open('verification_report.md', 'w') as f:

    f.write(report)

6.7 附件

  • 观测数据：cleaned_observation_data.csv

  • 初始预测数据：predicted_pm25.csv

  • 改进后的预测数据：new_predicted_pm25.csv

  • 气象数据：met_data.csv

  • 调整后的气象数据：adjusted_met_data.csv

通过以上步骤，我们详细记录了模式验证与评估的过程，为后续研究提供了宝贵的数据和经验。

评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

kkchenjj

你的鼓励是我最大的动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值