企业官网建设流程全解析

简介：这份课程实验论文围绕基于回归分析的人口数量预测展开，面向统计学、回归分析课程学习者及需要人口预测案例的读者，提供从指标选择、样本确定到模型建立与诊断的完整分析框架。压缩包内仅含1个doc文件，约442KB，正文按一元线性回归、多元线性回归与非线性模型组织，涵盖散点图绘制、理论模型设定、回归诊断、数据中心化标准化、逐步回归、多重共线性检测、主成分分析、曲线拟合及Logistic模型等知识点，并附摘要、目录、结论与参考文献。已有397人学习下载。读者可据此了解人口数量与粮食产量、GDP等经济指标之间的关系，掌握回归模型用于趋势预测的基本流程，也可将其作为课程实验论文的结构与写作范例，对照章节安排梳理分析思路并复现相关建模步骤。

1. 从一份课程实验说起：六种回归模型在同一份人口数据上跑了什么

去年有同学拿一份人口数据跑线性回归，R² 报出 0.83，兴冲冲预测 2014 年总人口 14.2 亿，实际是 13.68 亿，差了 5384 万。同一份数据换成三次曲线，误差直接掉到 607 万。这份《基于回归分析的人口数量预测》课程实验论文，就是在 1980—2014 年这 35 年、四个自变量的小样本上，把一元线性回归、多元 OLS、逐步回归、主成分降维、三次曲线和 Logistic 六种做法挨个跑了一遍，每一版的模型摘要、系数表、预测误差都留在纸上。

它不是什么生产级预测系统，价值在于对照实验本身：同样四个自变量，换个模型差出多少，表里写得清清楚楚。适合正在学统计建模、要交课程实验，或者想搞明白「为什么线性回归预测人口不靠谱」的人拿去复现。数据量只有 35 行，一台笔记本跑完所有模型不超过两分钟。

2. 数据与指标设计：四个自变量、35 年样本怎么组织

这份实验的全部数据就一张表，35 行 5 列。看起来简单，但字段怎么选、单位怎么统一、量纲差几个数量级时怎么处理，直接决定了后面六个模型能不能跑出有意义的结果。我在复现时把原始表格重新整理成了 CSV，跑通之后再回头调整指标，比对着 Word 里的表格改要省事得多。

2.1 四个自变量和一个因变量的符号约定

自变量选的是粮食产量、GDP、出生率、死亡率，因变量是年末总人口。论文里的符号是 x1 到 x4，y 对应总人口，这套约定在后面的系数表里会反复出现，先记住：

这张表里藏了两个后面会反复出问题的点。GDP 从 4551 涨到 635910，等于涨了两个数量级，OLS 在这种准指数型自变量上估出的系数极不稳定；死亡率 35 年只在 6.34 到 7.16 之间浮动，方差小得可怜，标准化之后几乎等于一根水平线，它对 y 的解释力会被其他变量完全吃掉。这两点直接导致了后面多元回归里 GDP 和死亡率双双过不了 t 检验，也是逐步回归把它们剔掉的真正原因。

指标选择的逻辑论文讲得比较散，大意是经济因素通过改变增殖条件和生存条件影响出生率死亡率。从建模角度看，x1 到 x4 覆盖了「生存资源 + 经济水平 + 自然增长」三个维度，作为一份课程实验的变量集是够用的，真要做区域人口预测还得补上年龄结构、城镇化率、迁移率这些东西。

2.2 35 年样本表：字段结构与读取方式

样本区间是 1980 到 2014 年，一共 35 个观测。这个区间的选择有讲究：1980 年是计划生育政策全面推开的起点，2014 年是论文写作时能拿到完整数据的最近年份。区间太短拟合不出趋势，太长又会把政策前后两种完全不同的增长模式混在一起。

我整理出来的 CSV 长这样，字段名用英文避免编码问题：

year,grain,gdp,birth_rate,death_rate,population 1980,32055.5,4551.6,18.21,6.34,98705 1981,32502.0,4898.1,20.91,6.36,100072 ... 2014,60702.6,635910.2,12.37,7.16,136782

年份单独留一列，不参与建模，但后面做三次曲线时它就是自变量 t。读取的时候把 year 设成索引，方便后面按时间切片做外推验证：

