灰度预测算法实现

灰色预测的概念

灰色系统、白色系统和黑色系统

白的系统是指一个系统的内部特征是完全已知的，既系统信息是完全充分的
黑色系统是一个系统的内部信息对外界来说一无所知，只能通过它与外界的联系来加以预测研究
灰色介于白色和黑色之间，灰色系统的一部分信息是已知的，另外一部分信息是未知的，系统内各因素间有不确定的关系。

灰色预测法

灰色预测法是一种灰色预测系统的预测方法
灰色预测通过鉴别系统因素之间发展趋势的相异程度，即进行关联分析，并对原始数据进行生成处理来寻找系统变动的规律，生成有较强规律性的数据序列，然后建立相应的微分方程模型，从而预测事物未来发展趋势的状况。

灰色生成数列

灰色系统理论认为，尽管客观表象复杂，但总是有整体功能的，因此必然蕴含某种内在规律，关键是如何选择适当的方式挖掘和利用它。灰色系统是通过对原始数据的整理来寻求其变化规律的，这是一种就数据寻求数据的现实规律的途径，也就是灰色序列的生成。一切灰色序列都能通过某种生成，弱化其随机性，显示其规律性。数据生成的常用方式有累加生成，累减生成和加权累加生成。

累加生成（AGO）

累减生成（IAGO）

加权邻值生成

灰色模型GM（1，1）

灰色系统理论是基于关联空间，光滑离散函数等概念定义灰导数与灰微分方程，进而用离散数据列建立微分方程形式的动态模型，即灰色模型是利用离散随机数经过生成变为随机性被显著削弱而且交友规律的生成数，建立起微分方程形式的模型，这样便于对其变化过程进行研究和描述。

GM(1,1)的白化型

GM(1,1)灰色预测的步骤

数据检验和处理

建立GM(1,1)模型

检验预测值

灰度预测实现代码如下所示：

def GM11(x0):
    x1 = x0.cumsum() #1-AGO序列 生成累加数据集合
    z1 = (x1[:len(x1)-1] + x1[1:])/2.0 # 紧邻均值（means）生成序列 写的太棒了 简直无敌
    z1 = z1.reshape(len(z1),1)    # reshape函数，不改变元素的数据，但是传入的参数为int 或 tuple
    B = np.append(-z1,np.ones_like(z1),axis = 1)                # ones_like函数返回全是1的原类型数组。 append函数
    Yn = x0[1:].reshape(len(x0)-1,1)
    [[a],[b]] = np.dot(np.dot(np.linalg.inv(np.dot(B.T,B)),B.T),Yn)          #dot函数计算矩阵乘积, linalg.inv为求逆矩阵   (本质上计算参数)
    f = lambda k: (x0[0]-b/a)*np.exp(-a*(k-1))-(x0[0]-b/a)*np.exp(-a*(k-2))  # 计算还原值 exp函数计算指数
    delta = np.abs(x0-np.array([f(i) for i in range(1,len(x0)+1)]))
    C = delta.std()/x0.std()
    P = 1.0 * (np.abs(delta - delta.mean()) < 0.6745 * x0.std()).sum() / len(x0)
    return f, a, b, x0[0], C, P  # 返回灰色预测函数、a、b、首项、方差比、小残差概率

近期做了一个项目，通过财政收入影响因素分析及预测模型，利用灰度预测建立模型，实现数据预测。具体的灰度预测模型代码如下：

#-*- coding: utf-8 -*-
'''
title:灰度预测学习算法
备注：numpy和pandas 太强大了也太好用了，还要多学习和实践。
'''

import numpy as np
import pandas as pd

def GM11(x0):
    x1 = x0.cumsum() #1-AGO序列 生成累加数据集合
    z1 = (x1[:len(x1)-1] + x1[1:])/2.0 # 紧邻均值（means）生成序列 写的太棒了 简直无敌
    z1 = z1.reshape(len(z1),1)    # reshape函数，不改变元素的数据，但是传入的参数为int 或 tuple
    B = np.append(-z1,np.ones_like(z1),axis = 1)                # ones_like函数返回全是1的原类型数组。 append函数
    Yn = x0[1:].reshape(len(x0)-1,1)
    [[a],[b]] = np.dot(np.dot(np.linalg.inv(np.dot(B.T,B)),B.T),Yn)          #dot函数计算矩阵乘积, linalg.inv为求逆矩阵   (本质上计算参数)
    f = lambda k: (x0[0]-b/a)*np.exp(-a*(k-1))-(x0[0]-b/a)*np.exp(-a*(k-2))  # 计算还原值 exp函数计算指数
    delta = np.abs(x0-np.array([f(i) for i in range(1,len(x0)+1)]))
    C = delta.std()/x0.std()
    P = 1.0 * (np.abs(delta - delta.mean()) < 0.6745 * x0.std()).sum() / len(x0)
    return f, a, b, x0[0], C, P  # 返回灰色预测函数、a、b、首项、方差比、小残差概率


inputfile =  '/Users/mac/Desktop/data/data1.csv'
data = pd.read_csv(inputfile)
data.index = range(1994,2014)
data.loc[2014] = None    # loc函数选择指定的行和列进行操作，异常方便。
data.loc[2015] = None

l = ['x1','x2','x3','x4','x5','x7']
for i in l:
    f = GM11(data[i][range(1994, 2014)].as_matrix())[0]  # as_matrix 将列向量转化成行向量。
    data[i][2014] = f(len(data)-1)
    data[i][2015] = f(len(data))
    data[i] = data[i].round(2)

print data[l+['y']]