可以。
成都創(chuàng)新互聯(lián)基于分布式IDC數(shù)據(jù)中心構(gòu)建的平臺為眾多戶提供服務(wù)器托管雅安 四川大帶寬租用 成都機柜租用 成都服務(wù)器租用。
使用Python自帶的sum函數(shù)，sum函數(shù)是個內(nèi)置函數(shù)，可以求一個數(shù)字列表的和，并且可以帶初始值，如果不帶初始值的話，默認(rèn)是0。
首個參數(shù)為可迭代的列表，初始值默認(rèn)為0，也可以為其他值，比如說[]，空列表在Python里面，類型是動態(tài)類型，一種操作或接口，到底做何操作取決于對象本身比如說同樣是+，如果兩者都是數(shù)字1+1=2，如果兩者都是字符串，則'1'+'1'='11'所以如果這里的start本身為[]，則會執(zhí)行列表合并的操作。
1、什么是多元線性回歸模型？
當(dāng)y值的影響因素不唯一時,采用多元線性回歸模型。
y =y=β0+β1x1+β2x2+...+βnxn
例如商品的銷售額可能不電視廣告投入,收音機廣告投入,報紙廣告投入有關(guān)系,可以有 sales =β0+β1*TV+β2* radio+β3*newspaper.
2、使用pandas來讀取數(shù)據(jù)
pandas 是一個用于數(shù)據(jù)探索、數(shù)據(jù)分析和數(shù)據(jù)處理的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()
上面代碼的運行結(jié)果：
TV ?Radio ?Newspaper ?Sales
0 ?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
上面顯示的結(jié)果類似一個電子表格，這個結(jié)構(gòu)稱為Pandas的數(shù)據(jù)幀(data frame)，類型全稱：pandas.core.frame.DataFrame.
pandas的兩個主要數(shù)據(jù)結(jié)構(gòu)：Series和DataFrame：
Series類似于一維數(shù)組，它有一組數(shù)據(jù)以及一組與之相關(guān)的數(shù)據(jù)標(biāo)簽(即索引)組成。
DataFrame是一個表格型的數(shù)據(jù)結(jié)構(gòu)，它含有一組有序的列，每列可以是不同的值類型。DataFrame既有行索引也有列索引，它可以被看做由Series組成的字典。
[python]?view plain?copy
#?display?the?last?5?rows
data.tail()
只顯示結(jié)果的末尾5行
?TV ?Radio ?Newspaper ?Sales
195 ? 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的形狀,注意第一列的叫索引，和數(shù)據(jù)庫某個表中的第一列類似。
(200,4)?
3、分析數(shù)據(jù)
特征：
TV：對于一個給定市場中單一產(chǎn)品，用于電視上的廣告費用（以千為單位）
Radio：在廣播媒體上投資的廣告費用
Newspaper：用于報紙媒體的廣告費用
響應(yīng)：
Sales：對應(yīng)產(chǎn)品的銷量
在這個案例中，我們通過不同的廣告投入，預(yù)測產(chǎn)品銷量。因為響應(yīng)變量是一個連續(xù)的值，所以這個問題是一個回歸問題。數(shù)據(jù)集一共有200個觀測值，每一組觀測對應(yīng)一個市場的情況。
注意：這里推薦使用的是seaborn包。網(wǎng)上說這個包的數(shù)據(jù)可視化效果比較好看。其實seaborn也應(yīng)該屬于matplotlib的內(nèi)部包。只是需要再次的單獨安裝。
[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
返回的結(jié)果：
seaborn的pairplot函數(shù)繪制X的每一維度和對應(yīng)Y的散點圖。通過設(shè)置size和aspect參數(shù)來調(diào)節(jié)顯示的大小和比例?？梢詮膱D中看出，TV特征和銷量是有比較強的線性關(guān)系的，而Radio和Sales線性關(guān)系弱一些，Newspaper和Sales線性關(guān)系更弱。通過加入一個參數(shù)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()
結(jié)果顯示如下：
4、線性回歸模型
優(yōu)點：快速；沒有調(diào)節(jié)參數(shù)；可輕易解釋；可理解。
缺點：相比其他復(fù)雜一些的模型，其預(yù)測準(zhǔn)確率不是太高，因為它假設(shè)特征和響應(yīng)之間存在確定的線性關(guān)系，這種假設(shè)對于非線性的關(guān)系，線性回歸模型顯然不能很好的對這種數(shù)據(jù)建模。
線性模型表達(dá)式：?y=β0+β1x1+β2x2+...+βnxn?其中
y是響應(yīng)
β0是截距
β1是x1的系數(shù)，以此類推
在這個案例中：?y=β0+β1?TV+β2?Radio+...+βn?Newspaper
(1)、使用pandas來構(gòu)建X(特征向量)和y(標(biāo)簽列)
scikit-learn要求X是一個特征矩陣，y是一個NumPy向量。
pandas構(gòu)建在NumPy之上。
因此，X可以是pandas的DataFrame，y可以是pandas的Series，scikit-learn可以理解這種結(jié)構(gòu)。
[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
輸出結(jié)果如下：
TV ?Radio ?Newspaper
0 ?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()
輸出的結(jié)果如下：
0 ? ?22.1
1 ? ?10.4
2 ? ? 9.3
3 ? ?18.5
4 ? ?12.9
Name: Sales
（2）、構(gòu)建訓(xùn)練集與測試集
[html]?view plain?copy
pre?name="code"?class="python"span?style="font-size:14px;"##構(gòu)造訓(xùn)練集和測試集
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
輸出結(jié)果如下：
(150, 3)
(150,)
(50, 3)
(50,)
注：上面的結(jié)果是由train_test_spilit()得到的，但是我不知道為什么我的版本的sklearn包中居然報錯：
ImportError ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? Traceback (most recent call last)ipython-input-182-3eee51fcba5a in module() ? ? ?1 ###構(gòu)造訓(xùn)練集和測試集---- 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包。再一次調(diào)用時就沒有錯誤了。
2、自己寫函數(shù)來認(rèn)為的隨機構(gòu)造訓(xùn)練集和測試集。(這個代碼我會在最后附上。)
（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_
輸出的結(jié)果如下：
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
如何解釋各個特征對應(yīng)的系數(shù)的意義？
對于給定了Radio和Newspaper的廣告投入，如果在TV廣告上每多投入1個單位，對應(yīng)銷量將增加0.0466個單位。就是加入其它兩個媒體投入固定，在TV廣告上每增加1000美元（因為單位是1000美元），銷量將增加46.6（因為單位是1000）。但是大家注意這里的newspaper的系數(shù)居然是負(fù)數(shù)，所以我們可以考慮不使用newspaper這個特征。這是后話，后面會提到的。
（4）、預(yù)測
[python]?view plain?copy
y_pred?=?linreg.predict(X_test)
print?y_pred
[python]?view plain?copy
print?type(y_pred)
輸出結(jié)果如下：
[ 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) 評價測度
對于分類問題，評價測度是準(zhǔn)確率，但這種方法不適用于回歸問題。我們使用針對連續(xù)數(shù)值的評價測度(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預(yù)測的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
最后的結(jié)果如下：
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")?#顯示圖中的標(biāo)簽
plt.xlabel("the?number?of?sales")
plt.ylabel('value?of?sales')
plt.show()
顯示結(jié)果如下：（紅色的線是真實的值曲線，藍(lán)色的是預(yù)測值曲線）
直到這里整個的一次多元線性回歸的預(yù)測就結(jié)束了。
6、改進(jìn)特征的選擇
在之前展示的數(shù)據(jù)中，我們看到Newspaper和銷量之間的線性關(guān)系竟是負(fù)關(guān)系（不用驚訝，這是隨機特征抽樣的結(jié)果。換一批抽樣的數(shù)據(jù)就可能為正了），現(xiàn)在我們移除這個特征，看看線性回歸預(yù)測的結(jié)果的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
最后的到的系數(shù)與測度如下：
LinearRegression(copy_X=True, fit_intercept=True, normalize=False)
2.81843904823
[ 0.04588771 ?0.18721008]
RMSE by hand: 1.28208957507
然后再次使用ROC曲線來觀測曲線的整體情況。我們在將Newspaper這個特征移除之后，得到RMSE變小了，說明Newspaper特征可能不適合作為預(yù)測銷量的特征，于是，我們得到了新的模型。我們還可以通過不同的特征組合得到新的模型，看看最終的誤差是如何的。
備注：
之前我提到了這種錯誤：
注：上面的結(jié)果是由train_test_spilit()得到的，但是我不知道為什么我的版本的sklearn包中居然報錯：
ImportError ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? Traceback (most recent call last)ipython-input-182-3eee51fcba5a in module() ? ? ?1 ###構(gòu)造訓(xùn)練集和測試集---- 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包。再一次調(diào)用時就沒有錯誤了。
2、自己寫函數(shù)來認(rèn)為的隨機構(gòu)造訓(xùn)練集和測試集。(這個代碼我會在最后附上。)
這里我給出我自己寫的函數(shù)：
一、概觀
scipy中的optimize子包中提供了常用的最優(yōu)化算法函數(shù)實現(xiàn)。我們可以直接調(diào)用這些函數(shù)完成我們的優(yōu)化問題。optimize中函數(shù)最典型的特點就是能夠從函數(shù)名稱上看出是使用了什么算法。下面optimize包中函數(shù)的概覽：
1.非線性最優(yōu)化
fmin -- 簡單Nelder-Mead算法
fmin_powell -- 改進(jìn)型Powell法
fmin_bfgs -- 擬Newton法
fmin_cg -- 非線性共軛梯度法
fmin_ncg -- 線性搜索Newton共軛梯度法
leastsq -- 最小二乘
2.有約束的多元函數(shù)問題
fmin_l_bfgs_b ---使用L-BFGS-B算法
fmin_tnc ---梯度信息
fmin_cobyla ---線性逼近
fmin_slsqp ---序列最小二乘法
nnls ---解|| Ax - b ||_2 for x=0
3.全局優(yōu)化
anneal ---模擬退火算法
brute --強力法
4.標(biāo)量函數(shù)
fminbound
brent
golden
bracket
5.擬合
curve_fit-- 使用非線性最小二乘法擬合
6.標(biāo)量函數(shù)求根
brentq ---classic Brent (1973)
brenth ---A variation on the classic Brent（1980）ridder ---Ridder是提出這個算法的人名
bisect ---二分法
newton ---牛頓法
fixed_point
7.多維函數(shù)求根
fsolve ---通用
broyden1 ---Broyden’s first Jacobian approximation.
broyden2 ---Broyden’s second Jacobian approximationnewton_krylov ---Krylov approximation for inverse Jacobiananderson ---extended Anderson mixing
excitingmixing ---tuned diagonal Jacobian approximationlinearmixing ---scalar Jacobian approximationdiagbroyden ---diagonal Broyden Jacobian approximation8.實用函數(shù)
line_search ---找到滿足強Wolfe的alpha值
check_grad ---通過和前向有限差分逼近比較檢查梯度函數(shù)的正確性二、實戰(zhàn)非線性最優(yōu)化
fmin完整的調(diào)用形式是：
fmin(func, x0, args=(), xtol=0.0001, ftol=0.0001, maxiter=None, maxfun=None, full_output=0, disp=1, retall=0, callback=None)不過我們最常使用的就是前兩個參數(shù)。一個描述優(yōu)化問題的函數(shù)以及初值。后面的那些參數(shù)我們也很容易理解。如果您能用到，請自己研究。下面研究一個最簡單的問題，來感受這個函數(shù)的使用方法：f(x)=x**2-4*x+8，我們知道，這個函數(shù)的最小值是4，在x=2的時候取到。
from scipy.optimize import fmin #引入優(yōu)化包def myfunc(x):
return x**2-4*x+8 #定義函數(shù)
x0 = [1.3] #猜一個初值
xopt = fmin(myfunc, x0) #求解
print xopt #打印結(jié)果
運行之后，給出的結(jié)果是：
Optimization terminated successfully.
Current function value: 4.000000
Iterations: 16
Function evaluations: 32
[ 2.00001953]
程序準(zhǔn)確的計算得出了最小值，不過最小值點并不是嚴(yán)格的2，這應(yīng)該是由二進(jìn)制機器編碼誤差造成的。
除了fmin_ncg必須提供梯度信息外，其他幾個函數(shù)的調(diào)用大同小異，完全類似。我們不妨做一個對比：
from scipy.optimize import fmin,fmin_powell,fmin_bfgs,fmin_cgdef myfunc(x):
return x**2-4*x+8
x0 = [1.3]
xopt1 = fmin(myfunc, x0)
print xopt1
xopt2 = fmin_powell(myfunc, x0)
print xopt2
xopt3 = fmin_bfgs(myfunc, x0)
print xopt3
xopt4 = fmin_cg(myfunc,x0)
print xopt4
給出的結(jié)果是：
Optimization terminated successfully.
Current function value: 4.000000
Iterations: 16
Function evaluations: 32
[ 2.00001953]
Optimization terminated successfully.
Current function value: 4.000000
Iterations: 2
Function evaluations: 53
1.99999999997
Optimization terminated successfully.
Current function value: 4.000000
Iterations: 2
Function evaluations: 12
Gradient evaluations: 4
[ 2.00000001]
Optimization terminated successfully.
Current function value: 4.000000
Iterations: 2
Function evaluations: 15
Gradient evaluations: 5
[ 2.]
我們可以根據(jù)給出的消息直觀的判斷算法的執(zhí)行情況。每一種算法數(shù)學(xué)上的問題，請自己看書學(xué)習(xí)。個人感覺，如果不是純研究數(shù)學(xué)的工作，沒必要搞清楚那些推導(dǎo)以及定理云云。不過，必須了解每一種算法的優(yōu)劣以及能力所及。在使用的時候，不妨多種算法都使用一下，看看效果分別如何，同時，還可以互相印證算法失效的問題。
在from scipy.optimize import fmin之后，就可以使用help(fmin)來查看fmin的幫助信息了。幫助信息中沒有例子，但是給出了每一個參數(shù)的含義說明，這是調(diào)用函數(shù)時候的最有價值參考。
有源碼研究癖好的，或者當(dāng)你需要改進(jìn)這些已經(jīng)實現(xiàn)的算法的時候，可能需要查看optimize中的每種算法的源代碼。在這里：https:/ / github. com/scipy/scipy/blob/master/scipy/optimize/optimize.py聰明的你肯定發(fā)現(xiàn)了，順著這個鏈接往上一級、再往上一級，你會找到scipy的幾乎所有源碼！
重慶分公司
重慶分公司
