Numpy线性回归：散点图正常但拟合直线异常问题排查

我用numpy对全球各国性别不平等指数（x）和生育率（y）做线性回归，散点图显示正常，x均值、y均值和相关系数看起来是对的（161个数据点未手动验算），但截距（alpha）和斜率（beta）计算错误，预测的生育率范围在-13到19之间（实际数据为1到7），怀疑协方差与方差的计算存在问题。

原代码 #

import pandas as pd
import numpy as np
from matplotlib import pyplot as plt

data = pd.read_csv("GII Data.csv")
df = pd.DataFrame(data)
print(df, "\n")

xmean = np.mean(x)
ymean = np.mean(y)

x = data.GII 
y = data.fertility_rate

print(f'xmean = {xmean}')
print(f'ymean = {ymean}\n')

covariance = np.cov(x, y)
print(f'covariance = {covariance}\n')

variance = (np.var(x))
print(f'variance = {variance}\n')  

# Calculate beta and alpha
beta = covariance.sum() / variance.sum()
alpha = ymean - (beta * xmean)
print(f'alpha = {alpha}')
print(f'beta = {beta}\n')
correlation = np.corrcoef(x, y)
print(f'correlation coefficient = {correlation}\n')

ypred = alpha + beta * x 

plt.figure(figsize=(8, 4))
plt.plot(x, y, 'ro')   # scatter plot showing actual data
plt.title('Gender Inequality and Children per Woman: r = 0.78')
plt.xlabel('Gender Inequality Index (GII)')
plt.ylabel('Fertility Rate')
plt.plot(x, ypred)     # regression line
plt.show()

问题点分析 #

• 变量定义顺序错误：计算xmean和ymean时，x和y还未赋值，会导致取到旧变量值或直接报错，必须先定义x、y再计算均值。
• 协方差使用错误：np.cov(x,y)返回2×2的协方差矩阵，其中[0,1]位置是x与y的协方差，直接sum()会把矩阵所有元素相加，应该取交叉协方差项。
• 方差自由度不匹配：np.var(x)默认计算总体方差（除以n），而np.cov默认计算样本协方差（除以n-1），两者自由度不一致，会导致斜率计算偏差。
• 冗余的sum()调用：np.var(x)返回标量，variance.sum()完全多余，直接使用variance即可。

修正后的代码 #

import pandas as pd
import numpy as np
from matplotlib import pyplot as plt

data = pd.read_csv("GII Data.csv")
df = pd.DataFrame(data)
print(df, "\n")

# 先定义变量x、y，再计算均值
x = data.GII 
y = data.fertility_rate

xmean = np.mean(x)
ymean = np.mean(y)

print(f'xmean = {xmean}')
print(f'ymean = {ymean}\n')

# 获取x与y的样本协方差（取协方差矩阵的交叉项）
covariance = np.cov(x, y)[0, 1]
print(f'covariance = {covariance}\n')

# 计算x的样本方差，与协方差自由度保持一致
variance = np.var(x, ddof=1)
# 也可以用np.cov直接计算样本方差：variance = np.cov(x, x)[0,0]
print(f'variance = {variance}\n')  

# 正确计算斜率beta和截距alpha
beta = covariance / variance
alpha = ymean - (beta * xmean)
print(f'alpha = {alpha}')
print(f'beta = {beta}\n')

correlation = np.corrcoef(x, y)
print(f'correlation coefficient = {correlation}\n')

ypred = alpha + beta * x 

plt.figure(figsize=(8, 4))
plt.plot(x, y, 'ro')   # 散点图展示实际数据
plt.title('性别不平等指数与女性生育子女数: r = 0.78')
plt.xlabel('性别不平等指数（GII）')
plt.ylabel('生育率')
plt.plot(x, ypred)     # 回归线
plt.show()

内容的提问来源于stack exchange，提问作者BretteConnolly