麻瓜 | 数学建模日记 | 第三天

本文介绍了一位编程者使用Python实现灰色预测模型(GM(1,1))的过程,包括级比检查、预测函数和精度评估。通过对时间序列数据(例如2008年至2018年的某指标数据)的建模,展示了如何预测未来值并计算相对误差、预测准确率和级比偏差。最后,通过可视化展示预测结果。
摘要由CSDN通过智能技术生成

python 灰色模型预测

建模体验

搞代码的一天,还是十分迷糊

不过幸好已确定以后用 python 来进行处理接下来的编程问题

麻瓜还可以!

灰色预测方法

附上未完成的代码

import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
def step_ratio(x0):
    n = len(x0)
    ratio = [x0[i]/x0[i+1] for i in range (len(x0)-1)]
    print (f"级比:{ratio}")
    min_la,max_la = min(ratio),max(ratio)
    thred_la = [np.exp(-2/(n+2)),np.exp(2/(n+2))]
    if min_la <thred_la[0] or max_la > thred_la[-1]:
        print ("级比超过灰色模型的范围")
    else:
        print ("级比满足要求,可用GM(1,1)模型")
    return ratio,thred_la

def predict(x0):
    n = len(x0)
    x1 = np.cumsum(x0)
    z = np.zeros(n-1)
    for i in range (n-1):
        z[i] = 0.5*(x1[i]+x1[i+1])
    B = [-z,[1]*(n-1)]
    Y = x0[1:]
    u = np.dot(np.linalg.inv(np.dot(B,np.transpose(B))),np.dot(B,Y))
    x1_solve = np.zeros(n)
    x0_solve = np.zeros(n)
    x1_solve[0] = x0_solve[0] = x0[0]
    for i in range (1,n):
        x1_solve[i] = (x0[0]-u[1]/u[0])*np.exp(-u[0]*i)+u[1]/u[0]
    for i in range (1,n):
        x0_solve[i] = x1_solve[i]-x1_solve[i-1]
    return x0_solve,x1_solve,u

def accuracy(x0,x0_sovle,ratio,u):
    epsilon = x0-x0_solve
    delta = abs(epsilon/x0)
    print (f"相对误差:{delta}")
    #Q = np.mean(delta)
    #C = np.std(epsilon)/np.std(x0)
    S1 = np.std(x0)
    S1_new = S1*0.6745
    temp_P = epsilon[abs(epsilon-np.mean(epsilon))<S1_new]
    P = len(temp_P)/len(x0)
    print (f"预测准确率:{P*100}%")
    ratio_solve = [x0_solve[i]/x0_solve[i+1] for i in range (len(x0_solve)-1)]
    rho = [1-(1-0.5*u[0]/u[1])/(1+0.5)*u[0]/u[1]*(ratio[i]/ratio_solve[i])for i in range (len(ratio))]
    print (f"级比偏差:{rho}")
    return epsilon,delta,rho,P

if __name__ == '__main__':
    data = pd.DataFrame(data = {"year":[2008,2009,2010,2011,2012,2013,2014,2015,2016,2017,2018],\
                                "eqL":[418.65,435.08,460.11,463.36,467.40,470.00,472.66,475.55,477.50,477.70,483.00]})
    x0 = np.array(data.iloc[:,1])
    ratio,thred_la = step_ratio(x0)
    x0_solve,x1_solve,u = predict(x0)
    epsilon,delta,rho,P = accuracy(x0,x0_solve,ratio,u)
    plt.title("Predict Graph")
    plt.xlabel("Years")
    plt.ylabel("Result")
    plt.plot(data.year,data.eqL)
    plt.legend()
    plt.show()


评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值