python3 线性回归验证方法
作者:SunWuKong_Hadoop 发布时间:2022-12-03 10:49:09
标签:python3,线性,回归,验证
如下所示:
#-*- coding: utf-8 -*-
import pandas as pd
import numpy as np
from patsy.highlevel import dmatrices
#2.7里面是from patsy import dmatrices
from statsmodels.stats.outliers_influence import variance_inflation_factor
import statsmodels.api as sm
import scipy.stats as stats
from sklearn.metrics import mean_squared_error
import seaborn as sns
import matplotlib.pyplot as plt
import matplotlib.mlab as mlab
import matplotlib
#数据获取
ccpp = pd.read_excel('CCPP.xlsx')
ccpp.describe()
#绘制各变量之间的散点图
sns.pairplot(ccpp)
plt.show()
#发电量(PE)与自变量之间的相关系数
a = ccpp.corrwith(ccpp.PE)
print(a)
#将因变量PE,自变量AT,V,AP和截距项(值为1的1维数值)以数据框的形式组合起来
y,x = dmatrices('PE~AT+V+AP',data = ccpp,return_type = 'dataframe')
#构造空的数据框
vif = pd.DataFrame()
vif[""VIF Factor""] = [variance_inflation_factor(x.values,i) for i in range(x.shape[1])]
vif[""features""] = x.columns
print (vif)
#构建PE与AT,V和AP之间的线性模型
fit = sm.formula.ols('PE~AT+V+AP',data=ccpp).fit()
b = fit.summary()
# print(b)
#计算模型的RMSE值
pred = fit.predict()
c = np.sqrt(mean_squared_error(ccpp.PE,pred))
print(c)
#离群点检验
outliers = fit.get_influence()
#高杠杆值点(帽子矩阵)
leverage = outliers.hat_matrix_diag
#dffits值
dffits = outliers.dffits[0]
#学生化残差
resid_stu = outliers.resid_studentized_external
#cook距离
cook = outliers.cooks_distance[0]
#covratio值
covratio = outliers.cov_ratio
#将上面的几种异常值检验统计量与原始数据集合并
contat1 = pd.concat([pd.Series(leverage,name = 'leverage'),pd.Series(dffits,name ='dffits'),
pd.Series(resid_stu,name = 'resid_stu'),pd.Series(cook,name = 'cook'),
pd.Series(covratio,name ='covratio'),],axis = 1)
ccpp_outliers = pd.concat([ccpp,contat1],axis = 1)
d = ccpp_outliers.head()
print(d)
#计算异常值数量的比例
outliers_ratio = sum(np.where((np.abs(ccpp_outliers.resid_stu)>2),1,0))/ccpp_outliers.shape[0]
e = outliers_ratio
print(e)
#删除异常值
ccpp_outliers = ccpp_outliers.loc[np.abs(ccpp_outliers.resid_stu)<=2,]
#重新建模
fit2 = sm.formula.ols('PE~AT+V+AP',data = ccpp_outliers).fit()
f = fit2.summary()
# print(f)
pred2 = fit2.predict()
g = np.sqrt(mean_squared_error(ccpp_outliers.PE,pred2))
print(g)
#
#残差的正态性检验(直方图法)
resid = fit2.resid
#中文和负号的正常显示
# plt.rcParams['font.sans=serif'] = ['Microsoft YaHei']
plt.rcParams['font.sans-serif'] = ['SimHei']
# plt.rcParams['font.sans=serif'] = 'sans-serif'
plt.rcParams['axes.unicode_minus'] = False
plt.hist(resid,bins = 100,normed = True,color = 'steelblue',edgecolor = 'k')
#设置坐标轴标签和标题
plt.title('残差直方图')
plt.ylabel('密度值')
#生成正态曲线的数据
x1 = np.linspace(resid.min(),resid.max(),1000)
normal = mlab.normpdf(x1,resid.mean(),resid.std())
#绘制正态分布曲线
plt.plot(x1,normal,'r-',linewidth = 2,label = '正态分布曲线')
#生成核密度曲线的数据
kde = mlab.GaussianKDE(resid)
x2 = np.linspace(resid.min(),resid.max(),1000)
#绘制核密度曲线
plt.plot(x2,kde(x2),'k-',linewidth = 2,label = '核密度曲线')
#去除图形顶部边界和右边界的刻度
plt.tick_params(top = 'off',right = 'off')
#显示图例
plt.legend(loc='best')
#显示图形
plt.show()
#生成的正态曲线的数据
pp_qq_plot = sm.ProbPlot(resid)
pp_qq_plot.ppplot(line = '45')
plt.title('P-P图')
pp_qq_plot.qqplot(line = 'q')
plt.title('Q-Q图')
plt.show()
#残差的正态性检验(非参数法)
standard_resid = (resid-np.mean(resid))/np.std(resid)
g = stats.kstest(standard_resid,'norm')
print(g)
# 总结:由于shapiro正态性检验对样本量的需求是5000以内,而本次数据集样本量有9000多,故选择k-s来完成正态性检验。
# 从k-s检验的p值来看,拒绝了残差服从正态分布的假设,即认为残差并不满足正态性假设这个前提。
# 如果残差不服从正态分布的话,建议对Y变量进行box-cox变换处理。
# 由于fit2模型的残差并没有特别明显的偏态(偏度为0.058,接近于0),故这里就不对Y进行变换。
#
# import scipy.stats as stats
# #找到box-cox变换的Lambda系数
# lamd = stats.boxcox_normmax(vif.y,method = 'mle')
# #对y进行变换
# vif['trans_y'] = stats.boxcox(vif.y,lamd)
# #建模
# fit3 = sm.formula.ols('y~x1+x2...',data = vif).fit()
# fit3.summary()
来源:https://blog.csdn.net/SunWuKong_Hadoop/article/details/80254848
0
投稿
猜你喜欢
- HTMLParser是python用来解析html的模块。它可以分析出html里面的标签、数据等等,是一种处理html的简便途径。HTMLP
- Python“json.decoder.JSONDecodeError: Expecting value: line 1
- 数据库的使用过程中由于程序方面的问题有时候会碰到重复数据,重复数据导致了数据库部分设置不能正确设置……方法一以下为引用的内容:declare
- 最近经常遇到xmlHttp组件的问题, 今天终于有空将解决过程中的经验整理下来了! ^o^一、在运用xmlhttp组件编
- 本文实例讲述了Python基于Matplotlib库简单绘制折线图的方法。分享给大家供大家参考,具体如下:Matplotlib画折线图,有一
- Software as a service 软件即服务,21世纪开始兴起的一种完全创新的软件应用模式。客户通过互联网向厂商定购所需的应用软件
- 1.我们可以为每一个实例对象增加方法。也就是说我们在每次使用‘类'之外的方法时候,都需要创建一次。 function D
- 最近在做一个电信的MIS系统,由于数据库的数据量庞大(最少也有500万),发现了一个sql server 长时间占有内存的现象。当你查询数据
- 占位符说明1、%c,格式化字符及其ASCII码2、%s,格式化字符串3、%d,格式化整数4、%u,格式化无符号整数5、%o,格式化无符号八进
- Python 输出 "Hello, World!",英文没有问题,但是如果你输出中文字符"你好,世界"
- 代码如下:USE TestDB declare @conversation uniqueidentifier while exists (s
- 第一步:创建转向控制页面创建网站默认的首页文件(通常为"index.asp"或"default.asp&quo
- 笛卡尔坐标系对于平面坐标系,任一射线OP与x轴夹角θ的范围,可以取[0,2π)或者(-π,&
- ASP生成柱型体,折线图,饼图源代码。一:纯ASP代码生成图表函数2——折线图;二:纯ASP代码生成图表函数1——柱状图 ;三:纯
- 错误代码错误消息说 明ASP0100Out of memory内存不足(不能分配要求的内存)ASP0101Unexpected
- zip文件是我们经常使用的打包格式之一,python解压和压缩zip效率非凡。 python解压zip文档:#/usr/bin/python
- 注意,在改变数值之前锁定应用,确保一段时间里只有一个客户执行该语句。<SCRIPT LANGUAGE="VBScr
- 问题你写了一个装饰器作用在某个函数上,但是这个函数的重要的元信息比如名字、文档字符串、注解和参数签名都丢失了。解决方案任何时候你定义装饰器的
- PDOStatement::errorInfoPDOStatement::errorInfo — 获取跟上一次语句句柄操作相关的扩展错误信息
- 先谈一下面包屑的由来:很久很久以前,在大森林的边上住着一个贫穷的樵夫,他妻子和两个孩子与他相依为命。他的儿子名叫汉赛尔,女儿名叫格莱特。后来