2020年发生的新冠肺炎疫情灾难给各国人民的生活造成了诸多方面的影响。
有的同学是收入上的，有的同学是感情上的，有的同学是心理上的，还有的同学是体重上的。
本文基于中国2020年3月之前的疫情数据，建立时间序列RNN模型，对中国的新冠肺炎疫情结束时间进行预测。

  1. import os
  2. import datetime
  3. import importlib 
  4. import torchkeras
  6. #打印时间
  7. def printbar():
  8.     nowtime = datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S')
  9.     print("\n"+"=========="*8 + "%s"%nowtime)
  11. #mac系统上pytorch和matplotlib在jupyter中同时跑需要更改环境变量
  12. os.environ["KMP_DUPLICATE_LIB_OK"]="TRUE"

准备数据

  1. import numpy as np
  2. import pandas as pd 
  3. import matplotlib.pyplot as plt
  1. %matplotlib inline
  2. %config InlineBackend.figure_format = 'svg'
  4. df = pd.read_csv("../data/covid-19.csv",sep = "\t")
  5. df.plot(x = "date",y = ["confirmed_num","cured_num","dead_num"],figsize=(10,6))
  6. plt.xticks(rotation=60);

image.png

  1. dfdata = df.set_index("date")
  2. dfdiff = dfdata.diff(periods=1).dropna()
  3. dfdiff = dfdiff.reset_index("date")
  5. dfdiff.plot(x = "date",y = ["confirmed_num","cured_num","dead_num"],figsize=(10,6))
  6. plt.xticks(rotation=60)
  7. dfdiff = dfdiff.drop("date",axis = 1).astype("float32")

image.png

  1. dfdiff.head()

image.png

下面我们通过继承torch.utils.data.Dataset实现自定义时间序列数据集。
torch.utils.data.Dataset是一个抽象类，用户想要加载自定义的数据只需要继承这个类，并且覆写其中的两个方法即可：

  • __len__:实现len(dataset)返回整个数据集的大小。
  • __getitem__:用来获取一些索引的数据，使dataset[i]返回数据集中第i个样本。

不覆写这两个方法会直接返回错误。

  1. import torch 
  2. from torch import nn 
  3. from torch.utils.data import Dataset,DataLoader,TensorDataset
  6. #用某日前8天窗口数据作为输入预测该日数据
  7. WINDOW_SIZE = 8
  9. class Covid19Dataset(Dataset):
  11.     def __len__(self):
  12.         return len(dfdiff) - WINDOW_SIZE
  14.     def __getitem__(self,i):
  15.         x = dfdiff.loc[i:i+WINDOW_SIZE-1,:]
  16.         feature = torch.tensor(x.values)
  17.         y = dfdiff.loc[i+WINDOW_SIZE,:]
  18.         label = torch.tensor(y.values)
  19.         return (feature,label)
  21. ds_train = Covid19Dataset()
  23. #数据较小，可以将全部训练数据放入到一个batch中，提升性能
  24. dl_train = DataLoader(ds_train,batch_size = 38)

定义模型

使用Pytorch通常有三种方式构建模型：

  1. 使用nn.Sequential按层顺序构建模型；
  2. 继承nn.Module基类构建自定义模型；
  3. 继承nn.Module基类构建模型并辅助应用模型容器进行封装。

此处选择第二种方式构建模型。由于接下来使用类形式的训练循环，我们进一步将模型封装成torchkeras中的Model类来获得类似Keras中高阶模型接口的功能。Model类实际上继承自nn.Module类。

  1. import torch
  2. from torch import nn 
  3. import importlib 
  4. import torchkeras 
  6. torch.random.seed()
  8. class Block(nn.Module):
  9.     def __init__(self):
  10.         super(Block,self).__init__()
  12.     def forward(self,x,x_input):
  13.         x_out = torch.max((1+x)*x_input[:,-1,:],torch.tensor(0.0))
  14.         return x_out
  16. class Net(nn.Module):
  17.     def __init__(self):
  18.         super(Net, self).__init__()
  19.         # 3层lstm
  20.         self.lstm = nn.LSTM(input_size = 3,hidden_size = 3,num_layers = 5,batch_first = True)
  21.         self.linear = nn.Linear(3,3)
  22.         self.block = Block()
  24.     def forward(self,x_input):
  25.         x = self.lstm(x_input)[0][:,-1,:]
  26.         x = self.linear(x)
  27.         y = self.block(x,x_input)
  28.         return y
  30. net = Net()
  31. model = torchkeras.Model(net)
  32. print(model)
  34. model.summary(input_shape=(8,3),input_dtype = torch.FloatTensor)

