概念
在統(tǒng)計(jì)學(xué)中，回歸分析（regression analysis)指的是確定兩種或兩種以上變量間相互依賴的定量關(guān)系的一種統(tǒng)計(jì)分析方法?；貧w分析按照涉及的變量的多少，分為一元回歸和多元回歸分析；按照因變量的多少，可分為簡單回歸分析和多重回歸分析；按照自變量和因變量之間的關(guān)系類型，可分為線性回歸分析和非線性回歸分析。
在大數(shù)據(jù)分析中，回歸分析是一種預(yù)測性的建模技術(shù)，它研究的是因變量（目標(biāo)）和自變量（預(yù)測器）之間的關(guān)系。這種技術(shù)通常用于預(yù)測分析，時(shí)間序列模型以及發(fā)現(xiàn)變量之間的因果關(guān)系。
回歸分析的步驟
- 確定回歸方程的解釋變量和被解釋變量。
- 確定回歸模型，建立回歸方程。
- 對回歸方程進(jìn)行各種驗(yàn)證。
- 利用回歸方程進(jìn)行預(yù)測。
一元線性回歸
一元線性回歸模型
如果只有一個(gè)自變量X，而且因變量Y和自變量X之間的數(shù)量變化關(guān)系呈近似線性關(guān)系，就可以建立一元線性回歸方程，由自變量X的值來預(yù)測因變量Y的值，這就是一元線性回歸預(yù)測。
如果因變量Y和自變量X之間呈線性相關(guān)，那就是說，對于自變量X的某一值 xi ，因變量Y對應(yīng)的取值 yi 不是唯一確定的，而是有很多的可能取值，它們分布在一條直線的上下，這是因?yàn)閅還受除自變量以外的其他因素的影響。這些因素的影響大小和方向都是不確定的，通常用一個(gè)隨機(jī)變量(記為ε)來表示，也稱為隨機(jī)擾動(dòng)項(xiàng)。于是，Y和X之間的依存關(guān)系可表示為：
α —— 截距
β —— 斜率
ε —— 誤差項(xiàng)，反映除x和y的線性關(guān)系之外的隨機(jī)因素對y的影響
一元線性回歸方程
描述因變量y的期望值如何依賴于自變量x的方程稱為回歸方程：
- 如果回歸方程的參數(shù)已知，對于一個(gè)給定的x值，利用回歸方程就能計(jì)算出y的期望值。
- 用樣本統(tǒng)計(jì)量代替回歸方程中的位置參數(shù)，就得到估計(jì)回歸方程，簡稱回歸直線。
參數(shù)的最小二乘法估計(jì)
對于回歸直線，關(guān)鍵在于求解參數(shù)，常用高斯提出的最小二乘法，它是使因變量的觀察值y與估計(jì)值之間的離差平方和達(dá)到最小來求解。
展開可得：
我們求Q的極小值點(diǎn)，對其求偏導(dǎo)(β0和β1)：
求解得：
利用回歸直線進(jìn)行估計(jì)和預(yù)測
- 點(diǎn)估計(jì)：利用估計(jì)的回歸方程，對于x的某一個(gè)特定得值，求出y得一個(gè)估計(jì)值。
- 區(qū)間估計(jì)：利用估計(jì)得回歸方程，對于x的一個(gè)特定值，求出y的一個(gè)估計(jì)值的區(qū)間。
估計(jì)標(biāo)準(zhǔn)誤差的計(jì)算
為了度量回歸方程的可靠性，通常計(jì)算估計(jì)的標(biāo)準(zhǔn)誤差。它度量觀察值回繞著回歸直線的變化程度或分散程度。
估計(jì)平均誤差：
- 自由度為：n-2（由兩個(gè)不能改變參數(shù)）。
- 估計(jì)標(biāo)準(zhǔn)誤差越大，則數(shù)據(jù)點(diǎn)圍繞回歸直線的分散程度越大，回歸方程大的代表性越小。
- 估計(jì)標(biāo)準(zhǔn)誤差越小，則數(shù)據(jù)點(diǎn)圍繞回歸直線的分散程度越小，回歸方程大的代表性越大，即可靠性高。
置信區(qū)間估計(jì)
在1 — α置信水平下預(yù)測區(qū)間
影響區(qū)間寬度的因素
- 置信水平（1 - α), 區(qū)間寬度隨置信水平的增大而增大。
- 數(shù)據(jù)的離散程度Se，區(qū)間寬度隨離散程度的增大而增大。
- 樣本容量，區(qū)間寬度隨樣本容量的增大而減小。
- X0與X均值之間的差異，隨差異程度的增大而增大。
回歸直線的擬合優(yōu)度
回歸直線與各觀測點(diǎn)的接近程度稱為回歸直線對數(shù)據(jù)的擬合優(yōu)度。
-
總平方和（SST）：
-
回歸平方和（SSR）：
-
殘差平方和（SSE）：
總平方和可以分解為回歸平方和、殘差平方和兩部分：SST=SSR+SSE
- SST反映因變量的n個(gè)觀察值與其均值的總離差。
- SSR反映了y的總變差中，由于x與y的線性關(guān)系引起的y的變化部分。
- SSE反映了除了x對y的線性影響之外的其他因素對y變差的作用，是不能由回歸直線來解釋的y的變差部分。
判定系數(shù)
回歸平方和占總平方的比例，用R^2表示，其值在0到1之間。
- R^2 == 0：說明y的變化與x無關(guān)，x完全無助于解釋y的變差
- R^2 == 1：說明殘差平方和為0，擬合是完全的，y的變化只與x有關(guān)（理想情況下）
顯著性檢驗(yàn)
顯著性檢驗(yàn)的主要目的是根據(jù)所建立的估計(jì)方程用自變量x來估計(jì)或預(yù)測因變量y的取值。當(dāng)建立了估計(jì)方程后，還不能馬上進(jìn)行估計(jì)或預(yù)測，因?yàn)樵摴烙?jì)方程是根據(jù)樣本數(shù)據(jù)得到的，它是否真實(shí)的反映了變量x和y之間的關(guān)系，需要通過檢驗(yàn)后才能證實(shí)。
根據(jù)樣本數(shù)據(jù)擬合回歸方程時(shí)，實(shí)際上就已經(jīng)假定變量x與y之間存在著線性關(guān)系，并假定誤差項(xiàng)是一個(gè)服從正態(tài)分布的隨機(jī)變量，且具有相同的方差。但這些假設(shè)是否成立需要檢驗(yàn)，顯著性檢驗(yàn)包括兩方面：
- 線性關(guān)系檢驗(yàn)
- 回歸系數(shù)檢驗(yàn)
線性關(guān)系檢驗(yàn)
線性關(guān)系檢驗(yàn)是檢驗(yàn)自變量x和因變量y之間的線性關(guān)系是否顯著或者說，它們之間能否用一個(gè)線性模型來表示。
講均方回歸（MSR）同均方殘差（MSE）加以比較，應(yīng)用F檢驗(yàn)來分析二者之間的差別是否顯著。
- 均方回歸：回歸平方和SSR除以相應(yīng)的自由度（自變量個(gè)數(shù)K）
- 均方殘差：殘差平方和SSE除以相應(yīng)的自由度（n-k-1）
H0：β1=0所有回歸系數(shù)與0無顯著差異，y與全體x的線性關(guān)系不顯著。
計(jì)算檢驗(yàn)統(tǒng)計(jì)量F：
回歸系數(shù)檢驗(yàn)
回歸系數(shù)顯著性檢驗(yàn)的目的是通過檢驗(yàn)回歸系數(shù)β的值與0是否有顯著性差異，來判斷Y與X之間是否有顯著性線性關(guān)系。若β=0，則總體回歸方程中不含X項(xiàng)（即Y不隨X變動(dòng)而變動(dòng)），因此變量Y與X之間不存在線性關(guān)系；若β≠0，說明變量Y與X之間存在顯著的線性關(guān)系。
計(jì)算檢驗(yàn)的統(tǒng)計(jì)量：
兩個(gè)檢驗(yàn)的區(qū)別
線性關(guān)系的檢驗(yàn)是檢驗(yàn)自變量與因變量是否可以用線性來表達(dá)，而回歸系數(shù)的檢驗(yàn)是對樣本數(shù)據(jù)計(jì)算的回歸系數(shù)檢驗(yàn)總體中回歸系數(shù)是否為0。
- 在一元線性回歸中，自變量只有一個(gè)，線性關(guān)系檢驗(yàn)與回歸系數(shù)檢驗(yàn)是等價(jià)的。
- 多元回歸分析中，這兩種檢驗(yàn)的意義是不同的。線性關(guān)系檢驗(yàn)只能用來檢驗(yàn)總體回歸關(guān)系的顯著性，而回歸系數(shù)檢驗(yàn)可以對各個(gè)回歸系數(shù)分別進(jìn)行檢驗(yàn)。
多元線性回歸
經(jīng)常遇到某一現(xiàn)象的發(fā)展變化取決于幾個(gè)影響因素的情況，也就是一個(gè)因變量和幾個(gè)自變量有依存關(guān)系的情況，這時(shí)需用多元線性回歸分析。
- 多元線性回歸分析預(yù)測法，是指通過對兩及其以上的自變量與一個(gè)因變量的相關(guān)性分析，建立預(yù)測模型進(jìn)行預(yù)測和控制的方法。
多元線性回歸預(yù)測模型一般式為：
調(diào)整的多重判定系數(shù)
用樣本容量n和自變量的個(gè)數(shù)k去修正R^2得到：
避免增加自變量而高估R^2
曲線回歸分析
直線關(guān)系是兩變量間最簡單的一種關(guān)系，曲線回歸分析的基本任務(wù)是通過兩個(gè)相關(guān)變量x與y的實(shí)際觀測數(shù)據(jù)建立曲線回歸方程，以揭示x與y間的曲線聯(lián)系的形式。
曲線分析最困難和首要的工作是確定自變量與因變量的曲線關(guān)系的類型，曲線回歸分析的基本過程：
- 將x或y進(jìn)行變量轉(zhuǎn)換。
- 對新變量進(jìn)行直線回歸分析、建立直線回歸方程并進(jìn)行顯著性檢驗(yàn)和區(qū)間估計(jì)。
- 將新變量還原為原變量，由新變量的直線回歸方程和置信區(qū)間得出原變量的曲線回歸方程和置信區(qū)間。
由于曲線回歸模型種類繁多，所以沒有通用的回歸方程可直接使用。但是對于某些特殊的回歸模型，可以通過變量代換、取對數(shù)等方法對其線性化，然后使用標(biāo)準(zhǔn)方程求解參數(shù)，再將參數(shù)帶回原方程就是所求。
例如：
多重共線性
回歸模型中兩個(gè)或兩個(gè)以上的自變量彼此相關(guān)的現(xiàn)象
多重共線性帶來的問題：
- 回歸系數(shù)估計(jì)值的不穩(wěn)定增強(qiáng)
- 回歸系數(shù)假設(shè)檢驗(yàn)的結(jié)果不顯著等
多重共線性檢驗(yàn)的主要方法
- 容忍度
- 方差膨脹因子（VIF）
容忍度
Ri是解釋變量xi與方程中其他解釋變量間的復(fù)相關(guān)系數(shù)；
容忍度在0~1之間，越接近于0，表示多重共線性越強(qiáng)。
方差膨脹因子
方差膨脹因子是容忍度的倒數(shù)：
VIFi越大，特別是>=10，說明解釋變量xi與方程中其他解釋變量之間有嚴(yán)重的多重共線性；
VIFi越接近于1，表明解釋變量xi和其他解釋變量之間的多重共線性越弱。
Python工具包介紹
Statsmodels
用來做統(tǒng)計(jì)分析
Statsmodels 是一個(gè) Python 模塊，它提供用于估計(jì)許多不同統(tǒng)計(jì)模型以及進(jìn)行統(tǒng)計(jì)測試和統(tǒng)計(jì)數(shù)據(jù)探索的類和函數(shù)。 每個(gè)估算器都有一個(gè)廣泛的結(jié)果統(tǒng)計(jì)列表。 對照現(xiàn)有的統(tǒng)計(jì)數(shù)據(jù)包對結(jié)果進(jìn)行測試，以確保它們是正確的。 該軟件包是在開源修改 BSD（3 條款）許可下發(fā)布的。
一元線性回歸
# 導(dǎo)入功能包，產(chǎn)生數(shù)據(jù)
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
nsample = 20 # 20個(gè)樣本點(diǎn)
x = np.linspace(0,10,nsample) # 從0到20之間選擇20個(gè)數(shù)(等距)
print(x)
[ 0. 0.52631579 1.05263158 1.57894737 2.10526316 2.63157895
3.15789474 3.68421053 4.21052632 4.73684211 5.26315789 5.78947368
6.31578947 6.84210526 7.36842105 7.89473684 8.42105263 8.94736842
9.47368421 10. ]
# 加一列1，為了與常數(shù)項(xiàng)進(jìn)行組合
X = sm.add_constant(x)
print(X)
[[ 1. 0. ]
[ 1. 0.52631579]
[ 1. 1.05263158]
[ 1. 1.57894737]
[ 1. 2.10526316]
[ 1. 2.63157895]
[ 1. 3.15789474]
[ 1. 3.68421053]
[ 1. 4.21052632]
[ 1. 4.73684211]
[ 1. 5.26315789]
[ 1. 5.78947368]
[ 1. 6.31578947]
[ 1. 6.84210526]
[ 1. 7.36842105]
[ 1. 7.89473684]
[ 1. 8.42105263]
[ 1. 8.94736842]
[ 1. 9.47368421]
[ 1. 10. ]]
# 產(chǎn)生β0與β1（真實(shí)情況）
beta = np.array([2,5])
print(beta)
[2 5]
# 誤差項(xiàng)(設(shè)置抖動(dòng)，為了讓數(shù)據(jù)更真實(shí))
e = np.random.normal(size=nsample)
print(e)
[ 0.82944189 1.67229293 0.39780164 1.24295698 1.475433 0.8642918
-0.34345566 1.03677286 0.26148177 -1.03382748 -0.09721028 0.66122185
-0.54803173 0.40733307 0.32082923 1.38139669 -0.81902716 0.45718195
1.31683565 1.03592956]
# 實(shí)際值y
y = np.dot(X,beta) + e
print(y)
[ 2.82944189 6.30387188 7.66095953 11.13769382 14.00174879 16.02218654
17.44601802 21.45782549 23.31411335 24.65038305 28.2185792 31.60859027
33.03091564 36.61785938 39.1629345 42.8550809 43.286236 47.19402405
50.6852567 53.03592956]
# 最小二乘法
model = sm.OLS(y,X)
# 擬合數(shù)據(jù)
res = model.fit()
# 計(jì)算的回歸系數(shù)β0和β1
print(res.params)
[2.77235974 4.95072454]
# 產(chǎn)生全部結(jié)果
res.summary()
# 擬合的估計(jì)值
y_ = res.fittedvalues
print(y_)
[ 2.77235974 5.37800424 7.98364873 10.58929322 13.19493772 15.80058221
18.4062267 21.0118712 23.61751569 26.22316018 28.82880467 31.43444917
34.04009366 36.64573815 39.25138265 41.85702714 44.46267163 47.06831613
49.67396062 52.27960511]
fig, ax = plt.subplots(figsize=(8,6))
ax.plot(x,y,'o',label='data') # 原始數(shù)據(jù)
ax.plot(x,y_,'r--',label='test') # 擬合數(shù)據(jù)
ax.legend(loc='best')
plt.show()
高階回歸
# Y=5+2X+3X^2
nsample=50
x=np.linspace(0,10,nsample)
X=np.column_stack((x,x**2)) # 添加原始數(shù)據(jù)
X=sm.add_constant(X)
# 則X的第一列為1，第二列為x，第三列為x^2
# 構(gòu)造beta
beta=np.array([5,2,3]) # 真實(shí)值
e=np.random.normal(size=nsample) # 誤差
y=np.dot(X,beta)+e # 產(chǎn)生y真實(shí)值
model=sm.OLS(y,X) # 最小二乘法
results=model.fit()
results.params # 計(jì)算最終系數(shù)值
[4.93081682, 2.22120911, 2.97636437]
# 產(chǎn)生全部結(jié)果
results.summary()
# 畫圖展示結(jié)果
y_fitted=results.fittedvalues
fit,ax=plt.subplots(figsize=(8,6))
ax.plot(x,y,'o',label='data')
ax.plot(x,y_fitted,'r--.',label='OLS')
ax.legend(loc='best')
plt.show()
分類變量
假設(shè)分類變量有3個(gè)取值（a，b，c），比如考試成績有3個(gè)等級(jí)。那么a就是（1，0，0），b（0，1，0），c（0，0，1），這個(gè)時(shí)候就需要3個(gè)系數(shù)β0，β1，β2，也就是β0x0+β1x1+β2x2。
# 構(gòu)造0列表
nsample=50
groups=np.zeros(nsample,int)
groups
array([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
0, 0, 0, 0, 0, 0])
# 轉(zhuǎn)化為矩陣
groups[20:40]=1
groups[40:]=2
dummy=sm.categorical(groups,drop=True) #轉(zhuǎn)化類型
dummy
# Y=5+2X+3Z1+6*Z2+9*Z3.
x = np.linspace(0,20,nsample)
X=np.column_stack((x,dummy))
X=sm.add_constant(X)
beta=[5,2,3,6,9]
e=np.random.normal(size=nsample)
y=np.dot(X,beta)+e
result=sm.OLS(y,X).fit()
result.summary()
# 畫圖
fig, ax=plt.subplots(figsize=(8,6))
ax.plot(x,y,'o',label="data")
ax.plot(x,result.fittedvalues,'r--',label="OLS")
ax.legend(loc='best')
plt.show()
Scikit-learn
用來做機(jī)器學(xué)習(xí)
Scikit-learn（以前稱為scikits.learn，也稱為sklearn）是針對Python 編程語言的免費(fèi)軟件機(jī)器學(xué)習(xí)庫 。它具有各種分類，回歸和聚類算法，包括支持向量機(jī)，隨機(jī)森林，梯度提升，k均值和DBSCAN，并且旨在與Python數(shù)值科學(xué)庫NumPy和SciPy聯(lián)合使用。
實(shí)戰(zhàn)：汽車價(jià)格預(yù)測
數(shù)據(jù)字典
數(shù)據(jù)來源： https://archive.ics.uci.edu/ml/datasets/Automobile
數(shù)據(jù)讀取與分析
# 加載功能包
import numpy as np
import pandas as pd
import datetime
# 數(shù)據(jù)可視化和缺失值
import matplotlib.pyplot as plt
import seaborn as sns # 高級(jí)可視化
import missingno as msno # 對缺失值進(jìn)行可視化展示
%matplotlib inline
# 統(tǒng)計(jì)數(shù)據(jù)
from statsmodels.distributions.empirical_distribution import ECDF
from sklearn.metrics import mean_squared_error, r2_score
# 機(jī)器學(xué)習(xí)
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import Lasso, LassoCV
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.ensemble import RandomForestRegressor
seed=123 # 隨機(jī)種子(控制隨機(jī)因素)
# 讀取數(shù)據(jù)集
data = pd.read_csv("CarPrice_Assignment.csv", na_values = '?')
data.columns
Index([‘car_ID’, ‘symboling’, ‘CarName’, ‘fueltype’, ‘a(chǎn)spiration’,
‘doornumber’, ‘carbody’, ‘drivewheel’, ‘enginelocation’, ‘wheelbase’,
‘carlength’, ‘carwidth’, ‘carheight’, ‘curbweight’, ‘enginetype’,
‘cylindernumber’, ‘enginesize’, ‘fuelsystem’, ‘boreratio’, ‘stroke’,
‘compressionratio’, ‘horsepower’, ‘peakrpm’, ‘citympg’, ‘highwaympg’,
‘price’],
dtype=‘object’)
# 查看數(shù)據(jù)類型
data.dtypes
car_ID int64
symboling int64
CarName object
fueltype object
aspiration object
doornumber object
carbody object
drivewheel object
enginelocation object
wheelbase float64
carlength float64
carwidth float64
carheight float64
curbweight int64
enginetype object
cylindernumber object
enginesize int64
fuelsystem object
boreratio float64
stroke float64
compressionratio float64
horsepower int64
peakrpm int64
citympg int64
highwaympg int64
price float64
dtype: object
# 查看數(shù)據(jù)本身
print("IN total:",data.shape)
data.head(5)
205條數(shù)據(jù)，26個(gè)特征
# 顯示CSV中的數(shù)值信息
data.describe()
缺失值處理（NaN）
當(dāng)然本案例是沒有缺失值的，若有需要用missingno功能包進(jìn)行可視化展示，要么去掉缺失值，要么用平均值或者眾數(shù)、中位數(shù)替代。
# 查看缺失值
sns.set(style = "ticks") # 指定風(fēng)格
msno.matrix(data)
若有缺失值：
對于缺失值較為嚴(yán)重的數(shù)據(jù)
# 展示缺失值數(shù)據(jù)
data[pd.isnull(data['normalized-losses'])].head()
# 進(jìn)行數(shù)據(jù)分布分析，然后選擇填充數(shù)據(jù)
sns.set(style = "ticks")
plt.figure(figsize = (12, 5))
c = '#366DE8'
# ECDF
plt.subplot(121)
cdf = ECDF(data['normalized-losses']) #查看連續(xù)分布，累計(jì)結(jié)果
plt.plot(cdf.x, cdf.y, label = "statmodels", color = c);
plt.xlabel('normalized losses'); plt.ylabel('ECDF');
# overall distribution
plt.subplot(122)
plt.hist(data['normalized-losses'].dropna(),
bins = int(np.sqrt(len(data['normalized-losses']))),
color = c);
# 按照不同的風(fēng)險(xiǎn)等級(jí)分成不同的組，然后觀察每組缺失值情況
data.groupby([symboling])["normalized-losses"].describe()
# 使用最相關(guān)數(shù)值進(jìn)行填充
data = data.dropna(subset = ['price', 'bore', 'stroke', 'peak-rpm', 'horsepower', 'num-of-doors']) #對于缺失值少的幾列，直接刪掉缺失值
data['normalized-losses'] = data.groupby('symboling')['normalized-losses'].transform(lambda x: x.fillna(x.mean())) #填充缺失值
print('In total:', data.shape)
data.head()
特征相關(guān)性
cormatrix = data.corr()
cormatrix
cormatrix *= np.tri(*cormatrix.values.shape, k=-1).T #返回函數(shù)的上三角矩陣，把對角線上的置0，讓他們不是最高的。
cormatrix
#某一指標(biāo)與其他指標(biāo)的關(guān)系
cormatrix = cormatrix.stack()
cormatrix
# 按照從大到小排序
cormatrix = cormatrix.reindex(cormatrix.abs().sort_values(ascending=False).index).reset_index()
cormatrix
# 選擇相關(guān)系數(shù)大的一個(gè)影響因素即可
# 并選擇本身性質(zhì)關(guān)系的一個(gè)即可（如長寬高—>體積）
data['volume'] = data.carlength * data.carwidth * data.carheight
data.drop(['carwidth', 'carlength', 'carheight',
'curbweight', 'citympg'],
axis = 1, # 1 for columns
inplace = True)
# 新變量
data.columns
Index([‘car_ID’, ‘symboling’, ‘CarName’, ‘fueltype’, ‘a(chǎn)spiration’,
‘doornumber’, ‘carbody’, ‘drivewheel’, ‘enginelocation’, ‘wheelbase’,
‘enginetype’, ‘cylindernumber’, ‘enginesize’, ‘fuelsystem’, ‘boreratio’,
‘stroke’, ‘compressionratio’, ‘horsepower’, ‘peakrpm’, ‘highwaympg’,
‘price’, ‘volume’],
dtype=‘object’)
# 畫出相關(guān)性的熱度圖
# 計(jì)算相關(guān)矩陣
corr_all = data.corr()
# 為上三角形生成蒙版
mask = np.zeros_like(corr_all, dtype = np.bool)
mask[np.triu_indices_from(mask)] = True
# 設(shè)置 matplotlib 圖
f, ax = plt.subplots(figsize = (11, 9))
# 使用蒙版和正確的縱橫比繪制熱圖
sns.heatmap(corr_all, mask = mask,
square = True, linewidths = .5, ax = ax, cmap = "BuPu")
plt.show()
# 統(tǒng)計(jì)每一個(gè)指標(biāo)的情況
sns.pairplot(data, hue = 'fueltype', palette = 'plasma')
# 返回的是當(dāng)前對于不同變量來說的各種圖
# 按照hue的指標(biāo)，把數(shù)據(jù)分為幾種
# 某幾個(gè)指標(biāo)之間的關(guān)系
sns.lmplot('price', 'horsepower', data,
hue = 'fueltype', col = 'fueltype', row = 'doornumber',
palette = 'plasma',
fit_reg = True);
#參數(shù)：1指定橫軸；2指定縱軸；傳進(jìn)數(shù)據(jù)data；hue按照哪個(gè)變化（參數(shù)的個(gè)數(shù)）；col按列區(qū)分；row按行區(qū)分
預(yù)處理
對連續(xù)值進(jìn)行標(biāo)量化
# 目標(biāo)和特點(diǎn)
target = data.price
regressors = [x for x in data.columns if x not in ['price']]
features = data.loc[:, regressors]
num = ['symboling', 'volume', 'horsepower', 'wheelbase',
'boreratio', 'stroke','compressionratio', 'peakrpm']
# 縮放數(shù)據(jù)
standard_scaler = StandardScaler() # 進(jìn)行特征標(biāo)準(zhǔn)化
features[num] = standard_scaler.fit_transform(features[num])
# 展示結(jié)果
features.head()
# 對分類屬性進(jìn)行one-bot編碼
# 分類變量
classes = ['CarName', 'fueltype', 'aspiration', 'doornumber',
'carbody', 'drivewheel', 'enginelocation',
'enginetype', 'cylindernumber', 'fuelsystem']
# 僅使用 continios vars 創(chuàng)建新數(shù)據(jù)集
dummies = pd.get_dummies(features[classes])
features = features.join(dummies).drop(classes, axis = 1)
# 新數(shù)據(jù)集
print('In total:', features.shape)
features.head()
205條數(shù)據(jù)，196個(gè)特征
# 劃分?jǐn)?shù)據(jù)集，按照30%劃分（訓(xùn)練集和測試集）
X_train, X_test, y_train, y_test = train_test_split(features, target, test_size = 0.3,random_state = seed)
print("Train", X_train.shape, "and test", X_test.shape)
Train (143, 196) and test (62, 196)
Lasso回歸
多加了一個(gè)絕對值想來懲罰過大的系數(shù)：alphas=0那就是最小二乘了。
# 對數(shù)刻度：以 2 為底的對數(shù)
# 高值將更多變量歸零
alphas = 2. ** np.arange(2, 20) #指定alphas的范圍
scores = np.empty_like(alphas)
# 測試集
for i, a in enumerate(alphas):
lasso = Lasso(random_state = seed) # 指定模型
lasso.set_params(alpha = a)
lasso.fit(X_train, y_train)
scores[i] = lasso.score(X_test, y_test)
# 交叉驗(yàn)證cross validation
lassocv = LassoCV(cv = 10, random_state = seed)
lassocv.fit(features, target)
lassocv_score = lassocv.score(features, target)
lassocv_alpha = lassocv.alpha_
plt.figure(figsize = (10, 4))
plt.plot(alphas, scores, '-ko')
plt.axhline(lassocv_score)
plt.xlabel(r'$\alpha$')
plt.ylabel('CV Score')
plt.xscale('log', basex = 2)
sns.despine(offset = 15)
print('CV results:', lassocv_score, lassocv_alpha)
# 特征重要性分析
# lassocv 系數(shù)
coefs = pd.Series(lassocv.coef_, index = features.columns)
# 打印出選擇/消除特征的數(shù)量
print("Lasso picked " + str(sum(coefs != 0)) + " features and eliminated the other " + \
str(sum(coefs == 0)) + " features.")
# 展示前5個(gè)和后5個(gè)
coefs = pd.concat([coefs.sort_values().head(5), coefs.sort_values().tail(5)])
plt.figure(figsize = (10, 4))
coefs.plot(kind = "barh")
plt.title("Coefficients in the Lasso Model")
plt.show()
# 選擇上面計(jì)算出來的最好的Alphas，進(jìn)行預(yù)測
model_l1 = LassoCV(alphas = alphas, cv = 10, random_state = seed).fit(X_train, y_train)
y_pred_l1 = model_l1.predict(X_test)
model_l1.score(X_test, y_test)
0.8156778549898953
# residual plot 殘差圖——畫圖表示實(shí)際值和預(yù)測值之間的差異
plt.rcParams['figure.figsize'] = (6.0, 6.0)
preds = pd.DataFrame({"preds": model_l1.predict(X_train), "true": y_train})
preds["residuals"] = preds["true"] - preds["preds"]
preds.plot(x = "preds", y = "residuals", kind = "scatter", color = "c")
文章來源：http://www.zghlxwxcb.cn/news/detail-464410.html
#計(jì)算指標(biāo)：MSE和R2
def MSE(y_true,y_pred):
mse = mean_squared_error(y_true, y_pred)
print('MSE: %2.3f' % mse)
return mse
def R2(y_true,y_pred):
r2 = r2_score(y_true, y_pred)
print('R2: %2.3f' % r2)
return r2
MSE(y_test, y_pred_l1); R2(y_test, y_pred_l1);
MSE: 6033979.103
R2: 0.816文章來源地址http://www.zghlxwxcb.cn/news/detail-464410.html
# 結(jié)果預(yù)測
d = {'true' : list(y_test),
'predicted' : pd.Series(y_pred_l1)
}
pd.DataFrame(d).head()
到了這里，關(guān)于統(tǒng)計(jì)分析——回歸分析的文章就介紹完了。如果您還想了解更多內(nèi)容，請?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！
