灰色预测模型：原理、Python实现与销量预测实战)
做数据分析这几年，最头疼的往往不是算法不够高级，而是手头数据实在太少。刚接到一个客户项目时，历史销量可能只有几个月的样本，别说对照周期，连基本统计显著性都凑不出来。这时候硬上ARIMA、回归这类大样本模型，结果通常惨不忍睹。灰色预测就是在这种“小样本、贫信息”场景下出现的工具，而GM(1,1)又是其中最基础、最普及的模型。可一旦你真拿GM(1,1)去跑几组线性增长或近似线性的数据，就会发现它的模拟值和真实值总隔着一层，预测曲线也偏“平”。DNGM(1,1)这个改进模型的思路，就是绕开GM(1,1)对微分方程离散化时产生的那层误差，直接建立在离散域上带线性趋势项的递推结构，让模型在保留指数趋势的同时，也能兼顾线性增长。这篇文章就把DNGM(1,1)的原理、推导、手算案例、Python代码和避坑经验一次性讲清楚。无论你是在做销量预测、能源负荷预测，还是参加数学建模竞赛，只要适合小样本预测的场景，这套方法都值得放进工具箱。
1. 为什么需要 DNGM(1,1)：从小样本预测的痛点说起
1.1 灰色预测到底解决了什么问题
灰色系统理论的核心思想是“部分信息已知、部分信息未知”。我们拿到的数据往往只是系统输出的一条轨迹，内部结构、外部干扰都不清楚，属于典型的贫信息场景。统计模型依赖大量样本去估计分布，机器学习模型依赖大量样本去拟合映射，而灰色预测走的是另一条路：通过对原始序列做累加生成，把随机波动一步步“攒”起来，让原本毛毛躁躁的曲线变成一条平滑、有规律的递增曲线，然后再对这个规律建模。
这里最关键的一步就是累加生成。你手里的数据可能是这样：
120, 132, 141, 150, 158, 165, 171, 178, 184, 190
看着有点规律，但噪声也不小。做一次累加之后变成：
120, 252, 393, 543, 701, 866, 1037, 1215, 1399, 1589
这条曲线就非常光滑了，近似一条带微弧的直线。灰色模型就是在累加序列上建方程，得出结果后再通过累减还原回原始序列。这个过程本质上是在“用累加压制噪声”，而不是像传统统计那样“用分布描述噪声”。
灰色预测能火这么多年，原因很现实：它对样本量要求极低，通常4到10个点就能建模；计算量小，Excel都能手算；短期预测精度在同类小样本方法里相当能打。但传统GM(1,1)有一个绕不开的毛病——它对“近似齐次指数增长”的数据效果好，对线性增长、非指数趋势的数据就明显吃力。DNGM(1,1)就是为了解决这个毛病来的。
1.2 传统 GM(1,1) 的数学骨架和四大局限
先把GM(1,1)的建模流程梳理一遍。设原始非负序列为：
[ x^{(0)} = (x^{(0)}(1), x^{(0)}(2), \dots, x^{(0)}(n)) ]
一次累加生成（1-AGO）得到：
[ x^{(1)}(k)=\sum_{i=1}^{k}x^{(0)}(i) ]
然后定义背景值：
[ z^{(1)}(k)=0.5x^{(1)}(k)+0.5x^{(1)}(k-1) ]
GM(1,1)对应的白化微分方程是：
[ \frac{dx^{(1)}}{dt}+a x^{(1)}=b ]
参数 (a,b) 通过最小二乘估计得到，时间响应式为：
[ \hat{x}^{(1)}(k+1)=\left(x^{(0)}(1)-\frac{b}{a}\right)e^{-ak}+\frac{b}{a} ]
最后累减还原得到预测值。
这套流程看起来顺理成章，但实际用起来有四个明显的坑：
第一，微分方程在离散数据上天然存在近似误差。方程里那个导数本来是需要连续信息的，可我们只有离散点。GM(1,1)用梯形法构造背景值，本质上是用“前后两点均值”去近似曲线下的面积。数据变化越剧烈，这个近似误差越大。
第二，白化方程的解是指数形式的 (e^{-ak})。这意味着GM(1,1)默认数据服从指数增长规律。碰到近似线性的序列，模型会硬生生把直线掰成指数曲线来拟合，结果就是模拟值系统性偏离真实值，预测值往往偏大或偏小。
第三，初值被强行固定为第一个数据点 (x^{(0)}(1))。从最小二乘的角度看，第一个点不一定是最优的估计起点。尤其当第一个点本身是异常值时，整个预测曲线都会被带偏。
第四，对数据要求高，级比检验不通过时精度骤降。灰色模型虽然小样本能用，但不是任何小样本都能用。数据必须满足一定的级比条件，否则累加序列并不呈现准指数规律，模型硬建模就是自欺欺人。
这四个坑，正是DNGM(1,1)要逐一击破的对象。
2. DNGM(1,1) 的建模原理：一条差分方程如何改写预测逻辑
2.1 从“微分方程近似”到“离散递推”的思路转变
DNGM(1,1)模型全称通常写作 Discrete New Grey Model，也就是“新离散灰色模型”。它的核心改动，是把传统GM(1,1)那条连续微分方程，替换成一条直接在离散域上成立的差分递推方程：
[ x^{(1)}(k+1)=\beta_1 x^{(1)}(k)+\beta_2 k+\beta_3 ]
其中 (k=1,2,\dots,n-1)。
这条式子看起来平平无奇，但信息量很大。(\beta_1 x^{(1)}(k)) 是自回归项，负责捕捉序列自身的惯性增长；(\beta_2 k) 是线性趋势项，负责捕捉随时间匀速增加的部分；(\beta_3) 是常数项，相当于基线水平。相比GM(1,1)只有一个指数项，DNGM把“指数规律”和“线性规律”同时放进了模型里。
这里顺便说一句，不同文献里这类模型的叫法并不完全统一，有的叫NDGM(1,1)，有的叫DGM(1,1)的扩展形式，有的直接叫新离散灰色模型。命名虽然不同，但结构本质是一样的：在离散累积序列上加入线性时间项，再通过最小二乘估计参数。DNGM(1,1)就可以看作是这一类模型中的一种简洁表达。
为什么一定要从“微分方程离散化”改成“直接离散建模”？我个人的理解是，连续模型离散化之后，误差来源有两个：一是微分方程本身对真实系统的近似，二是离散化数值方法带来的近似。DNGM直接丢弃了第一层“微分方程近似”，直接从离散数据到差分方程，少了一层误差传递。这就好比你要测量一段弯曲山路，GM是先用微积分算出理论长度，再用梯形近似；DNGM则是直接沿着山路一步步走，每步都落在真实数据上。
2.2 参数估计：最小二乘与矩阵表达
DNGM(1,1)的参数估计并不复杂。把递推式展开，对 (k=1) 到 (n-1) 可以写出：
[ \begin{cases} x^{(1)}(2)=\beta_1 x^{(1)}(1)+\beta_2 \cdot 1+\beta_3 \ x^{(1)}(3)=\beta_1 x^{(1)}(2)+\beta_2 \cdot 2+\beta_3 \ \vdots \ x^{(1)}(n)=\beta_1 x^{(1)}(n-1)+\beta_2 \cdot (n-1)+\beta_3 \end{cases} ]
写成矩阵形式就是：
[ Y=\Phi \beta ]
其中：
[ Y=\begin{bmatrix} x^{(1)}(2) \ x^{(1)}(3) \ \vdots \ x^{(1)}(n) \end{bmatrix}, \quad \Phi=\begin{bmatrix} x^{(1)}(1) & 1 & 1 \ x^{(1)}(2) & 2 & 1 \ \vdots & \vdots & \vdots \ x^{(1)}(n-1) & n-1 & 1 \end{bmatrix}, \quad \beta=\begin{bmatrix} \beta_1 \ \beta_2 \ \beta_3 \end{bmatrix} ]
最小二乘解为：
[ \hat{\beta}=(\Phi^T\Phi)^{-1}\Phi^T Y ]
这里有个实操细节值得多说一句。实际计算时，我不建议直接写代码求矩阵逆，哪怕数据量小。因为当 (\Phi^T\Phi) 接近奇异时，直接求逆的数值误差会被放大好几倍。更好的做法是用svd或最小二乘求解器，比如Python里numpy.linalg.lstsq，或者MATLAB里的\运算。后面代码部分我会再演示。
2.3 为什么它能突破传统局限
回到第一节讲的那四个坑，我们逐条看DNGM的回应。
针对背景值近似误差：DNGM压根不构造背景值，也就不存在梯形近似的误差。模型直接从离散点之间的递推关系出发，每个方程都是用真实累加数据写出来的。
针对指数假设过于狭窄：DNGM在递推式里加入了 (\beta_2 k) 这一项。当数据线性增长时，这一项能直接捕获线性趋势，不再需要靠指数去“硬掰”。当数据是指数增长时，(\beta_1) 项会逼近一个合适的比例系数，模型同样能适应。也就是说，DNGM是“指数+线性”混合结构，覆盖面比GM宽得多。
针对初值固定问题：虽然递推的第一步仍然用 (x^{(1)}(1)) 作为启动值，但后续所有模拟值都由递推方程自己演化，并不要求第一个原始点必须精确落在拟合曲线上。这一点和GM的时间响应式不同，GM的解是理论上的指数曲线，第一个点天然被当作曲线上的点；DNGM的递推解是一条由数据“喂”出来的轨迹，更像让数据自己说话。
针对级比条件：DNGM的适用面更宽，对近似线性序列、带趋势项的序列都有不错的拟合效果。当然这并不意味着它可以违反灰色预测的基本前提，数据如果毫无规律，什么灰色模型都救不了。
3. 完整案例手算：10周销量数据的 DNGM(1,1) 建模过程
3.1 数据准备与级比检验
用一个实际业务场景来演示。某线上店铺最近10周销售额（单位：万元）如下：
| 周次 k | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 |
|---|---|---|---|---|---|---|---|---|---|---|
| 销售额 x(0)(k) | 120 | 132 | 141 | 150 | 158 | 165 | 171 | 178 | 184 | 190 |
先用级比检验确认数据适不适合做灰色预测。级比定义为：
[ \lambda(k)=\frac{x^{(0)}(k-1)}{x^{(0)}(k)}, \quad k=2,3,\dots,n ]
对于n=10的数据，常用的可容覆盖区间是：
[ (e^{-2/(n+1)},\ e^{2/(n+1)}) = (0.8338,\ 1.1994) ]
我们这批数据的级比大约为：
1.100, 1.068, 1.064, 1.053, 1.044, 1.036, 1.041, 1.034, 1.033
全部落在(0.8338, 1.1994)区间内，说明序列满足准指数规律，可以直接建模。
这里补充一个经验：如果某个级比稍微超出区间，不必急着推翻重来，可以先尝试平移处理，再重新检验。平移量 (c) 一般取 (|x_{\min}|+1)，把整个序列抬高到正数范围，灰色模型对平移后的序列依然有效，预测结果再平移回去即可。
3.2 累加生成与矩阵求解
对原始序列做一次累加：
| k | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 |
|---|---|---|---|---|---|---|---|---|---|---|
| x(1)(k) | 120 | 252 | 393 | 543 | 701 | 866 | 1037 | 1215 | 1399 | 1589 |
接下来构造矩阵。对 (k=1) 到 (n-1=9)，有：
[ Y=\begin{bmatrix} 252 \ 393 \ 543 \ 701 \ 866 \ 1037 \ 1215 \ 1399 \ 1589 \end{bmatrix} ]
[ \Phi=\begin{bmatrix} 120 & 1 & 1 \ 252 & 2 & 1 \ 393 & 3 & 1 \ 543 & 4 & 1 \ 701 & 5 & 1 \ 866 & 6 & 1 \ 1037 & 7 & 1 \ 1215 & 8 & 1 \ 1399 & 9 & 1 \end{bmatrix} ]
对 (\Phi) 和 (Y) 做最小二乘解，得到：
[ \beta_1 \approx 0.8968, \quad \beta_2 \approx 23.72, \quad \beta_3 \approx 119.43 ]
注意看这里的 (\beta_2) 相当大，说明这组数据里线性趋势占了主导，传统GM(1,1)只有指数项，自然会吃亏。(\beta_1) 略小于1，保证了累加序列的递推不会失控，整体结构是稳的。
递推式写出来就是：
[ \hat{x}^{(1)}(k+1)=0.8968,\hat{x}^{(1)}(k)+23.72,k+119.43 ]
以 (x^{(1)}(1)=120) 为启动值，依次递推得到累加序列的模拟值，再做累减还原，就得到原始序列的拟合值。
3.3 模拟结果与误差分析
用上面的递推式，得到的结果如下：
| 周次 k | 真实值 | DNGM拟合值 | 相对误差 |
|---|---|---|---|
| 1 | 120 | 120.00 | 0.00% |
| 2 | 132 | 130.77 | 0.93% |
| 3 | 141 | 140.89 | 0.08% |
| 4 | 150 | 150.20 | 0.13% |
| 5 | 158 | 158.39 | 0.25% |
| 6 | 165 | 165.60 | 0.36% |
| 7 | 171 | 172.35 | 0.79% |
| 8 | 178 | 178.24 | 0.14% |
| 9 | 184 | 183.65 | 0.19% |
| 10 | 190 | 188.32 | 0.88% |
平均相对误差大约0.38%，这个精度在短期预测里已经相当不错。我拿同样数据跑了一遍标准GM(1,1)，平均相对误差大概在1.5%到2%左右，虽然也不差，但明显比DNGM高出一个量级。特别是第6到第10期，GM的模拟误差往往越往后越大，而DNGM的误差始终在小范围内波动。
外推预测未来3周，也就是第11到第13周的销售额，递推会继续往k=10、11、12的方向推进，得到：
| 周次 | 预测值 |
|---|---|
| 11 | 约192.45 |
| 12 | 约196.40 |
| 13 | 约200.36 |
这里需要提醒一句，外推步数越多，不确定性越大。灰色模型的短期预测优势明显，长期预测只能做趋势参考，不建议直接拿来当最终决策依据。
4. Python 复现与代码实现
4.1 核心函数实现（可直接复制）
先给一份完整的Python实现。代码不搞花哨，核心就是构造矩阵、最小二乘求解、递推模拟外推三步。
import numpy as np def dngm_forecast(x0, steps=5): """ DNGM(1,1) 模型 参数: x0 : 原始数据序列，list或np.ndarray steps: 需要外推预测的步数 返回: x0_hat: 拟合值 + 预测值 beta : 模型参数 [beta1, beta2, beta3] """ x0 = np.asarray(x0, dtype=float).reshape(-1) n = len(x0) # 1. 一次累加生成 x1 = np.cumsum(x0) # 2. 构造矩阵 # 递推式: x1(k+1) = beta1 * x1(k) + beta2 * k + beta3, k=1,2,...,n-1 Y = x1[1:].reshape(-1, 1) k_index = np.arange(1, n).reshape(-1, 1) Phi = np.hstack([x1[:-1].reshape(-1, 1), k_index, np.ones((n - 1, 1))]) # 3. 最小二乘求解 beta, _, _, _ = np.linalg.lstsq(Phi, Y, rcond=None) beta = beta.flatten() # 4. 递推模拟与外推 m = n + steps x1_hat = np.zeros(m) x1_hat[0] = x0[0] for k in range(1, m): x1_hat[k] = beta[0] * x1_hat[k - 1] + beta[1] * k + beta[2] # 5. 累减还原 x0_hat = np.empty(m) x0_hat[0] = x0[0] x0_hat[1:] = np.diff(x1_hat) return x0_hat, beta这段代码的核心逻辑就一句话：先累加，用最小二乘拟合递推系数，再用递推式从已知推未知。如果你打算自己写一遍，我建议重点关注矩阵构造那一段，很容易把 (k) 的起始位置搞错。递推式里的 (k) 是从1开始，而不是从0开始，所以Phi矩阵的第二列是arange(1, n)，不是arange(n-1)。
4.2 调用示例与结果说明
仍然用10周销量数据来演示调用：
x0 = [120, 132, 141, 150, 158, 165, 171, 178, 184, 190] pred, beta = dngm_forecast(x0, steps=3) print('参数 beta:', beta) print('拟合+预测结果:') for i, v in enumerate(pred, 1): print(f'第{i}期: {v:.2f}')输出结果和第三节手算的一致，第1到第10期是拟合值，第11到第13期是外推预测值。实际项目里我通常会把数据、拟合值、真实值、残差拼成一个DataFrame输出，方便直接画图对比。
还有个小技巧：外推预测时，如果业务上每周都有新数据进来，建议滚动更新。比如这周拿到第11周真实值，就把第2到第11周的数据作为训练集重新建模，预测第12周，而不是继续沿用上一次的模型曲线。这种做法在预测领域叫“新陈代谢”，能让模型始终贴着最新趋势走，比一次性外推远好几步要稳得多。
5. 常见问题与避坑指南
5.1 数据预处理：灰色预测最容易翻车的地方
灰色预测不是万能药，它有几个硬前提。第一个是数据必须非负。如果序列里有0或负数，累加序列会出现负值，级比检验也很可能失效，模型输出会变得很奇怪。解决办法是平移，统一加一个常数 (c=\max(|x_{\min}|+1, 0))，建完模再平移回去。
第二个是数据不能太少。最少3个点才能估计出3个参数，但3个点的模型基本没有泛化能力。我的经验是5到15个点比较合适，少于5个点不如直接用移动平均或指数平滑，多于15个点时传统统计模型或机器学习模型逐渐变得可行，灰色预测的优势就不明显了。
第三个是级比检验不通过时的处理。常见手段包括取对数变换、开方变换、平移变换。平移是最温和的，不改变数据形态；取对数是把指数趋势掰成线性趋势，如果掰完之后级比通过了，说明数据本质上还是适合灰色建模的。
第四个是警惕数据里的突变点。灰色模型对异常值非常敏感，因为累加会把异常点的误差传递到后续所有值。建模前最好先画个折线图看看有没有明显离群点，有的话先做平滑或修正，否则预测结果会被带偏。
5.2 精度检验指标与等级判断
模型建完之后，不能只看拟合曲线顺不顺眼，要有一套客观的精度检验标准。常见的指标有两个：平均相对误差MAPE和后验差比C。
平均相对误差的计算很简单：
[ MAPE = \frac{1}{n}\sum_{k=1}^{n}\left|\frac{x^{(0)}(k)-\hat{x}^{(0)}(k)}{x^{(0)}(k)}\right|\times 100% ]
后验差比C的计算稍微复杂一点。设残差序列为 (\varepsilon(k)=x^{(0)}(k)-\hat{x}^{(0)}(k))，那么：
[ C=\frac{S_2}{S_1} ]
其中 (S_1) 是原始序列的标准差，(S_2) 是残差序列的标准差。C越小说明残差波动相对原始波动越小，模型越稳定。实际判断标准参考下面的表：
| 精度等级 | MAPE | 后验差比C |
|---|---|---|
| 优 | <1% | <0.35 |
| 合格 | 1%~5% | 0.35~0.50 |
| 勉强 | 5%~10% | 0.50~0.65 |
| 不合格 | >10% | >0.65 |
注意，MAPE和后验差比C是两套独立的评价体系，最好同时看。有的模型MAPE很低，但残差存在明显的趋势性波动，说明还有信息没提取干净，这时候可以考虑残差修正。
5.3 模型效果还能怎么往上提
DNGM(1,1)本身已经比传统GM(1,1)强，但实际项目中我经常还会叠加下面三个技巧。
第一个是残差修正。把模型的残差序列再拿去做一次灰色预测，然后叠加回原预测结果。这个方法对系统性偏差的修正效果很明显，尤其当残差里还有线性趋势的时候。不过残差建模的输入数据往往更“毛糙”，样本量又更小，所以修正模型也不宜建得太复杂。
第二个是滚动的新陈代谢预测。前面已经提过，每次拿到新真实值后，把最旧的数据丢弃、加入新数据重新建模。这样做的好处是模型永远跟随最新趋势，不会因为某个历史极端值被一直影响。
第三个是组合预测。把DNGM(1,1)当作基学习器，再用一个机器学习模型去拟合它的残差，或者和ARIMA、指数平滑的结果做加权融合。灰色模型擅长提取整体趋势，统计模型或机器学习模型擅长捕捉局部波动，两者结合常常能拿下一个相当不错的预测精度。这类组合方式在数学建模竞赛里也比较讨巧，因为评委喜欢看到“为什么单一模型不够，组合为什么更好”的完整逻辑。
最后再分享一个我的个人习惯。现在拿到一组小样本数据，我不会直接套模型，而是先把数据点画出来看一眼增长形态。如果接近指数型，传统GM(1,1)足够；如果接近线性或带有趋势扰动，我第一个跑的就是DNGM(1,1)；如果数据本身波动很大、看不出规律，那就先用平移或平滑处理，再做级比检验决定要不要走灰色这条路。工具从来不嫌多，关键是每把工具用在合适的地方。DNGM(1,1)的出现不是为了取代GM(1,1)，而是给预测工具箱里增加了一把更趁手的改锥——当数据呈现出线性特征时，它比传统模型更贴合真实规律，这就是它真正的价值所在。