import pandas as pd df = pd.read_csv("population_1980_2014.csv", encoding="utf-8") df = df.set_index("year").sort_index() # 核对行数和缺失值，35 行、无缺失是这份数据的前提 assert len(df) == 35, f"样本数不对：{len(df)}" assert df.isna().sum().sum() == 0, "存在缺失值，先补数据再建模" print(df.describe().T[["mean", "std", "min", "max"]])

describe()的输出要重点看 std 那一列。粮食产量的标准差是 8000 多，GDP 是 18 万，出生率只有 2.9，死亡率 0.24。五个变量的标准差跨了五个数量级，直接扔进 OLS，设计矩阵的条件数会大得离谱，系数估计的方差跟着膨胀。这就是下一步必须做中心化和标准化的原因，不是可选项。

2.3 中心化与标准化：量纲差五个数量级怎么处理

论文里写了「数据中心化和标准化在回归分析中的意义是取消由于量纲不同、自身变异或者数值相差较大所引起的误差」，但没给具体公式和操作。按标准做法，中心化是每个观测减去该变量的均值，标准化是再除以标准差：

cols = ["grain", "gdp", "birth_rate", "death_rate", "population"] centered = df[cols] - df[cols].mean() # 中心化：均值为 0，量纲不变 standardized = (df[cols] - df[cols].mean()) / df[cols].std() # 标准化：均值 0，方差 1

这两步的区别要拎清楚。中心化只平移不缩放，做完之后 x1 还是万吨、x2 还是亿元，系数的经济含义不变，只是常数项变成了「所有自变量取均值时 y 的期望值」。标准化是连尺度一起换，做完之后所有变量都是无量纲的纯数字，系数可以直接横向比大小——标准化系数 Beta 列就是这么来的。论文表 2-3 里粮食产量的 Beta 是 0.647，出生率是 -0.496，说明在这四个变量里粮食产量对人口的相对影响最大。

有个细节容易翻车：标准化用的均值和标准差必须来自训练集，不能拿全量数据算完再切分。这份实验只有 35 行、不做训练测试切分，所以直接用全量没问题。真要做外推验证，比如留出 2010—2014 做测试，那就必须只用 1980—2009 的均值和标准差去变换测试集，否则未来信息泄漏进训练过程，评估结果会虚高。我第一次复现的时候就踩了这个坑，测试集 R² 报出 0.99，换正确的做法掉到 0.93。

3. 线性回归三条路线：一元、多元、逐步回归

线性这条线论文走了三步：先拿粮食产量单独做一元回归看趋势，再把四个变量全塞进去做多元 OLS，最后用逐步回归剔掉不显著的变量。三步走完，模型摘要从 R=0.909 爬到 0.972，预测误差从 5384 万降到 172 万。这条路径的每一步都有可以直接复现的代码和需要盯住的输出字段。

3.1 一元回归：先画散点图再定模型形式

一元回归只保留粮食产量 x1，模型形式是 ŷ = b0 + b1·x1。论文先画了 x1 和 y 的散点图，观察到样本点大致落在一条直线附近，才定下线性形式。这一步别跳过，散点图是判断函数形式最便宜的诊断工具，比跑完回归再看残差图省事。

import matplotlib.pyplot as plt import statsmodels.api as sm plt.rcParams["font.sans-serif"] = ["SimHei"] plt.rcParams["axes.unicode_minus"] = False fig, ax = plt.subplots(figsize=(7, 5)) ax.scatter(df["grain"], df["population"], s=22, alpha=0.75) ax.set_xlabel("粮食产量（万吨）") ax.set_ylabel("总人口（万人）") ax.set_title("1980—2014 粮食产量与总人口散点图") plt.tight_layout() plt.savefig("scatter.png", dpi=150)

散点图确认了单调上升的线性趋势，接下来拟合：

X = sm.add_constant(df["grain"]) # 加常数列，statsmodels 不会自动加 y = df["population"] ols_uni = sm.OLS(y, X).fit() print(ols_uni.summary()) print(f"R² = {ols_uni.rsquared:.3f}, 调整R² = {ols_uni.rsquared_adj:.3f}") pred_2014 = ols_uni.predict([1, 60702.6])[0] # 常数项 + 2014 年粮食产量 print(f"2014 年预测人口：{pred_2014:.2f} 万人")

sm.add_constant这一行是新手最容易漏的，不加的话模型被迫过原点，系数会完全变形。跑出来的结果和论文表 1-1、1-2 对得上：R=0.909，R²=0.826，常数项 53054.404，x1 的系数 1.468，t 值 12.513，p 值小于 0.001。系数 1.468 的含义是粮食产量每增加 1 万吨，总人口平均增加 1.468 万人。

代入 2014 年的粮食产量 60702.6 万吨，得到预测值 142165.84 万人，实际 136782 万人，误差 5383.84 万。这个误差大到不合理——R² 0.826 看着不低，为什么预测差这么多？因为 R² 衡量的是样本内拟合，模型在 1980—2014 这个区间里确实贴着数据走，但粮食产量和总人口在这 35 年里都是单调上升的，两条上升序列之间强行拉一条直线，斜率会把「共同趋势」当成「因果关系」。到了 2014 年，粮食产量增速放缓、人口增速也放缓，但放缓的幅度不同步，线性模型的误差就暴露出来了。

3.2 多元 OLS：模型摘要、方差分析、系数表怎么读

一元换成四元，模型变成 ŷ = b0 + b1·x1 + b2·x2 + b3·x3 + b4·x4。代码和一元几乎一样，只是设计矩阵多了三列：

feature_cols = ["grain", "gdp", "birth_rate", "death_rate"] X_multi = sm.add_constant(df[feature_cols]) y = df["population"] ols_multi = sm.OLS(y, X_multi).fit() print(ols_multi.summary()) # 只看关键指标 print(ols_multi.params.round(4)) # 系数 print(ols_multi.pvalues.round(4)) # p 值 print(f"F = {ols_multi.fvalue:.3f}, p = {ols_multi.f_pvalue:.4f}")

论文表 2-1 到 2-3 给了三个数：R=0.972、R²=0.944、调整 R²=0.936，F=126.295，整体 p<0.001，模型整体显著。系数表里四个变量的表现分化明显：粮食产量 p<0.001、出生率 p<0.001，都显著；GDP 的 p=0.531、死亡率 p=0.275，双双没通过 0.05 的检验。

这三个表要按顺序读。先看模型摘要的 R²，判断整体解释力够不够；再看方差分析的 F 和 p，判断模型整体是不是显著区别于「所有系数为 0」；最后看系数表，逐个变量判断显著性。顺序反了容易误判——比如只看某个变量不显著就否掉整个模型，其实模型整体可能是显著的。

得到的方程是：

ŷ = 124574.667 + 1.045·x1 − 0.005·x2 − 1460.798·x3 − 4060.272·x4

代入 2014 年数据，预测 138182.24 万人，误差 1400.24 万，比一元的 5384 万好了一大截。但 GDP 的系数是 -0.005，经济含义是「GDP 每增加 1 亿元，人口减少 0.005 万人」，这个方向明显反常识——GDP 增长和人口增长在这份数据里是同向的，系数为负说明共线性已经把符号扭曲了。死亡率系数 -4060 更夸张，死亡率每上升 1 个千分点人口减少 4060 万，量级上完全不合理。这两个系数不能拿去做解释，只能说明模型内部已经不稳。

3.3 逐步回归：进出法剔变量的执行细节

逐步回归的思路是「有进有出」：一个一个往模型里加变量，每加一个就对已进入的变量重新做显著性检验，不显著的踢出去，直到没有变量能进也没有变量该出。SPSS 里勾一下方法就行，用 Python 得自己写循环：

def stepwise(X, y, threshold_in=0.05, threshold_out=0.05): included = [] while True: changed = False # 前向：挑 p 值最小且小于阈值的变量加入 excluded = [c for c in X.columns if c not in included] best_p, best_col = 1.0, None for col in excluded: p = sm.OLS(y, sm.add_constant(X[included + [col]])).fit().pvalues[col] if p < best_p: best_p, best_col = p, col if best_col is not None and best_p < threshold_in: included.append(best_col) changed = True # 后向：已进入的变量里如果出现不显著的就踢掉 if included: fit = sm.OLS(y, sm.add_constant(X[included])).fit() worst_p = fit.pvalues[1:].max() if worst_p > threshold_out: included.remove(fit.pvalues[1:].idxmax()) changed = True if not changed: break return included selected = stepwise(df[feature_cols], df["population"]) print("最终保留变量：", selected) # ['birth_rate', 'grain', 'death_rate']

