应用回归分析-一元线性回归——应用回归分析中一元线性回归的python代码

2301_78456108

于 2024-06-16 19:46:23 发布

阅读量467

点赞数 9

文章标签：回归线性回归 python

本文链接：https://blog.csdn.net/2301_78456108/article/details/139724638

版权

一元线性回归模型的数学形式

只考虑两个变量的关系，回归模型可表示为

y=β0+β1x+ε

其中，

y 是因变量（dependent variable），
x 是自变量（independent variable），
β0 是截距（intercept），
β1 是斜率（slope），
ε 是误差项（error term）。

满足的基本假设

线性关系假设：y=β0+β1x+ε
独立性假设：E(εi)=0 且Cov(εi,εj)=0 对所有i不等于j成立
正态性假设：ε∼N(0,σ2)
等方差性假设：Var(εi)=σ2 对所有 i成立
无多重共线性：y=β0+β1x1+β2x2+...+βkxk+ε 中的自变量 xk 之间不存在高度相关性。
无影响异常值和杠杆点：异常值和杠杆点对结果的影响可以通过评估离群值统计来检测。

参数β0，β1的估计

根据最小二乘法所得到的结果代码如下：

为方便代码编写，我们以取β0 =c ， β1=m。


import pandas as pd

df=pd.read_csv("数据示例.txt",sep="\t",header=None)
print(df)
x=df.loc[:,1]
y=df.loc[:,0]
n=len(x)

# 计算均值
mean_x = sum(x)/len(x)
mean_y = sum(y)/len(y)

# 计算斜率 m 和截距 c
Lxy = sum([xi*yi for xi, yi in zip(x, y)]) - n * mean_x * mean_y
Lxx = sum([xi**2 for xi in x]) - n * mean_x**2
m = Lxy/Lxx
c = mean_y - m * mean_x

输出结果并绘画出观测到以及回归直线的图形

import matplotlib.pyplot as plt

# 输出结果
print("斜率 m:", m)
print("截距 c:", c)
plt.rcParams['font.sans-serif'] = ['SimHei']  # 设置中文字体为黑体
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号

# 绘制数据和拟合直线
plt.plot(x, y, 'o', label='原始数据', markersize=10)
plt.plot(x, m*x + c, 'r', label='回归直线')
plt.legend()
plt.show()

回归方程的显著性检验

t检验

#回归方程的显著性检验T
#H0:m=0,y与x之间没有真正的线性关系
#H1：m!=0,y与x之间有真正的线性关系


residuals = [yi - (c + m * xi) for xi, yi in zip(x, y)]
v = sum([residual**2 for residual in residuals]) / (n - 2)#回归方差
SE_slope = np.sqrt(v / Lxx)#m的标准差

# 计算t统计量
t_statistic = m / SE_slope
# 在自由度为 n-2 的 t 分布下查找p值
p_value = stats.t.sf(np.abs(t_statistic), n-2) * 2

# 根据p值判断显著性
alpha = 0.05
if p_value < alpha:
    print("回归方程显著")
else:
    print("回归方程不显著")

print("回归系数的t统计量:", t_statistic)
print("回归系数的p值:", p_value)

F检验

SST=SSR+SSE

构造F统计量为F=SSR/p/SSE/（n-p-1）其中p=1


#回归方程的显著性检验F
residuals =[yi - (c + m * xi) for xi, yi in zip(x, y)]
SSE = sum([residual**2 for residual in residuals])
residuals =[(yi - mean_y) for yi in y]
SST=sum([residual**2 for residual in residuals])
residuals =[(c + m * xi)-mean_y for xi in x]
SSR=sum([residual**2 for residual in residuals])
#print(SST,SSE+SSR)
F=(SSR/1)/(SSE/(n-2))
print("F值",F)
p_value = stats.f.sf(F, 1, n-2)
print("F值对应的P值", p_value)

回归系数的区间估计以及因变量的区间预测

vei为残差的方差，h为杠杆值


#残差的方差
a=[(xi-mean_x)**2 for xi in x]
vei=[(1-1/n-s/Lxx)*v for s in a]
h=[0]*n
for i in range(n):
    h[i]=1/n+(x[i]-mean_x)**2/Lxx
print(h)
# 获取自由度为1的t值
t_value = stats.t.ppf(1-0.05/2, n-2)
print(t_value)
#回归系数的区间估计
print("B1的1-a置信区间为：",(m-t_value*(v/Lxx)**0.5,m+t_value*(v/Lxx)**0.5))
print("B0的1-a置信区间为：",(c-t_value*(v*(1/n+(mean_x*mean_x)/Lxx))**0.5,c+t_value*(v*(1/n+(mean_x*mean_x)/Lxx))**0.5))
x0=3.5
y0=m*x0+c
print("y0的点估计值:",m*x0+c)
print("单个新值：",(y0-t_value*((1+h[0])*v)**0.5,y0+t_value*((1+h[0])*v)**0.5))
print("平均值E（y0）：",(y0-t_value*(h[0]*v)**0.5,y0+t_value*(h[0]*v)**0.5))

