1、什麽是多元線性回歸模型?
當y值的影響因素不唯壹時,采用多元線性回歸模型。
y =y=β0+β1x1+β2x2+...+βnxn?
例如商品的銷售額可能不電視廣告投入,收音機廣告投入,報紙廣告投入有關系,可以有 sales =β0+β1*TV+β2* radio+β3*newspaper.
2、使用pandas來讀取數據
pandas 是壹個用於數據探索、數據分析和數據處理的python庫
[python]?view plain?copy
import?pandas?as?pd?
[html]?view plain?copy
<pre?name="code"?class="python">#?read?csv?file?directly?from?a?URL?and?save?the?results
data?=?pd.read_csv('/home/lulei/Advertising.csv')?
#?display?the?first?5?rows?
data.head()?
上面代碼的運行結果:
TV ?Radio ?Newspaper ?Sales0 ?230.1 ? 37.8 ? 69.2 ? 22.1
1 ? 44.5 ? 39.3 ? 45.1 ? 10.4
2 ? 17.2 ? 45.9 ? 69.3 9.3
3 ?151.5 ? 41.3 ? 58.5 ? 18.5
4 ?180.8 ? 10.8 ? 58.4 ? 12.9
上面顯示的結果類似壹個電子表格,這個結構稱為Pandas的數據幀(data frame),類型全稱:pandas.core.frame.DataFrame.
pandas的兩個主要數據結構:Series和DataFrame:
Series類似於壹維數組,它有壹組數據以及壹組與之相關的數據標簽(即索引)組成。
DataFrame是壹個表格型的數據結構,它含有壹組有序的列,每列可以是不同的值類型。DataFrame既有行索引也有列索引,它可以被看做由Series組成的字典。
[python]?view plain?copy
#?display?the?last?5?rows?
data.tail()?
只顯示結果的末尾5行
TV ?Radio ?Newspaper ?Sales195 ? 38.2 3.7 ? 13.8 7.6
196 ? 94.2 4.9 8.1 9.7
197 ?177.0 9.3 6.4 ? 12.8
198 ?283.6 ? 42.0 ? 66.2 ? 25.5
199 ?232.1 8.6 8.7 ? 13.4
[html]?view plain?copy
#?check?the?shape?of?the?DataFrame(rows,?colums)?
data.shape?
查看DataFrame的形狀,註意第壹列的叫索引,和數據庫某個表中的第壹列類似。
(200,4)?
3、分析數據
特征:
TV:對於壹個給定市場中單壹產品,用於電視上的廣告費用(以千為單位)
Radio:在廣播媒體上投資的廣告費用
Newspaper:用於報紙媒體的廣告費用
響應:
Sales:對應產品的銷量
在這個案例中,我們通過不同的廣告投入,預測產品銷量。因為響應變量是壹個連續的值,所以這個問題是壹個回歸問題。數據集壹***有200個觀測值,每壹組觀測對應壹個市場的情況。
註意:這裏推薦使用的是seaborn包。網上說這個包的數據可視化效果比較好看。其實seaborn也應該屬於matplotlib的內部包。只是需要再次的單獨安裝。
[python]?view plain?copy
import?seaborn?as?sns?
import?matplotlib.pyplot?as?plt
#?visualize?the?relationship?between?the?features?and?the?response?using?scatterplots?
sns.pairplot(data,?x_vars=['TV','Radio','Newspaper'],?y_vars='Sales',?size=7,?aspect=0.8)?
plt.show()#註意必須加上這壹句,否則無法顯示。?
[html]?view plain?copy
這裏選擇TV、Radio、Newspaper?作為特征,Sales作為觀測值?
[html]?view plain?copy
返回的結果:?
seaborn的pairplot函數繪制X的每壹維度和對應Y的散點圖。通過設置size和aspect參數來調節顯示的大小和比例。可以從圖中看出,TV特征和銷量是有比較強的線性關系的,而Radio和Sales線性關系弱壹些,Newspaper和Sales線性關系更弱。通過加入壹個參數kind='reg',seaborn可以添加壹條最佳擬合直線和95%的置信帶。
[python]?view plain?copy
sns.pairplot(data,?x_vars=['TV','Radio','Newspaper'],?y_vars='Sales',?size=7,?aspect=0.8,?kind='reg')?
plt.show()?
結果顯示如下:
4、線性回歸模型
優點:快速;沒有調節參數;可輕易解釋;可理解。
缺點:相比其他復雜壹些的模型,其預測準確率不是太高,因為它假設特征和響應之間存在確定的線性關系,這種假設對於非線性的關系,線性回歸模型顯然不能很好的對這種數據建模。
線性模型表達式:?y=β0+β1x1+β2x2+...+βnxn?其中
y是響應
β0是截距
β1是x1的系數,以此類推
在這個案例中:?y=β0+β1?TV+β2?Radio+...+βn?Newspaper
(1)、使用pandas來構建X(特征向量)和y(標簽列)
scikit-learn要求X是壹個特征矩陣,y是壹個NumPy向量。
pandas構建在NumPy之上。
因此,X可以是pandas的DataFrame,y可以是pandas的Series,scikit-learn可以理解這種結構。
[python]?view plain?copy
#create?a?python?list?of?feature?names?
feature_cols?=?['TV',?'Radio',?'Newspaper']?
#?use?the?list?to?select?a?subset?of?the?original?DataFrame?
X?=?data[feature_cols]?
#?equivalent?command?to?do?this?in?one?line?
X?=?data[['TV',?'Radio',?'Newspaper']]?
#?print?the?first?5?rows?
print?X.head()?
#?check?the?type?and?shape?of?X?
print?type(X)?
print?X.shape?
輸出結果如下:
?TV ?Radio ?Newspaper0 ?230.1 ? 37.8 ? 69.2
1 ? 44.5 ? 39.3 ? 45.1
2 ? 17.2 ? 45.9 ? 69.3
3 ?151.5 ? 41.3 ? 58.5
4 ?180.8 ? 10.8 ? 58.4
<class 'pandas.core.frame.DataFrame'>
(200, 3)
[python]?view plain?copy
#?select?a?Series?from?the?DataFrame?
y?=?data['Sales']?
#?equivalent?command?that?works?if?there?are?no?spaces?in?the?column?name?
y?=?data.Sales?
#?print?the?first?5?values?
print?y.head()?
輸出的結果如下:
0 22.11 10.4
2 9.3
3 18.5
4 12.9
Name: Sales
(2)、構建訓練集與測試集
[html]?view plain?copy
<pre?name="code"?class="python"><span?style="font-size:14px;">##構造訓練集和測試集?
from?sklearn.cross_validation?import?train_test_split?#這裏是引用了交叉驗證?
X_train,X_test,?y_train,?y_test?=?train_test_split(X,?y,?random_state=1)?
#default split is 75% for training and 25% for testing
[html]?view plain?copy
print?X_train.shape?
print?y_train.shape?
print?X_test.shape?
print?y_test.shape?
輸出結果如下:
(150, 3)(150,)
(50, 3)
(50,)
註:上面的結果是由train_test_spilit()得到的,但是我不知道為什麽我的版本的sklearn包中居然報錯:
ImportError ? Traceback (most recent call last)<ipython-input-182-3eee51fcba5a> in <module>() ?1 ###構造訓練集和測試集----> 2 from sklearn.cross_validation import train_test_split ?3 #import sklearn.cross_validation ?4 X_train,X_test, y_train, y_test = train_test_split(X, y, random_state=1) ?5 # default split is 75% for training and 25% for testingImportError: cannot import name train_test_split處理方法:1、我後來重新安裝sklearn包。再壹次調用時就沒有錯誤了。
2、自己寫函數來認為的隨機構造訓練集和測試集。(這個代碼我會在最後附上。)
(3)sklearn的線性回歸
[html]?view plain?copy
from?sklearn.linear_model?import?LinearRegression?
linreg?=?LinearRegression()?
model=linreg.fit(X_train,?y_train)?
print?model?
print?linreg.intercept_?
print?linreg.coef_?
輸出的結果如下:
LinearRegression(copy_X=True, fit_intercept=True, normalize=False)2.66816623043
[ 0.04641001 ?0.19272538 -0.00349015]
[html]?view plain?copy
#?pair?the?feature?names?with?the?coefficients?
zip(feature_cols,?linreg.coef_)?
輸出如下:
[('TV', 0.046410010869663267),('Radio', 0.19272538367491721),
('Newspaper', -0.0034901506098328305)]
y=2.668+0.0464?TV+0.192?Radio-0.00349?Newspaper
如何解釋各個特征對應的系數的意義?
對於給定了Radio和Newspaper的廣告投入,如果在TV廣告上每多投入1個單位,對應銷量將增加0.0466個單位。就是加入其它兩個媒體投入固定,在TV廣告上每增加1000美元(因為單位是1000美元),銷量將增加46.6(因為單位是1000)。但是大家註意這裏的newspaper的系數居然是負數,所以我們可以考慮不使用newspaper這個特征。這是後話,後面會提到的。
(4)、預測
[python]?view plain?copy
y_pred?=?linreg.predict(X_test)?
print?y_pred?
[python]?view plain?copy
print?type(y_pred)?
輸出結果如下:
[ 14.58678373 ? 7.92397999 ?16.9497993 ? 19.35791038 ? 7.36360284
7.35359269 ?16.08342325 ? 9.16533046 ?20.35507374 ?12.63160058 22.83356472 ? 9.66291461 ? 4.18055603 ?13.70368584 ?11.4533557 4.16940565 ?10.31271413 ?23.06786868 ?17.80464565 ?14.53070132 15.19656684 ?14.22969609 ? 7.54691167 ?13.47210324 ?15.00625898 19.28532444 ?20.7319878 ? 19.70408833 ?18.21640853 ? 8.50112687 9.8493781 9.51425763 ? 9.73270043 ?18.13782015 ?15.41731544 5.07416787 ?12.20575251 ?14.05507493 ?10.6699926 7.16006245 11.80728836 ?24.79748121 ?10.40809168 ?24.05228404 ?18.44737314 20.80572631 ? 9.45424805 ?17.00481708 ? 5.78634105 ? 5.10594849]<type 'numpy.ndarray'>
5、回歸問題的評價測度
(1) 評價測度
對於分類問題,評價測度是準確率,但這種方法不適用於回歸問題。我們使用針對連續數值的評價測度(evaluation metrics)。
這裏介紹3種常用的針對線性回歸的測度。
1)平均絕對誤差(Mean Absolute Error, MAE)
(2)均方誤差(Mean Squared Error, MSE)
(3)均方根誤差(Root Mean Squared Error, RMSE)
這裏我使用RMES。
[python]?view plain?copy
<pre?name="code"?class="python">#計算Sales預測的RMSE?
print?type(y_pred),type(y_test)?
print?len(y_pred),len(y_test)?
print?y_pred.shape,y_test.shape?
from?sklearn?import?metrics?
import?numpy?as?np?
sum_mean=0?
for?i?in?range(len(y_pred)):?
sum_mean+=(y_pred[i]-y_test.values[i])**2?
sum_erro=np.sqrt(sum_mean/50)?
#?calculate?RMSE?by?hand?
print?"RMSE?by?hand:",sum_erro?
最後的結果如下:
<type 'numpy.ndarray'> <class 'pandas.core.series.Series'>50 50
(50,) (50,)
RMSE by hand: 1.42998147691
(2)做ROC曲線
[python]?view plain?copy
import?matplotlib.pyplot?as?plt?
plt.figure()?
plt.plot(range(len(y_pred)),y_pred,'b',label="predict")?
plt.plot(range(len(y_pred)),y_test,'r',label="test")?
plt.legend(loc="upper?right")?#顯示圖中的標簽?
plt.xlabel("the?number?of?sales")?
plt.ylabel('value?of?sales')?
plt.show()?
顯示結果如下:(紅色的線是真實的值曲線,藍色的是預測值曲線)
直到這裏整個的壹次多元線性回歸的預測就結束了。
6、改進特征的選擇
在之前展示的數據中,我們看到Newspaper和銷量之間的線性關系竟是負關系(不用驚訝,這是隨機特征抽樣的結果。換壹批抽樣的數據就可能為正了),現在我們移除這個特征,看看線性回歸預測的結果的RMSE如何?
依然使用我上面的代碼,但只需修改下面代碼中的壹句即可:
[python]?view plain?copy
#create?a?python?list?of?feature?names?
feature_cols?=?['TV',?'Radio',?'Newspaper']?
#?use?the?list?to?select?a?subset?of?the?original?DataFrame?
X?=?data[feature_cols]?
#?equivalent?command?to?do?this?in?one?line?
#X?=?data[['TV',?'Radio',?'Newspaper']]#只需修改這裏即可<pre?name="code"?class="python"?style="font-size:?15px;?line-height:?35px;">X?=?data[['TV',?'Radio']]?#去掉newspaper其他的代碼不變?
# print the first 5 rowsprint X.head()# check the type and shape of Xprint type(X)print X.shape
最後的到的系數與測度如下:
LinearRegression(copy_X=True, fit_intercept=True, normalize=False)
2.81843904823[ 0.04588771 ?0.18721008]
RMSE by hand: 1.28208957507
然後再次使用ROC曲線來觀測曲線的整體情況。我們在將Newspaper這個特征移除之後,得到RMSE變小了,說明Newspaper特征可能不適合作為預測銷量的特征,於是,我們得到了新的模型。我們還可以通過不同的特征組合得到新的模型,看看最終的誤差是如何的。
備註:
之前我提到了這種錯誤:
註:上面的結果是由train_test_spilit()得到的,但是我不知道為什麽我的版本的sklearn包中居然報錯:
ImportError ? Traceback (most recent call last)<ipython-input-182-3eee51fcba5a> in <module>() ?1 ###構造訓練集和測試集----> 2 from sklearn.cross_validation import train_test_split ?3 #import sklearn.cross_validation ?4 X_train,X_test, y_train, y_test = train_test_split(X, y, random_state=1) ?5 # default split is 75% for training and 25% for testingImportError: cannot import name train_test_split處理方法:1、我後來重新安裝sklearn包。再壹次調用時就沒有錯誤了。
2、自己寫函數來認為的隨機構造訓練集和測試集。(這個代碼我會在最後附上。)
這裏我給出我自己寫的函數: