灰色神经网络预测模型：PGM(1,1)与贝叶斯正则化实现TFP高精度预测

简介：一篇发表于《西南师范大学学报（自然科学版）》的学术论文，聚焦经济增长中全要素生产率（TFP）的预测问题，面向经济学研究者、量化建模爱好者及数据科学从业人员。文章将PGM(1,1)灰色模型与贝叶斯正则化神经网络相结合，构建组合预测模型，并以中国全要素生产率为实例，验证其在非线性、不确定性经济数据下的预测效果，为经济数据建模提供了新思路。压缩包内共1个文件，文件类型为PDF，大小约260KB，可直接阅读或打印。目前已有106人浏览学习。读者可从文中了解传统TFP测算方法的局限、贝叶斯正则化神经网络的核心原理及模型构建流程，适合用于论文选题、课题研究或量化建模参考。

1. 灰色神经网络预测模型：TFP预测不是算不准，而是单模型都输在了非线性上

全要素生产率预测的难点不在数据量，而在数据的非线性。2001—2010年中国全要素生产率数据介于0.58—0.66之间，带有明显的波动特征，索洛残差法、隐性变量法这类传统测算思路建立在线性假设上，拿去预测基本翻车。这篇《经济增长中全要素生产率的灰色神经网络预测模型》发表于《西南师范大学学报（自然科学版）》2016年第5期，给出了一条能复现的组合路线：先用PGM(1,1)模型对原始序列做一阶累加，压缩随机性，再用贝叶斯正则化神经网络拟合非线性残差，实测2009年和2010年预测相对误差为1.78%和0.37%，比单纯BP网络稳得多。适合做经济时序预测、灰色系统与深度学习交叉研究的读者，尤其适合时间序列预测模型的小样本场景。下面从选型逻辑、公式推导到复现参数逐层拆开。

2. 为什么是PGM(1,1)加贝叶斯正则化：单模型各自的坑与互补逻辑

2.1 传统TFP测算与预测方法的三个死穴

论文开篇把以往TFP研究的问题说得很透：过去大量成果集中在测算，而不是预测。索洛残差法把“余值”当作技术进步的贡献，但这个余值同时包含了资源配置效率、政策扰动、统计误差等所有说不清的因素，过于宽泛；隐性变量法把TFP视为从C-D生产函数残差中分离出的独立状态变量，仍然被新古典经济学的理论基础和规模报酬不变假设框住；前沿生产函数法依赖产出缺口的估算，而产出缺口本身的估算就有偏差。这三种方法在理论上各有定位，但共同的问题是：都假定经济系统存在稳定的线性结构。

回到实际数据，TFP受知识发展、规模经济、资源配置、政策变化共同影响，这几个因素之间的交互是非线性的，而且随着时间动态变化。用线性框架去外推，结果是误差随预测步长迅速放大。更现实的问题是，每一次测算都需要做大幅度的经济数据调研统计，人力物力成本高。论文的思路是换个方向：不追求把TFP拆解得多精细，而是建立一个与TFP数据特点吻合的预测模型，通过历史观测值直接推断未来趋势。这样既节约测算成本，又能为政策制定提供参考。这个转向在2016年的研究里不算主流，放到今天看就是典型的“用数据建模代替机理建模”，和现在流行的prophet时序预测模型、xgboost回归预测模型做时序外推时的思路本质相同。

2.2 PGM(1,1)与GM(1,1)的差别：背景值精确化是关键

灰色模型适合小样本、贫信息场景，TFP数据恰好符合这个特征：年度数据，样本量通常只有十年左右。GM(1,1)是最常见的灰色预测模型，但它的一个固有问题是背景值构造存在误差。传统GM(1,1)用一阶累加后相邻两点的平均值作为背景值，这个近似在序列变化平缓时够用，但当序列在某一段出现拐点或增速突变时，梯形面积近似会偏离实际积分值，导致发展系数a和灰作用量b的估计有偏。

PGM(1,1)的改进点就是对模型的权值p做精确求解。论文引用了李树峰、陈志丹关于pGM(1,1)模型权值p精确求解的结果，通过重新构造数据阵求出新的灰参数列，消除传统模型固有偏差。通俗地说，GM(1,1)的背景值公式固定取相邻两点算术平均，PGM则把这个平均的权重当作未知参数来估计，让背景值更贴合实际数据形态。这个差异在小样本下尤其明显，因为样本少，任何系统性偏差都无法通过数据量摊平。