代码完整版如下

from scipy import stats
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
df=pd.read_csv("数据示例.txt",sep="\t",header=None)
print(df)
x=df.loc[:,1]
y=df.loc[:,0]
n=len(x)

# 计算均值
mean_x = sum(x)/len(x)
mean_y = sum(y)/len(y)

# 计算斜率 m 和截距 c
Lxy = sum([xi*yi for xi, yi in zip(x, y)]) - n * mean_x * mean_y
Lxx = sum([xi**2 for xi in x]) - n * mean_x**2
m = Lxy/Lxx
c = mean_y - m * mean_x

# 输出结果
print("斜率 m:", m)
print("截距 c:", c)
plt.rcParams['font.sans-serif'] = ['SimHei']  # 设置中文字体为黑体
plt.rcParams['axes.unicode_minus'] = False  # 用来正常显示负号

# 绘制数据和拟合直线
plt.plot(x, y, 'o', label='原始数据', markersize=10)
plt.plot(x, m*x + c, 'r', label='回归直线')
plt.legend()
plt.show()

#回归方程的显著性检验T
#H0:m=0,y与x之间没有真正的线性关系
#H1：m!=0,y与x之间有真正的线性关系

residuals = [yi - (c + m * xi) for xi, yi in zip(x, y)]
v = sum([residual**2 for residual in residuals]) / (n - 2)#回归方差
SE_slope = np.sqrt(v / Lxx)#m的标准差

# 计算t统计量
t_statistic = m / SE_slope
# 在自由度为 n-2 的 t 分布下查找p值
p_value = stats.t.sf(np.abs(t_statistic), n-2) * 2

# 根据p值判断显著性
alpha = 0.05
if p_value < alpha:
    print("回归方程显著")
else:
    print("回归方程不显著")

print("回归系数的t统计量:", t_statistic)
print("回归系数的p值:", p_value)

#回归方程的显著性检验F
residuals =[yi - (c + m * xi) for xi, yi in zip(x, y)]
SSE = sum([residual**2 for residual in residuals])
residuals =[(yi - mean_y) for yi in y]
SST=sum([residual**2 for residual in residuals])
residuals =[(c + m * xi)-mean_y for xi in x]
SSR=sum([residual**2 for residual in residuals])
#print(SST,SSE+SSR)
F=(SSR/1)/(SSE/(n-2))
print("F值",F)
p_value = stats.f.sf(F, 1, n-2)
print("F值对应的P值", p_value)

#相关系数的显著性检验
Lyy=SST
r=Lxy/(Lxx*Lyy)**(1/2)
print(r)
print(m*(Lxx/SST)**0.5)
t=(n-2)**0.5*r/(1-r**2)**0.5
p_value = stats.t.sf(np.abs(t), n-2) * 2
print(p_value)

#残差的方差
a=[(xi-mean_x)**2 for xi in x]
vei=[(1-1/n-s/Lxx)*v for s in a]
h=[0]*n
for i in range(n):
    h[i]=1/n+(x[i]-mean_x)**2/Lxx
print(h)
# 获取自由度为1的t值
t_value = stats.t.ppf(1-0.05/2, n-2)
print(t_value)
#回归系数的区间估计
print("B1的1-a置信区间为：",(m-t_value*(v/Lxx)**0.5,m+t_value*(v/Lxx)**0.5))
print("B0的1-a置信区间为：",(c-t_value*(v*(1/n+(mean_x*mean_x)/Lxx))**0.5,c+t_value*(v*(1/n+(mean_x*mean_x)/Lxx))**0.5))
x0=3.5
y0=m*x0+c
print("y0的点估计值:",m*x0+c)
print("单个新值：",(y0-t_value*((1+h[0])*v)**0.5,y0+t_value*((1+h[0])*v)**0.5))
print("平均值E（y0）：",(y0-t_value*(h[0]*v)**0.5,y0+t_value*(h[0]*v)**0.5))

数据格式如下

数据示例.txt
31.30   4.30
24.00   2.10
17.30   1.10
43.20   6.10
36.40   4.80
26.10   3.80

2301_78456108

关注

9
点赞
踩
6

收藏

觉得还不错? 一键收藏
1
评论
应用回归分析-一元线性回归——应用回归分析中一元线性回归的python代码

应用回归分析中一元线性回归的python代码
复制链接

扫一扫