2024/12/1周报

摘要

本周使用LSTM，Transformer等基础模型对伶俐污水处理厂的数据进行预测，预测出水水质的一些参数，效果还可以，未来将用其他厂更多更合理的数据进行预测尝试。此外，还对基于人工智能技术在污水处理中的应用数学模型和变量进行系统性表达，以此全面描述污水处理过程的多维效率，并在实际运行中应用智能优化工具提升处理效果和效率。

Abstract

This week, LSTM, Transformer and other basic models are used to predict the data of Lingli sewage treatment plant, and some parameters of effluent water quality are predicted, and the effect is okay, and more and more reasonable data from other plants will be used to predict in the future. In addition, the mathematical models and variables based on artificial intelligence technology in wastewater treatment are systematically expressed, so as to comprehensively describe the multi-dimensional efficiency of the sewage treatment process, and intelligent optimization tools are applied to improve the treatment effect and efficiency in actual operation.

COD含氧量预测

LSTM

以下首先是使用LSTM模型对出水水质COD含氧量进行预测，该代码实现了一个基于优化LSTM的时间序列预测模型，用于预测"出水水质COD801"。数据预处理包括异常值截断、对数变换、时间特征提取、移动平均与滞后特征生成，并使用标准化处理。模型通过多层LSTM、层归一化和Dropout增强泛化能力，采用MSE损失函数和Adam优化器，动态调整学习率进行训练。最终，通过测试集验证模型性能，评估指标包括RMSE、MAE、MAPE和R²，并通过可视化实际值与预测值对比分析预测效果：

import numpy as np
import pandas as pd
import torch
import torch.nn as nn
from sklearn.preprocessing import MinMaxScaler
from sklearn.model_selection import train_test_split
from torch.utils.data import DataLoader, Dataset
from sklearn.metrics import mean_squared_error, mean_absolute_error, r2_score
import matplotlib.pyplot as plt

# 设置随机种子
def set_random_seeds(seed):
    import random
    random.seed(seed)
    np.random.seed(seed)
    torch.manual_seed(seed)
    torch.cuda.manual_seed_all(seed)
    torch.backends.cudnn.deterministic = True
    torch.backends.cudnn.benchmark = False

# 设置随机种子
seed = 42
set_random_seeds(seed)

# 检查 GPU 是否可用
device = torch.device("cuda" if torch.cuda.is_available() else "cpu")
print(f"Using device: {device}")

# 加载数据
file_path = "C:/Users/52253/Desktop/伶俐园区/cod/COD23.1.xlsx"
data = pd.read_excel(file_path, sheet_name='Sheet1')

# 数据预处理
data['时间'] = pd.to_datetime(data['时间'])
data.set_index('时间', inplace=True)

# 处理目标变量异常值（如大于 95% 分位数的值进行截断）
cod801_upper_limit = data['出水水质COD801'].quantile(0.95)
data['出水水质COD801'] = np.clip(data['出水水质COD801'], None, cod801_upper_limit)

# 对目标变量进行对数变换，减少分布偏态
data['出水水质COD801'] = np.log1p(data['出水水质COD801'])

# 添加时间特征
data['hour'] = data.index.hour
data['day_of_week'] = data.index.dayofweek
data['month'] = data.index.month
data['hour_sin'] = np.sin(2 * np.pi * data['hour'] / 24)
data['hour_cos'] = np.cos(2 * np.pi * data['hour'] / 24)

# 添加移动平均特征
data['rolling_mean_3'] = data['出水水质COD801'].rolling(window=3).mean()
data['rolling_mean_7'] = data['出水水质COD801'].rolling(window=7).mean()

# 为 pH 值添加移动平均特征
data['ph201_rolling_mean_3'] = data['进水PH201'].rolling(window=3).mean()
data['ph801_rolling_mean_3'] = data['出水PH801'].rolling(window=3).mean()

# 添加滞后特征
for lag in range(1, 4):
    data[f'lag_{lag}'] = data['出水水质COD801'].shift(lag)
    data[f'ph201_lag_{lag}'] = data['进水PH201'].shift(lag)
    data[f'ph801_lag_{lag}'] = data['出水PH801'].shift(lag)

# 删除包含 NaN 的行（由于滞后特征引入的 NaN）
data.dropna(inplace=True)

