托馬茲·卓巴斯的《資料分析實戰》，2018年6月出版，本系列為讀書筆記。主要是為了系統整理，加深記憶。 第6章涵蓋了許多回歸模型，有線性的，也有非線性的。我們還會複習随機森林和支援向量機，它們可用來解決分類或回歸問題。
python資料分析個人學習讀書筆記-目錄索引
第6章涵蓋了許多回歸模型，有線性的，也有非線性的。我們還會複習随機森林和支援向量機，它們可用來解決分類或回歸問題。
本章會介紹一些技術，以預測發電廠生産的電量。你将學習以下主題：
·識别并解決資料中的多重共線性
·建構線性回歸模型，預測發電廠生産的電量
·使用OLS預測生産的電量
·使用CART估算發電廠生産的電量
·将kNN模型用于回歸問題
·将随機森林模型用于回歸分析
·使用SVM預測發電廠生産的電量
·訓練神經網絡，預測發電廠生産的電量
6.1導論
現實世界中一個常見的問題就是預測數量，或者用更抽象的說法，是在自變量和因變量之間找到聯系。本章中，我們專注于預測發電廠生産的電量。
本章使用的資料集來自于美國聯邦能源資訊管理局（U.S.Energy Information Adminis-tration）。從網站上下載下傳2014年的資料：http://www.eia.gov/electricity/data/eia923/xls/f923_2014.zip。
下載下傳頁面：https://www.eia.gov/electricity/data/eia923/ 右邊選擇ZIP下載下傳即可。
我們僅使用EIA923_Schedules_2_3_4_5_M_12_2014_Final_Revision.xlsx檔案中Generation and Fuel Data這一sheet的資料。我們将預測淨産量（兆瓦時）。所有資料都是類别變量（州名或燃料類型），我們決定将它們寫死。
最終，我們的資料集僅有4494行記錄，是整個資料集的一個子集。我們僅選擇少數州2014年超出100MWh的發電廠。我們也隻選取使用特定燃料類型的發電廠——Aggregate Fuel Code（AER）：coal（COL），distillate（DFO），hydroelectric conventional（HYC），biogenic municipal solid waste and landfill gas（MLG），natural gas（NG），nuclear（NUC），solar（SUN）以及wind（WND）。另外，隻選取特定發動機的發電廠：combined-cycle combustion turbine part（CT），combustion turbine（GT），hydraulic turbine（HY），internal combustion（IC），photovoltaic（PV），steam turbine（ST）以及wind turbine（WT）。
發電廠的電量資料方差很大，并且高度偏斜，中位數是18219MWh，平均數是448213MWh；25%的發電廠（2017家）在2014年的電量少于3496MWh，而美國人均用電量是12kWh左右（2010年的資料，http://energyalmanac.ca.gov/electricity/us_per_capita_electricity-2010.html），幾乎相當于291人的小鎮，如下圖所示：
大部分發電廠使用天然氣（超過1300家），将近750家使用水力或光生發電，接近580家使用風力發電。這樣，我們得到了我們資料集中使用可再生能源的發電廠比例：
NG天然氣，HYC水力，SUN太陽能、光生，WND風力
接近2500家發電廠使用蒸汽渦輪生産電力，緊接着的是水輪機和瓦斯渦輪：
ST蒸汽渦輪，HY水輪，GT瓦斯輪機
我們資料集中加利福尼亞的發電廠最多，超過1200家。德克薩斯以超過500家的數量位居第二，排第三的是馬薩諸塞：
6.2識别并解決資料中的多重共線性
多重共線性指的是這樣的場景，一個（或多個）自變量可以表示為其他自變量的線性組合。
例如，假設我們用州中的人口、家庭數和發電廠數來預測用電量。這種情況下，可以想見，州中人口越多，家庭數應該也越多，也就是說，家庭數可以表示為州中人口的某些（接近）線性關系。
現在，如果我們要基于多重共線的資料估算模型，很有可能會有一個（甚至所有的共線變量）是不顯著的。将這個變量移除（并隻保留和因變量最相關的變量，這樣，保留了大部分對方差的貢獻），并不會降低模型的解釋能力。
準備:需裝好pandas、NumPy和Scikit。我們使用Scikit降低次元。
要判斷我們的資料是否共線，我們需要檢視自變量協方差矩陣的特征值；如果資料共線，對協方差矩陣的分析有助于辨識變量（regression_multicollinearity.py檔案）：
1 # this is needed to load helper from the parent folder
2 import sys
3 sys.path.append('..')
4
5 # the rest of the imports
6 import helper as hlp
7 import pandas as pd
8 import numpy as np
9 import sklearn.decomposition as dc
10
11 def reduce_PCA(x, n):
12 '''
13 Reduce the dimensions using Principal Component
14 Analysis
15 '''
16 # create the PCA object
17 pca = dc.PCA(n_components=n, whiten=True)
18
19 # learn the principal components from all the features
20 return pca.fit(x)
21
22 # the file name of the dataset
23 r_filename = '../../Data/Chapter06/power_plant_dataset.csv'
24
25 # read the data
26 csv_read = pd.read_csv(r_filename)
27
28 x = csv_read[csv_read.columns[:-1]].copy()
29 y = csv_read[csv_read.columns[-1]]
30
31 # produce correlation matrix for the independent variables
32 #生成自變量協方差矩陣
33 corr = x.corr()
34
35 # and check the eigenvectors and eigenvalues of the matrix
36 #檢查矩陣的特征向量與特征值
37 w, v = np.linalg.eig(corr)
38 print('Eigenvalues: ', w)
39
40 # values that are close to 0 indicate multicollinearity
41 #接近0的值意味着多重共載
42 s = np.nonzero(w < 0.01)
43 # inspect which variables are collinear
44 #找出共線變量
45 print('Indices of eigenvalues close to 0:', s[0])
46
47 all_columns = []
48 for i in s[0]:
49 print('\nIndex: {0}. '.format(i))
50
51 t = np.nonzero(abs(v[:,i]) > 0.33)
52 all_columns += list(t[0]) + [i]
53 print('Collinear: ', t[0]) 原理：首先，我們和往常一樣載入必需的子產品并讀入資料集。然後，我們将資料集拆分成自變量x（我們為原始的資料集建立一個副本）和因變量y。
接下來，我們生成所有自變量之間的協方差矩陣，使用NumPy在linalg子產品中的eig（...）方法找到特征值和特征向量。eig（...）方法要求輸入一個方陣。
隻有方陣才有特征向量和特征值。非方形的矩陣有奇異值（https://www.math.washington.edu/~greenbau/Math_554/Course_Notes/ch1.5.pdf）。
接近于0的特征值意味着（多重）共線性。要找到這些值（我們的例子中，<0.01）的索引，我們先建立一個真值向量：每個元素不是False（0）就是True（1）。簡單一句w<0.001就可以完成；這個操作的結果就是到原始向量w的等長真值向量，這個向量中小于0.001的值标為True。然後我們可以使用NumPy的.nonzero（...）方法傳回非零元素的索引清單（我們的例子中，傳回的就是True元素的索引）。你很有可能在不止一個位置獲得特征值接近0的元素。對于我們的資料集，你可能得到類似的結果：
/*
Indices of eigenvalues close to 0: [28 29 30 31 32]
*/ 我們周遊索引，找到相應的特征向量。特征向量中顯著大于0的元素會幫我們找到共線的變量。由于特征向量可以有負值，我們傳回絕對值>0.33的元素索引：
/*
Index: 28.
Collinear: [0 1 3]
Index: 29.
Collinear: [ 9 11 12 13]
Index: 30.
Collinear: [15]
Index: 31.
Collinear: [ 2 10]
Index: 32.
Collinear: [ 4 14]
*/ 注意這些變量的重複。我們将這些值存到all_columns清單。等我們周遊所有接近0的特征值，我們來看看哪些變量是共線的：
for i in np.unique(all_columns):
print('Variable {0}: {1}'.format(i, x.columns[i])) 我們周遊all_columns清單中去重後的值（使用NumPy的.unique（...）方法），輸出x中所有列相應的變量名：
/*
Variable 0: fuel_aer_NG
Variable 1: fuel_aer_DFO
Variable 2: fuel_aer_HYC
Variable 3: fuel_aer_SUN
Variable 4: fuel_aer_WND
Variable 9: mover_GT
Variable 10: mover_HY
Variable 11: mover_IC
Variable 12: mover_PV
Variable 13: mover_ST
Variable 14: mover_WT
Variable 15: state_CA
Variable 28: state_OH
Variable 29: state_GA
Variable 30: state_WA
Variable 31: total_fuel_cons
Variable 32: total_fuel_cons_mmbtu
*/ 現在我們能清楚看出來資料集中的共線性了：某種發動機當然與某種燃料相關。比如，風力發動機用的是風，光能發動機用的是太陽能，水力發電靠的是水。另外，加利福尼亞州也出現了；加利福尼亞是擁有最多風力渦輪機的州（第二位是德克薩斯，http://www.awea.org/resources/statefactsheets.aspx?itemnumber=890）。
（除了删減模型中的變量而外）解決多重共線性的一個方法就是對資料集降維。我們可以通過PCA做到這一點（參考本書5.3節）。
我們重用reduce_PCA（...）方法：
1 # and reduce the data keeping only 5 principal components
2 n_components = 5
3 z = reduce_PCA(x, n=n_components)
4 pc = z.transform(x)
5
6 # how much variance each component explains?
7 #每個成分對方差貢獻多少大？
8 print('\nVariance explained by each principal component: ',
9 z.explained_variance_ratio_)
10
11 # and total variance accounted for
12 #一共貢獻多少方差
13 print('Total variance explained: ',
14 np.sum(z.explained_variance_ratio_)) 我們希望得到5個PC（principal components，主成分）。reduce_PCA（...）方法估算模型，.transform（...）方法轉換我們的資料，傳回5個主成分。我們可以檢視每個主成分對整體方差的貢獻：
第一個主成分大約貢獻了14.5%，第二個約13%，第三個約11%。剩下兩個的貢獻加起來大約16%，是以一共是：
/*
Variance explained by each principal component: [0.14446578 0.13030196 0.11030824 0.0935897 0.06694611]
Total variance explained: 0.5456117974906106
*/ 我們将主成分和其他變量一起存到一個檔案中：
1 # append the reduced dimensions to the dataset
2 for i in range(0, n_components):
3 col_name = 'p_{0}'.format(i)
4 x[col_name] = pd.Series(pc[:, i])
5
6 x[csv_read.columns[-1]] = y
7 csv_read = x
8
9 # output to file
10 w_filename = '../../Data/Chapter06/power_plant_dataset_pc.csv'
11 with open(w_filename, 'w',newline='') as output:
12 output.write(csv_read.to_csv(index=False)) 首先，我們将主成分附加到x資料集後面；我們希望将因變量放在最後。然後，我們将資料集輸出到power_plant_dataset_pc.csv供後續使用。
6.3建構線性回歸模型
線性回歸無疑是最易于建構的模型。當你知道因變量和自變量之間是線性關系時就可以選擇。
準備：需裝好pandas、NumPy和Scikit。
用Scikit估算回歸模型很簡單（regression_linear.py檔案）：
1 # this is needed to load helper from the parent folder
2 import sys
3 sys.path.append('..')
4
5 # the rest of the imports
6 import helper as hlp
7 import pandas as pd
8 import numpy as np
9 import sklearn.linear_model as lm
10
11 @hlp.timeit
12 def regression_linear(x,y):
13 '''
14 Estimate a linear regression估算線性回歸
15 '''
16 # create the regressor object建立回歸對象
17 linear = lm.LinearRegression(fit_intercept=True,
18 normalize=True, copy_X=True, n_jobs=-1)
19
20 # estimate the model估算模型
21 linear.fit(x,y)
22
23 # return the object
24 return linear
25
26 # the file name of the dataset
27 r_filename = '../../Data/Chapter06/power_plant_dataset_pc.csv'
28
29 # read the data
30 csv_read = pd.read_csv(r_filename)
31
32 # select the names of columns
33 dependent = csv_read.columns[-1]
34 independent_reduced = [
35 col
36 for col
37 in csv_read.columns
38 if col.startswith('p')
39 ]
40
41 independent = [
42 col
43 for col
44 in csv_read.columns
45 if col not in independent_reduced
46 and col not in dependent
47 ]
48
49 # split into independent and dependent features
50 #拆成自變量和因變量特征
51 x = csv_read[independent]
52 x_red = csv_read[independent_reduced]
53 y = csv_read[dependent]
54
55
56 # estimate the model using all variables (without PC)
57 #使用所有變量估算模型
58 regressor = regression_linear(x,y)
59
60 # print model summary
61 print('\nR^2: {0}'.format(regressor.score(x,y)))
62 coeff = [(nm, coeff)
63 for nm, coeff
64 in zip(x.columns, regressor.coef_)]
65 intercept = regressor.intercept_
66 print('Coefficients: ', coeff)
67 print('Intercept', intercept)
68 print('Total number of variables: ',
69 len(coeff) + 1) 原理：按照慣例，我們先導入必需的子產品并讀入資料集。然後将資料拆成子集。
首先，我們取出代表因變量的列名（資料集的最後一列）、主成分（independent_reduced）以及自變量。選取主成分時，由于所有的主成分列都以p_開頭（參見本書6.2節中的“更多”部分），我們使用Python内建的.startswith（...）方法。選取自變量時，我們取出所有既不是independent_reduced也不是dependent的列。然後我們從資料集中提取這些列。
我們調用regression_linear（...）方法估算模型。這個方法接受兩個參數：自變量x和因變量y。
regression_linear（...）方法内部，我們先用.LinearRegression（...）方法建立模型對象。我們指定模型要估算常數（fit_intercept=True），正規化資料，并複制自變量。
Python（預設情況下）傳參是傳引用（不像C或C++那樣是傳值），是以方法改變了傳入的資料，這一可能性總是存在的。http://www.python-course.eu/passing_arguments.php。
我們也将并行跑的任務數設為計算機的核數。盡管對于線性回歸來說，這一般不是個問題。
有了模型對象，我們應用（.fit（...）方法）資料，即估算并傳回模型。
線性回歸模型的性能可以R2(R平方)的形式表示。這個名額可以了解為模型對方差的貢獻的度量。
在本書6.4節，我們會介紹如何手動計算R2。
LinearRegression對象提供了.score（...）方法。我們的模型有如下的結果：
/*
R^2: 0.9965751657989851
*/ 我們使用清單表達式建立了coeff清單，以将所有變量的清單及相應的系數列印出來。我們先（用zip（...）方法）将x的列和回歸器的coeff_屬性壓縮在一起；這建立了一個元組構成的清單，我們将周遊清單中的元組以建立coeff清單：
注意州的系數之間的差距有多大。另外，所有電力相關的變量（燃料和發動機）都是負的。這展現了州對我們的模型的整體表現有一個最小的影響，特别是與截距相比時：
/*
Intercept -11531009013.842264
*/ 我們已由前一技巧知道，燃料總量是相關變量，是以我們可以安心地将其移除。
别忘了資料集中有多重共線性，我們隻用主成分估算模型：
1 # estimate the model
2 regressor_nm = regression_linear(x_no_state,y)
3
4 # print model summary
5 print('\nR^2: {0}'.format(regressor_nm.score(x_no_state,y)))
6 coeff = [(nm, coeff)
7 for nm, coeff
8 in zip(x_no_state.columns, regressor_nm.coef_)]
9 intercept = regressor_nm.intercept_
10 print('Coefficients: ', coeff)
11 print('Intercept', intercept)
12 print('Total number of variables: ',
13 len(coeff) + 1) 這個模型的R2結果更差：
/*
R^2: 0.12740548278756392
*/ 由于州的影響不大，我們決定試試不考慮這個因素。我們也不考慮total_fuel_consumption：
1 # removing the state variables and keeping only fuel and state
2 columns = [col for col in independent if 'state' not in col and col != 'total_fuel_cons']
3 x_no_state = x[columns] R2有下降，但不顯著：
/*
R^2: 0.9957289744500172
*/ 不過，現在我們模型的自變量由35個（包括截距）變成了18個。有個精神要把握住，類似的表現下，你應該選解釋變量更少的模型。
要找到最佳直線（并且你隻有一個解釋變量），你可以使用Seaborn。盡管結果并不直覺，但我們可以使用主成分來展示。
首先，由于我們将主成分存在列中，我們将它們放在棧頂，這樣資料集中隻有三列——PC是主成分，x是主成分的值，y是正規化的發電量：
1 # stack up the principal components
2 #将主成分放在棧頂
3 pc_stack = pd.DataFrame()
4
5 # stack up the principal components
6 #将主成分放在棧頂
7 for col in x_red.columns:
8 series = pd.DataFrame()
9 series['x'] = x_red[col]
10 series['y'] = y
11 series['PC'] = col
12 pc_stack = pc_stack.append(series) 我們先建立一個空的pandas DataFrame。然後，我們周遊x_red中的列（隻有主成分）。DataFrame的内部序列将主成分的值放在x，将因變量的值放在y，将主成分的名字放在PC。然後我們将序列附加到pc_stack。
現在繪制圖線：
1 # Show the results of a linear regression within each
2 # principal component
3 sns.lmplot(x='x', y='y', col='PC', hue='PC', data=pc_stack,
4 col_wrap=2, size=5)
5
6 pl.savefig('../../Data/Chapter06/charts/regression_linear.png',
7 dpi=300) .lmplot（...）方法将線性回歸用到我們的資料上。col參數指定了圖表的分組，我們将為5個主成分生成圖表。hue參數将為每個主成分更換顔色。col_wrap=2意味着我們在一行中放兩個圖表，每個都是5英寸高。
結果如下：
Tips:
/*
D:\tools\Python37\lib\site-packages\seaborn\regression.py:546: UserWarning: The `size` paramter has been renamed to `height`; please update your code.
warnings.warn(msg, UserWarning)
*/ 解決方案：size=5改為height=5
/*
sns.lmplot(x='x', y='y', col='fuel', hue='fuel',
data=fuel_stack, col_wrap=2, height=5)
*/ 我們也檢視兩個燃料變量，看是否存線上性關系：
# select only the fel consumption
fuel_cons = ['total_fuel_cons','total_fuel_cons_mmbtu']
x = csv_read[fuel_cons] 下面的圖顯示了，total_fuel_cons_mmbtu存線上性關系，total_fuel_cons也存在接近線性的關系：
/*
MMBTU，million British Thermal Units代表百萬英熱機關，百萬英制熱機關。
1mmBtu=2.52X10^8cal(卡) =2.52X10^5 kcal(千卡)
1桶原油=5.8 MMBTU
*/ 燃料消耗總量不考慮存儲量；隻考慮使用的資源總體積（或重量）。但是，一磅的鈾或钚比起一磅的煤，産生的能量要多得多。是以total_fuel_cons在x的低值一端有很大的偏離。當我們考慮特定的能量（或對液體燃料來說的能量密度），total_fuel_cons_mmbtu，我們在能量的消耗和産出之間看到了線性關系。
你可以從圖中看出來，這條線的趨勢<1，代表着系統的損失。注意，有些燃料是可再生的，風力或太陽能發電時我們幾乎不用承擔輸入相關（燃料）的成本。
6.4使用OLS預測生産的電量
OLS（Ordinary Least Squares，最小二乘法）也是一個線性模型。使用最小二乘法實際上也估算了一個線性回歸。然而，即使自變量與因變量之間不是線性關系，隻要參數之間是線性關系，OLS就可以估算出一個模型。
準備:需裝好pandas和Statsmodels。
按照慣例，為模型估算包一層函數（regression_ols.py檔案）：
1 import statsmodels.api as sm
2
3 @hlp.timeit
4 def regression_ols(x,y):
5 '''
6 Estimate a linear regression
7 '''
8 # add a constant to the data
9 x = sm.add_constant(x)
10
11 # create the model object
12 model = sm.OLS(y, x)
13
14 # and return the fit model
15 return model.fit()
16
17 # the file name of the dataset
18 r_filename = '../../Data/Chapter06/power_plant_dataset_pc.csv'
19
20 # read the data
21 csv_read = pd.read_csv(r_filename)
22
23 # select the names of columns
24 dependent = csv_read.columns[-1]
25 independent_reduced = [
26 col
27 for col
28 in csv_read.columns
29 if col.startswith('p')
30 ]
31
32 independent = [
33 col
34 for col
35 in csv_read.columns
36 if col not in independent_reduced
37 and col not in dependent
38 ]
39
40 # split into independent and dependent features
41 x = csv_read[independent]
42 y = csv_read[dependent]
43
44 # estimate the model using all variables (without PC)
45 regressor = regression_ols(x,y)
46 print(regressor.summary()) 原理：首先，我們按照慣例讀入必需的資料并選取自變量和因變量。
然後我們估算模型。先往資料中加一個常量；這僅僅是加一列1。然後用.OLS（...）方法建立回歸器模型；這個方法與Scikit中的對應方法不同，接受自變量和因變量。最後，我們應用并傳回模型。接着，列印出模型總結（這個在Scikit的線性回歸中沒有實作）。我們簡化了表以節省空間：
/*
D:\tools\Python37\lib\site-packages\numpy\core\fromnumeric.py:2495: FutureWarning: Method .ptp is deprecated and will be removed in a future version.
Use numpy.ptp instead.
return ptp(axis=axis, out=out, **kwargs)
*/ 解決方案：無
可見，燃料類型、發動機和州名這幾個變量在統計上是不顯著的：p-value大于0.05。這些變量可以放心地從模型中移除，不會損失太多精度。注意，數字2告訴我們資料集中有很嚴重的多重共線性問題——這一點我們已經知道了。
R2以及total_fuel_cons與total_fuel_cons_mmbtu統計上的顯著性證明了我們之前的發現：僅靠這兩個變量就決定了發電廠的發電量。使用什麼燃料并沒有多大差别：同樣的輸入能量會得到同樣的電量，不管是什麼燃料。其實發電廠本質上就是個能量轉換器：它将一種形式的能量轉化為另一種。問題來了——為什麼使用化石燃料，而不是可再生能源？
我們看下如果僅僅使用這兩個變量，模型是不是同樣有效：
1 # remove insignificant variables
2 #移除不顯著的變量
3 significant = ['total_fuel_cons', 'total_fuel_cons_mmbtu']
4 x_red = x[significant]
5
6 # estimate the model with limited number of variables
7 # 用僅有的變量估算模型
8 regressor = regression_ols(x_red,y)
9 print(regressor.summary()) 模型的表現并沒有變差：
如你所見，R2從0.997降到了0.996，差異可忽略不計。注意這個常量，盡管統計上不顯著，對最後的結果幾乎沒有影響。實際上，模型中值得保留的變量隻有total_fuel_cons_mmbtu，因為它與發電量有最強的（線性）關系。
我們也可以使用MLPY估算OLS模型（regression_ols_alternative.py檔案）：
1 import mlpy as ml
2
3 @hlp.timeit
4 def regression_linear(x,y):
5 '''
6 Estimate a linear regression
7 '''
8 # create the model object
9 ols = ml.OLS()
10
11 # estimate the model
12 ols.learn(x, y)
13
14 # and return the fit model
15 return ols
16
17 # the file name of the dataset
18 r_filename = '../../Data/Chapter06/power_plant_dataset_pc.csv'
19
20 # read the data
21 csv_read = pd.read_csv(r_filename)
22
23 # remove insignificant variables
24 significant = ['total_fuel_cons_mmbtu']
25 x = csv_read[significant]
26
27 # x = np.array(csv_read[independent_reduced[0]]).reshape(-1,1)
28 y = csv_read[csv_read.columns[-1]]
29
30 # estimate the model using all variables (without PC)
31 regressor = regression_linear(x,y)
32
33 # predict the output
34 predicted = regressor.pred(x)
35
36 # and calculate the R^2
37 score = hlp.get_score(y, predicted)
38 print('R2: ', score) MLPY的.OLS（）方法隻能建構簡單的線性模型。我們使用最顯著的變量，total_fuel_cons_mmbtu。由于MLPY不提供關于模型表現的資訊，我們決定手寫一個評分函數（helper.py檔案）：
1 def get_score(y, predicted):
2 '''
3 Method to calculate R^2
4 '''
5 # calculate the mean of actuals
6 mean_y = y.mean()
7
8 # calculate the total sum of squares and residual
9 # sum of squares
10 sum_of_square_total = np.sum((y - mean_y)**2)
11 sum_of_square_resid = np.sum((y - predicted)**2)
12
13 return 1 - sum_of_square_resid / sum_of_square_total 這個方法遵循了決定系數的标準定義：用1減去殘差平方和與總平方和的比值。總平方和是因變量中每個觀測值與平均值的差的平方和。殘差平方和是因變量每個實際觀測值與模型輸出值的差的平方和。
get_score（...）方法要求兩個輸入參數：因變量的實際觀測值和模型預測的值。
對于本模型，我們算出下述值：
/*
The method regression_linear took 0.00 sec to run.
R2: 0.9951670803634637
*/ 參考:通路這個網址，看看OLS的另一種做法：http://www.datarobot.com/blog/ordinary-least-squares-in-python/。
6.5使用CART估算發電廠生産的電量
之前，我們用決策樹給銀行電話分類過（參考本書3.6節）。對于回歸問題來說，分類樹和回歸樹是一碼事。
準備：需裝好pandas和Scikit。
用Scikit估算CART很簡單（regression_cart.py檔案）：
1 import sklearn.tree as sk
2
3 @hlp.timeit
4 def regression_cart(x,y):
5 '''
6 Estimate a CART regressor
7 '''
8 # create the regressor object
9 cart = sk.DecisionTreeRegressor(min_samples_split=80,
10 max_features="auto", random_state=66666,
11 max_depth=5)
12
13 # estimate the model
14 cart.fit(x,y)
15
16 # return the object
17 return cart 原理：按照慣例，我們先載入資料，提取出因變量y和自變量x。同樣以主成分的形式将自變量存到x_red。
估算模型時，和之前所有用Scikit估算的模型一樣，我們建立模型對象。.DecisionTree Regressor（...）方法接受一系列參數。min_samples_split是節點中為了實作拆分至少要有的樣本數。max_features是樹可以使用的最多變量數；auto參數是說，如果需要的話，模型可以使用所有的自變量。random_state決定了模型的初始随機狀态，max_depth決定了樹最多能有多少層。
然後應用資料，傳回估算的模型。我們看看模型表現如何：
# print out the results
print('R2: ', regressor.score(x,y)) 結果看上去和之前差不多：
/*
The method regression_cart took 0.01 sec to run.
R2: 0.9595455136505057
The method regression_cart took 0.01 sec to run.
R: 0.7035525409649943
0. p_0: 0.17273577812298352
1. p_1: 0.08688236875805344
2. p_2: 0.0
3. p_3: 0.5429677683710962
4. p_4: 0.19741408474786684
*/ 決策樹的美在于，模型根本不會使用不顯著的自變量。我們看看哪些是顯著的：
1 for counter, (nm, label) \
2 in enumerate(
3 zip(x.columns, regressor.feature_importances_)
4 ):
5 print("{0}. {1}: {2}".format(counter, nm,label)) /*
0. fuel_aer_NG: 0.0
1. fuel_aer_DFO: 0.0
2. fuel_aer_HYC: 0.0
3. fuel_aer_SUN: 0.0
4. fuel_aer_WND: 0.0
5. fuel_aer_COL: 0.0
6. fuel_aer_MLG: 0.0
7. fuel_aer_NUC: 0.0
8. mover_CT: 0.0
9. mover_GT: 0.0
10. mover_HY: 0.0
11. mover_IC: 0.0
12. mover_PV: 0.0
13. mover_ST: 0.0
14. mover_WT: 0.0
15. state_CA: 0.0
16. state_TX: 0.0
17. state_NY: 0.0
18. state_FL: 0.0
19. state_MN: 0.0
20. state_MI: 0.0
21. state_NC: 0.0
22. state_PA: 0.0
23. state_MA: 0.0
24. state_WI: 0.0
25. state_NJ: 0.0
26. state_IA: 0.0
27. state_IL: 0.0
28. state_OH: 0.0
29. state_GA: 0.0
30. state_WA: 0.0
31. total_fuel_cons: 0.0
32. total_fuel_cons_mmbtu: 1.0
*/ 我們周遊regressor（.DecisionTreeRegressor（...））的feature_importances_屬性，檢視使用zip（...）的相應名字。下面是簡略的清單：
可以看到，模型隻認為total_fuel_cons_mmbtu是顯著的，具有100%的預測能力；其他變量都對模型的表現沒有影響。
可以通過最後的決策樹确認這一點：
# and export to a .dot file
sk.export_graphviz(regressor,
out_file='../../Data/Chapter06/CART/tree.dot') 要将.dot檔案轉換為易讀格式，可以在指令行使用dot指令。切到路徑下（假設你在Codes/Chapter6檔案夾）：
要為我們的樹生成一個PDF檔案，隻需執行下面這行指令：
/*
D:\PythonLearn\資料分析實戰-托馬茲·卓巴斯\tools\GraphViz\release\bin\dot -Tpdf D:\Java2018\practicalDataAnalysis\Data\Chapter06\CART\tree.dot -o D:\Java2018\practicalDataAnalysis\Data\Chapter06\CART\tree.pdf
D:\PythonLearn\資料分析實戰-托馬茲·卓巴斯\tools\GraphViz\release\bin\dot -Tpdf D:\Java2018\practicalDataAnalysis\Data\Chapter06\CART\tree_red.dot -o D:\Java2018\practicalDataAnalysis\Data\Chapter06\CART\tree_red.pdf
*/ 得到的樹應該隻有一個決定變量：
X[32]，看一下列名清單，果然就是total_fuel_cons_mmbtu啊。
現在看看如果模型隻有主成分，表現如何：
1 # estimate the model using Principal Components only
2 regressor_red = regression_cart(x_red,y)
3
4 # print out the results
5 print('R: ', regressor_red.score(x_red,y)) R2的值也可接受：
/*
R^2: 0.7035525409649943
*/ 驚奇的是，（如果你還記得每個主成分對方差的貢獻）第三個主成分是最重要的，而第二個是完全沒用的：
/*
0. p_0: 0.17273577812298352
1. p_1: 0.08688236875805344
2. p_2: 0.0
3. p_3: 0.5429677683710962
4. p_4: 0.19741408474786684
*/ （通常）使用主成分的問題在于無法直接看懂結果，要得出真實意義就有點繞（https://onlinecourses.science.psu.edu/stat505/node/54）。
得到的樹比前一棵更複雜：
參考:這裡對這裡對CART的介紹挺好：http://www.stat.wisc.edu/~loh/treeprogs/guide/wires11.pdf。
6.6将kNN模型用于回歸問題
第3章中用到的k鄰近模型，盡管多用于解決分類問題，但其實也可以用于解決回歸問題。本技巧中來談談怎麼用。
用Scikit估算模型很簡單（regression_knn.py檔案）：
1 import sklearn.neighbors as nb
2 import sklearn.cross_validation as cv
3
4 @hlp.timeit
5 def regression_kNN(x,y):
6 '''
7 Build the kNN classifier
8 '''
9 # create the classifier object
10 knn = nb.KNeighborsRegressor(n_neighbors=80,
11 algorithm='kd_tree', n_jobs=-1)
12
13 # fit the data
14 knn.fit(x,y)
15
16 #return the classifier
17 return knn 原理：首先，我們讀入資料，并拆成因變量y和自變量x_sig；我們隻選取之前發現是顯著的變量total_fuel_cons和total_fuel_cons_mmbtu。我們也會測試隻使用主成分建構的模型；這種情況下，我們使用x_principal。
我們使用Scikit的.KNeighborsRegressor（...）方法建構kNN回歸模型。我們将n_neighbors設為80，指定kd_tree作為使用的算法。
檢視kd-樹是如何形成的：http://www.alglib.net/other/nearestneighbors.php。
指定并行任務數為-1，模型會使用與處理器核數相同的任務數。
估算的模型R2分值為0.94。
但是問題來了：模型的表現對資料來說有多敏感？我們用Scikit的cross_validation子產品來測一下這個：
1 import sklearn.cross_validation as cv
2 # test the sensitivity of R2
3 scores = cv.cross_val_score(regressor, x_sig, y, cv=100)
4 print('Expected R2: {0:.2f} (+/- {1:.2f})'\
5 .format(scores.mean(), scores.std()**2)) 首先，我們加載子產品。然後，我們使用.cross_val_score（...）方法測試模型。第一個參數是要測試的模型，第二個參數是自變量的集合，第三個參數是因變量的集合。最後一個參數cv，指定了折疊的數量：原始資料集随機劃分成的分塊數。通常情況下，交叉驗證隻用到一個分塊，剩下的用作訓練資料集。
我們的結果：
/*
The method regression_kNN took 0.00 sec to run.
R^2: 0.9433911310360819
Expected R^2: 0.91 (+/- 0.03)
*/ /*
ModuleNotFoundError: No module named 'sklearn.cross_validation'
*/ 解決方案：改為從 sklearn.model_selection
/*
import sklearn.model_selection as cv
*/ 我們運作得到的R^2分數是0.94。不過，交叉驗證得到的期望值是0.91，标準差是0.03（是以我們運作得到的結果在這個範圍之内）。
我們測試隻使用主成分的模型：
1 # estimate the model using Principal Components only
2 #隻使用主成分估算模型
3 regressor_principal = regression_kNN(x_principal,y)
4
5 print('R^2: ', regressor_principal.score(x_principal,y)) 得到的結果你也許會意外（特别是如果你學到的說法是R^2的取值範圍在0到1之間）：
/*
The method regression_kNN took 0.01 sec to run.
R^2: 0.4602300520714968
Expected R^2: -18.55 (+/- 14427.43)
*/ 是以，我們運作得到的R^2分數是0.46。然而，期望值是負的（标準差極高）。
我們怎麼會得到負值？考慮這個例子：我們有一組資料，差不多遵循這個關系y=x^2。現在，我們假設有些模型（隻能估算線性關系）給我們一個遵循y_pred=3*x函數的預測。
這樣的預測，R^2的值就在-0.4左右。如果你分析下R^2的計算過程，你就會明白某些情況下R^2的值會變成負的。
參考StackOverflow上這個回答：http://stats.stackexchange.com/questions/12900/when-is-r-squared-negative。
6.7将随機森林模型用于回歸分析
類似于決策樹，随機森林也可以用于解決回歸問題。我們之前用它歸類電話（參考本書3.7節）。這裡，我們使用随機森林預測電廠的發電量。
随機森林是一組模型。這個例子的代碼來自第3章（regression_randomForest.py檔案）：
1 import sys
2 sys.path.append('..')
3
4 # the rest of the imports
5 import helper as hlp
6 import pandas as pd
7 import numpy as np
8 import sklearn.ensemble as en
9 import sklearn.model_selection as cv
10
11 @hlp.timeit
12 def regression_rf(x,y):
13 '''
14 Estimate a random forest regressor
15 估算随機森林回歸器
16 '''
17 # create the regressor object
18 random_forest = en.RandomForestRegressor(
19 min_samples_split=80, random_state=666,
20 max_depth=5, n_estimators=10)
21
22 # estimate the model
23 random_forest.fit(x,y)
24
25 # return the object
26 return random_forest 原理：按照慣例，我們讀入資料并拆成自變量和因變量。後面會用于估算模型。
.RandomForestRegressor（...）方法建立了模型對象。和決策樹一樣，min_samples_split決定拆分時節點最少的樣本數。random_state決定了模型的初始随機狀态，max_depth決定了樹最多能有多少層。最後一個參數n_estimators，告訴模型森林中估算多少樹。
估算後，我們列印出R^2評分，測試它對輸入的敏感度，以清單形式列印出變量及其重要性：
1 # print out the results
2 print('R^2: ', regressor_red.score(x_red,y))
3
4 # test the sensitivity of R2 測試r^2敏感度
5 scores = cv.cross_val_score(regressor_red, x_red, y, cv=100)
6 print('Expected R^2: {0:.2f} (+/- {1:.2f})'\
7 .format(scores.mean(), scores.std()**2))
8
9 # print features importance
10 for counter, (nm, label) \
11 in enumerate(
12 zip(x_red.columns, regressor_red.feature_importances_)
13 ):
14 print("{0}. {1}: {2}".format(counter, nm,label)) Tips:
/*
Traceback (most recent call last):
File "D:\Java2018\practicalDataAnalysis\Codes\Chapter06\regression_randomForest.py", line 79, in <module>
x_red = csv_read[features[0]]
KeyError: "None of [Int64Index([32], dtype='int64')] are in the [columns]"
*/ 解決方案：
改為：x_red = csv_read.iloc[:,features[0]]
代碼的結果如下（已簡化）：
/*
The method regression_rf took 0.06 sec to run.
R^2: 0.9704592485243615
Expected R^2: 0.82 (+/- 0.21)
0. fuel_aer_NG: 0.0
1. fuel_aer_DFO: 0.0
2. fuel_aer_HYC: 0.0
......
8. mover_CT: 0.0
9. mover_GT: 0.0
10. mover_HY: 0.0
11. mover_IC: 0.0
12. mover_PV: 0.0
13. mover_ST: 1.9909417055476134e-05
14. mover_WT: 0.0
15. state_CA: 0.0
16. state_TX: 0.0
......
28. state_OH: 0.0
29. state_GA: 0.0
30. state_WA: 0.0
31. total_fuel_cons: 0.0
32. total_fuel_cons_mmbtu: 0.9999800905829446
*/ 又一次，我們（果然）看到total_fuel_cons_mmbtu起決定作用。mover_ST也出現了，但它的重要性幾乎等于沒有，是以我們可以移除它。
隻用total_fuel_cons_mmbtu構模組化型的結果如下：
/*
The method regression_rf took 0.03 sec to run.
R^2: 0.9704326205779759
Expected R^2: 0.81 (+/- 0.22)
0. total_fuel_cons_mmbtu: 1.0
*/ 如你所見，R2和敏感度的差别可以忽略不計。
6.8使用SVM預測發電廠生産的電量
SVM（Support Vector Machine，支援向量機）已流行多年。模型使用多種核技巧，自變量和因變量之間最複雜的關系也可以模組化。
本技巧中，使用人工生成的資料，我們展示SVM的威力。
準備：需裝好pandas、NumPy、Scikit和Matplotlib。
本技巧中，我們用4種核測試用于回歸的SVM（regression_svm.py檔案）：
1 import sys
2 sys.path.append('..')
3
4 # the rest of the imports
5 import helper as hlp
6 import pandas as pd
7 import numpy as np
8 import sklearn.svm as sv
9 import matplotlib.pyplot as plt
10
11 @hlp.timeit
12 def regression_svm(x, y, **kw_params):
13 '''
14 Estimate a SVM regressor
15 '''
16 # create the regressor object
17 svm = sv.SVR(**kw_params)
18
19 # estimate the model
20 svm.fit(x,y)
21
22 # return the object
23 return svm
24
25 # simulated dataset
26 x = np.arange(-2, 2, 0.004)
27 errors = np.random.normal(0, 0.5, size=len(x))
28 y = 0.8 * x**4 - 2 * x**2 + errors
29
30 # reshape the x array so its in a column form
31 x_reg = x.reshape(-1, 1)
32
33 models_to_test = [
34 {'kernel': 'linear'},
35 {'kernel': 'poly','gamma': 0.5, 'C': 0.5, 'degree': 4},
36 {'kernel': 'poly','gamma': 0.5, 'C': 0.5, 'degree': 6},
37 {'kernel': 'rbf','gamma': 0.5, 'C': 0.5}
38 ] 原理：首先，我們導入所有必需子產品。然後，我們模拟資料。我們先用NumPy的.arange（...）方法生成x向量；資料集從-2到2，步長為0.004。然後引入一些誤差：.random.normal（...）方法傳回和x長度相等的樣本數，x是從一個均值為0标準差為0.5的正态分布中取出的。我們的因變量這樣指定：y=0.8x4-2x2+errors。
為了估算模型，我們要将x轉為列的形式；我們使用.reshape（...）方法重塑x。
models_to_test清單包括了4個關鍵詞參數，用于regression_svm（...）方法；我們估算一個線性核SVM，兩個多項式核SVM（一個度數為4，一個度數為6），以及一個RBF核SVM。
regression_svm（...）方法和所有模型估算方法的邏輯一緻：先建立模型對象，然後估算模型，最後傳回估算的模型。
本技巧中，我們生成一個圖表，以具象化模型的預測。
建立圖，調整框：
1 # create chart frame
2 plt.figure(figsize=(len(models_to_test) * 2 + 3, 9.5))
3 plt.subplots_adjust(left=.05, right=.95,
4 bottom=.05, top=.96, wspace=.1, hspace=.15) 建立指定大小的圖，調整子圖：為每個子圖調整left和right，bottom和top邊距，以及子圖間的空間，wspace和hspace。
現在周遊要建構的模型，加到圖中：
1 for i, model in enumerate(models_to_test):
2 # estimate the model
3 regressor = regression_svm(x_reg, y, **model)
4
5 # score
6 score = regressor.score(x_reg, y)
7
8 # plot the chart
9 plt.subplot(2, 2, i + 1)
10 if model['kernel'] == 'poly':
11 plt.title('Kernel: {0}, deg: {1}'\
12 .format(model['kernel'], model['degree']))
13 else:
14 plt.title('Kernel: {0}'.format(model['kernel']))
15 plt.ylim([-4, 8])
16 plt.scatter(x, y)
17 plt.plot(x, regressor.predict(x_reg), color='r')
18 plt.text(.9, .9, ('R^2: {0:.2f}'.format(score)),
19 transform=plt.gca().transAxes, size=15,
20 horizontalalignment='right') 先用regression_svm（...）方法估算模型。估算好後，計算R^2分值。
然後用.subplot（...）方法指定（2乘2的網格上）子圖的位置，将目前的循環次數i加1；圖表中，子圖從1開始計數，而清單索引從0開始。如果是多項式核，圖表的标題中應展現度數；否則給出核的名稱就足夠了。
然後，我們用原始資料x和y建立散布圖。之後，我們用回歸器的預測值繪圖；為了醒目，我們使用紅色。在各個子圖中，y軸的範圍保持在-4到8。
然後，我們在右上角放上R^2的值；指定坐标為（.9，.9）并使用轉換可以做到這一點。字型設為15，右對齊。
最後儲存到檔案：
plt.savefig('../../data/Chapter06/charts/regression_svm.png',
dpi= 300) 從得到的圖中可以看出，RBF核對于高度非線性資料表現有多好。我們又一次看到了線性核的負的R2值。可以注意到，多項式核的表現也并沒有足夠好：
更多：使用RBF核估算SVM的主要難點在于決定為C（誤差懲罰參數）和gamma（控制對特征向量差異的敏感度的RBF核參數）參數選擇何值。
最後，我們可以使用格點搜尋法，一個計算密集型的做法。不過，我們将使用Optunity，一個多種環境下各種優化器的庫。
關于Optunity可以參考：http://optunity.readthedocs.org/en/latest/。
Anaconda不直接支援Optunity，我們需要自行安裝。
獨立安裝 optunity
>pip install optunity
/*
Looking in indexes: https://pypi.tuna.tsinghua.edu.cn/simple
Collecting optunity
Downloading https://pypi.tuna.tsinghua.edu.cn/packages/32/4d/d49876a49e105b56755eb5ba06a4848ee8010f7ff9e0f11a13aefed12063/Optunity-1.1.1.tar.gz (4.6MB)
.......
Successfully built optunity
Installing collected packages: optunity
Successfully installed optunity-1.1.1
FINISHED
*/ 既然有了Optunity，我們來看看代碼（regression_svm_alternative.py檔案）。這裡我跳過了各種引入；完整的引入清單可以檢視檔案：
1 import optunity as opt
2 import optunity.metrics as mt
3
4 # simulated dataset
5 x = np.arange(-2, 2, 0.002)
6 errors = np.random.normal(0, 0.5, size=len(x))
7 y = 0.8 * x**4 - 2 * x**2 + errors
8
9 # reshape the x array so its in a column form
10 x_reg = x.reshape(-1, 1)
11
12 @opt.cross_validated(x=x_reg, y=y, num_folds=10, num_iter=5)
13 def regression_svm(
14 x_train, y_train, x_test, y_test, logC, logGamma):
15 '''
16 Estimate a SVM regressor
17 '''
18 # create the regressor object
19 svm = sv.SVR(kernel='rbf',
20 C=0.1 * logC, gamma=0.1 * logGamma)
21
22 # estimate the model
23 svm.fit(x_train,y_train)
24
25 # decision function
26 decision_values = svm.decision_function(x_test)
27
28 # return the object
29 return mt.roc_auc(y_test, decision_values) 首先，類似前一個例子，我們建立模拟資料。
然後用一個Optunity中的方法裝飾regression_svm（...）方法。.cross_validated（...）方法接受一系列參數：x是自變量，y是因變量。
num_folds是用于交叉驗證的塊數。我們将1/10的樣本用于交叉驗證，9/10用于訓練。num_iter是交叉驗證的循環次數。
注意regression_svm（...）方法的定義變了。cross_validated（...）方法先注入訓練的x和y樣本，然後是測試的x和y，以及指定的其他參數：我們的例子中，logC和logGamma。對于後兩個參數我們要為其尋找最優值。
方法從估算模型開始：先用指定的參數建立模型，應用資料。decision_function（...）方法傳回一個清單，放的是x中元素到超平面的距離。我們傳回的是在ROC曲線之下的區域。（ROC參考本書3.2節）
使用Optunity提供的maximize（...）方法找到C和gamma的最佳值。第一個參數是我們要最大化的函數，即regression_svm（...）。num_evals指定了允許的函數評估的次數。最後兩個指定了logC和logGamma搜尋空間的邊界。
優化之後，我們就可以用找到的值估算最終的模型并儲存圖檔：
1 # and the values are...
2 print('The optimal values are: C - {0:.2f}, gamma - {1:.2f}'\
3 .format(0.1 * hps['logC'], 0.1 * hps['logGamma']))
4
5 # estimate the model with optimal values
6 regressor = sv.SVR(kernel='rbf',
7 C=0.1 * hps['logC'],
8 gamma=0.1 * hps['logGamma'])\
9 .fit(x_reg, y)
10
11 # predict the output
12 predicted = regressor.predict(x_reg)
13
14 # and calculate the R^2
15 score = hlp.get_score(y, predicted)
16 print('R2: ', score)
17
18 # plot the chart
19 plt.scatter(x, y)
20 plt.plot(x, predicted, color='r')
21 plt.title('Kernel: RBF')
22 plt.ylim([-4, 8])
23 plt.text(.9, .9, ('R^2: {0:.2f}'.format(score)),
24 transform=plt.gca().transAxes, size=15,
25 horizontalalignment='right')
26
27 plt.savefig(
28 '../../data/Chapter06/charts/regression_svm_alt.png',
29 dpi=300
30 ) 大部分代碼之前已經解釋過了，這裡不再重複。
得到的值應該是這樣：
/*
D:\tools\Python37\lib\site-packages\optunity\metrics.py:204: RuntimeWarning: invalid value encountered in true_divide
return float(TP) / (TP + FN)
...
D:\tools\Python37\lib\site-packages\optunity\metrics.py:204: RuntimeWarning: invalid value encountered in true_divide
return float(TP) / (TP + FN)
The optimal values are: C - 0.38, gamma - 1.44
R^2: 0.8662968996734814
*/ 和前一個例子相比，R^2值上升了；現在回歸更比對資料了：
參考：強烈推薦熟悉一下Optunity：http://optunity.readthedocs.org/en/latest/。
6.9訓練神經網絡，預測發電廠生産的電量
神經網絡是強有力的回歸器，可以應用到任何資料上。然而，本例中，我們已經知道資料遵循線性關系了，是以就使用一個簡單的模型來應用資料。
準備：需裝好pandas和PyBrain。
類似第3章中的方式，我們指定fitANN（...）方法（regression_ann.py檔案）。有些引入語句已略去：
1 import sys
2 sys.path.append('..')
3
4 # the rest of the imports
5 import helper as hlp
6 import pandas as pd
7 import pybrain.structure as st
8 import pybrain.supervised.trainers as tr
9 import pybrain.tools.shortcuts as pb
10
11 @hlp.timeit
12 def fitANN(data):
13 '''
14 Build a neural network regressor
15 '''
16 # determine the number of inputs and outputs
17 inputs_cnt = data['input'].shape[1]
18 target_cnt = data['target'].shape[1]
19
20 # create the regressor object
21 ann = pb.buildNetwork(inputs_cnt,
22 inputs_cnt * 3,
23 target_cnt,
24 hiddenclass=st.TanhLayer,
25 outclass=st.LinearLayer,
26 bias=True
27 )
28
29 # create the trainer object
30 trainer = tr.BackpropTrainer(ann, data,
31 verbose=True, batchlearning=False)
32
33 # and train the network
34 trainer.trainUntilConvergence(maxEpochs=50, verbose=True,
35 continueEpochs=2, validationProportion=0.25)
36
37 # and return the regressor
38 return ann 原理：我們先讀入資料，拆分資料，并建立ANN訓練資料：
1 # split the data into training and testing
2 train_x, train_y, \
3 test_x, test_y, \
4 labels = hlp.split_data(csv_read,
5 y='net_generation_MWh', x=['total_fuel_cons_mmbtu'])
6
7 # create the ANN training and testing datasets
8 training = hlp.prepareANNDataset((train_x, train_y),
9 prob='regression')
10 testing = hlp.prepareANNDataset((test_x, test_y),
11 prob='regression') 自變量隻使用total_fuel_cons_mmbtu。注意.prepareANNDataset（...）方法使用了新參數。
prob（代表problem，問題）指定了要用資料集解決一個回歸問題：對于分類問題，我們需要有類别名額及其到1的補充abs（item-1），參考本書3.8節。
現在可以應用模型了：
1 # train the model
2 regressor = fitANN(training)
3
4 # predict the output from the unseen data
5 predicted = regressor.activateOnDataset(testing)
6
7 # and calculate the R^2
8 score = hlp.get_score(test_y, predicted[:, 0])
9 print('R^2: ', score) fitANN（...）方法先建立神經網絡。我們的神經網絡有一個輸入神經元，三個（帶tanh激活函數的）隐藏層神經元，一個（帶線性激活函數的）輸出神經元。
和第3章一樣使用線上訓練（batchlearning=False）。我們指定maxEpochs為50，即使到時候網絡沒有達到最小的誤差，我們也終止訓練。我們将四分之一的資料用于驗證。
訓練之後，我們通過預測發電量來測試網絡，計算R^2：
/*
The method fitANN took 68.70 sec to run.
R^2: 0.9843443717345242
*/ 可以看出，我們神經網絡的表現不輸其他模型。
參考
/* Traceback (most recent call last):
File "D:\Java2018\practicalDataAnalysis\Codes\Chapter06\regression_ann.py", line 8, in <module>
import pybrain.structure as st
File "D:\tools\Python37\lib\site-packages\pybrain\__init__.py", line 1, in <module>
from structure.__init__ import *
ModuleNotFoundError: No module named 'structure'
*/ 解決方案（均失敗）：
/*
pip install git+https://github.com/pybrain/pybrain.git
This seems to be an issue with the library when used from Python 3 and installed using pip. The latest version listed is 0.3.3, but I'm still getting v0.3.0 when installing.
Try installing numpy and scipy in Anaconda, then (after activating the conda environment) install directly from the git repo:
pip install git+https://github.com/pybrain/[email protected]
嘗試方法1：pip install git+https://github.com/pybrain/[email protected]
Running command git clone -q https://github.com/pybrain/pybrain.git 'd:\Temp\pip-req-build-0soc1yjx'
Looking in indexes: https://pypi.tuna.tsinghua.edu.cn/simple
Collecting git+https://github.com/pybrain/[email protected]
Cloning https://github.com/pybrain/pybrain.git (to revision 0.3.3) to d:\temp\pip-req-build-0soc1yjx
Running command git checkout -q 882cf66a89b6b6bb8112f41d200c55461464aa06
Requirement already satisfied: scipy in d:\tools\python37\lib\site-packages (from PyBrain==0.3.1) (1.3.3)
Requirement already satisfied: numpy>=1.13.3 in d:\tools\python37\lib\site-packages (from scipy->PyBrain==0.3.1) (1.17.4)
Building wheels for collected packages: PyBrain
Building wheel for PyBrain (setup.py): started
Building wheel for PyBrain (setup.py): finished with status 'done'
Created wheel for PyBrain: filename=PyBrain-0.3.1-cp37-none-any.whl size=472111 sha256=07f80767b5065f40565e595826df759d2dc66740d3fb6db709d7343cbdc2a8c0
Stored in directory: d:\Temp\pip-ephem-wheel-cache-h4bvm49q\wheels\a3\1c\d1\9cfc0e65ef0ed5559fd4f2819e4423e9fa4cfe02ff3a5c3b3e
Successfully built PyBrain
Installing collected packages: PyBrain
Found existing installation: PyBrain 0.3
Uninstalling PyBrain-0.3:
Successfully uninstalled PyBrain-0.3
Successfully installed PyBrain-0.3.1
FINISHED
嘗試方法2：
D:\tools\Python37\Scripts>pip3 install https://github.com/pybrain/pybrain/archive/0.3.3.zip
Looking in indexes: https://pypi.tuna.tsinghua.edu.cn/simple
Collecting https://github.com/pybrain/pybrain/archive/0.3.3.zip
Downloading https://github.com/pybrain/pybrain/archive/0.3.3.zip
/ 1.9MB 163kB/s
Requirement already satisfied: scipy in d:\tools\python37\lib\site-packages (from PyBrain==0.3.1) (1.3.3)
Requirement already satisfied: numpy>=1.13.3 in d:\tools\python37\lib\site-packages (from scipy->PyBrain==0.3.1) (1.17.4)
Building wheels for collected packages: PyBrain
Building wheel for PyBrain (setup.py) ... done
Created wheel for PyBrain: filename=PyBrain-0.3.1-cp37-none-any.whl size=468236 sha256=ce9591da503451919d71cd322324930d3d2259c7fdbe7c3041a8309bcf1ff9b2
Stored in directory: d:\Temp\pip-ephem-wheel-cache-ouwcf150\wheels\0b\04\38\2f174aa3c578350870947ca6ab12e0eb89aef3478c9610eb0a
Successfully built PyBrain
Installing collected packages: PyBrain
Successfully installed PyBrain-0.3.1
D:\tools\Python37\Scripts>
*/ 最終反複檢查：
pip install https://github.com/pybrain/pybrain/archive/0.3.3.zip
安裝0.3.3，其實還是0.3.1
/*
File "D:\tools\Python37\lib\site-packages\pybrain\tools\functions.py", line 4, in <module>
from scipy.linalg import inv, det, svd, logm, expm2
ImportError: cannot import name 'expm2' from 'scipy.linalg' (D:\tools\Python37\lib\site-packages\scipy\linalg\__init__.py)
*/ 解決方案補充：
根據提示錯誤是“D:\tools\Python37\lib\site-packages\pybrain\tools\functions.py”檔案中的import scipy.linalg import inv, det, svd, logm, expm2出錯，查找到一些英文資料說明SciPy中expm2和expm3魯棒性很差，官方建議不再使用這兩個，而是使用expm，是以，在前的所提到的路徑中找到functions.py檔案，找到語句“from scipy.linalg import inv, det, svd, logm, expm2”，将expm2更改成expm，儲存，退出，再次嘗試，引用成功。
/* Total error: 0.0234644034879
Total error: 0.00448033176755
Total error: 0.00205475261566
Total error: 0.000943766244883
Total error: 0.000441922091808
Total error: 0.00023137206216
Total error: 0.00014716265305
Total error: 0.000113546124273
Total error: 0.000100591874178
Total error: 9.48731641822e-05
Total error: 9.17187548879e-05
Total error: 9.00545474786e-05
Total error: 8.8754314954e-05
Total error: 8.78622211966e-05
Total error: 8.64148298491e-05
Total error: 8.52580031627e-05
Total error: 8.45496751335e-05
Total error: 8.31773350186e-05
Total error: 8.22848381841e-05
Total error: 8.13736136075e-05
Total error: 8.01448506931e-05
Total error: 7.91777696701e-05
Total error: 7.82493713543e-05
Total error: 7.75330633562e-05
Total error: 7.63823501292e-05
Total error: 7.55286037883e-05
Total error: 7.44680561538e-05
Total error: 7.36680418929e-05
Total error: 7.27012120487e-05
Total error: 7.18980979088e-05
Total error: 7.09005144214e-05
Total error: 6.99157006989e-05
Total error: 6.9222957236e-05
Total error: 6.82360145737e-05
Total error: 6.71353795523e-05
Total error: 6.63631396761e-05
Total error: 6.56685823709e-05
Total error: 6.53158010387e-05
Total error: 6.4518956904e-05
Total error: 6.35146379129e-05
Total error: 6.28304190101e-05
Total error: 6.22125569291e-05
Total error: 6.13769678787e-05
Total error: 6.09232354752e-05
Total error: 6.0155987076e-05
Total error: 5.94907225251e-05
Total error: 5.88562745717e-05
Total error: 5.80850744217e-05
Total error: 5.75905974532e-05
Total error: 5.65098746457e-05
Total error: 5.60159718614e-05
('train-errors:', '[0.023464 , 0.00448 , 0.002055 , 0.000944 , 0.000442 , 0.000231 , 0.000147 , 0.000114 , 0.000101 , 9.5e-05 , 9.2e-05 , 9e-05 , 8.9e-05 , 8.8e-05 , 8.6e-05 , 8.5e-05 , 8.5e-05 , 8.3e-05 , 8.2e-05 , 8.1e-05 , 8e-05 , 7.9e-05 , 7.8e-05 , 7.8e-05 , 7.6e-05 , 7.6e-05 , 7.4e-05 , 7.4e-05 , 7.3e-05 , 7.2e-05 , 7.1e-05 , 7e-05 , 6.9e-05 , 6.8e-05 , 6.7e-05 , 6.6e-05 , 6.6e-05 , 6.5e-05 , 6.5e-05 , 6.4e-05 , 6.3e-05 , 6.2e-05 , 6.1e-05 , 6.1e-05 , 6e-05 , 5.9e-05 , 5.9e-05 , 5.8e-05 , 5.8e-05 , 5.7e-05 , 5.6e-05 ]')
('valid-errors:', '[2.98916 , 0.006323 , 0.003087 , 0.001464 , 0.000735 , 0.00041 , 0.000263 , 0.000197 , 0.000168 , 0.000154 , 0.000147 , 0.000144 , 0.00014 , 0.000143 , 0.000136 , 0.000132 , 0.000131 , 0.000131 , 0.000136 , 0.000127 , 0.000126 , 0.000124 , 0.000124 , 0.000122 , 0.000123 , 0.000121 , 0.00012 , 0.000117 , 0.000116 , 0.000118 , 0.000114 , 0.000113 , 0.000112 , 0.000117 , 0.00011 , 0.000108 , 0.000107 , 0.000107 , 0.000105 , 0.000105 , 0.000105 , 0.000103 , 0.000101 , 0.000102 , 0.000103 , 9.8e-05 , 9.7e-05 , 9.6e-05 , 9.5e-05 , 9.6e-05 , 9.7e-05 , 9.4e-05 ]')
The method fitANN took 68.70 sec to run.
R^2: 0.9843443717345242
*/ 檢視使用PyBrain訓練神經網絡的例子：http://fastml.com/pybrain-a-simple-neural-networks-library-in-python/。
第6章完。
杭州這個季節的窗外細雨打在窗上，感覺有點陰冷。東方既白，還要早起找警察叔叔。
随書源碼官方下載下傳：
http://www.hzcourse.com/web/refbook/detail/7821/92
邀月注：本文版權由邀月和部落格園共同所有，轉載請注明出處。
助人等于自助! [email protected]