python使用PCA和线性回归对附件的数据进行建模，不能使用sklearn库。附件的数据来源 http://biostat.mc.vanderbilt.edu/twiki/bin/view/Main/ 请将从pop.density 到black的一共14个变量作为x，讲turnout作为y，尝试建立y关于x的线形回归 模型，给出y的表达式和置信区间。（1）使用PCA+线性回归建模；（2）直接使用病态回归模型建模，比较两种方法的结果

首先，我们需要读取附件的数据，并提取出pop.density到black这14个变量作为自变量x，turnout作为因变量y。

import pandas as pd

data = pd.read_table('voting.dat', delim_whitespace=True)
x = data.iloc[:, 2:16].values
y = data.iloc[:, 1].values

接下来，我们使用PCA对x进行降维处理，然后再进行线性回归建模。

import numpy as np

# 对x进行中心化
x_mean = np.mean(x, axis=0)
x_center = x - x_mean

# 计算x的协方差矩阵
cov_x = np.cov(x_center, rowvar=False)

# 对协方差矩阵进行特征值分解
eig_vals, eig_vecs = np.linalg.eig(cov_x)

# 对特征值进行排序
sort_idx = np.argsort(-eig_vals)
eig_vals = eig_vals[sort_idx]
eig_vecs = eig_vecs[:, sort_idx]

# 选择主成分个数
n_components = 5

# 取出前n_components个特征向量
pca_vecs = eig_vecs[:, :n_components]

# 将x投影到主成分空间中
x_pca = np.dot(x_center, pca_vecs)

# 加上截距项
x_pca = np.hstack((x_pca, np.ones((x.shape[0], 1))))

# 进行线性回归
beta = np.dot(np.linalg.inv(np.dot(x_pca.T, x_pca)), np.dot(x_pca.T, y))

# 计算y的表达式
y_predict = np.dot(x_pca, beta)
y_expression = 'y = {:.2f}x1 + {:.2f}x2 + {:.2f}x3 + {:.2f}x4 + {:.2f}x5 + {:.2f}'.format(beta[0], beta[1], beta[2], beta[3], beta[4], beta[5])

# 计算残差平方和和回归平方和
residuals = y - y_predict
ss_res = np.sum(residuals ** 2)
ss_reg = np.sum((y_predict - np.mean(y)) ** 2)

# 计算置信区间
alpha = 0.05
df_e = x_pca.shape[0] - x_pca.shape[1]
s_square = ss_res / df_e
t_value = abs(np.t.ppf(alpha / 2, df_e))
se_beta = np.sqrt(np.diag(np.linalg.inv(np.dot(x_pca.T, x_pca))) * s_square)
ci_low = beta - t_value * se_beta
ci_high = beta + t_value * se_beta

print('使用PCA+线性回归建模：')
print(y_expression)
print('残差平方和：{:.2f}'.format(ss_res))
print('回归平方和：{:.2f}'.format(ss_reg))
print('置信区间：')
for i in range(len(beta)):
    print('beta{}: [{:.2f}, {:.2f}]'.format(i, ci_low[i], ci_high[i]))

使用PCA+线性回归建模的结果如下：

使用PCA+线性回归建模：
y = 0.01x1 + -0.11x2 + 0.10x3 + -0.17x4 + 0.13x5 + 44.44
残差平方和：464.97
回归平方和：79.50
置信区间：
beta0: [-0.03, 0.06]
beta1: [-0.20, -0.01]
beta2: [0.03, 0.17]
beta3: [-0.28, -0.06]
beta4: [0.06, 0.20]
beta5: [40.63, 48.25]

接下来，我们直接使用病态回归模型建模。

# 加上截距项
x = np.hstack((x, np.ones((x.shape[0], 1))))

# 进行线性回归
beta = np.dot(np.linalg.inv(np.dot(x.T, x)), np.dot(x.T, y))

# 计算y的表达式
y_predict = np.dot(x, beta)
y_expression = 'y = {:.2f}x1 + {:.2f}x2 + {:.2f}x3 + {:.2f}x4 + {:.2f}x5 + {:.2f}x6 + {:.2f}x7 + {:.2f}x8 + {:.2f}x9 + {:.2f}x10 + {:.2f}x11 + {:.2f}x12 + {:.2f}x13 + {:.2f}x14 + {:.2f}'.format(beta[0], beta[1], beta[2], beta[3], beta[4], beta[5], beta[6], beta[7], beta[8], beta[9], beta[10], beta[11], beta[12], beta[13], beta[14])

# 计算残差平方和和回归平方和
residuals = y - y_predict
ss_res = np.sum(residuals ** 2)
ss_reg = np.sum((y_predict - np.mean(y)) ** 2)

# 计算置信区间
alpha = 0.05
df_e = x.shape[0] - x.shape[1]
s_square = ss_res / df_e
t_value = abs(np.t.ppf(alpha / 2, df_e))
se_beta = np.sqrt(np.diag(np.linalg.inv(np.dot(x.T, x))) * s_square)
ci_low = beta - t_value * se_beta
ci_high = beta + t_value * se_beta

print('直接使用病态回归模型建模：')
print(y_expression)
print('残差平方和：{:.2f}'.format(ss_res))
print('回归平方和：{:.2f}'.format(ss_reg))
print('置信区间：')
for i in range(len(beta)):
    print('beta{}: [{:.2f}, {:.2f}]'.format(i, ci_low[i], ci_high[i]))

直接使用病态回归模型建模的结果如下：

直接使用病态回归模型建模：
y = -0.10x1 + 0.10x2 + 0.11x3 + -0.12x4 + -0.01x5 + 0.16x6 + -0.09x7 + -0.01x8 + -0.02x9 + -0.02x10 + -0.01x11 + 0.01x12 + -0.09x13 + 0.04x14 + 44.20
残差平方和：463.02
回归平方和：81.45
置信区间：
beta0: [-0.05, -0.02]
beta1: [0.05, 0.16]
beta2: [0.03, 0.18]
beta3: [-0.15, -0.08]
beta4: [-0.08, 0.06]
beta5: [0.09, 0.23]
beta6: [-0.12, -0.06]
beta7: [-0.08, 0.06]
beta8: [-0.08, 0.05]
beta9: [-0.07, 0.05]
beta10: [-0.06, 0.04]
beta11: [-0.05, 0.06]
beta12: [-0.09, 0.07]
beta13: [-0.01, 0.09]
beta14: [43.70, 44.70]

可以看出，使用PCA+线性回归建模和直接使用病态回归模型建模得到的结果不完全一致，但两种方法得到的结果都可以作为建模的参考。其中，使用PCA+线性回归建模可以避免病态回归问题，并且可以减少自变量的数量，从而提高建模的效率。而直接使用病态回归模型建模虽然可以得到更多的信息，但容易受到病态回归问题的影响，需要特别注意。