# 提取特征与目标
features = ['进水水质COD201', '溶氧仪DOT301', '溶氧仪DOT302', '出水总磷TP801', '出水氨氮NHT801', '进水氨氮',
            'hour_sin', 'hour_cos', 'day_of_week', 'month', 'rolling_mean_3', 'rolling_mean_7',
            'lag_1', 'lag_2', 'lag_3', '进水PH201', '出水PH801', 'ph201_rolling_mean_3',
            'ph801_rolling_mean_3', 'ph201_lag_1', 'ph201_lag_2', 'ph201_lag_3',
            'ph801_lag_1', 'ph801_lag_2', 'ph801_lag_3']
target = ['出水水质COD801']

# 数据标准化
scaler_x = MinMaxScaler()
scaler_y = MinMaxScaler()

X = scaler_x.fit_transform(data[features])
y = scaler_y.fit_transform(data[target])

# 划分训练集和测试集
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, shuffle=False, random_state=seed)

# 定义数据集类
class TimeSeriesDataset(Dataset):
    def __init__(self, X, y, seq_length):
        self.X = X
        self.y = y
        self.seq_length = seq_length

    def __len__(self):
        return len(self.X) - self.seq_length

    def __getitem__(self, index):
        return (
            torch.tensor(self.X[index:index+self.seq_length], dtype=torch.float32),
            torch.tensor(self.y[index+self.seq_length], dtype=torch.float32)
        )

# 序列长度
seq_length = 20

train_dataset = TimeSeriesDataset(X_train, y_train, seq_length)
test_dataset = TimeSeriesDataset(X_test, y_test, seq_length)

# 固定数据加载顺序
def worker_init_fn(worker_id):
    np.random.seed(seed + worker_id)

train_loader = DataLoader(train_dataset, batch_size=32, shuffle=True, worker_init_fn=worker_init_fn)
test_loader = DataLoader(test_dataset, batch_size=32, shuffle=False, worker_init_fn=worker_init_fn)

# 定义优化后的 LSTM 模型
class OptimizedLSTM(nn.Module):
    def __init__(self, input_size, hidden_size, output_size, num_layers=1, dropout=0.3):
        super(OptimizedLSTM, self).__init__()
        self.hidden_size = hidden_size
        self.num_layers = num_layers

        self.lstm = nn.LSTM(input_size, hidden_size, num_layers, batch_first=True, dropout=dropout)
        self.layer_norm = nn.LayerNorm(hidden_size)  # 添加层归一化
        self.fc = nn.Linear(hidden_size, output_size)

    def forward(self, x):
        h0 = torch.zeros(self.num_layers, x.size(0), self.hidden_size).to(x.device)
        c0 = torch.zeros(self.num_layers, x.size(0), self.hidden_size).to(x.device)
        out, _ = self.lstm(x, (h0, c0))
        out = self.layer_norm(out[:, -1, :])  # 添加归一化
        out = self.fc(out)
        return out

# 模型参数
input_size = len(features)
hidden_size = 128
output_size = 1
num_layers = 3
learning_rate = 0.0005
num_epochs = 150
dropout = 0.3

# 初始化模型、损失函数和优化器
model = OptimizedLSTM(input_size, hidden_size, output_size, num_layers, dropout).to(device)
criterion = nn.MSELoss()
optimizer = torch.optim.Adam(model.parameters(), lr=learning_rate)

# 学习率调度器
scheduler = torch.optim.lr_scheduler.StepLR(optimizer, step_size=30, gamma=0.5)

# 训练模型
model.train()
for epoch in range(num_epochs):
    train_loss = 0.0
    for X_batch, y_batch in train_loader:
        X_batch, y_batch = X_batch.to(device), y_batch.to(device)
        optimizer.zero_grad()
        outputs = model(X_batch)
        loss = criterion(outputs, y_batch)
        loss.backward()
        optimizer.step()
        train_loss += loss.item()
    scheduler.step()  # 更新学习率

    print(f'Epoch {epoch+1}/{num_epochs}, Loss: {train_loss/len(train_loader):.4f}')

# 测试模型
model.eval()
predictions = []
actuals = []
with torch.no_grad():
    for X_batch, y_batch in test_loader:
        X_batch, y_batch = X_batch.to(device), y_batch.to(device)
        outputs = model(X_batch)
        predictions.append(outputs.cpu().numpy())
        actuals.append(y_batch.cpu().numpy())

# 数据逆归一化
predictions = np.concatenate(predictions, axis=0)
actuals = np.concatenate(actuals, axis=0)

