文章最前： 我是Octopus，這個名字來源于我的中文名–章魚；我熱愛編程、熱愛算法、熱愛開源。所有源碼在我的個人github
；這博客是記錄我學習的點點滴滴，如果您對 Python、Java、AI、算法有興趣，可以關注我的動態，一起學習，共同進步。
2020年發生的新冠肺炎疫情災難給各國人民的生活造成了諸多方面的影響。
有的同學是收入上的，有的同學是感情上的，有的同學是心理上的，還有的同學是體重上的。
本文基于中國2020年3月之前的疫情數據，建立時間序列RNN模型，對中國的新冠肺炎疫情結束時間進行預測。
import torch
print("torch.__version__ = ", torch.__version__)torch.__version__ = 2.0.1
公眾號 算法美食屋 回復關鍵詞：pytorch， 獲取本項目源碼和所用數據集百度云盤下載鏈接。
import os#mac系統上pytorch和matplotlib在jupyter中同時跑需要更改環境變量
os.environ["KMP_DUPLICATE_LIB_OK"]="TRUE"
一，準備數據
本文的數據集取自tushare，獲取該數據集的方法參考了以下文章。
《https://zhuanlan.zhihu.com/p/109556102》
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt%matplotlib inline
%config InlineBackend.figure_format = 'svg'df = pd.read_csv("./eat_pytorch_datasets/covid-19.csv",sep = "\t")
df.plot(x = "date",y = ["confirmed_num","cured_num","dead_num"],figsize=(10,6))
plt.xticks(rotation=60);
dfdata = df.set_index("date")
dfdiff = dfdata.diff(periods=1).dropna()
dfdiff = dfdiff.reset_index("date")dfdiff.plot(x = "date",y = ["confirmed_num","cured_num","dead_num"],figsize=(10,6))
plt.xticks(rotation=60)
dfdiff = dfdiff.drop("date",axis = 1).astype("float32")
dfdiff.head()
| confirmed_num | cured_num | dead_num | |
|---|---|---|---|
| 0 | 457.0 | 4.0 | 16.0 |
| 1 | 688.0 | 11.0 | 15.0 |
| 2 | 769.0 | 2.0 | 24.0 |
| 3 | 1771.0 | 9.0 | 26.0 |
| 4 | 1459.0 | 43.0 | 26.0 |
下面我們通過繼承torch.utils.data.Dataset實現自定義時間序列數據集。
torch.utils.data.Dataset是一個抽象類，用戶想要加載自定義的數據只需要繼承這個類，并且覆寫其中的兩個方法即可：
__len__:實現len(dataset)返回整個數據集的大小。__getitem__:用來獲取一些索引的數據，使dataset[i]返回數據集中第i個樣本。
不覆寫這兩個方法會直接返回錯誤。
import torch
from torch import nn
from torch.utils.data import Dataset,DataLoader,TensorDataset#用某日前8天窗口數據作為輸入預測該日數據
WINDOW_SIZE = 8class Covid19Dataset(Dataset):def __len__(self):return len(dfdiff) - WINDOW_SIZEdef __getitem__(self,i):x = dfdiff.loc[i:i+WINDOW_SIZE-1,:]feature = torch.tensor(x.values)y = dfdiff.loc[i+WINDOW_SIZE,:]label = torch.tensor(y.values)return (feature,label)ds_train = Covid19Dataset()#數據較小，可以將全部訓練數據放入到一個batch中，提升性能
dl_train = DataLoader(ds_train,batch_size = 38)for features,labels in dl_train:break #dl_train同時作為驗證集
dl_val = dl_train二，定義模型
使用Pytorch通常有三種方式構建模型：使用nn.Sequential按層順序構建模型，繼承nn.Module基類構建自定義模型，繼承nn.Module基類構建模型并輔助應用模型容器進行封裝。
此處選擇第二種方式構建模型。
import torch
from torch import nn
import importlib
import torchkeras torch.random.seed()class Block(nn.Module):def __init__(self):super(Block,self).__init__()def forward(self,x,x_input):x_out = torch.max((1+x)*x_input[:,-1,:],torch.tensor(0.0))return x_outclass Net(nn.Module):def __init__(self):super(Net, self).__init__()# 3層lstmself.lstm = nn.LSTM(input_size = 3,hidden_size = 3,num_layers = 5,batch_first = True)self.linear = nn.Linear(3,3)self.block = Block()def forward(self,x_input):x = self.lstm(x_input)[0][:,-1,:]x = self.linear(x)y = self.block(x,x_input)return ynet = Net()
print(net)Net((lstm): LSTM(3, 3, num_layers=5, batch_first=True)(linear): Linear(in_features=3, out_features=3, bias=True)(block): Block()
)
Net((lstm): LSTM(3, 3, num_layers=5, batch_first=True)(linear): Linear(in_features=3, out_features=3, bias=True)(block): Block()
)
from torchkeras import summary
summary(net,input_data=features);
--------------------------------------------------------------------------
Layer (type) Output Shape Param #
==========================================================================
LSTM-1 [-1, 8, 3] 480
Linear-2 [-1, 3] 12
Block-3 [-1, 3] 0
==========================================================================
Total params: 492
Trainable params: 492
Non-trainable params: 0
--------------------------------------------------------------------------
Input size (MB): 0.000069
Forward/backward pass size (MB): 0.000229
Params size (MB): 0.001877
Estimated Total Size (MB): 0.002174
--------------------------------------------------------------------------
三，訓練模型
訓練Pytorch通常需要用戶編寫自定義訓練循環，訓練循環的代碼風格因人而異。
有3類典型的訓練循環代碼風格：腳本形式訓練循環，函數形式訓練循環，類形式訓練循環。
此處我們通過引入torchkeras庫中的KerasModel工具來訓練模型，無需編寫自定義循環。
torchkeras詳情: https://github.com/lyhue1991/torchkeras
注：循環神經網絡調試較為困難，需要設置多個不同的學習率多次嘗試，以取得較好的效果。
from torchmetrics.regression import MeanAbsolutePercentageErrordef mspe(y_pred,y_true):err_percent = (y_true - y_pred)**2/(torch.max(y_true**2,torch.tensor(1e-7)))return torch.mean(err_percent)net = Net()
loss_fn = mspe
metric_dict = {"mape":MeanAbsolutePercentageError()}optimizer = torch.optim.Adam(net.parameters(), lr=0.03)
lr_scheduler = torch.optim.lr_scheduler.StepLR(optimizer, step_size=10, gamma=0.0001)
from torchkeras import KerasModel
model = KerasModel(net,loss_fn = loss_fn,metrics_dict= metric_dict,optimizer = optimizer,lr_scheduler = lr_scheduler) dfhistory = model.fit(train_data=dl_train,val_data=dl_val,epochs=100,ckpt_path='checkpoint',patience=10,monitor='val_loss',mode='min',callbacks=None,plot=True,cpu=True)[0;31m<<<<<< 🐌 cpu is used >>>>>>[0m
18.00% [18/100] [00:02<00:10]
████████████████████100.00% [1/1] [val_loss=0.4363, val_mape=0.5570]
████████████████████100.00% [1/1] [val_loss=0.4363, val_mape=0.5570]
[0;31m<<<<<< val_loss without improvement in 10 epoch,early stopping >>>>>>
[0m
四，評估模型
評估模型一般要設置驗證集或者測試集，由于此例數據較少，我們僅僅可視化損失函數在訓練集上的迭代情況。
model.evaluate(dl_val)100%|█████████████████████████████████| 1/1 [00:00<00:00, 63.91it/s, val_loss=0.384, val_mape=0.505]{'val_loss': 0.38373321294784546, 'val_mape': 0.5048269033432007}
五，使用模型
此處我們使用模型預測疫情結束時間，即 新增確診病例為0 的時間。
#使用dfresult記錄現有數據以及此后預測的疫情數據
dfresult = dfdiff[["confirmed_num","cured_num","dead_num"]].copy()
dfresult.tail()
| confirmed_num | cured_num | dead_num | |
|---|---|---|---|
| 41 | 143.0 | 1681.0 | 30.0 |
| 42 | 99.0 | 1678.0 | 28.0 |
| 43 | 44.0 | 1661.0 | 27.0 |
| 44 | 40.0 | 1535.0 | 22.0 |
| 45 | 19.0 | 1297.0 | 17.0 |
#預測此后1000天的新增走勢,將其結果添加到dfresult中
for i in range(1000):arr_input = torch.unsqueeze(torch.from_numpy(dfresult.values[-38:,:]),axis=0)arr_predict = model.forward(arr_input)dfpredict = pd.DataFrame(torch.floor(arr_predict).data.numpy(),columns = dfresult.columns)dfresult = pd.concat([dfresult,dfpredict],ignore_index=True)
dfresult.query("confirmed_num==0").head()# 第50天開始新增確診降為0，第45天對應3月10日，也就是5天后，即預計3月15日新增確診降為0
# 注：該預測偏樂觀| confirmed_num | cured_num | dead_num | |
|---|---|---|---|
| 50 | 0.0 | 999.0 | 0.0 |
| 51 | 0.0 | 948.0 | 0.0 |
| 52 | 0.0 | 900.0 | 0.0 |
| 53 | 0.0 | 854.0 | 0.0 |
| 54 | 0.0 | 810.0 | 0.0 |
dfresult.query("cured_num==0").head()
# 第137天開始新增治愈降為0，第45天對應3月10日，也就是大概3個月后，即6月12日左右全部治愈。
# 注: 該預測偏悲觀，并且存在問題，如果將每天新增治愈人數加起來，將超過累計確診人數。
| confirmed_num | cured_num | dead_num | |
|---|---|---|---|
| 137 | 0.0 | 0.0 | 0.0 |
| 138 | 0.0 | 0.0 | 0.0 |
| 139 | 0.0 | 0.0 | 0.0 |
| 140 | 0.0 | 0.0 | 0.0 |
| 141 | 0.0 | 0.0 | 0.0 |
六，保存模型
模型權重保存在了model.ckpt_path路徑。
print(model.ckpt_path)
checkpoint
model.load_ckpt('checkpoint') #可以加載權重