我一般会额外做一个动作：在进入PGM之前先对原始序列做一次平稳性检查。如果数据有明显趋势突变，比如政策调整导致TFP突然跳升，一阶累加后仍然可能出现异常斜率段，这时候需要考虑对原始数据进行分段建模或先做对数化处理。论文里的中国TFP数据走势相对平稳，直接建模没有问题，但换到行业级数据时这一条要格外留意。

2.3 传统BP网络的过拟合困境与贝叶斯正则化的自动剪枝

传统BP神经网络因为万能逼近定理，理论上可以拟合任意非线性函数，这是它在洪水预测、焊缝外观预测、城市生活垃圾预测等领域被反复验证过的能力。但BP的一个臭名昭著的问题是过拟合：当网络节点数多于样本信息承载量时，训练集误差可以降到很低，测试集误差却大得离谱。TFP预测场景里，训练样本可能只有六七个，而一个三层网络稍微设宽一点就有几十个权重参数，参数数量远超样本数，过拟合几乎是必然的。

贝叶斯正则化神经网络的思路是给目标函数加一个惩罚项。把训练的目标从“最小化均方误差”改为“最小化均方误差与网络复杂度的加权和”，并且让这个权重由数据自己决定。Mackay在1992年的工作把这个过程放在贝叶斯框架下，将网络连接权值看作随机向量，假设训练集和连接权值的先验概率服从高斯分布，然后通过后验概率最大化求解最优的正则化系数。训练过程中网络能算出一个关键量——有效参数个数γ，它反映网络真正在使用的参数数量。如果网络结构庞大，但数据只需要少量参数就能解释，γ会自动收缩；如果数据确实复杂，γ会保持在一个较高的水平。

这就是3-25-1这个看似夸张的结构（25个隐含节点）在样本量很小的情况下依然能工作的原因。放到实际对比来看，用同样数据跑随机森林回归预测模型或标准BP，前者在非线性外推上缺乏有效外插能力，后者在无正则化时训练误差可以降到接近零，预测结果却完全漂移。贝叶斯正则化把“网络规模”和“样本规模”的匹配问题变成了一个自动优化问题，而不是人工试错问题。这一点对任何做小样本机器学习风险预测模型的人来说都是最重要的启示：样本量不够时，与其拼命降节点数，不如引入机制让网络自己控制有效复杂度。

3. PGM(1,1)建模全过程：从一阶累加到灰参数求解的完整计算链

3.1 一阶累加与背景值构造：把随机波动压成平滑曲线

PGM(1,1)的第一步是对原始序列做一阶累加生成。设原始TFP观测序列为X(0) = (x(0)(1), x(0)(2), …, x(0)(n))，一阶累加后得到X(1)，其中x(1)(k) = Σx(0)(m)，m从1到k。这一步的作用本质上是积分操作：原始年度数据里夹杂的政策冲击、统计误差等随机扰动，在累加过程中互相抵消，序列变成单调递增形式。累加后的序列规律性显著增强，这正是灰色系统理论的核心假设——系统行为数据中蕴含的信息可以通过累加生成来显性化。

第二步是计算背景值。PGM(1,1)的背景值采用相邻累加值的算术平均：

z(1)(k) = (x(1)(k) + x(1)(k+1)) / 2，k = 1, 2, …, n−1。

这个z序列替代原始的累加序列进入参数估计，目的是让微分方程的白化模型与离散数据之间有更合理的对应关系。

提示：如果直接拿原始序列做GM(1,1)建模而不做累加，模型会退化成指数平滑的某种变体，完全失去灰色模型处理贫信息系统的优势。累加是灰建模的根基，不能跳过。

论文的实例中，以2001年中国TFP数据为初始值，用7维序列长度建立PGM(1,1)模型。n=7这个取值不是拍脑袋定的：样本区间是2001—2010年共10个数据，按滑动预测的设定需要保留最后两年做验证，可用训练数据是前8年，再减去初始值占用的位置，7维是信息量与数据量的折中点。如果数据更长，建议做5维到10维的敏感性扫描，观察参数a的稳定性。

3.2 灰参数求解的四步矩阵运算：最小二乘估计的完整套路

灰参数列a = [a, b]^T的求解是建模的核心。计算流程如下表：

第四步的预测模型表达式为：

x̂(1)(k+1) = [x(0)(1) − b/a]·e^(−ak) + b/a，k = 1, 2, …。

得到累加序列的预测值后，再通过相邻累减还原出原始尺度的预测值。整个流程都是标准线性代数运算，可以在Excel里实现，也可以用任何编程语言几十行代码完成。

