帝王谷资源网

Python进行统计建模
脚本专栏 2024/12/25 佚名
3 2
帝王谷资源网 Design By www.wdxyy.com

前言

大家好，在之前的文章中我们已经讲解了很多Python数据处理的方法比如读取数据、缺失值处理、数据降维等，也介绍了一些数据可视化的方法如Matplotlib、pyecharts等，那么在掌握了这些基础技能之后，要进行更深入的分析就需要掌握一些常用的建模方法，本文将讲解如何利用Python进行统计分析。和之前的文章类似，本文只讲如何用代码实现，不做理论推导与过多的结果解释(事实上常用的模型可以很轻松的查到完美的推导与解析)。因此读者需要掌握一些基本的统计模型比如回归模型、时间序列等。

Statsmodels简介

在Python 中统计建模分析最常用的就是Statsmodels模块。Statsmodels是一个主要用来进行统计计算与统计建模的Python库。主要有以下功能：

  • 探索性分析：包含列联表、链式方程多重插补等探索性数据分析方法以及与统计模型结果的可视化图表，例如拟合图、箱线图、相关图、时间序列图等
  • 回归模型：线性回归模型、非线性回归模型、广义线性模型、线性混合效应模型等
  • 其他功能：方差分析、时间序列分析等模型的参数估计与估计参数的假设检验等

安装 brew install Statsmodels
文档 github.com/statsmodels/statsmodels

线性回归模型：普通最小二乘估计

线性模型有普通最小二乘(OLS)广义最小二乘(GLS)、加权最小二乘(WLS)等，Statsmodels对线性模型有较好的支持，来看个最简单的例子：普通最小二乘(OLS)

首先导入相关包

%matplotlib inline
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
from statsmodels.sandbox.regression.predstd import wls_prediction_std
np.random.seed(9876789)

然后创建数据，先设置样本量为100

nsample = 100 #样本数量

然后设置x1和x2，x1是0到10等差排列，x2是x1的平方

x = np.linspace(0, 10, 100)
X = np.column_stack((x, x**2))

再设置beta、误差项与响应变量y

beta = np.array([1, 0.1, 10])
e = np.random.normal(size=nsample)
X = sm.add_constant(X)
y = np.dot(X, beta) + e

接着建立回归模型

model = sm.OLS(y, X) 
results = model.fit()
print(results.summary())

查看模型结果

Python进行统计建模

是不是和R语言输出的结果形式很接近？回归系数值、P-value、R-squared等评估回归模型的参数值全部都有，还可以使用dir(results)获得全部变量的值并调取出来

print('Parameters: ', results.params)
print('R2: ', results.rsquared)

那么回归模型的就是y=1.3423-0.0402x1+10.0103x2，当然这个模型可以继续优化那么就交给读者完成。接下来我们来绘制一下样本点与回归曲线

y_fitted = results.fittedvalues
fig, 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')

Python进行统计建模

时间序列：ARMA

关于时间序列的模型有很多，我们选择ARMA模型示例，首先导入相关包并生成数据

%matplotlib inline
import numpy as np
import statsmodels.api as sm
import pandas as pd
from statsmodels.tsa.arima_process import arma_generate_sample
np.random.seed(12345)

arparams = np.array([.75, -.25])
maparams = np.array([.65, .35])

arparams = np.r_[1, -arparams]
maparams = np.r_[1, maparams]
nobs = 250
y = arma_generate_sample(arparams, maparams, nobs)

接着，我们可以添加一些日期信息。对于本例，我们将使用pandas时间序列并建立模型

dates = sm.tsa.datetools.dates_from_range('1980m1', length=nobs)
y = pd.Series(y, index=dates)
arma_mod = sm.tsa.ARMA(y, order=(2,2))
arma_res = arma_mod.fit(trend='nc', disp=-1)

Python进行统计建模

最后再做一下预测

import matplotlib.pyplot as plt
fig, ax = plt.subplots(figsize=(10,8))
fig = arma_res.plot_predict(start='1999-06-30', end='2001-05-31', ax=ax)
legend = ax.legend(loc='upper left')

Python进行统计建模

回归诊断：估计回归模型

首先导入相关包

