main.py/GM11.py

# -*- coding: utf-8 -*-
import numpy as np


def GM11(x0):
    """
    灰色预测GM(1,1)模型
    参数:
        x0: 原始序列（numpy数组或列表）
    返回:
        f: 预测函数
        a, b: 模型参数
        x0_0: 首项
        C: 方差比
        P: 小残差概率
    """
    x0 = np.array(x0).flatten()
    x1 = x0.cumsum()  # 1-AGO序列

    # 紧邻均值生成序列
    z1 = (x1[:len(x1) - 1] + x1[1:]) / 2.0
    z1 = z1.reshape((len(z1), 1))

    # 构建矩阵B
    B = np.append(-z1, np.ones_like(z1), axis=1)
    Yn = x0[1:].reshape((len(x0) - 1, 1))

    # 最小二乘法计算参数 [a, b]
    params = np.dot(np.dot(np.linalg.inv(np.dot(B.T, B)), B.T), Yn)
    a, b = params[0][0], params[1][0]

    # 预测函数
    def f(k):
        return (x0[0] - b / a) * np.exp(-a * (k - 1)) - (x0[0] - b / a) * np.exp(-a * (k - 2))

    # 模型检验
    y_hat = np.array([f(i) for i in range(1, len(x0) + 1)])
    delta = np.abs(x0 - y_hat)

    # 方差比C
    C = delta.std() / x0.std() if x0.std() != 0 else np.inf

    # 小残差概率P
    e = np.abs(delta - delta.mean())
    S0 = 0.6745 * x0.std()
    P = (e < S0).sum() / len(e)

    return f, a, b, x0[0], C, P