《利用Python進行數據分析》13章（中二)建模庫介紹

來自專欄 Python程序員7 人贊了文章

前文傳送門：

【翻譯】《利用Python進行數據分析·第2版》第1章 準備工作

【翻譯】《利用Python進行數據分析·第2版》第2章（上）Python語法基礎，IPython和Jupyter

【翻譯】《利用Python進行數據分析·第2版》第2章（中）Python語法基礎，IPython和Jupyter

【翻譯】《利用Python進行數據分析·第2版》第2章（下）Python語法基礎，IPython和Jupyter

【翻譯】《利用Python進行數據分析·第2版》第3章（上）Python的數據結構、函數和文件

【翻譯】《利用Python進行數據分析·第2版》第3章（中）Python的數據結構、函數和文件

【翻譯】《利用Python進行數據分析·第2版》第3章（下）Python的數據結構、函數和文件

【翻譯】《利用Python進行數據分析·第2版》第4章（上）NumPy基礎：數組和矢量計算

【翻譯】《利用Python進行數據分析·第2版》第4章（中）NumPy基礎：數組和矢量計算

【翻譯】《利用Python進行數據分析·第2版》第4章（下）NumPy基礎：數組和矢量計算

【翻譯】《利用Python進行數據分析·第2版》第5章（上）pandas入門

【翻譯】《利用Python進行數據分析·第2版》第5章（中）pandas入門

【翻譯】《利用Python進行數據分析·第2版》第5章（下）pandas入門

【翻譯】《利用Python進行數據分析·第2版》第6章（上） 數據載入、存儲與文件格式

【翻譯】《利用Python進行數據分析·第2版》第6章（中） 數據載入、存儲與文件格式

【翻譯】《利用Python進行數據分析·第2版》第6章（下） 數據載入、存儲與文件格式

【翻譯】《利用Python進行數據分析·第2版》第7章（上）數據清洗和準備

【翻譯】《利用Python進行數據分析·第2版》第7章（中） 數據清洗和準備

【翻譯】《利用Python進行數據分析·第2版》第7章（下） 數據清洗和準備

【翻譯】《利用Python進行數據分析·第2版》第8章（上） 數據規整：聚合、合併和重塑

【翻譯】《利用Python進行數據分析·第2版》第8章（中） 數據規整：聚合、合併和重塑

【翻譯】《利用Python進行數據分析·第2版》第8章（下） 數據規整：聚合、合併和重塑

【翻譯】《利用Python進行數據分析·第2版》第9章（上） 繪圖和可視化

【翻譯】《利用Python進行數據分析·第2版》第9章（中） 繪圖和可視化

【翻譯】《利用Python進行數據分析·第2版》第9章（下） 繪圖和可視化

【翻譯】《利用Python進行數據分析·第2版》第10章（上） 數據聚合與分組運算

【翻譯】《利用Python進行數據分析·第2版》第10章（中） 數據聚合與分組運算

【翻譯】《利用Python進行數據分析·第2版》第10章（下） 數據聚合與分組運算

【翻譯】《利用Python進行數據分析·第2版》第11章（上） 時間序列

【翻譯】《利用Python進行數據分析·第2版》第11章（中） 時間序列

【翻譯】《利用Python進行數據分析·第2版》第11章（中二） 時間序列

【翻譯】《利用Python進行數據分析·第2版》第11章（下） 時間序列

【翻譯】《利用Python進行數據分析·第2版》第12章（上） pandas高級應用

【翻譯】《利用Python進行數據分析·第2版》第12章（中） pandas高級應用

【翻譯】《利用Python進行數據分析·第2版》第12章（下） pandas高級應用

【翻譯】《利用Python進行數據分析·第2版》第13章（上） Python建模庫介紹

【翻譯】《利用Python進行數據分析·第2版》第13章（中） Python建模庫介紹

13.3 statsmodels介紹

statsmodels是Python進行擬合多種統計模型、進行統計試驗和數據探索可視化的庫。Statsmodels包含許多經典的統計方法，但沒有貝葉斯方法和機器學習模型。

statsmodels包含的模型有：

• 線性模型，廣義線性模型和健壯線性模型
• 線性混合效應模型
• 方差（ANOVA）方法分析
• 時間序列過程和狀態空間模型
• 廣義矩估計

下面，我會使用一些基本的statsmodels工具，探索Patsy公式和pandasDataFrame對象如何使用模型介面。

估計線性模型

statsmodels有多種線性回歸模型，包括從基本（比如普通最小二乘）到複雜（比如迭代加權最小二乘法）的。

statsmodels的線性模型有兩種不同的介面：基於數組，和基於公式。它們可以通過API模塊引入：

import statsmodels.api as smimport statsmodels.formula.api as smf

為了展示它們的使用方法，我們從一些隨機數據生成一個線性模型：

def dnorm(mean, variance, size=1): if isinstance(size, int): size = size, return mean + np.sqrt(variance) * np.random.randn(*size)# For reproducibilitynp.random.seed(12345)N = 100X = np.c_[dnorm(0, 0.4, size=N), dnorm(0, 0.6, size=N), dnorm(0, 0.2, size=N)]eps = dnorm(0, 0.1, size=N)beta = [0.1, 0.3, 0.5]y = np.dot(X, beta) + eps

