简介：《基于回归分析的人口数量预测》是一份面向统计学、应用统计及相关专业学生的课程实验论文文档，适合作为回归分析课程设计、实验报告撰写与人口预测建模练习的参考材料。压缩包内共1个doc文件，大小约442KB，内容以论文正文形式组织，完整呈现从指标选择、样本确定到模型建立与诊断的分析流程。文档围绕人口数量与粮食产量、经济指标之间的关系展开，涵盖一元线性回归的散点图绘制、理论模型设定与回归诊断，多元线性回归的数据中心化标准化、逐步回归、多重共线性检测及主成分分析，并延伸至曲线回归与Logistic模型等非线性方法。已有396人学习，可帮助读者理解回归建模的完整步骤、论文目录结构与结果表达方式，并对照完成课程实验或开展类似预测分析。
1. 人口预测的悖论：R² 最高的模型，往往最先崩
同一张 35 行的年度数据表，做同一件事——基于回归分析的人口数量预测，五个模型给出的 2014 年预测值从 136174 万到 142180 万，最离谱的那个差了 5400 万。更反直觉的是：把 R² 写到 1.000 的三次曲线并不是误差最小的那个里最可信的，而参数置信区间宽到正负四十亿的 Logistic 模型，软件照样吐出了一个像模像样的"预测值"。这说明一个模型能不能用，跟它小数点后留几位没关系，只跟三件事有关：自变量有没有真的携带时间信息、参数是不是可识别、误差有没有被数据本身的噪声盖过去。后面的章节就沿着这三条线，把一份课程实验里的五个模型拆开看。
2. 一元线性回归：53054 + 1.468x 这条直线是怎么算出来的，又错在哪
2.1 最小二乘的几何：残差平方和最小到底在最小化什么
把 35 个年份摊在平面上，横轴是粮食产量，纵轴是年末人口总量，你会看到一团向右上倾斜的点云。一元线性回归要做的事，是找一条直线，让所有点到这条直线的竖直距离的平方和最小。写成矩阵形式，就是解正规方程 $(X^TX)\beta = X^Ty$，其中 $X$ 的第一列全是 1（截距项），第二列是粮食产量。这个方程有唯一解的充要条件是 $X^TX$ 可逆，也就是自变量之间不能存在精确的线性关系。
这份实验里拟合出来的是 $\hat{y} = 53054.404 + 1.468x$。斜率 1.468 的字面意思是：粮食产量每增加 1 万吨，人口总量增加 1.468 万人。这个数字在业务上其实很荒谬——它等于在说粮食产量翻一倍，人口也差不多翻一倍。原因不难找：1980 到 2014 年，粮食产量从 32055 万吨涨到 60702 万吨，人口从 98705 万涨到 136782 万，两条序列都在单调上升，把它们放在一起做回归，本质上是在用一个"时间"的替身去解释另一个"时间"的替身。这就是所谓的伪回归味道：R² 很高，但斜率没有因果含义，只有统计含义。
理解这一点之后，后面所有模型的成败就都有了解释框架：谁离"时间"更近，谁的拟合就更好看；谁离"机制"更近，谁的系数才谈得上可解释。
2.2 系数表逐行读：标准误、t 值、R² 各自回答什么问题
实验结果里的模型摘要显示 $R = 0.909$、$R^2 = 0.826$、调整后 $R^2 = 0.821$、回归标准误 4985.997。系数表里，常数项 53054.404、标准误 5500.013、t 值 9.646；粮食产量系数 1.468、标准误 0.117、标准化系数 Beta 0.909、t 值 12.513。
这三组数字回答的是三个完全不同的问题，混着看就会翻车。R² 回答"模型解释了因变量多少比例的变异"，0.826 意味着人口总量的波动里有 82.6% 能被粮食产量的线性变化"跟着解释"掉。t 值回答"这个系数是不是显著不为零"，它等于系数除以标准误，12.513 对应 p 值远小于 0.05，说明斜率确实不是抽风抽出来的。而回归标准误 4986 回答的是"平均而言，真实值离回归线有多远"——单位是万人，也就是平均每个年份偏了近 5000 万人。
这三个指标里，最容易被误读的就是 R²。它衡量的是相对变异，不是绝对误差。当因变量本身量级很大（人口从 9.8 亿涨到 13.7 亿），哪怕残差有 5000 万，R² 依然能轻松上 0.8。所以看到"拟合程度很高"这种描述时，第一反应应该是去翻回归标准误，而不是点头。
另外，标准化系数 Beta 0.909 和相关系数 R 0.909 数值相同，这不是巧合。一元回归里，Beta 就等于自变量与因变量的皮尔逊相关系数，这也解释了为什么一元模型的 R² 恰好是相关系数的平方。
2.3 5383 万人的预测误差：为什么它没超出预测区间
拿 2014 年的粮食产量 60702.6 万吨代入方程，得到 $\hat{y} = 142165.8355$ 万人，而实际值是 136782 万人，误差 5383.8 万，相对误差 3.9%。很多人看到这里就下结论"模型不行"，但真正的诊断要做的是把点预测换成区间预测。
预测区间的公式是 $\hat{y}0 \pm t{\alpha/2,,n-2}\cdot s\sqrt{1 + \frac{1}{n} + \frac{(x_0-\bar{x})^2}{S_{xx}}}$。这里 $s=4986$，$n=35$，$x_0=60702.6$ 是所有样本里最大的粮食产量，$\bar{x}$ 大约在 46000 附近，粮食产量的样本标准差约 7500，于是 $S_{xx}\approx 34\times 7500^2 \approx 1.9\times10^9$，杠杆项 $\frac{(x_0-\bar{x})^2}{S_{xx}} \approx \frac{2.16\times10^8}{1.9\times10^9}\approx 0.114$。代进去，$\sqrt{1+0.0286+0.114}\approx 1.07$，半宽约 $2.04\times4986\times1.07 \approx 10900$ 万人。
也就是说，95% 的预测区间是 142166 ± 10900，粗略在 13.1 亿到 15.3 亿之间。真实值 136782 稳稳落在区间里面。这个结论很关键：模型并没有"算错"，它只是"没用"——区间宽到可以装下小半个亚洲，任何落在这个区间里的预测都谈不上有价值。而如果只算均值响应的置信区间（把根号里的 1 去掉，只留 $\frac{1}{n}+h$），半宽大约是 3800 万，明显窄一截。这两个区间的差别，正是"预测一个新个体"和"估计回归线在该点的高度"的差别，很多人做预测时贴的是置信区间，报出去的却是预测误差，这就是口径混淆。
用代码把这件事一次算清楚：
import numpy as np import pandas as pd import statsmodels.api as sm # 读入 1980—2014 年年度数据，列：year, grain(万吨), gdp(亿元), birth, death, pop(万人) df = pd.read_csv("pop_1980_2014.csv").sort_values("year").reset_index(drop=True) y = df["pop"] X = sm.add_constant(df[["grain"]]) # 必须补截距列，否则 OLS 会强制回归线过原点 fit1 = sm.OLS(y, X).fit() print(fit1.params) # const=53054.404, grain=1.468 print(fit1.bse) # 标准误 5500.013 / 0.117 print(fit1.rsquared) # 0.826 print(np.sqrt(fit1.mse_resid)) # 回归标准误 4985.997，单位万人 # 2014 年外推：同时给出均值置信区间和单点预测区间 X0 = sm.add_constant(pd.DataFrame({"grain": [60702.6]})) frame = fit1.get_prediction(X0).summary_frame(alpha=0.05) print(frame[["mean", "mean_ci_lower", "mean_ci_upper", "obs_ci_lower", "obs_ci_upper"]])sm.add_constant的作用是给设计矩阵补一列 1，对应截距项；漏掉它，模型会被迫经过原点，斜率会被严重扭曲。mse_resid是残差平方和除以自由度 $n-2$，开方即回归标准误。get_prediction().summary_frame()一次返回两组区间：mean_ci_*是回归线高度的置信区间，obs_ci_*是单个观测的预测区间，后者永远更宽，因为多了观测本身那层随机误差。做人口这类外推预测时，报obs_ci才算诚实。
3. 多元回归与多重共线性：GDP 的系数为什么会变成负数
3.1 中心化标准化改的不是 R²，是 X'X 的条件数
实验里把粮食产量、GDP、出生率、死亡率四个自变量一起放进模型之前，先做了一步中心化和标准化，理由是"取消量纲不同、数值相差较大所引起的误差"。这个说法方向对，但机制说得含糊。
真实的机制是：OLS 对自变量的线性变换是不变的。你把 GDP 从亿元换成万元，系数会跟着缩放，但 t 值、p 值、R²、残差全都一模一样。所以标准化绝对不会"提高拟合优度"。它真正改变的是 $X^TX$ 的条件数。原始数据里，GDP 在 4551 到 636463 之间（注意这里单位标的是亿元，但出生率那一列标着 %、数值却是 12.37 这种量级，显然是千分率，表头单位本身就是个坑），粮食产量在 3 万到 6 万之间，两者数量级差了近两个量级。$X^TX$ 的对角元因此相差好几个数量级，矩阵接近奇异，直接求解时浮点误差会被放大，系数可能在最后几位上抖。
标准化之后，每个变量的方差都归一，$X^TX$ 变成相关矩阵，对角元全是 1，条件数大幅下降。这跟模型好不好没关系，跟"解方程解得准不准"有关系。判断一个建模流程是否专业，看它是否在拟合前检查条件数、在拟合后检查 VIF，比看它用了多少个变量靠谱得多。
3.2 特征值、条件指数、VIF：三套指标测的是同一件事
四个自变量一起进模型，结果确实好看了：$R = 0.972$、$R^2 = 0.944$、调整后 $R^2 = 0.936$、回归标准误降到 2967.568，F 检验 126.295、p<0.001。但系数表露馅了：GDP 的系数是 -0.005，标准化 Beta -0.078，t 值 -0.635，p 值 0.531；死亡率系数 -4060.272，t 值 -1.112，p 值 0.275。两个变量都不显著，而且 GDP 的符号是负的——经济总量越大，人口越少，这在常识上说不通。
共线性诊断给出了三组互相印证的证据。第一组是特征值：四个维度里，有三个特征值接近 0（0.028、0.068、0.003），只剩一个 1.972 撑着；第二组是条件指数：1.000、8.447、6.557、33.798、7.506、36.436、86.828、186.63，其中有三个超过 10，一个接近 190；第三组是 VIF，剔除 GDP 后剩下的三个变量 VIF 分别是 2.502、2.502、4.146，全部低于 10 这条常用阈值。
这三组指标的数学来源其实是同一个。VIF 的定义是 $1/(1-R_j^2)$，把第 j 个自变量对其余自变量做回归，得到的 $R_j^2$ 越接近 1，VIF 就越大，容忍度（VIF 的倒数）就越接近 0。条件指数则是把标准化后的设计矩阵做奇异值分解，用最大奇异值除以第 j 个奇异值，本质上衡量的是 $X^TX$ 在某个方向上的"塌陷"程度。而特征值就是奇异值的平方。所以三者测的是同一件事：某个自变量能不能被其余自变量的线性组合几乎完全表示出来。
共线性带来的直接后果，是系数估计的方差膨胀。方差的膨胀因子就是 VIF，VIF=4 意味着系数方差被放大了 4 倍，标准误变成 2 倍。标准误一变大，t 值就变小，本来显著的变量会变得不显著；更糟的是，在强共线下，系数的点估计本身会剧烈震荡，符号翻转是常见现象。GDP 系数为负，就是这种震荡的产物，而不是什么"经济增长抑制人口"的经济学发现。
用代码把 VIF 算出来只需要几行：
from statsmodels.stats.outliers_influence import variance_inflation_factor Xm = sm.add_constant(df[["grain", "gdp", "birth", "death"]]) vif = pd.Series( [variance_inflation_factor(Xm.values, i) for i in range(Xm.shape[1])], index=Xm.columns, ) print(vif) # 看 grain / gdp / birth / death 四行，const 行没有意义variance_inflation_factor每次只算一列，它内部做的事就是把该列对其余列回归再取 $1/(1-R^2)$。要提醒的是，这个函数对常数项也照算不误，const 那一行的数值没有解释力，直接忽略。另外，VIF 必须在同一批变量上比较才有意义，删一个变量，其余所有变量的 VIF 都会变。
3.3 逐步回归的"有进有出"：剔除 GDP 后误差缩到 1724 万
既然 GDP 和死亡率不显著，实验里用了逐步回归来筛变量。逐步回归的机制是"有进有出"：一个一个往模型里加变量，每加一个就对已在模型里的变量重新做显著性检验，不显著的踢出去，反复迭代直到没有变量能进、也没有变量该出。
三步之后稳定下来的模型只留了三个变量：常数 134678.878、粮食产量 1.012、出生率 -1366.202、死亡率 -5681.557，三个 t 值分别是 9.412、7.184、-6.396，全部显著；VIF 分别是 2.868、4.146、1.808，全部低于 10。方程写成：
$$\hat{y} = 134678.878 + 1.012x_{\text{粮食}} - 1366.202x_{\text{出生率}} - 5681.557x_{\text{死亡率}}$$
代入 2014 年数据得到预测值 138506.21662 万人。这里有个必须点出来的细节：原稿写"与真实值相差 172.21662"，但 138506.21662 减去 136782 等于 1724.21662，是 1724 万而不是 172 万。这个小数点错位在后文的模型对比里会直接把结论带偏，因为它会让读者以为逐步回归的精度是三次曲线（误差 607 万）的三分之一，而实际上逐步回归的样本外误差是 1724 万，比全模型（1400 万）还差。
这个反常结果恰恰是最值得琢磨的地方。逐步回归在样本内更"干净"：变量全部显著、共线性全部消除、方程更简洁。但样本外误差反而变大了。原因是逐步回归优化的目标是样本内的信息准则（F 检验、AIC 之类），剔除变量降低的是方差、抬高的是偏差，在这份数据上偏差抬高的幅度超过了方差下降的收益。换句话说，逐步回归买到的是"解释得清楚"，不是"猜得准"。把这两件事混为一谈，是回归建模里最常见的后悔药。
4. 主成分回归：94.18% 的方差贡献率买到了什么
4.1 特征值就是方差：主成分如何把 4 个共线变量摊成 4 根正交轴
面对共线性，除了删变量还有一条路：把自变量本身换掉。主成分分析做的事，是对标准化后的设计矩阵做特征分解，得到四根互相正交的新坐标轴，每根轴上承载的方差大小就是对应的特征值。
实验里的结果是：第一主成分特征值 3.040，解释 76.004% 的方差；第二主成分 0.727，解释 18.176%；前两个累计 94.18%；第三 0.152（3.793%），第四 0.081（2.028%）。四个特征值加起来正好是 4——因为标准化后每个变量方差为 1，总方差等于变量个数。
这组数字透露了两个信息。好消息是前两个主成分就吃掉了 94.18% 的信息，降到二维几乎不损失什么；坏消息是，第一主成分独占 76%，意味着四个自变量之间存在一根非常强的公共轴，几乎所有变量都在这根轴上同向变化。这根轴是什么？从业务上看，它几乎就是"时间"——粮食产量、GDP、人口都在随年份单调上升，出生率和死亡率在下降，它们被同一股力量推着走。主成分分析没有消除共线性，它只是把共线性重新参数化了：原来四个变量彼此相关，现在变成第一主成分承载公共趋势、第二主成分承载剩余对比。
代码层面这一步很直白：
from sklearn.preprocessing import StandardScaler from sklearn.decomposition import PCA Xz = StandardScaler().fit_transform(df[["grain", "gdp", "birth", "death"]]) pca = PCA().fit(Xz) print(pca.explained_variance_) # 特征值 3.040 / 0.727 / 0.152 / 0.081 print(pca.explained_variance_ratio_) # 76.0% / 18.2% / 3.8% / 2.0% print(pca.components_[:2]) # 前两个主成分的载荷向量 Z = pca.transform(Xz)[:, :2] # 只保留前两维得分必须先StandardScaler再做 PCA，顺序反了第一主成分会被量级最大的那一列垄断，这里就会变成 GDP 的独角戏。explained_variance_是协方差矩阵的特征值，components_的每一行是一根主轴的载荷向量（也就是特征向量），transform得到的是样本在这些新轴上的坐标。降维就是取前 k 列。
4.2 得分回归再还原：为什么系数又变回 1.468
主成分回归的完整流程有两步变换：先把原始自变量投影成主成分得分，用得分对人口做回归，再把系数还原回原始变量的尺度上。实验里最终给出的还原结果是：常数 53075.711、粮食产量系数 1.468、GDP 系数 $2.494\times10^{-5}$。
看到 1.468 这个数，应该立刻警觉——它和一元回归里那个 1.468 几乎一模一样。这不是巧合。第一主成分承载了 76% 的方差且几乎与粮食产量同向，第二主成分承载 18% 的方差，主要刻画的是 GDP 相对于其他变量的偏离。用前两个主成分回归时，人口对第一主成分的依赖极强，对第二主成分几乎没反应，还原回去之后，GDP 的系数就被压到了 $10^{-5}$ 这个量级——等于宣告"控制住粮食产量之后，GDP 对人口没有额外解释力"。
这个结论本身是合理的，和逐步回归剔除 GDP 的判断一致。但代价也很明显：预测值变成了 142180.52727 万人，误差 5398.5 万，比一元回归的 5383.8 万还差一点。降维降掉的那 5.82% 方差里，恰好包含了模型真正需要的那点信号。
原稿把误差变大归因于"主成分含义模糊、解释性差"，这个说法只对了一半。解释性差确实是 PCA 的固有代价，但预测变差的直接原因是：主成分是按方差大小排序的，不是按与因变量的相关性排序的。排在前面的主成分可能方差很大，对 y 却没什么用；排在后面的主成分方差小，反而可能携带关键的判别信息。PCA 是无监督方法，它压根不知道 y 的存在。
4.3 降维的代价：无监督排序不等于对 y 有用
把 PCA 和逐步回归放在一起对比，能看清两条降维路线的本质区别。逐步回归是有监督的，它每一步都盯着 y，靠显著性检验决定去留；主成分分析是无监督的，它只看 X 自身的方差结构，不关心 y。前者容易过拟合样本、后者容易丢掉有用信息，两种风险方向相反。
这里还有一个更隐蔽的问题：主成分回归的系数还原公式是 $\beta = V_k \cdot \gamma$，其中 $V_k$ 是保留的 k 个载荷向量组成的矩阵，$\gamma$ 是得分回归的系数。还原出来的 $\beta$ 是有偏估计——它被约束在 $V_k$ 张成的 k 维子空间里，永远不可能等于真实的最优系数。偏差大小取决于被丢弃的那些主成分方向上，真实的 $\beta$ 分量有多大。如果真实系数恰好落在被丢弃的方向上，误差会非常大。
所以主成分回归适合什么场景？适合自变量高度共线、且你确信主要变异方向就是主要信号方向的场合，比如光谱数据、图像特征。人口预测这种"所有指标都被时间推着走"的场景里，第一主成分就是时间本身，降维等于把问题退化成趋势外推，那还不如直接对时间做回归来得干脆。这也就顺理成章地把讨论推到最后一个模型。
5. 排查手册：五个让回归预测悄悄失真的陷阱
5.1 五条踩坑记录：现象、原因、解决
第一个坑：系数符号和业务常识反着来。现象是 GDP 的系数出现 -0.005 且 p 值 0.531。原因不是经济学规律，而是共线性导致系数方差膨胀、点估计剧烈震荡，符号翻转是标准症状；同时 GDP 与粮食产量、出生率高度同步，控制住后者的前提下 GDP 的偏效应已经接近于零，剩下的只是噪声。解决办法是先算 VIF 和条件指数，超过阈值就做变量筛选、岭回归或者主成分回归，不要硬解释符号。
第二个坑：样本内 R² 漂亮，样本外一塌糊涂。现象是三次曲线 R² 达到 1.000、回归标准误只有 365，看起来无可挑剔。原因在于样本内拟合优度衡量的是"模型对已见数据的贴合程度"，参数越多越容易贴合；而预测关心的是"对未见数据的泛化能力"，两者不是一回事。解决办法是做滚动起点回测：用 1980—2005 拟合、预测 2006—2010，再把窗口往后滑，看样本外误差的均值和波动，而不是看全样本 R²。
第三个坑：单个系数不显著，却当成通过了检验。现象是三次曲线里 $t^3$ 项的系数 -0.070、标准误 0.077、t 值 -0.908、p 值 0.371，明显不显著，但结论里写的是"每个不同次幂的 t 都经过了显著性检验"。原因是把整体 F 检验的显著性当成了每个系数的显著性，而多项式各次幂之间天然高度相关（$t$、$t^2$、$t^3$ 在 1 到 35 的区间上相关系数极高），单项 t 检验本来就容易失灵。解决办法是对多项式先做正交化（用正交多项式基），或者干脆用中心化后的 $(t-\bar{t})$ 各次幂，再检查每一项的边际贡献。
第四个坑：非线性最小二乘的参数根本不可识别。现象是 Logistic 拟合给出 $a = 1.217\times10^5$、$b = 1.176\times10^5$、$k = 13.21$，而 b 和 k 的 95% 置信区间分别是 $(-4.741\times10^9,\ 4.742\times10^9)$ 和 $(-4.032\times10^4,\ 4.035\times10^4)$，宽度跨越九个数量级。原因是增长曲线还处在上升段，拐点没出现，上渐近线 $a$ 和增长率 $k$ 之间存在强 trade-off：把 $a$ 调大、$k$ 调小，前半段曲线几乎不变。数据里没有能钉住 $a$ 的信息，优化器只能在一条长长的谷底里乱走。更致命的是，拟出来的上渐近线是 121700 万人，而样本后半段人口已经涨到 136782 万，等于天花板在样本内就被顶穿了。解决办法是把 $a$ 固定成外部给定的上限（比如按生育率情景推算），只估 $b$ 和 $k$；或者换成分段增长模型、灰色模型，别硬拟合一个 S 形。
第五个坑：数据的底噪比模型之间的差异还大。现象是五个模型的误差从 607 万到 5399 万，看起来差异巨大，但普查十年一次，年度人口数据多数是靠抽样调查和推算得到的，本身就有百万量级的误差。原因在于把"模型差异"和"数据噪声"混在一起比较。解决办法是先估计数据本身的误差水平（用相邻年份差分、普查年与推算年的差值来估），如果模型差异没有明显超过这个底噪，就不要宣称"某模型显著更优"。
5.2 上线前的三个硬检查：残差、滚动回测、参数扰动
第一项是残差图。把残差对拟合值、对时间各画一张，看有没有喇叭口（异方差）、有没有 U 形（漏了二次项）、有没有连续同号（自相关）。人口这种时序数据，残差几乎必然正自相关，Durbin-Watson 统计量通常会明显小于 2，这意味着 OLS 的标准误被低估、t 值被高估，显著性结论要打折看。
第二项是滚动回测。固定一个起点，每次往后扩一年做拟合、预测下一年，累计几十次之后看误差的分布而不是单点误差。这一步能直接暴露"模型在样本尾部是否失稳"，而尾部恰恰是外推时最依赖的区域。
第三项是参数扰动。对数据做 bootstrap 重采样，把模型重估几百次，看系数和预测值的分布宽度。像 Logistic 那种置信区间宽到九个数量的情况，扰动一次就会暴露，比看拟合摘要快得多。
# 滚动起点回测：每次用截至当年的数据拟合，预测下一年 errors = [] for end in range(2000, 2015): tr = df[df.year <= end] te = df[df.year == end + 1] m = sm.OLS(tr["pop"], sm.add_constant(tr[["grain"]])).fit() pred = m.predict(sm.add_constant(te[["grain"]]))[0] errors.append(pred - te["pop"].values[0]) print(np.mean(np.abs(errors)), np.std(errors)) # 看平均 <p> <a href="https://download.csdn.net/download/2201_75683007/87851862" style="color:#ec7500;font-size:14px;"> 本文还有配套的精品资源，点击获取 </a> <img alt="menu-r.4af5f7ec.gif" src="https://csdnimg.cn/release/wenkucmsfe/public/img/menu-r.4af5f7ec.gif" style="width:16px;margin-left:4px;vertical-align:text-bottom;cursor:text;"> </p>
