本文目錄一覽：
python 的LinearRegression包，怎麼導出回歸模型公式？
線性回歸是機器學習算法中最簡單的算法之一，它是監督學習的一種算法，主要思想是在給定訓練集上學習得到一個線性函數，在損失函數的約束下，求解相關係數，最終在測試集上測試模型的回歸效果。
也就是說 LinearRegression 模型會構造一個線性回歸公式
y’ = w^T x + b
，其中 w 和 x 均為向量，w 就是係數，截距是 b，得分是根據真實的 y 值和預測值 y’ 計算得到的。
python回歸模型保存
1、首先需要使用公式將回歸結果計算出來。
2、其次選擇回歸。
3、最後將其另存為，另存為到word中就可以保存了。Python由荷蘭數學和計算機科學研究學會的吉多范羅蘇姆於1990年代初設計，作為一門叫做ABC語言的替代品。Python提供了高效的高級數據結構，還能簡單有效地面向對象編程。
python多元線性回歸怎麼計算
1、什麼是多元線性回歸模型？
當y值的影響因素不唯一時,採用多元線性回歸模型。
y =y=β0+β1×1+β2×2+…+β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 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
上面顯示的結果類似一個電子表格，這個結構稱為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 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的形狀,注意第一列的叫索引，和數據庫某個表中的第一列類似。
(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+β1×1+β2×2+…+β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 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()
輸出的結果如下：
0 22.1
1 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、自己寫函數來認為的隨機構造訓練集和測試集。(這個代碼我會在最後附上。)
這裡我給出我自己寫的函數：
怎麼看python中邏輯回歸輸出的解釋
以下為python代碼，由於訓練數據比較少，這邊使用了批處理梯度下降法，沒有使用增量梯度下降法。
##author:lijiayan##data:2016/10/27
##name:logReg.pyfrom numpy import *import matplotlib.pyplot as pltdef loadData(filename):
data = loadtxt(filename)
m,n = data.shape print ‘the number of examples:’,m print ‘the number of features:’,n-1 x = data[:,0:n-1]
y = data[:,n-1:n] return x,y#the sigmoid functiondef sigmoid(z): return 1.0 / (1 + exp(-z))#the cost functiondef costfunction(y,h):
y = array(y)
h = array(h)
J = sum(y*log(h))+sum((1-y)*log(1-h)) return J# the batch gradient descent algrithmdef gradescent(x,y):
m,n = shape(x) #m: number of training example; n: number of features x = c_[ones(m),x] #add x0 x = mat(x) # to matrix y = mat(y)
a = 0.0000025 # learning rate maxcycle = 4000 theta = zeros((n+1,1)) #initial theta J = [] for i in range(maxcycle):
h = sigmoid(x*theta)
theta = theta + a * (x.T)*(y-h)
cost = costfunction(y,h)
J.append(cost)
plt.plot(J)
plt.show() return theta,cost#the stochastic gradient descent (m should be large,if you want the result is good)def stocGraddescent(x,y):
m,n = shape(x) #m: number of training example; n: number of features x = c_[ones(m),x] #add x0 x = mat(x) # to matrix y = mat(y)
a = 0.01 # learning rate theta = ones((n+1,1)) #initial theta J = [] for i in range(m):
h = sigmoid(x[i]*theta)
theta = theta + a * x[i].transpose()*(y[i]-h)
cost = costfunction(y,h)
J.append(cost)
plt.plot(J)
plt.show() return theta,cost#plot the decision boundarydef plotbestfit(x,y,theta):
plt.plot(x[:,0:1][where(y==1)],x[:,1:2][where(y==1)],’ro’)
plt.plot(x[:,0:1][where(y!=1)],x[:,1:2][where(y!=1)],’bx’)
x1= arange(-4,4,0.1)
x2 =(-float(theta[0])-float(theta[1])*x1) /float(theta[2])
plt.plot(x1,x2)
plt.xlabel(‘x1’)
plt.ylabel((‘x2’))
plt.show()def classifyVector(inX,theta):
prob = sigmoid((inX*theta).sum(1)) return where(prob = 0.5, 1, 0)def accuracy(x, y, theta):
m = shape(y)[0]
x = c_[ones(m),x]
y_p = classifyVector(x,theta)
accuracy = sum(y_p==y)/float(m) return accuracy
調用上面代碼：
from logReg import *
x,y = loadData(“horseColicTraining.txt”)
theta,cost = gradescent(x,y)print ‘J:’,cost
ac_train = accuracy(x, y, theta)print ‘accuracy of the training examples:’, ac_train
x_test,y_test = loadData(‘horseColicTest.txt’)
ac_test = accuracy(x_test, y_test, theta)print ‘accuracy of the test examples:’, ac_test
學習速率=0.0000025，迭代次數=4000時的結果：
似然函數走勢（J = sum(y*log(h))+sum((1-y)*log(1-h))），似然函數是求最大值，一般是要穩定了才算最好。
下圖為計算結果，可以看到訓練集的準確率為73%，測試集的準確率為78%。
這個時候，我去看了一下數據集，發現沒個特徵的數量級不一致，於是我想到要進行歸一化處理：
歸一化處理句修改列loadData(filename)函數：
def loadData(filename):
data = loadtxt(filename)
m,n = data.shape print ‘the number of examples:’,m print ‘the number of features:’,n-1 x = data[:,0:n-1]
max = x.max(0)
min = x.min(0)
x = (x – min)/((max-min)*1.0) #scaling y = data[:,n-1:n] return x,y
在沒有歸一化的時候，我的學習速率取了0.0000025（加大就會震蕩，因為有些特徵的值很大，學習速率取的稍大，波動就很大），由於學習速率小，迭代了4000次也沒有完全穩定。現在當把特徵歸一化後（所有特徵的值都在0~1之間），這樣學習速率可以加大，迭代次數就可以大大減少，以下是學習速率=0.005，迭代次數=500的結果：
此時的訓練集的準確率為72%，測試集的準確率為73%
從上面這個例子，我們可以看到對特徵進行歸一化操作的重要性。
原創文章，作者：小藍，如若轉載，請註明出處：https://www.506064.com/zh-hk/n/291670.html
微信掃一掃 
支付寶掃一掃 