這裡，我使用了「真實」模型和可知參數beta。此時，dnorm可用來生成正太分布數據，帶有特定均值和方差。現在有：

In [66]: X[:5]Out[66]: array([[-0.1295, -1.2128, 0.5042], [ 0.3029, -0.4357, -0.2542], [-0.3285, -0.0253, 0.1384], [-0.3515, -0.7196, -0.2582], [ 1.2433, -0.3738, -0.5226]])In [67]: y[:5]Out[67]: array([ 0.4279, -0.6735, -0.0909, -0.4895,-0.1289])

像之前Patsy看到的，線性模型通常要擬合一個截距。sm.add_constant函數可以添加一個截距的列到現存的矩陣：

In [68]: X_model = sm.add_constant(X)In [69]: X_model[:5]Out[69]: array([[ 1. , -0.1295, -1.2128, 0.5042], [ 1. , 0.3029, -0.4357, -0.2542], [ 1. , -0.3285, -0.0253, 0.1384], [ 1. , -0.3515, -0.7196, -0.2582], [ 1. , 1.2433, -0.3738, -0.5226]])

sm.OLS類可以擬合一個普通最小二乘回歸：

In [70]: model = sm.OLS(y, X)

這個模型的fit方法返回了一個回歸結果對象，它包含估計的模型參數和其它內容：

In [71]: results = model.fit()In [72]: results.paramsOut[72]: array([ 0.1783, 0.223 , 0.501 ])

對結果使用summary方法可以列印模型的詳細診斷結果：

In [73]: print(results.summary())OLS Regression Results ==============================================================================Dep. Variable: y R-squared: 0.430Model: OLS Adj. R-squared: 0.413Method: Least Squares F-statistic: 24.42Date: Mon, 25 Sep 2017 Prob (F-statistic): 7.44e-12Time: 14:06:15 Log-Likelihood: -34.305No. Observations: 100 AIC: 74.61Df Residuals: 97 BIC: 82.42Df Model: 3 Covariance Type: nonrobust ============================================================================== coef std err t P>|t| [0.025 0.975]------------------------------------------------------------------------------x1 0.1783 0.053 3.364 0.001 0.073 0.283x2 0.2230 0.046 4.818 0.000 0.131 0.315x3 0.5010 0.080 6.237 0.000 0.342 0.660==============================================================================Omnibus: 4.662 Durbin-Watson: 2.201Prob(Omnibus): 0.097 Jarque-Bera (JB): 4.098Skew: 0.481 Prob(JB): 0.129Kurtosis: 3.243 Cond. No.1.74==============================================================================Warnings:[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

這裡的參數名為原始的名字x1, x2等等。假設所有的模型參數都在一個DataFrame中：

In [74]: data = pd.DataFrame(X, columns=[col0, col1, col2])In [75]: data[y] = yIn [76]: data[:5]Out[76]: col0 col1 col2 y0 -0.129468 -1.212753 0.504225 0.4278631 0.302910 -0.435742 -0.254180 -0.6734802 -0.328522 -0.025302 0.138351 -0.0908783 -0.351475 -0.719605 -0.258215 -0.4894944 1.243269 -0.373799 -0.522629 -0.128941

現在，我們使用statsmodels的公式API和Patsy的公式字元串：

In [77]: results = smf.ols(y ~ col0 + col1 + col2, data=data).fit()In [78]: results.paramsOut[78]: Intercept 0.033559col0 0.176149col1 0.224826col2 0.514808dtype: float64In [79]: results.tvaluesOut[79]: Intercept 0.952188col0 3.319754col1 4.850730col2 6.303971dtype: float64

觀察下statsmodels是如何返回Series結果的，附帶有DataFrame的列名。當使用公式和pandas對象時，我們不需要使用add_constant。

給出一個樣本外數據，你可以根據估計的模型參數計算預測值：

In [80]: results.predict(data[:5])Out[80]: 0 -0.0023271 -0.1419042 0.0412263 -0.3230704 -0.100535dtype: float64

statsmodels的線性模型結果還有其它的分析、診斷和可視化工具。除了普通最小二乘模型，還有其它的線性模型。

估計時間序列過程

statsmodels的另一模型類是進行時間序列分析，包括自回歸過程、卡爾曼濾波和其它態空間模型，和多元自回歸模型。

用自回歸結構和雜訊來模擬一些時間序列數據：

init_x = 4import randomvalues = [init_x, init_x]N = 1000b0 = 0.8b1 = -0.4noise = dnorm(0, 0.1, N)for i in range(N): new_x = values[-1] * b0 + values[-2] * b1 + noise[i] values.append(new_x)

這個數據有AR(2)結構（兩個延遲），參數是0.8和-0.4。當你和AR模型，你可能不知道滯後項的個數，因此可以用較多的滯後量來擬合這個模型：

In [82]: MAXLAGS = 5In [83]: model = sm.tsa.AR(values)In [84]: results = model.fit(MAXLAGS)

結果中的估計參數首先是截距，其次是前兩個參數的估計值：

In [85]: results.paramsOut[85]: array([-0.0062, 0.7845, -0.4085, -0.0136, 0.015 , 0.0143])

更多的細節以及如何解釋結果超出了本書的範圍，可以通過statsmodels文檔學習更多。