predictions = scaler_y.inverse_transform(predictions)
actuals = scaler_y.inverse_transform(actuals)

# 计算评估指标
rmse = np.sqrt(mean_squared_error(actuals, predictions))
mae = mean_absolute_error(actuals, predictions)
mape = np.mean(np.abs((actuals - predictions) / actuals)) * 100
r2 = r2_score(actuals, predictions)

# 打印评估指标
print(f"Evaluation Metrics:")
print(f"Root Mean Squared Error (RMSE): {rmse:.4f}")
print(f"Mean Absolute Error (MAE): {mae:.4f}")
print(f"Mean Absolute Percentage Error (MAPE): {mape:.2f}%")
print(f"R^2 Score: {r2:.4f}")

# 可视化实际值与预测值对比
plt.figure(figsize=(12, 6))
plt.plot(actuals, label='Actual', alpha=0.7)
plt.plot(predictions, label='Predicted', alpha=0.7)
plt.legend()
plt.title('Optimized LSTM Prediction of 出水水质COD801')
plt.show()

GRU

以下使用将LSTM替换成GRU：

class OptimizedGRU(nn.Module):
    def __init__(self, input_size, hidden_size, output_size, num_layers=1, dropout=0.3):
        super(OptimizedGRU, self).__init__()
        self.hidden_size = hidden_size
        self.num_layers = num_layers
        self.gru = nn.GRU(input_size, hidden_size, num_layers, batch_first=True, dropout=dropout)
        self.layer_norm = nn.LayerNorm(hidden_size)
        self.fc = nn.Linear(hidden_size, output_size)

    def forward(self, x):
        h0 = torch.zeros(self.num_layers, x.size(0), self.hidden_size).to(x.device)
        out, _ = self.gru(x, h0)
        out = self.layer_norm(out[:, -1, :])
        out = self.fc(out)
        return out

其他模型

还尝试使用了如下一些模型：

class CNNTransformerModel(nn.Module):
    def __init__(self, input_size, seq_length, d_model=128, n_heads=4, e_layers=3, d_ff=256, dropout=0.3):
        super(CNNTransformerModel, self).__init__()
        self.seq_length = seq_length
        self.conv1 = nn.Conv1d(in_channels=input_size, out_channels=64, kernel_size=3, padding=1)
        self.conv2 = nn.Conv1d(in_channels=64, out_channels=d_model, kernel_size=3, padding=1)
        self.relu = nn.ReLU()
        self.encoder = nn.TransformerEncoder(
            nn.TransformerEncoderLayer(d_model=d_model, nhead=n_heads, dim_feedforward=d_ff, dropout=dropout),
            num_layers=e_layers
        )
        self.projector = nn.Linear(d_model, 1)
        self.dropout = nn.Dropout(dropout)

# 定义 CNN-LSTM-自注意力模型
class CLATT_SelfAttention(nn.Module):
    def __init__(self, input_size, hidden_size, output_size, cnn_kernel_size=3, lstm_layers=1):
        super(CLATT_SelfAttention, self).__init__()
        self.cnn = nn.Sequential(
            nn.Conv1d(input_size, hidden_size, kernel_size=cnn_kernel_size, padding=1),
            nn.ReLU(),
            nn.BatchNorm1d(hidden_size)
        )
        self.lstm = nn.LSTM(hidden_size, hidden_size, num_layers=lstm_layers, batch_first=True)
        self.attention = SelfAttention(hidden_size)
        self.fc = nn.Linear(hidden_size, output_size)

    def forward(self, x):
        x = x.permute(0, 2, 1)
        cnn_out = self.cnn(x)
        cnn_out = cnn_out.permute(0, 2, 1)
        lstm_out, _ = self.lstm(cnn_out)
        attn_out, attn_weights = self.attention(lstm_out)
        context = attn_out.mean(dim=1)
        output = self.fc(context)
        return output

实验结果

以下是测试的一些结果：

文献阅读

题目

时间序列预测中的深度学习技术：在汇率中的应用

引言

在特定的汇率下，预测金融序列在决策、投资组合配置、风险管理和衍生品定价问题中起着关键作用。近年来，机器学习模型由于其处理非线性关系的能力而被应用于时间序列预测。这项研究评估了神经网络模型，如MLP，LSTM，GRU和拟议的AE-BiGRU-Swish预测汇率序列。使用深度学习技术调整了这些模型的参数，比较了神经网络模型预测汇率序列的性能与传统的方法，如ARIMA和ARFIMA模型的时间范围为1，7，14，21和30天的提前预测。结果表明，混合神经网络模型优于其他方法来预测汇率。

