仔細斟酌了幾秒鍾,我終於決定寫下這篇博客,或許應該取一個高大上的名字,比如“Linear Regression to predict PM2.5 using Gradient Descent”,不過回頭看了看慘淡的排名…
言歸正傳,本文僅適合看見該題目無從下手的小白進行參考,用的算法也是入門爆款的線性回歸的梯度下降,也是我用py擼的第一個梯度下降,下面切入正題。
題目鏈接:https://www.kaggle.com/c/ml2019spring-hw1/overview
train.csv因為中文亂碼所以我手動把亂碼(標題)改成了英文,然后把代表地區的那一列刪掉了,因為預測的都是這個地區,這一列沒啥用。
首先,看到該題目,我有很多設想,題目說:給定一個月前20天每個小時的PM2.5的值(當然還有別的參數),然后在一個月最后十天里,每十小時為一組,給出前9個小時,預測第10個小時的PM2.5值。
那么我首先就設想,那我把每個小時都單獨拎出來進行回歸,比如:1月1號的1點,1月2號的1點,1月3號的一點等等,可是問題來了,即便我以這種方式把每一天的每一個小時的PM2.5單獨求出一個回歸方程,那我如何預測后一天的PM2.5呢?-_-|| 然后這個方法就被我pass掉了。而且事實上,我通過繪圖分析,這個方法似乎也不太好。
下圖就是我分析的每天0點的的PM2.5值,可以看出每天同一時刻的PM2.5值差別都挺大的。
然后我接着想,PM2.5是每個小時測的,那么理論上每兩個相鄰間隔下的PM2.5應該差別不會太大,那可不可以把每個小時的PM2.5都做到一個回歸里,即不考慮哪天是哪天,而只以小時進行考慮,盡管這種方式跟上面的思路存在同樣的問題(即我即便將該回歸方程擬合了出來,我又該怎么預測下一時刻的PM2.5值呢?)我還是對前幾天的每個小時的PM2.5的圖像進行了分析
代碼如下,我分析了前10天每天的每個小時的PM2.5的變化趨勢
限於篇幅,這里只粘貼前3天的數據,可以看出每天PM2.5的波動還是有一定規律的
盡管每天的PM2.5值的變化有規律,但是最后提供的預測數據是以10個小時為一個單位的,而且上圖中也可以看出,相鄰小時內PM2.5的值變化可能會很大,比如上午時分上報高峰期時,PM2.5激增。
說了這么多,也沒想出什么能用的思路來,那就先直接來最簡單粗暴的,我同樣把所有數據以10小時為單位分割,用前9小時來預測第10小時,不過為了增加數據量,我采取了一點小策略,下面打個比方:
假設一共測試了12個小時,那么我的數據是:1~ 10小時,2~ 11小時,3~ 12小時,如此如此這般這般,也算是盡可能增加數據量了(但是我並不確定這樣更好,或許這就是使結果偏離正確答案的一個重大因素,因為它盡可能將這個思路的效果發揮到了極致,但這個思路或許(事實上也是)本身就有問題)
說這個思路有問題,是因為接下來我采用的就是非常常規的線性回歸了,hypothesis就是下圖所示了,但這玩意模擬出來的是直線啊,盡管它有10維,但該是直線還是直線,所以它本身就存在一些問題,不過誰叫它最簡單呢,現在最重要的是先做出來嘛,做出來了再講后續優化。
下面以代碼為主了。
首先導入一些工具庫
# initialize
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
%matplotlib inline
import seaborn as sns
然后導入數據,最下面一行打印出測試數據的維數,結果是(240,23),測試數據量確實挺少的
train_data_path = './train.csv'
test_data_path = './test.csv'
train_data = pd.read_csv(train_data_path)
#這里是把不需要的值給去掉,只留下數據
train_data = train_data[(train_data['kind'] == 'PM2.5')].iloc[:, 3:].astype(int)
train_data.shape
上面的train_data經過處理后結構如下:
下面就是初始化一些值,trainX和trainY就是上面講的思路將每10組數據都給抽取出來,而且因為用上面的思路給擴充了訓練集的大小,因此此時trainX一共有5000多行了
trainX = []
trainY = []
wholeVector = []
# 該循環將訓練集中所有數據全部存到一個一維的列表wholeVector中以方便后續操作
for i in range(train_data.shape[0]):
temp = train_data.iloc[i].values
for j in range(len(temp)):
wholeVector.append(temp[j])
for i in range(len(wholeVector) - 9):
tempx = wholeVector[i: i + 9]
tempy = wholeVector[i + 9]
trainX.append(tempx)
trainY.append(tempy)
# 測試集
trainX = np.array(trainX)
# 給trainX的第一列前加入一列全1的列,方便后續矩陣操作(不理解的可以去看看吳恩達的視頻)
trainX = np.c_[np.ones(trainX.shape[0]), trainX]
trainY = np.array(trainY)
定義計算損失值的函數,並將損失值保存下來,用來判斷算法執行情況,也可以判斷learning rate選擇的好不好
def getLossValue(X, Y, theta, loss_history):
''' 計算損失值,並將損失值存入loss_history列表中 '''
tempvec = Y - (theta.T.dot(X.T)).T
loss = tempvec.T.dot(tempvec)
loss_history.append(loss)
return loss
下面就正式開始進行模型的訓練了,畢竟很簡單,關鍵代碼就兩三行
# 學習率
learning_rate = 0.00000001
# 參數
theta = np.zeros(10)
# 損失函數的歷史取值
loss_history = []
X = trainX
Y = trainY
for i in range(100000):
# 梯度下降前先計算並保存一下損失值,以便之后審查
getLossValue(trainX, trainY, theta, loss_history)
# 更新參數
tempV = (Y - (theta.T.dot(X.T)).T).T.dot(X)
temp2 = -2 * tempV
theta = theta - learning_rate * temp2
print(loss_history)
然后這里最后是打印了loss_history的,即所有的損失值記錄
可以看到一開始的損失值高達400余萬,還是蠻嚇人的,而最后收斂到也有20多萬
這么看畢竟不直觀,因此我把損失值畫了個圖,可以看出收斂的還挺快的,看起來似乎迭代一兩千次就差不多收斂了,而后迭代十萬次百萬次得到的結果也不會好太多,因此在這個思路上這個方法似乎已經達到極限了。
損失值從400多萬降到20多萬,雖然還是蠻高的,不過也是降低了不少了,那就准備提交一下試試。
# 開始預測
test_dataa = pd.read_csv(test_data_path)
test_data = test_dataa[test_dataa['AMB_TEMP'] == 'PM2.5'].iloc[:, 2:].astype(int)
test_result = []
for i in range(test_data.shape[0]):
t = test_data.iloc[i].values
t = np.r_[1, t]
temp = theta.T.dot(t)
if temp < 0:
temp = 0
test_result.append(temp)
test_result中發現了一個負數,當然只比0小了一點點,這是不允許出現的,所以我直接讓它等於0了,別的看起來還挺正常的,沒有損失值那么嚇人。
然后導出數據(PS:因為對pandas不熟悉所以最終還是決定寫了個循環):
df = pd.DataFrame([])
for i in range(len(test_result)):
tp = {'id': 'id_' + str(i), 'value': test_result[i]}
tp = pd.DataFrame(tp, index=[i])
df = df.append(tp)
print(df)
output = pd.DataFrame(df)
output.to_csv('submission.csv', index=False)
然后趕緊在kaggle提交一下(注:上傳不了文檔的需要科學上網哦):
這個得分啥水平呢?
一打開排行榜就是Private Leaderboard,這個得分在這個榜上似乎是相當不錯了(一陣竊喜)
咦,旁邊還有個Public Leaderboard,點開看看,然后一直滑到最下面
倒數第七啊哈哈哈…
這個比賽因為已經結束了所以上不了榜,咱也不知道這個成績該上哪個榜,不過我估摸着,還是Public這個吧…
總之,到這兒也算是把代碼擼出來了,這個算法也算是為了梯度下降而去下降,並沒有從多個方面進行考慮,因此,我在最后還想提出一些優化思想。
1、 是否可以借鑒隨機森林的思想(只是說的高大上一點,其實就是取多組值進行預測然后求平均值)?
2、 這條線不應該是直線,而應該是曲線,所以是否可以考慮引入X的平方或者更高次方?
3、 本算法中為了增加測試數據而增加的測試數據是否將本算法本就不太准的准星帶的更偏?
4、 最后的預測值能否與時間掛鈎,即對於每一個小時都設計一個梯度下降,比如每天的10點,都用當天1~ 9點的數據進行訓練,然后需要預測的時間如果是10點就采用這個模型
5、 本題中除PM2.5外還提供了許多別的參數,能否將另外一些參數也考慮在內?
以后或許會有更新…
結語:搞機器學習,數學很重要!而且還得熟練掌握基礎的數據分析技巧…