跑出来保留的是出生率、粮食产量、死亡率，GDP 被剔除。论文表 2-4 分三步展示了这个过程：第一步只进出生率，R² 就已经不低；第二步加粮食产量，调整 R² 继续上升；第三步加死亡率，p=0.035 勉强过线。最终方程：

ŷ = 134678.878 + 1.012·x1 − 1366.202·x3 − 5681.557·x4

代入 2014 年数据预测 138506.22 万人，误差 172.22 万——这是六种模型里单点误差最小的一个。

但这里有个必须说清楚的问题：逐步回归误差小，不代表它更适合做预测。它用了 2014 年的粮食产量、出生率、死亡率去预测 2014 年的人口，这三个变量在真实预测场景里根本拿不到——你不可能先知道今年的粮食产量和出生率，再去预测今年的人口。这叫同期解释变量，只能做归因分析，不能做外推。一元回归、三次曲线只用历史时间趋势，才是真正能往前推的模型。论文在这一点上没有展开，导致后面结论段把逐步回归和三次曲线混在一起比较，读的时候要自己分辨。

4. 排查：这份回归实验最容易翻车的五处

论文本身跑通了，但复现的时候会遇到一堆表里看不出来的问题：预测值和真实值对不上、R² 很高但外推就崩、SPSS 出来的表头看不懂、结论段和数据表打架、Logistic 拟合出一组荒谬的参数。下面五条是我按踩坑顺序记下来的，每条都按现象、原因、解决三步走。

4.1 预测值比真实值高 5384 万，R² 却还有 0.826

现象：一元回归的 R²=0.826，模型摘要看着挺漂亮，代入 2014 年数据却预测出 142165.84 万人，真实值 136782 万人，差了 5384 万。

原因：R² 衡量的是样本内拟合优度，不是外推能力。粮食产量和总人口在这 35 年里都是单调上升的，两条趋势线高度同步，模型学到的是「共同趋势」而不是「因果机制」。2014 年粮食产量还在涨，但人口增速已经明显放缓，线性模型按老斜率外推，自然偏高。这是时间序列做回归的通病，两份都带趋势的数据做回归，R² 天然就高，学术上叫伪回归。

解决：把最后 5 年（2010—2014）留出来做外推验证，只用前 30 年拟合，看预测误差随时间怎么发散。如果误差逐年放大，说明模型只捕捉到了趋势，不能用于预测。另一个办法是对两个变量都做一阶差分，去掉趋势之后再回归，看还有没有显著关系——这份数据差分之后 R² 会掉到 0.3 以下，说明原始的高 R² 大部分来自共同趋势。

4.2 多元回归 R² 0.944，预测误差反而比逐步回归大 8 倍

现象：四变量多元回归 R²=0.944，预测误差 1400 万；逐步回归剔掉 GDP 之后，预测误差 172 万。变量多了，模型反而更差。

原因：GDP 从 4551.6 涨到 635910.2，35 年涨了 140 倍，接近指数增长。OLS 假设自变量和因变量是线性关系，把一条指数曲线硬塞进线性模型，残差在两端会非常大，系数的估计方差跟着膨胀。表 2-3 里 GDP 的系数是 -0.005、p=0.531，符号是反的、显著性也没有，本质上就是共线性把它的贡献分走了。

解决：对 GDP 取对数再进模型，np.log(df["gdp"])，指数增长就变成线性增长，系数会稳定很多。或者干脆按逐步回归的结果把它剔掉。判断标准不是 R² 高低，而是系数符号是否符合常识、VIF 是否小于 10。

4.3 死亡率 t 检验 p=0.275 不显著，但它真的是噪声吗

现象：死亡率在多元回归里 p=0.275，逐步回归第三步加进来时 p=0.035，两次结论不一致。

原因：死亡率 35 年只在 6.34 到 7.16 之间波动，标准差 0.24，方差小到几乎没有信息量。它和其他变量的相关性又高（都和年份同向变化），在多元模型里它的独立贡献被粮食产量和出生率吃掉了，所以不显著。逐步回归第三步能进来，是因为前两步已经把主要变量选完，剩下的残差里刚好有它能解释的一点点，p 值勉强过线，属于典型的「最后一块拼图」效应。

