RNN時間序列預測空氣質量

摘要： 關於RNN可以參考這篇文章博文地址 背景 時間序列是一個比較場景的問題， 獲取資料 從UIC下載學習庫下載地址下載其中的csv檔案 圖片.png 其各個...

關於RNN可以參考這篇文章博文地址

背景

時間序列是一個比較場景的問題，

獲取資料

從UIC下載學習庫下載地址下載其中的csv檔案

其各個資料特徵如下

處理資料

先匯入需要用到的模組

這裡用到了keras的API關於keras的使用可以看官方教程點選前往

from tensorflow import keras
from tensorflow.keras import layers
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt

匯入資料

data = pd.read_csv('D:\machine\project/testbeijinpm2.5/PRSA_data_2010.1.1-2014.12.31.csv')

檢視下列資訊看是否匯入成功

data.columns

檢視資料情況可以用data.head在jupyter會顯示前5行

data.head()

pm2.5就是我們需要的label通過下面操作檢視是否有空的資料

data['pm2.5'].isna().sum()

可以看到有2067行資料包含空的pm2.5

對於這些空的資料，如果採取刪除的動作，那樣會導致資料不連續，可能對於訓練的結果產生影響，這裡採用的處理方式是講空值處理成它最近的點的值，這樣可以保證資料的連續性和前後關聯性，有不至於有太大誤差

執行如下操作處理

data=data.iloc[24:].fillna(method='ffill')

可以觀察到時間的年月日小時分佈在不同列

它們應該分佈在同一列屬於同一個特徵

這裡匯入python的時間模組進行處理

import datetime

將資料 中的year，month，day，hour整合成一條資料賦值給data['time']

#時間組合
data['time']=data.apply(lambdax: datetime.datetime(year=x['year'],
month=x['month'],
day=x['day'],
hour=x['hour']),
axis=1
)

刪除原本的year，month，day，no列

data.drop(columns=['year','month','day','hour','No'],inplace=True)

檢視下資料的情況

設定成按照時間索引

#設定成按時間索引
data=data.set_index('time')

從資料中可以看到cbwd不是一些數值標量，這需要對它進行轉換，先檢視cbwd有哪些內容

#找到cbwd有哪些內容
data.cbwd.unique()

將這幾個標識編碼

#把風向編碼化
data=data.join(pd.get_dummies(data.cbwd))

編碼化完成後刪除cbwd列

del data['cbwd']

查看下錶格，可以看到4個標識已經編碼化了

data.tail()

可以下面程式碼檢視最後1000次pm2.5變化情況

data['pm2.5'][-1000:].plot()

最後1000次溫度變化情況

data['TEMP'][-1000:].plot()

可以通過檢視列數檢視特徵數量

#特徵數量
len(data.columns)

時間序列預測是觀測過去多久的值，來預測未來多久的值，因此這裡定義兩個變數seq_length代表過去5×24小時的資料，來預測未來24小時的資料。

seq_length = 5*24
delay = 24

然後我們就可以採用資料，以seq_length和delay來取樣

#取樣的資料
data_ = []
for t in range(len(data) - seq_length - delay):
data_.append(data.iloc[t:t + seq_length + delay ])

檢視下采樣後的資料形狀

data_[0].shape

可以看到144個樣例，11個特徵

將取樣的資料轉換成np的array

data_ = np.array([df.values for df in data_])
data_.shape

洗牌資料

np.random.shuffle(data_)

特徵切片

#特徵切片
x = data_[:, :5*24, :]
# label
y = data_[:,-1, 0]

檢視下x,y的資料形狀

x.shape
y.shape

劃分80%資料用於訓練集20%用於測試集合

split_b = data_.shape[0]*0.8
split_b = int(split_b)
train_x = x[ :split_b]
train_y = y[ :split_b]
test_x = x[split_b: ]
test_y = y[split_b: ]

資料歸一化，讓不同特徵規範統一範圍，有效利用神經網路的訓練避免引起誤差。

#計算均值
#資料歸一化，
mean = test_x.mean(axis = 0)
std = train_x.std(axis = 0)
train_x = (train_x - mean)/std
test_x = (test_x - mean) / std

搭建神經網路

全連結神經網路方式

batch_size = 128 #每次訓練128個數據

model = keras.Sequential()
model.add(layers.Flatten(input_shape=(train_x.shape[1:])))
model.add(layers.Dense(32,activation='relu'))
model.add(layers.Dense(1))

model.compile(
optimizer='adam',#優化演算法 
loss = 'mse',#均方差 
metrics= ['mae'] #平均絕對誤差
)

history = model.fit(train_x, train_y,
batch_size =batch_size,
epochs = 50,
validation_data=(test_x,test_y)
)

#畫圖看情況
import matplotlib.pyplot as plt
%matplotlib inline

plt.plot(history.epoch,history.history['mean_absolute_error'],c='r') #訓練資料
plt.plot(history.epoch,history.history['val_mean_absolute_error'],c='g') #測試資料

LSTM方式

在訓練過程中可以看到速度明顯要比全連線的方式要好慢上許多，而訓練的效果也比全連結網路的要要很多。

train_x.shape#第一維度大小，第二維度觀測次數，第三維度每次觀測的特徵

這隻有一層的LSTM

model = keras.Sequential()
model.add(layers.LSTM(32,input_shape=(120,11)))
model.add(layers.Dense(1))

編譯選擇模型的優化方式

model.compile(
optimizer='adam',#優化演算法 
loss = 'mse',#均方差 
metrics= ['mae'] #平均絕對誤差
)

history = model.fit(train_x, train_y,
batch_size =batch_size,
epochs = 200,
validation_data=(test_x,test_y)
)

plt.plot(history.epoch,history.history.get('mean_absolute_error'),'y',label='Training loss')
plt.plot(history.epoch,history.history.get('val_mean_absolute_error'),'b',label= 'Test loss')
plt.legend()

LSTM層優化

新增多層優化擬合能力，需要新增

return_sequences，最後一層不需要

可以觀察到訓練速度更加緩慢，但訓練效果卻是跨越式的進展。

model = keras.Sequential()
model.add(layers.LSTM(32,input_shape=(120,11),return_sequences=True))
model.add(layers.LSTM(32,return_sequences=True))
model.add(layers.LSTM(32,return_sequences=True))
model.add(layers.LSTM(32))
model.add(layers.Dense(1))

#在訓練過程中降低學習速率，加快訓練速度同時增加訓練的效果
#連續的3個val_LOSS沒有降低則降低學習速率0.3,最小的學習速率為0.000001
lr_reduce = keras.callbacks.ReduceLROnPlateau('val_loss',patience=3,factor=0.3,min_lr=0.00001)

model.compile(
optimizer='adam',#優化演算法 
loss = 'mse',#均方差 
metrics= ['mae'] #平均絕對誤差
)

history = model.fit(train_x, train_y,
batch_size =batch_size,
epochs = 200,
callbacks = [lr_reduce],
validation_data=(test_x,test_y)
)

plt.plot(history.epoch,history.history.get('mean_absolute_error'),'y',label='Training loss')
plt.plot(history.epoch,history.history.get('val_mean_absolute_error'),'b',label= 'Test loss')
plt.legend()

可以很明顯的看到效果的提升