statsmodels︱python常规统计模型库

之前看sklearn线性模型没有R方,F检验,回归系数T检验等指标,于是看到了statsmodels这个库,看着该库输出的结果真是够怀念的。。



1 安装

pip install statsmodels

不过有可能会报错:

ImportError: cannot import name 'factorial' from 'scipy.misc' 
(E:\Anaconda3.7\lib\site-packages\scipy\misc\__init__.py)

是跟scipy版本不匹配,笔者是删掉之前的pip uninstall statsmodels,再重新安装了一下就好了:

pip install --pre statsmodels -i https://pypi.tuna.tsinghua.edu.cn/simple

2 相关模型介绍

相关文档可见:https://www.statsmodels.org/stable/examples/index.html

在这里插入图片描述

包含的模型有:

2.1 线性模型

在这里插入图片描述

2.2 离散选择模型(Discrete Choice Model, DCM)

在这里插入图片描述

参考:离散选择模型(Discrete Choice Model, DCM)简介——之一

离散选择模型(Discrete Choice Model, DCM)在经济学领域和社会学领域都有广泛的应用。
例如,消费者在购买汽车的时候通常会比较几个不同的品牌,如福特、本田、大众,等等。
如果将消费者选择福特汽车记为Y=1,选择本田汽车记为Y=2,选择大众汽车记为Y=3;那么在研究消费者选择何种汽车品牌的时候,由于因变量不是一个连续的变量(Y=1, 2, 3),传统的线性回归模型就有一定的局限(见DCM系列文章第2篇)。
再比如,在交通安全研究领域,通常将交通事故的严重程度划分为3大类:

  • (1)仅财产损失(Property Damage Only, PDO),
  • (2)受伤(Injury),
  • (3)死亡(Fatality);
    在研究各类因素(如道路坡度、弯道曲率等、车龄、光照、天气条件等)对事故严重程度的影响的时候,由于因变量(事故严重程度)是一个离散变量(仅3个选项),使用离散选择模型可以提供一个有效的建模途径。
    在这里插入图片描述

2.3 非参数统计

在这里插入图片描述

2.4 广义线性模型 - Generalized Linear Models

在这里插入图片描述

2.5 稳健回归——Robust Regression

在这里插入图片描述

2.6 广义估计方程

在这里插入图片描述

2.7 方差分析

在这里插入图片描述

2.8 时间序列分析——Time Series Analysis

在这里插入图片描述

2.9 空间计量必备:状态空间模型——State space models

在这里插入图片描述

在这里插入图片描述

2.10 多元统计模型——因子/主成分分析

在这里插入图片描述


3 相关模型demo

3.1 线性回归模型

可参考:https://www.statsmodels.org/stable/examples/notebooks/generated/ols.html

#  线性模型
import statsmodels.api as sm
import numpy as np
x = np.linspace(0,10,100)
y = 3*x + np.random.randn()+ 10
# Fit and summarize OLS model
X = sm.add_constant(x)
mod = sm.OLS(y,X)
result = mod.fit()
print('Parameters: ', result .params)
print('Standard errors: ', result .bse)
print('Predicted values: ', result .predict())
print(result.summary())


# 预测数据
print(result.predict(X[:5]))

输出结果超级熟悉。

  • result.params是回归系数
  • result.summary()把模型相关系数都打印出来
    其中,预测的时候,如果不给入参数result.predict(),则默认是X

在这里插入图片描述

在这里插入图片描述

3.2 广义线性模型——GLM

参考:https://www.statsmodels.org/stable/examples/notebooks/generated/glm.html

import statsmodels.formula.api as smf
star98 = sm.datasets.star98.load_pandas().data
formula = 'SUCCESS ~ LOWINC + PERASIAN + PERBLACK + PERHISP + PCTCHRT + \
           PCTYRRND + PERMINTE*AVYRSEXP*AVSALK + PERSPENK*PTRATIO*PCTAF'
dta = star98[['NABOVE', 'NBELOW', 'LOWINC', 'PERASIAN', 'PERBLACK', 'PERHISP',
              'PCTCHRT', 'PCTYRRND', 'PERMINTE', 'AVYRSEXP', 'AVSALK',
              'PERSPENK', 'PTRATIO', 'PCTAF']].copy()
endog = dta['NABOVE'] / (dta['NABOVE'] + dta.pop('NBELOW'))
del dta['NABOVE']
dta['SUCCESS'] = endog

mod1 = smf.glm(formula=formula, data=dta, family=sm.families.Binomial()).fit()
mod1.summary()
mod1.predict(dta)

formula是常规的公式,其中所有X/Y数据都放在一个dataframe之中。
在这里插入图片描述

print('Total number of trials:',  data.endog[0].sum())
print('Parameters: ', res.params)
print('T-values: ', res.tvalues)

在这里插入图片描述

包括了回归系数,T检验值

3.3 稳健回归

参考:https://www.statsmodels.org/stable/examples/notebooks/generated/robust_models_0.html

nsample = 50
x1 = np.linspace(0, 20, nsample)
X = np.column_stack((x1, (x1-5)**2))
X = sm.add_constant(X)
sig = 0.3   # smaller error variance makes OLS<->RLM contrast bigger
beta = [5, 0.5, -0.0]
y_true2 = np.dot(X, beta)
y2 = y_true2 + sig*1. * np.random.normal(size=nsample)
y2[[39,41,43,45,48]] -= 5   # add some outliers (10% of nsample)

X2 = X[:,[0,1]]
res2 = sm.OLS(y2, X2).fit()
print(res2.params)
print(res2.bse)


resrlm2 = sm.RLM(y2, X2).fit()
print(resrlm2.params)
print(resrlm2.bse)
print(resrlm2.summary())

在这里插入图片描述


4 其他

4.1 模型结果如何CSV导出?

可以通过as_csv()将模型导出

resrlm2 = sm.RLM(y, x).fit()
resrlm2.summary()
with open( 'model_rlm.csv', 'w') as fh:
    fh.write(resrlm2.summary().as_csv())

不过导出的格式比较奇怪:
在这里插入图片描述

4.2 画模型图以及保存

import statsmodels.api as sm
import numpy as np
import matplotlib.pyplot as plt

# 准备数据
x = np.linspace(0,10,100)
y = 3*x + np.random.randn()+ 10

# Fit and summarize OLS model
res = sm.OLS(y,x).fit()
print(res.params)
print(res.summary())
# 稳健回归
resrlm = sm.RLM(y, x).fit()

# 画图
fig, ax = plt.subplots(figsize=(8,6))
ax.plot(x, y, 'o', label="truey ")
ax.plot(x, res.predict(), 'o', label="ols")  # res2.predict(X2) == res2.predict()
ax.plot(x, resrlm.predict(), 'b-', label="rlm")# resrlm2.predict(X2) == resrlm2.predict()
legend = ax.legend(loc="best")

# 图保存
plt.savefig( 'image.jpg')

4.3 快速获取模型输出参数:P检验、F检验、P统计量

def get_model_param(res2,name = 'all'):
    model_param_dict = {'name':name,  # 模型的名字
                        'rsquared':res2.rsquared, # R方
                        'fvalue':res2.fvalue, # F值,整个模型
                        'f_pvalue':res2.f_pvalue, # P值,整个模型
                        'params':res2.params[0], # 回归系数
                        'pvalues':res2.pvalues[0], # 回归系数 P检验 0.000
                        'tvalues':res2.tvalues[0]} # 回归系数 T检验 276.571
    return model_param_dict

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值