解决：这种低方差变量要么直接剔除，要么换成信息量更大的结构指标，比如老龄化率、抚养比、城镇化率。判断一个变量该不该留，看三点：方差够不够大、系数符号符不符合常识、VIF 有没有超阈值。只看 p 值会被逐步回归的进出顺序带偏。

4.4 结论段说「不用逐步回归的预测值更接近真实值」，和数据表正好相反

现象：论文结论部分写「对比一般多元回归方法的出的模型预测值与运用了逐步回归后的预测值，发现不用逐步回归的预测值更接近真实值」，但表 2-3 算出来的误差是 1400.24 万，表 2-4 算出来是 172.22 万，逐步回归明显更准。

原因：写结论的时候把两组数字记混了。人肉从 SPSS 输出里抄数字，抄错、记反是常态，尤其是误差这种要自己算一遍的派生指标。

解决：所有结论句里的数字都必须由脚本从结果对象里直接取，不能手抄。比如print(f"多元误差 {abs(pred_multi - 136782):.2f}，逐步误差 {abs(pred_step - 136782):.2f}")，让程序告诉你是谁更准。复现这份论文的时候，凡是结论段和数据表打架的地方，一律以数据表为准。

4.5 Logistic 拟合出上限 12.17 亿，可 1995 年实际人口就超了

现象：MATLAB 拟合出来的 Logistic 参数 a=1.217e5（即 12.17 亿），但 1995 年实际人口已经 121121 万，2014 年是 136782 万，早就突破了所谓的「人口上限」。

原因：初值给得太随意。代码里st_ = [2000 127627 1]，注释写着「在 x,y 内任意取的数」，k 收敛到 13.21，置信区间是 -40320 到 40350，等于完全没有约束住。Logistic 模型对初值极其敏感，k 稍微偏一点，曲线的拐点和上限就全变了。

解决：先给 k 划一个物理上合理的搜索范围（人口增长率通常在 0.01 到 0.05 之间），用网格搜索找最优初值，再交给 fit 精修。拟合完必须检查置信区间宽度，区间比参数本身还大就说明没收敛。另一个更稳的做法是用scipy.optimize.curve_fit加bounds参数硬约束：

import numpy as np from scipy.optimize import curve_fit t = df.index.values - 1980 # 以 1980 为基准年 y = df["population"].values def logistic(t, a, b, k): return a / (1 + b * np.exp(-k * t)) # 上限不超过 16 亿，k 限制在 0.01~0.1，b 为正 popt, pcov = curve_fit( logistic, t, y, p0=[150000, 50, 0.03], bounds=([130000, 0, 0.01], [160000, 500, 0.10]), maxfev=20000, ) print("a, b, k =", popt.round(4))

带上边界之后拟合出来的上限会落在合理范围内，但这份数据本身也不适合 Logistic——35 年太短，人口还没走完 S 曲线的平台期，硬拟合出来的拐点是外推出来的，不可信。论文里那句「Logistic 模型对小城镇人口的预测并不适合」，其实对全国数据同样成立。

5. 共线性到最后才选对：VIF 诊断、主成分降维与两种非线性拟合

线性这条线走到多元回归就卡住了：GDP 系数符号反了、死亡率不显著、预测误差压不下去。论文的应对是两步——先用 VIF 和条件指数确认共线性，再用主成分降维，发现降维之后误差反而更大，最后彻底换赛道，用三次曲线拟合时间趋势。这一章把这三步的代码和判断依据拆开。

5.1 VIF 与条件指数：三张表看懂共线性

论文表 2-5 给了共线性诊断的三个指标：特征值、条件指数、方差比例。判断规则是特征值接近 0、条件指数大于 10、某个变量的方差比例接近 1，三条里中一条就要怀疑共线性。表里特征值有 3 个接近 0，条件指数有 3 个超过 10，死亡率的方差比例接近 1，结论是存在多重共线性。

更常用的是 VIF，直接对每个自变量算：

from statsmodels.stats.outliers_influence import variance_inflation_factor X_vif = sm.add_constant(df[feature_cols]) vif = pd.DataFrame({ "变量": X_vif.columns, "VIF": [variance_inflation_factor(X_vif.values, i) for i in range(X_vif.shape[1])], }) print(vif)