本文中，作者提出了使用非线性模型，特别是非递归多层感知器神经网络，长短期记忆，门控递归单元网络和AutoEncoder双向门控递归单元多步提前（1，7，14，21和30天）预测汇率，试图验证这些模型是否可以捕捉汇率序列中的非线性关系。
为了评估非线性模型的性能，作者还调整了自回归积分移动平均（ARIMA）和自回归分数积分移动平均（ARFIMA）模型，并根据MAE，MAPE，MSE，RMSE，ACCURACY和CE等误差指标比较结果。

模型与方法

多层感知器神经网络MLP

MLP是静态的，模型输出表示如下：

其中，φ（·）是激活函数， ∑是输入的线性聚合函数，由权重W加权。
MLP训练需要确定通过误差反向传播算法（BP）获得的权重向量w。
BP算法基于误差校正规则，这意味着应用迭代的基于梯度的方法，该方法寻求最小化估计输出y和期望输出y。信息处理有两个步骤：向前和向后。前向步骤通过网络传播输入，生成输出（传播），后向步骤包括基于误差 计算成本函数的梯度，通过以下公式获得：

BP算法使用基于向下梯度的搜索技术：

最后用下式更新W：

其中η ∈ [0，1]是学习率或自适应步长， ∇(·)是梯度算子。

LSTM

LSTM单元由输入信号Xi和输出yi-1反馈，LSTM单元在所有阶段共享相同的参数。

LSTM使用误差反向传播时间算法（BPTT）进行估计，该方法旨在最小化模型提供的输出与期望输出之间的成本函数：

BPTT计算成本函数关于向量w的导数：

使用以下公式更新W：

GRU

门控递归单元（GRU）是LSTM单元的简化版本，以相同的方式工作。
GRU是一个动态模型，它提供递归映射并存储先前状态的信息，以影响当前时刻的预测。GRU与其他网络的不同之处在于，它将细胞的整个状态暴露给网络中的其他单元，而LSTM细胞控制其内存内容的暴露。GRU单元有两个因子：更新因子zt，即遗忘门和输入门在单个门上的连接点，以及复位因子rt。

AE-BiGRU-Swish神经网络

本文提出混合模型AE-BiGRU-Swish，结合了监督和无监督学习。该模型与其他模型一起进行了跨不同时间范围（长达30天）的多步预测测试。自动编码器（AE）是一种利用无监督学习的神经网络模型。它由两个主要阶段组成：编码器和解码器。在编码器中，输入数据被转换为较低的维度，而解码器的目标是基于编码层提供的压缩表示来重新创建原始输入。在训练期间，目标是最大化原始输入数据的重建并最小化重建误差（等式9）。编码器（CODEC）包括两个连续的双向GRU（BiGRU），然后是具有Swish激活函数的dense层。另一方面，解码器（D）由dense层和两个BiGRU层组成。解码器的输出是高光谱向量，表示为D(Ξ(xi))，具有与原始向量相同的维度。

所提出的AE-BiGRU-Swish具有两个明显的特征：AE结构提供的非线性低维表示能力和使用Swish激活函数训练的BiGRU结构进行的系列预测。

实验

实验过程

本研究旨在评估MLP，LSTM，GRU和AE-BiGRU-Swish神经网络模型在预测雷亚尔/美元和欧元/美元汇率的1，7，14，21和30天的时间范围内的表现。我们从yahoo finance中提取了这些系列，这些系列公开表示为汇率BRL=X（R$/US$）和EURUSD=X（Euro/US$）。这两个系列分别是2003年12月31日至2021年5月4日的雷亚尔/美元数据和2003年12月31日至2021年5月5日的欧元/美元数据的日平均值，分别为4525和4526次观测。
图1（a）和（B）说明了雷亚尔/美元和欧元/美元汇率的演变情况。

神经网络模型的参数和超参数如下表所示：

评估标准

平均绝对误差（MAE）、平均绝对百分比误差（MAPE）、均方误差（MSE）、均方根误差（RMSE）、方向精密度（ACCURACY）和Nash-Sutcliffe效率（CE），定义如下：

实验结果

比较了不同模型对30步前的预测结果如下表所示：

总结

下周将继续预测出水水质的其他参数，例如出水总氮、总磷、氨氮等，并扩展更多模型。