与传统GM(1,1)不同的是，PGM(1,1)在初步求解后会重新构造数据阵并求解新灰参数列，这一步的目的是修正背景值近似带来的系统偏差。具体做法是先把初值代入初步模型生成一轮拟合序列，用拟合序列重新计算背景值，再代入最小二乘求解一次，得到修正后的参数。本质上是一轮迭代精化，让参数估计不再依赖固定的算术平均权重。

3.3 参数敏感性：初值选取和数据长度对结果的影响

灰建模有两个容易被忽略的参数敏感点。第一个是初值x(0)(1)的选择，模型预测式直接包含x(0)(1)，而灰色模型对初值的响应不是线性的，初值偏差会通过指数项被放大。论文直接以2001年数据为初始值，这在数据噪声较小时没有问题；但如果初始年份恰好是经济异常波动年份，建议把初值也纳入优化，用后验数据反推最优初值。

第二个是建模长度n。n越长，最小二乘估计的样本量越大，参数更稳定，但代价是纳入更早的历史信息。如果经济结构发生过阶段性变化，早期信息反而会成为噪声。论文用7维序列在10年数据上做预测，是典型的保守选择。

还有一个工程细节：灰参数求解涉及矩阵求逆(B^T B)^(−1)，如果原始序列近乎常数，背景值序列之间的差异极小，B^T B会接近奇异，求逆结果数值不稳定。遇到这种情况，可以给B^T B加一个极小的对角扰动项（比如1e-8），或者用岭回归替代普通最小二乘。这个坑在论文里没有展开，但复现时很容易碰到，尤其是数据经过平滑处理后。

4. 贝叶斯正则化神经网络与组合流程：3-25-1结构背后的可复现细节

4.1 目标函数重构与有效参数γ：贝叶斯正则化在解决什么问题

贝叶斯正则化神经网络与标准BP的本质区别在目标函数。标准BP的训练目标是最小化均方误差，贝叶斯正则化把目标函数改写为：

F = β·ED + α·EW

其中ED是训练集均方误差，EW是网络权值平方和的惩罚项，α和β是正则化系数。当β远大于α时，传统无正则化训练占主导，网络倾向把训练误差压到极小，结果是过拟合；当α远大于β时，惩罚项占主导，网络会收缩权值规模，但整体训练误差会变大，可能出现欠拟合。训练的目标就是找到最优的α和β，让模型在拟合能力和复杂度之间取得平衡。

Mackay在贝叶斯框架下将网络权值视为随机向量，假设权值先验分布服从高斯分布，通过在整个权值空间上学习获取后验条件概率。后验概率最大化可以得到α和β的最优解，解的形式涉及网络有效参数个数γ：

γ = N − 2α·tr(H)^(−1)

其中N是网络总权值个数，H = β·∇²ED + α·∇²EW是目标函数的Hessian矩阵，tr(H)^(−1)是Hessian矩阵逆的迹。γ的物理含义非常直观：它表示网络真正被数据支撑的参数个数。一个1000参数的网络如果γ只有30，说明网络实际只用30个自由度在拟合，其余参数被正则化压扁了。

Hessian矩阵的精确计算在高维网络下计算量巨大，Foresee和Hagan在1997年提出用高斯牛顿法近似计算Hessian矩阵，这就是MATLAB神经网络工具箱中trainbr函数的实现基础。整个训练过程是自动的，α和β在每一轮迭代中根据当前网络状态重新估计，直到收敛。这意味着你不需要手动搜索正则化系数，这也是贝叶斯正则化相比L2正则化更省心的关键。

4.2 三步组合建模架构：从PGM拟合值到滑动窗口样本集

组合模型的三步流程如下：

步骤1：以TFP观测值序列中的x1为建模初始点，选择适当的序列维度n（论文用7），利用PGM(1,1)模型拟合观测值，得到拟合值序列。这一步的输出是质量更高的输入数据——随机性被压缩、规律性增强的拟合序列。

步骤2：确定预测周期数i（论文通过对比试验确定为3），以拟合值序列的前三个值为网络输入，以观测值序列的第四个值为网络输出，构造第一个训练样本；然后新陈代谢去掉最老的信息，补充新信息，以第二、三、四个拟合值为输入，第五个观测值为输出，构造第二个训练样本。依次滚动，得到完整训练样本集。

步骤3：应用训练好的最优网络，输入最近一期的拟合窗口值，预测第n+1期的TFP值。由于网络输出直接对应观测值，不需要再做累减还原，误差不会像灰色模型那样随累减步骤累积。