image.png

训练模型

训练Pytorch通常需要用户编写自定义训练循环，训练循环的代码风格因人而异。
有3类典型的训练循环代码风格：

  1. 脚本形式训练循环；
  2. 函数形式训练循环；
  3. 类形式训练循环。

此处介绍一种类形式的训练循环。我们仿照Keras定义了一个高阶的模型接口Model,实现 fit, validate，predict, summary 方法，相当于用户自定义高阶API。
注：循环神经网络调试较为困难，需要设置多个不同的学习率多次尝试，以取得较好的效果。

  1. def mspe(y_pred,y_true):
  2.     err_percent = (y_true - y_pred)**2/(torch.max(y_true**2,torch.tensor(1e-7)))
  3.     return torch.mean(err_percent)
  5. model.compile(loss_func = mspe,optimizer = torch.optim.Adagrad(model.parameters(),lr = 0.1))
  1. dfhistory = model.fit(100,dl_train,log_step_freq=10)

image.png

评估模型

评估模型一般要设置验证集或者测试集，由于此例数据较少，我们仅仅可视化损失函数在训练集上的迭代情况。

  1. %matplotlib inline
  2. %config InlineBackend.figure_format = 'svg'
  4. import matplotlib.pyplot as plt
  6. def plot_metric(dfhistory, metric):
  7.     train_metrics = dfhistory[metric]
  8.     epochs = range(1, len(train_metrics) + 1)
  9.     plt.plot(epochs, train_metrics, 'bo--')
  10.     plt.title('Training '+ metric)
  11.     plt.xlabel("Epochs")
  12.     plt.ylabel(metric)
  13.     plt.legend(["train_"+metric])
  14.     plt.show()
  16. plot_metric(dfhistory,"loss")

image.png

使用模型

此处我们使用模型预测疫情结束时间，即 新增确诊病例为0 的时间。

  1. #使用dfresult记录现有数据以及此后预测的疫情数据
  2. dfresult = dfdiff[["confirmed_num","cured_num","dead_num"]].copy()
  3. dfresult.tail()

image.png

  1. #预测此后500天的新增走势,将其结果添加到dfresult中
  2. for i in range(500):
  3.     arr_input = torch.unsqueeze(torch.from_numpy(dfresult.values[-38:,:]),axis=0)
  4.     arr_predict = model.forward(arr_input)
  6.     dfpredict = pd.DataFrame(torch.floor(arr_predict).data.numpy(),
  7.                 columns = dfresult.columns)
  8.     dfresult = dfresult.append(dfpredict,ignore_index=True)
  1. dfresult.query("confirmed_num==0").head()
  3. # 第50天开始新增确诊降为0，第45天对应3月10日，也就是5天后，即预计3月15日新增确诊降为0
  4. # 注：该预测偏乐观

image.png

  1. dfresult.query("cured_num==0").head()
  3. # 第365天开始新增治愈降为0，即大概1年后。
  4. # 注: 该预测偏悲观，并且存在问题，如果将每天新增治愈人数加起来，将超过累计确诊人数。

image.png

  1. dfresult.query("dead_num==0").head()
  3. # 第54天开始新增确诊降为0
  4. # 注：该预测偏乐观

image.png

保存模型

推荐使用保存参数方式保存Pytorch模型。

  1. print(model.net.state_dict().keys())

image.png

  1. # 保存模型参数
  3. torch.save(model.net.state_dict(), "model/1_4_model_parameter.pkl")
  5. net_clone = Net()
  6. net_clone.load_state_dict(torch.load("model/1_4_model_parameter.pkl"))
  7. model_clone = torchkeras.Model(net_clone)
  8. model_clone.compile(loss_func = mspe)
  10. # 评估模型
  11. model_clone.evaluate(dl_train)

image.png