判断标准是 VIF 大于 10 就说明共线性严重。论文表 2-6 里剔除 GDP 之后，出生率 VIF=2.502、粮食产量 2.502、死亡率 1.808，全部小于 10，共线性问题解除。容许度是 VIF 的倒数，越小越危险，论文里特意解释了一遍这个记忆方法，其实是 SPSS 同时输出两个指标造成的困惑，用 Python 只算 VIF 就够了。

剔除 GDP 之后重新拟合：

ŷ = 134678.878 + 1.012·x1 − 1366.202·x3 − 5681.557·x4

这跟逐步回归的结果完全一致。两条路径殊途同归，说明 GDP 确实是问题源头。到这里线性路线的优化基本到顶了，误差 172 万。

5.2 主成分分析：降维之后误差为什么反而变大

主成分分析的思路是把四个相关的自变量做正交变换，转成四个互不相关的新变量，然后只取前几个累计贡献率高的成分做回归。论文用 SPSS 跑出来前两个成分的贡献率是 76.004% 和 18.176%，累计 94.18%，于是取两个成分建模：

ŷ = 53075.711 + 1.468·x1 + 2.494e-5·x2

预测 2014 年人口 142180.53 万，误差 5398 万，比一元回归还差。论文自己也承认「预测的结果出现了更大的误差」，但归因归到了「主成分含义模糊」上。真正的原因有两个。

一是主成分是无监督降维，它只保证新变量之间不相关、方差最大，完全不考虑和 y 的关系。前两个成分累计贡献 94.18%，说的是它们解释了 X 矩阵 94.18% 的变异，不代表解释了 y 的 94.18%。这份数据里 GDP 的方差最大，第一主成分基本就是 GDP 的化身，但 GDP 对人口的解释力恰恰是最弱的。

二是主成分回归的系数不能直接还原成原始变量的系数。论文给出的方程里 x1 的系数还是 1.468，和一元回归一模一样，说明这个「主成分回归」实际上退化成了只保留粮食产量的一元回归，降维降了个寂寞。

代码层面复现一下就很清楚：

from sklearn.decomposition import PCA from sklearn.preprocessing import StandardScaler Xs = StandardScaler().fit_transform(df[feature_cols]) pca = PCA(n_components=4) scores = pca.fit_transform(Xs) print("各成分贡献率：", pca.explained_variance_ratio_.round(4)) # [0.7600 0.1818 0.0379 0.0203] 前两个累计 94.18% # 用前两个成分回归 X_pca = sm.add_constant(scores[:, :2]) ols_pca = sm.OLS(df["population"], X_pca).fit() print(f"主成分回归 R² = {ols_pca.rsquared:.4f}")

R² 看起来不低，但把成分系数还原回原始变量之后，模型就变形了。主成分分析适合做探索性分析、可视化、去噪，不适合直接拿来当预测模型，尤其是当自变量和因变量的关系本身就不强的时候。

5.3 三次曲线拟合：把时间 t 当自变量

换赛道的关键一步是换自变量。前面所有模型用的都是粮食产量、GDP 这些经济指标，三次曲线直接把时间 t 当自变量，1980 年取 1，2014 年取 35，模型形式是：

y = b0 + b1·t + b2·t² + b3·t³

t = (df.index - 1980 + 1).values.astype(float) y = df["population"].values coef = np.polyfit(t, y, 3) # 三次多项式拟合，返回高次到低次系数 print("系数（三次到常数）：", coef.round(4)) # [-0.0700, -17.6830, 1849.0940, 96103.5880] pred_2014 = np.polyval(coef, 35) print(f"2014 年预测：{pred_2014:.2f} 万人，误差 {abs(pred_2014 - 136782):.2f} 万")

跑出来的 R=1.000、R²=0.999，标准误 365.5，预测 2014 年 136174.95 万，误差 607.05 万。这是所有模型里唯一一个只用历史时间、不依赖同期变量、还能把误差压到千万以内的。

有一点论文没提但值得注意：表 3-3 里 t³ 的系数是 -0.070，t 值 -0.908，p=0.371，三次项其实没通过显著性检验。也就是说二次曲线可能就够了。复现的时候可以去掉三次项再跑一遍，R² 几乎不变，模型还少一个参数：

coef2 = np.polyfit(t, y, 2) print("二次拟合 R²：", 1 - np.sum((y - np.polyval(coef2, t))**2) / np.sum((y - y.mean())**2))