%matplotlib inline
from statsmodels.compat import lzip
import numpy as np
import pandas as pd
import statsmodels.formula.api as smf
import statsmodels.stats.api as sms
import matplotlib.pyplot as plt

然后加载数据

url = 'https://raw.githubusercontent.com/vincentarelbundock/Rdatasets/master/csv/HistData/Guerry.csv'
dat = pd.read_csv(url)

拟合模型

results = smf.ols('Lottery ~ Literacy + np.log(Pop1831)', data=dat).fit()

查看结果

print(results.summary())

Python进行统计建模

回归诊断：残差的正态性

Jarque-Bera test:

name = ['Jarque-Bera', 'Chi^2 two-tail prob.', 'Skew', 'Kurtosis']
test = sms.jarque_bera(results.resid)
lzip(name, test)
####结果
[('Jarque-Bera', 3.3936080248431666),
('Chi^2 two-tail prob.', 0.1832683123166337),
('Skew', -0.48658034311223375),
('Kurtosis', 3.003417757881633)]

Omni test:

name = ['Chi^2', 'Two-tail probability']
test = sms.omni_normtest(results.resid)
lzip(name, test)
####结果
[('Chi^2', 3.713437811597181), ('Two-tail probability', 0.15618424580304824)]

回归诊断：异方差

Breush-Pagan test:

name = ['Lagrange multiplier statistic', 'p-value',
    'f-value', 'f p-value']
test = sms.het_breuschpagan(results.resid, results.model.exog)
lzip(name, test)
###结果
[('Lagrange multiplier statistic', 4.893213374093957),
('p-value', 0.08658690502352209),
('f-value', 2.503715946256434),
('f p-value', 0.08794028782673029)]
Goldfeld-Quandt test

name = ['F statistic', 'p-value']
test = sms.het_goldfeldquandt(results.resid, results.model.exog)
lzip(name, test)
####结果
[('F statistic', 1.1002422436378152), ('p-value', 0.3820295068692507)]

回归诊断：多重共线性

检查多重共线性可以使用

np.linalg.cond(results.model.exog)

结果是702.1792145490062，说明存在较强多重共线性。

结束语

以上就是Statsmodels的基本功能介绍，如果熟悉R的读者会发现很多命令与R是类似的。最后想多说一句，全文没有出现太多模型的理论知识，因为这些模型的推导过程随便百度一搜都能得到十分详细的优质回答，因此在学会如何用计算机实现之后必须要回过头去理解模型里每一个参数是怎样得到，又有哪些含义才算真正搞定。

以上就是Python进行统计建模的详细内容，更多关于Python统计建模的资料请关注其它相关文章！

标签：
Python,统计,Python,建模,Python,统计建模
帝王谷资源网 Design By www.wdxyy.com
广告合作：本站广告合作请联系QQ：858582 申请时备注：广告合作（否则不回）
免责声明：本站文章均来自网站采集或用户投稿，网站不提供任何软件下载或自行开发的软件！ 如有用户或公司发现本站内容信息存在侵权行为，请邮件告知！ 858582#qq.com
帝王谷资源网 Design By www.wdxyy.com
www.wdxyy.com 帝王谷资源网
139,976互联网资源
144,792高清电影
21,817无损音乐
631,128技术资源

更新日志

友情链接

杰晶网络 DDR爱好者之家 桃源资源网 杰网资源 富贵资源网 南强小屋 铁雪资源网 幽灵资源网 万梅资源网 狼山资源网 白云岛资源网 昆仑资源网 相思资源网 明霞山资源网 内蒙古资源网 黑松山资源网 茶园资源网 饿虎岗资源网 大旗谷资源网 常春岛资源网 岱庙资源网 兴国资源网 快活林资源网 蝙蝠岛资源网 帝王谷资源网 白云城资源网 伏龙阁资源网 清风细雨楼 天枫庄资源网 圆月山庄资源网 无争山庄资源网 神水资源网 移花宫资源网 神剑山庄资源网 无为清净楼资源网 金钱帮资源网 丐帮资源网 华山资源网 极乐门资源网 小李飞刀资源网 凤求凰客栈 风云阁资源网 金狮镖局 鸳鸯亭资源网 千金楼资源网 更多链接