第9章--隨機森林項目實戰——氣溫預測（1/2）
第8章已經講解過隨機森林的基本原理，本章將從實戰的角度出發，借助Python工具包完成氣溫預測任務，其中涉及多個模塊，主要包含隨機森林建模、特征選擇、效率對比、參數調優等。這個例子實在太長了，分為三篇介紹。這是第一篇。
隨機森林建模：氣溫預測的任務目標就是使用一份天氣相關數據來預測某一天的最高溫度，屬於回歸任務，首先觀察一下數據集：
輸出結果中表頭的含義如下。
- year,moth,day,week：分別表示的具體的時間。
- temp_2：前天的最高溫度值。
- temp_1：昨天的最高溫度值。
- average：在歷史中，每年這一天的平均最高溫度值。
- actual：就是標簽值，當天的真實最高溫度。
- friend：這一列可能是湊熱鬧的，你的朋友猜測的可能值，不管它就好。
該項目實戰主要完成以下3項任務。
1.使用隨機森林算法完成基本建模任務：包括數據預處理、特征展示、完成建模並進行可視化展示分析。
2.分析數據樣本量與特征個數對結果的影響：在保證算法一致的前提下，增加數據樣本個數，觀察結果變化。重新考慮特征工程，引入新特征后，觀察結果走勢。
3.對隨機森林算法進行調參，找到最合適的參數：掌握機器學習中兩種經典調參方法，對當前模型選擇最合適的參數。
9.1.1特征可視化與預處理
拿到數據之后，一般都會看看數據的規模，做到心中有數：
print('數據維度:', features.shape) #數據維度: (348, 9)
輸出結果顯示該數據一共有348條記錄，每個樣本有9個特征。如果想進一步觀察各個指標的統計特性，可以用.describe()展示：
輸出結果展示了各個列的數量，如果有數據缺失，數量就會有所減少。由於各列的統計數量值都是348，所以表明數據集中並不存在缺失值，並且均值、標准差、最大值、最小值等指標都在這里顯示。
對於時間數據，也可以進行格式轉換，原因在於有些工具包在繪圖或者計算的過程中，用標准時間格式更方便：
1 # 處理時間數據 2 import datetime 3 4 # 分別得到年，月，日 5 years = features['year'] 6 months = features['month'] 7 days = features['day'] 8 9 # datetime格式 10 dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)] 11 dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in dates]
為了更直觀地觀察數據，最簡單有效的辦法就是畫圖展示，首先導入Matplotlib工具包，再選擇一個合適的風格（其實風格差異並不是很大）：
1 # 准備畫圖 2 import matplotlib.pyplot as plt 3 4 %matplotlib inline 5 6 # 指定默認風格 7 plt.style.use('fivethirtyeight')
開始布局，需要展示4項指標，分別為最高氣溫的標簽值、前天、昨天、朋友預測的氣溫最高值。既然是4個圖，不妨采用2×2的規模，這樣會更清晰，對每個圖指定好其圖題和坐標軸即可：
1 # 設置布局 2 fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(nrows=2, ncols=2, figsize = (10,10)) 3 fig.autofmt_xdate(rotation = 45) 4 5 # 標簽值 6 ax1.plot(dates, features['actual']) 7 ax1.set_xlabel(''); ax1.set_ylabel('Temperature'); ax1.set_title('Max Temp') 8 9 # 昨天 10 ax2.plot(dates, features['temp_1']) 11 ax2.set_xlabel(''); ax2.set_ylabel('Temperature'); ax2.set_title('Previous Max Temp') 12 13 # 前天 14 ax3.plot(dates, features['temp_2']) 15 ax3.set_xlabel('Date'); ax3.set_ylabel('Temperature'); ax3.set_title('Two Days Prior Max Temp') 16 17 # 我的逗逼朋友 18 ax4.plot(dates, features['friend']) 19 ax4.set_xlabel('Date'); ax4.set_ylabel('Temperature'); ax4.set_title('Friend Estimate') 20 21 plt.tight_layout(pad=2)
上述代碼可以生成圖9-1的輸出。
圖9-1 各項特征指標
由圖可見，各項指標看起來還算正常（由於是國外的天氣數據，在統計標准上有些區別）。接下來，考慮數據預處理的問題，原始數據中的week列並不是一些數值特征，而是表示星期幾的字符串，計算機並不認識這些數據，需要轉換一下。
圖9-2是常用的轉換方式，稱作one-hot encoding或者獨熱編碼，目的就是將屬性值轉換成數值。對應的特征中有幾個可選屬性值，就構造幾列新的特征，並將其中符合的位置標記為1，其他位置標記為0。
圖9-2 特征編碼
既可以用Sklearn工具包中現成的方法完成轉換，也可以用Pandas中的函數，綜合對比后覺得用Pandas中的.get_dummies()函數最容易：
1 # 獨熱編碼 2 features = pd.get_dummies(features) 3 features.head(5)
完成數據集中屬性值的預處理工作后，默認會把所有屬性值都轉換成獨熱編碼的格式，並且自動添加后綴，這樣看起來更清晰。
其實也可以按照自己的方式設置編碼特征的名字，在使用時，如果遇到一個不太熟悉的函數，想看一下其中的細節，一個更直接的方法，就是在Notebook中直接調用help工具來看一下它的API文檔，下面返回的就是get_dummies的細節介紹，也可以查閱在線文檔：
help(pd.get_dummies) Help on function get_dummies in module pandas.core.reshape.reshape: get_dummies(data, prefix=None, prefix_sep='_', dummy_na=False, columns=None, sparse=False, drop_first=False, dtype=None) -> 'DataFrame' Convert categorical variable into dummy/indicator variables. Parameters ---------- data : array-like, Series, or DataFrame Data of which to get dummy indicators. prefix : str, list of str, or dict of str, default None String to append DataFrame column names. Pass a list with length equal to the number of columns when calling get_dummies on a DataFrame. Alternatively, `prefix` can be a dictionary mapping column names to prefixes. prefix_sep : str, default '_' If appending prefix, separator/delimiter to use. Or pass a list or dictionary as with `prefix`. dummy_na : bool, default False Add a column to indicate NaNs, if False NaNs are ignored. columns : list-like, default None Column names in the DataFrame to be encoded. If `columns` is None then all the columns with `object` or `category` dtype will be converted. sparse : bool, default False Whether the dummy-encoded columns should be backed by a :class:`SparseArray` (True) or a regular NumPy array (False). drop_first : bool, default False Whether to get k-1 dummies out of k categorical levels by removing the first level. dtype : dtype, default np.uint8 Data type for new columns. Only a single dtype is allowed. .. versionadded:: 0.23.0 Returns ------- DataFrame Dummy-coded data. See Also -------- Series.str.get_dummies : Convert Series to dummy codes. Examples -------- >>> s = pd.Series(list('abca')) >>> pd.get_dummies(s) a b c 0 1 0 0 1 0 1 0 2 0 0 1 3 1 0 0 >>> s1 = ['a', 'b', np.nan] >>> pd.get_dummies(s1) a b 0 1 0 1 0 1 2 0 0 >>> pd.get_dummies(s1, dummy_na=True) a b NaN 0 1 0 0 1 0 1 0 2 0 0 1 >>> df = pd.DataFrame({'A': ['a', 'b', 'a'], 'B': ['b', 'a', 'c'], ... 'C': [1, 2, 3]}) >>> pd.get_dummies(df, prefix=['col1', 'col2']) C col1_a col1_b col2_a col2_b col2_c 0 1 1 0 0 1 0 1 2 0 1 1 0 0 2 3 1 0 0 0 1 >>> pd.get_dummies(pd.Series(list('abcaa'))) a b c 0 1 0 0 1 0 1 0 2 0 0 1 3 1 0 0 4 1 0 0 >>> pd.get_dummies(pd.Series(list('abcaa')), drop_first=True) b c 0 0 0 1 1 0 2 0 1 3 0 0 4 0 0 >>> pd.get_dummies(pd.Series(list('abc')), dtype=float) a b c 0 1.0 0.0 0.0 1 0.0 1.0 0.0 2 0.0 0.0 1.0
特征預處理完成之后，還要把數據重新組合一下，特征是特征，標簽是標簽，分別在原始數據集中提取一下：
print('Shape of features after one-hot encoding:', features.shape) #Shape of features after one-hot encoding: (348, 15)
1 # 數據與標簽 2 import numpy as np 3 4 # 標簽 5 labels = np.array(features['actual']) 6 7 # 在特征中去掉標簽 8 features= features.drop('actual', axis = 1) 9 10 # 名字單獨保存一下，以備后患 11 feature_list = list(features.columns) 12 13 # 轉換成合適的格式 14 features = np.array(features)
在訓練模型之前，需要先對數據集進行切分：
1 # 數據集切分 2 from sklearn.model_selection import train_test_split 3 4 train_features, test_features, train_labels, test_labels = train_test_split(features, labels, test_size = 0.25, 5 random_state = 42) 6 print('訓練集特征:', train_features.shape) 7 print('訓練集標簽:', train_labels.shape) 8 print('測試集特征:', test_features.shape) 9 print('測試集標簽:', test_labels.shape)
訓練集特征: (261, 14) 訓練集標簽: (261,) 測試集特征: (87, 14) 測試集標簽: (87,)
9.1.2隨機森林回歸模型
萬事俱備，開始建立隨機森林模型，首先導入工具包，先建立1000棵樹模型試試，其他參數暫用默認值，然后深入調參任務：
1 # 導入算法 2 from sklearn.ensemble import RandomForestRegressor 3 4 # 建模 5 rf = RandomForestRegressor(n_estimators= 1000, random_state=42) 6 7 # 訓練 8 rf.fit(train_features, train_labels)
RandomForestRegressor(bootstrap=True, ccp_alpha=0.0, criterion='mse', max_depth=None, max_features='auto', max_leaf_nodes=None, max_samples=None, min_impurity_decrease=0.0, min_impurity_split=None, min_samples_leaf=1, min_samples_split=2, min_weight_fraction_leaf=0.0, n_estimators=1000, n_jobs=None, oob_score=False, random_state=42, verbose=0, warm_start=False)
由於數據樣本量非常小，所以很快可以得到結果，這里選擇先用MAPE指標進行評估，也就是平均絕對百分誤差。
1 # 預測結果 2 predictions = rf.predict(test_features) 3 4 # 計算誤差 5 errors = abs(predictions - test_labels) 6 7 # mean absolute percentage error (MAPE) 8 mape = 100 * (errors / test_labels) 9 10 print ('MAPE:',np.mean(mape))
MAPE: 6.011244187972058
其實對於回歸任務，評估方法還是比較多的，下面列出幾種，都很容易實現，也可以選擇其他指標進行評估。
9.1.3樹模型可視化方法
得到隨機森林模型后，現在介紹怎么利用工具包對樹模型進行可視化展示，首先需要安裝Graphviz工具，其配置過程如下。
第① 步：下載安裝。
登錄網站https://graphviz.gitlab.io/_pages/Download/Download_windows.html
下載graphviz-2.38.msi，完成后雙擊這個msi文件，然后一直單擊next按鈕，即可安裝Graphviz軟件（注意：一定要記住安裝路徑，因為后面配置環境變量會用到路徑信息，系統默認的安裝路徑是C:\Program Files (x86)\Graphviz2.38）。
第②步：配置環境變量。
將Graphviz安裝目錄下的bin文件夾添加到Path環境變量中。
本例中：
D:\tools\GraphViz\bin
第③步：驗證安裝。
進入Windows命令行界面，輸入“dot–version”命令，然后按住Enter鍵，如果顯示Graphviz的相關版本信息，則說明安裝配置成功，
dot -version dot - graphviz version 2.38.0 (20140413.2041) libdir = "D:\tools\GraphViz\bin" Activated plugin library: gvplugin_dot_layout.dll Using layout: dot:dot_layout Activated plugin library: gvplugin_core.dll Using render: dot:core Using device: dot:dot:core The plugin configuration file: D:\tools\GraphViz\bin\config6 was successfully loaded. render : cairo dot fig gd gdiplus map pic pov ps svg tk vml vrml xdot layout : circo dot fdp neato nop nop1 nop2 osage patchwork sfdp twopi textlayout : textlayout device : bmp canon cmap cmapx cmapx_np dot emf emfplus eps fig gd gd2 gif gv imap imap_np ismap jpe jpeg jpg
metafile pdf pic plain plain-ext png pov ps ps2 svg svgz tif tiff tk vml vmlz vrml wbmp xdot xdot1.2 xdot1.4 loadimage : (lib) bmp eps gd gd2 gif jpe jpeg jpg png ps svg
最后還需安裝graphviz、pydot和pydotplus插件，在命令行中輸入相關命令即可，代碼如下：
1 pip3 install graphviz 2 pip3 install pydot2 3 pip3 install pydotplus 4 pip3 install pydot
上述工具包安裝完成之后，就可以繪制決策樹模型：
1 # 導入所需工具包 2 from sklearn.tree import export_graphviz 3 import pydot #pip install pydot 4 5 # 拿到其中的一棵樹 6 tree = rf.estimators_[5] 7 8 # 導出成dot文件 9 export_graphviz(tree, out_file = 'tree.dot', feature_names = feature_list, rounded = True, precision = 1) 10 11 # 繪圖 12 (graph, ) = pydot.graph_from_dot_file('tree.dot') 13 14 # 展示 15 graph.write_png('tree.png');
執行完上述代碼，會在指定的目錄下（如果只指定其名字，會在代碼所在路徑下）生成一個tree.png文件，這就是繪制好的一棵樹的模型，如圖9-8所示。樹模型看起來有點太大，觀察起來不太方便，可以使用參數限制決策樹的規模，還記得剪枝策略嗎？預剪枝方案在這里可以派上用場。
1 print('The depth of this tree is:', tree.tree_.max_depth) 2 #The depth of this tree is: 15
圖9-9對生成的樹模型中各項指標的含義進行了標識，看起來還是比較好理解，其中非葉子節點中包括4項指標：所選特征與切分點、評估結果、此節點樣本數量、節點預測結果（回歸中就是平均）。
圖9-9 樹模型可視化中各項指標含義
9.1.4特征重要性
講解隨機森林算法的時候，曾提到使用集成算法很容易得到其特征重要性，在sklearn工具包中也有現成的函數，調用起來非常容易：
1 # 得到特征重要性 2 importances = list(rf.feature_importances_) 3 4 # 轉換格式 5 feature_importances = [(feature, round(importance, 2)) for feature, importance in zip(feature_list, importances)] 6 7 # 排序 8 feature_importances = sorted(feature_importances, key = lambda x: x[1], reverse = True) 9 10 # 對應進行打印 11 [print('Variable: {:20} Importance: {}'.format(*pair)) for pair in feature_importances]
Variable: temp_1 Importance: 0.7 Variable: average Importance: 0.19 Variable: day Importance: 0.03 Variable: temp_2 Importance: 0.02 Variable: friend Importance: 0.02 Variable: month Importance: 0.01 Variable: year Importance: 0.0 Variable: week_Fri Importance: 0.0 Variable: week_Mon Importance: 0.0 Variable: week_Sat Importance: 0.0 Variable: week_Sun Importance: 0.0 Variable: week_Thurs Importance: 0.0 Variable: week_Tues Importance: 0.0 Variable: week_Wed Importance: 0.0
上述輸出結果分別打印了當前特征及其所對應的特征重要性，繪制成圖表分析起來更容易：
1 # 轉換成list格式 2 x_values = list(range(len(importances))) 3 4 # 繪圖 5 plt.bar(x_values, importances, orientation = 'vertical') 6 7 # x軸名字 8 plt.xticks(x_values, feature_list, rotation='vertical') 9 10 # 圖名 11 plt.ylabel('Importance'); plt.xlabel('Variable'); plt.title('Variable Importances');
上述代碼可以生成圖9-10的輸出，可以明顯發現，temp_1和average這兩個特征的重要性占據總體的絕大部分，其他特征的重要性看起來微乎其微。那么，只用最厲害的特征來建模，其效果會不會更好呢？其實並不能保證效果一定更好，但是速度肯定更快，先來看一下結果：
圖9-10 隨機森林特征重要性
1 # 選擇最重要的那兩個特征來試一試 2 rf_most_important = RandomForestRegressor(n_estimators= 1000, random_state=42) 3 4 # 拿到這倆特征 5 important_indices = [feature_list.index('temp_1'), feature_list.index('average')] 6 train_important = train_features[:, important_indices] 7 test_important = test_features[:, important_indices] 8 9 # 重新訓練模型 10 rf_most_important.fit(train_important, train_labels) 11 12 # 預測結果 13 predictions = rf_most_important.predict(test_important) 14 15 errors = abs(predictions - test_labels) 16 17 # 評估結果 18 19 mape = np.mean(100 * (errors / test_labels)) 20 21 print('mape:', mape)
mape: 6.229055723613811
從損失值上觀察，並沒有下降，反而上升了，說明其他特征還是有價值的，不能只憑特征重要性就否定部分特征數據，一切還要通過實驗進行判斷。
但是，當考慮時間效率的時候，就要好好斟酌一下是否應該剔除掉那些用處不大的特征以加快構建模型的速度。到目前為止，已經得到基本的隨機森林模型，並可以進行預測，下面來看看模型的預測值與真實值之間的差異：
1 # 日期數據 2 months = features[:, feature_list.index('month')] 3 days = features[:, feature_list.index('day')] 4 years = features[:, feature_list.index('year')] 5 6 # 轉換日期格式 7 dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)] 8 dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in dates] 9 10 # 創建一個表格來存日期和其對應的標簽數值 11 true_data = pd.DataFrame(data = {'date': dates, 'actual': labels}) 12 13 # 同理，再創建一個來存日期和其對應的模型預測值 14 months = test_features[:, feature_list.index('month')] 15 days = test_features[:, feature_list.index('day')] 16 years = test_features[:, feature_list.index('year')] 17 18 test_dates = [str(int(year)) + '-' + str(int(month)) + '-' + str(int(day)) for year, month, day in zip(years, months, days)] 19 20 test_dates = [datetime.datetime.strptime(date, '%Y-%m-%d') for date in test_dates] 21 22 predictions_data = pd.DataFrame(data = {'date': test_dates, 'prediction': predictions}) 23 24 # 真實值 25 plt.plot(true_data['date'], true_data['actual'], 'b-', label = 'actual') 26 27 # 預測值 28 plt.plot(predictions_data['date'], predictions_data['prediction'], 'ro', label = 'prediction') 29 plt.xticks(rotation = '60'); 30 plt.legend() 31 32 # 圖名 33 plt.xlabel('Date'); plt.ylabel('Maximum Temperature (F)'); plt.title('Actual and Predicted Values');
通過上述輸出結果的走勢可以看出，模型已經基本能夠掌握天氣變化情況，接下來還需要深入數據，考慮以下幾個問題。
1.如果可利用的數據量增大，會對結果產生什么影響呢？
2.加入新的特征會改進模型效果嗎？此時的時間效率又會怎樣？
未完待續。
該書資源下載，請至異步社區：https://www.epubit.com