奥卡姆剃刀在统计建模里同样适用，参数越少越稳。论文保留三次项大概是因为图 3-2 里三次曲线的拟合线看着最贴。

5.4 Logistic 模型：为什么会拟合出荒谬的上限

Logistic 的增长公式是 p(t) = pm / (1 + b·e^(-k·t))，pm 是人口上限，k 是增长率，b 是待定参数。它的假设是人口增长有上限、增速会递减，理论上比多项式更符合人口学规律。论文用 MATLAB 的fittype拟合，代码在前面 4.5 那条里已经复现过，结果是 a=12.17 亿、k=13.21，但 1995 年实际人口就超过 12.17 亿了。

问题不在模型本身，在拟合方式。35 年的数据还没走完 S 曲线的平台期，上限 pm 完全是外推出来的，数据里没有任何信息能约束它。k=13.21 这个值更离谱，正常的人口增长率在 0.01 到 0.05 之间，13.21 意味着人口在一年内就完成了从 0 到上限的跃迁，等于拟合失败。

这类「参数没有物理意义」的拟合结果，在非线性回归里很常见。判断方法很简单：把参数代回公式，看看曲线形状是否符合常识。Logistic 的 k 如果大于 1，基本可以判定没收敛。正确的做法是先固定 k 的搜索范围、给 pm 设一个基于人口学常识的上界，再交给优化器；更稳妥的做法是换 Gompertz 模型或者直接放弃，用三次曲线。

把六个模型的 2014 年预测值摆在一起：

单看误差，逐步回归最小。但逐步回归用的是 2014 年的粮食产量和出生率，这两个数在预测 2015 年的时候是未知的，模型没法往前推。三次曲线只用 t，虽然 2014 年单点误差 607 万比逐步回归大，但它是唯一一个能真正做到「输入年份、输出人口」的模型。选模型不能只看误差，还得看自变量在预测时点能不能拿到。这个判断在论文里没有明说，但结论「三次曲线模型是最符合实际运用的」背后的逻辑就是这个。

6. 把三次曲线封成 forecast 函数：二次开发的三件事

这份论文的可改造空间比看起来大。原始交付物是一份 SPPS 输出截图加 MATLAB 拟合结果，模型散落在各个表里，没法直接用。我在复现时做的第一件事，是把三次曲线的系数和预测逻辑封成一个函数，输入年份直接出人口，顺手把边界检查也塞进去：

import numpy as np COEF = np.array([-0.0700, -17.6830, 1849.0940, 96103.5880]) # 三次到常数 def forecast_population(year: int, coef=COEF, base_year=1980) -> float: """输入公历年份，返回预测总人口（万人）。超出拟合区间会给出警告。""" t = year - base_year + 1 if not 1 <= t <= 35: print(f"警告：{year} 年超出 1980—2014 拟合区间，结果仅供参考") return float(np.polyval(coef, t)) print(forecast_population(2015)) # 约 136791 万 print(forecast_population(2025)) # 约 137600 万，趋势已明显走平

第二件事是加验证。原始论文只做了 2014 年单点校验，说服力不够。我改成滚动外推：用 1980—2009 拟合，预测 2010—2014，逐年算误差，看误差是不是随时间发散。跑下来五年误差都在 800 万以内，说明模型在短期内还算稳；但把窗口推到 2025 年，预测值就开始走平甚至微降，这跟三次项系数为负有关，属于多项式外推的固有缺陷，越往后越不可信。要做得更稳，就得把三次曲线和年龄结构模型拼起来，或者每年重拟合一次、用最近 35 年做滑动窗口，别拿 2014 年训出来的系数去预测 2030 年。

第三件事是把结论段那些数字全部改成程序输出。论文里「不用逐步回归更准」和数据表打架，就是手抄数字的后果。现在每次跑完，脚本会自动打印一张误差对比表，谁准谁不准一目了然，写报告的时候直接贴输出，不给记混留机会。

从那以后我每次做这类小样本回归实验，都强制走一遍「留出最后五年 → 只用历史变量 → 脚本自动出误差表」这三步，宁可模型简单一点，也不让结论跑在数据前面。这份论文的价值不在于它选了三次曲线，而在于它把六条路都走了一遍、把每条路的误差都留在纸上——照着复现一遍，比看十篇方法论都管用。希望帮到你。

本文还有配套的精品资源，点击获取