这里有一点需要特别说明：网络输入是PGM拟合值，输出是原始观测值。这个设计是组合模型的关键，它让神经网络学习的不是原始序列本身的形态，而是“灰色模型拟合值与真实值之间的残差规律”。灰色模型擅长捕捉趋势，神经网络擅长捕捉非线性映射，两者各管一段，互不干扰。如果直接用原始观测值进网络，相当于放弃了PGM的信息预处理能力，TFP数据的随机波动会让网络训练过程变成在噪声里捞信号。

4.3 网络结构、训练参数与MATLAB参考实现

论文最终确定的网络最优节点结构为3-25-1，训练误差设为0.01，预测周期数i=3。各关键参数的设定与理由整理如下表：

MATLAB参考实现的核心代码只有几行：

net = feedforwardnet(25); % 单隐含层，25个节点 net.trainFcn = 'trainbr'; % 贝叶斯正则化训练函数 net.trainParam.goal = 0.01; % 训练目标误差 [net, tr] = train(net, P_train, T_train); % P_train是PGM拟合值窗口 y_pred = sim(net, P_test); % 用最近窗口预测下一期

这里的feedforwardnet(25)创建单隐含层25节点网络，trainFcn设为trainbr就是启用贝叶斯正则化，trainParam.goal设为0.01对应论文的训练误差设定。P_train是滑动窗口构造的输入矩阵，每一列是一个3维窗口，T_train是对应的观测值。训练完成后用sim函数对新窗口做预测。

一个很容易翻车的细节是：trainbr在训练过程中会打印Mu（阻尼因子）和SSE（训练平方误差），SSE会停顿或反复，这是算法在更新α和β的正常表现，不需要按传统BP的习惯去手动调整学习率。如果你看到训练误差长时间停滞，检查的应该是窗口构造是否合理，而不是学习率。

5. 复现避坑清单：五个在实践中反复出现的错误与排查方法

5.1 训练误差降得很低，预测值却整体漂移

现象：网络在训练集上表现完美，误差远低于0.01，但预测出的2009年TFP值偏离真实值超过5%。

原因：这是典型的过拟合信号。训练误差过低本身就说明网络开始记忆训练样本的个体噪声，而不是学习规律。另一个常见原因是训练函数没设成trainbr，而是用了默认的trainlm（Levenberg-Marquardt），后者不具备自动正则化机制，在25个隐含节点的小样本场景下必然过拟合。

解决：把trainFcn强制设为trainbr，并把trainParam.goal保持在0.01左右，不要追求训练误差接近于零。训练完成后用tr结构体里的有效参数数量做诊断，25个隐含节点的网络如果γ小于10，说明网络自动剪枝工作正常。如果γ接近25甚至更大，说明正则化没有按预期生效，需要检查训练函数是否真的切换成功。

5.2 用原始TFP数据直接当网络输入，预测曲线锯齿状抖动

现象：跳过PGM步骤，直接把2001—2008年的原始TFP观测值作为窗口数据输入神经网络，预测结果在时间轴上出现剧烈的上下抖动。

原因：原始TFP数据本身包含随机扰动，神经网络在拟合扰动和拟合规律之间做了错误的取舍。PGM(1,1)在组合中的职责是压缩非线性、增强信息集成性，跳过它就等于放弃了灰色系统的预处理能力。

解决：严格执行组合流程，网络输入一律使用PGM(1,1)拟合值序列。建议对比一下两组实验的输入序列走势图：原始序列在2003—2005年有明显的波动，而PGM拟合值的曲线平滑得多，这就是神经网络真正需要的输入形态。如果两组曲线差异不大，说明数据本身的随机性不强，组合模型的优势就体现不出来——这本身也是一个有用的诊断结论。

5.3 训练阶段把真实观测值混进了输入窗口

现象：训练样本集的输入中有几列是直接用的真实TFP值而不是PGM拟合值，模型在训练集上表现尚可，在测试集上误差突然放大。

原因：数据泄漏。窗口输入必须严格来自PGM拟合值序列，输出是真实观测值。一旦把真实观测值放进输入，网络会学到“输入包含未来信息”这种假规律，虽然训练集上不会立刻暴露问题，但测试时输入窗口里没有未来信息，预测立刻失效。

解决：构造训练样本时区分两个序列。输入序列用拟合值，输出序列用真实观测值，在代码里用不同的变量名管理，样本构造完成后打印前两列的来源做抽查。更简单的一种做法是直接比较输入矩阵和输出矩阵的数值范围，如果两者高度重合，大概率发生了泄漏。这是所有小样本预测模型最容易踩的坑，不止这一篇论文。

5.4 隐含节点数随意改大，错误地期待更低的预测误差

现象：把隐含节点数从25改成50或100，训练误差确实下降了一些，但两个测试年份的预测误差反而上升了。

原因：节点数增加意味着网络参数数量增加，在贝叶斯正则化存在的情况下，网络会自动让大部分参数归零，但可搜索参数空间变大后，优化器找到全局最优的难度也增加了。更多节点不是免费的，它带来更高的训练不稳定性和局部最优风险。

解决：不要凭感觉扩大网络结构。论文中的25个节点是经过对比试验确定的，复现时先固定为25，再围绕18—35之间做扫描。每个结构至少训练3次取平均，因为trainbr每次初始权值不同，结果有随机波动。记录每次训练的有效参数个数γ，如果γ在某个节点数之后增长缓慢，说明继续加节点只是增加无效容量。

5.5 新数据预测时忘记保存归一化参数，结果整体偏移

现象：训练阶段预测误差正常，但换了一批新数据重新预测时，输出值整体比真实值低0.1左右。

原因：数据归一化参数不一致。MATLAB的mapminmax默认在每次调用时重新计算最小值和最大值，如果训练和测试分别调用了两次mapminmax，测试输入被映射到不同的区间，网络输出反归一化后自然偏移。

解决：用mapminmax('apply', X, ps)复用训练时得到的归一化结构体ps，反归一化用mapminmax('reverse', Y, ps)。参考代码如下：

[Pn_train, ps] = mapminmax(P_train, -1, 1); % 训练时保存归一化结构ps Pn_test = mapminmax('apply', P_test, ps); % 测试时复用ps，不重新计算 Y_pred = sim(net, Pn_test); Y_pred_inv = mapminmax('reverse', Y_pred, ps); % 反归一化回原始尺度

注意看代码里apply和reverse的用法：训练阶段输出的ps结构体记录了训练集的min和max，测试阶段apply计算输入映射用的是训练集的参数，反归一化同理。这是保证训练测试同分布的关键一步。所有回归类神经网络都存在这个问题，建议把“归一化参数保存到文件”作为每次训练的固定动作。

6. 让模型更可复现的进阶技巧：滚动验证与误差口径的统一

如果只按论文原流程跑一遍就收工，你大概率漏掉了对模型稳定性的验证。论文用2001—2008年数据预测2009和2010年两个值，属于单次留出验证，误差再小也只能代表一个切面。更稳的做法是滚动验证：先用2001—2007年训练、预测2008年；再用2001—2008年训练、预测2009年；最后用2001—2009年训练、预测2010年。每次预测都完整走一遍PGM(1,1)拟合、滑动窗口构造、网络训练三步，得到三个独立的预测误差。预测模型评估用滚动验证而不是单次留出，这是时间序列预测模型和普通回归预测模型在评估方式上最核心的差别。

评估指标建议统一用相对误差RE = |预测值−实际值| / 实际值 × 100%，论文给出的1.78%和0.37%就是这个口径。不要混用RMSE和MAE作为对比基准：小样本下RMSE对大误差敏感，一个离群预测会把整体指标拉得很差，而MAE没那么敏感，两个指标算出来的模型排名可能完全相反。跨模型比较（比如和随机森林回归预测模型、prophet时序预测模型对比）时，口径不统一会在汇报时引起无意义的争论。

还有一个值得做的二阶诊断：把PGM(1,1)的拟合残差画出来，观察残差序列是否有明显自相关。如果残差里还存在趋势片段，说明灰色模型部分没有把信息提取干净，可以尝试对残差再做一层灰色建模，或者把预测周期数i从3调到4，重新走一遍窗口构造。这类改进不需要动模型框架，改一个参数重跑的成本很低，但对最终误差的影响可能比换网络结构还大。

如果你打算直接跟着复现，下载论文PDF后重点核对三个地方：表2的PGM拟合值、表3的误差对比、以及第3节的组合模型步骤图——公式与代码之间的对应关系都在这些表里。从那以后，我每次做这类小样本时序预测，都强制走一遍固定检查清单：归一化参数是否复用、输入窗口是否全部来自拟合序列、trainbr是否生效、有效参数γ是否记录在案。这个五分钟能完成的检查，替我挡掉了好几次推倒重来的返工。希望帮到你。

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