优化的灰色模型

优化的是μ值和C值，若有错误希望大家多多指出

MATLAB 算法实现

clear;clc;
ysh =[0.851,0.66733333,0.603,0.779,0.857666667,0.374666667,0.669,0.813333333,0.252,0.536333333,0.190666667]; %输入的数据
n = length(ysh);
AGO = zeros(1,n);
AGO(1,1) = ysh(1,1);
for k = 2:n
    AGO(1,k) = AGO(1,k-1) + ysh(1,k);
end                                          %一次累加序列
for i = 1:n
    P(n) = ysh(n-1)/AGO(n);
    if P(n) <0.5
        continue
    else
        disp('不满足光滑性检验')
    end
end                                           %光滑性检验
for i = 1:n
    O(n) = AGO(n-1)/AGO(n);
    if 1<=O(n)<=1.5
        continue
    else
        disp('不满足指数规律检验要求')
    end
end                                            %指数规律检验
for i = 2:n
    S(n) = ysh(n)/ysh(n-1);
    if exp(-2/(n+1))<=S(n)<=exp(-2/(n+1))
        continue
    else
        disp('不满足级比检验要求')
    end
end                                             %级比检验
 
i = 0;
M = []; A = [];C = [];S = [];F = [];
for m = 0:0.001:1        %确定权重系数起始值、递增值和终止值
    i = i + 1;
    z = zeros(1,n-1);                            %计算背景值 
    for k = 1:n - 1
        z(1,k) = m*AGO(1,k+1) + (1-m)*AGO(1,k);   %邻值生成序列
    end
    B = zeros(n-1,2);
    Y = zeros(n-1,1);                             %计算数据矩阵 B 和数据向量 Y 
    for k = 1:n-1
        B(k,1) = -z(1,k);
        B(k,2) = 1;
        Y(k,1) = ysh(1,k+1);
    end
    N = zeros(2,1);
    N = inv(B'*B)*B'*Y;
    a = N(1,1);
    u = N(2,1);                                       %计算a，u值
    c1 = 0;c2 = 0;
    for k = 2:n
        c1 = c1 + ysh(1,k)*exp(-a*(k-1));
        c2 = c2 + exp(-2*a*(k-1));
    end
    c = ((ysh(1,1)-u/a)/(1-exp(a)) + c1)/((1-exp(a))^(-2)+c2);   %计算改进模型参数 C
    ycAGO = zeros(1,n);
    for k = 0:n-1
        ycAGO(1,k+1) = c*(1-exp(a))^(-1)*exp(-a*k)+u/a;          %计算预测值一次累加生成（1AGO）
    end
    yc = zeros(1,n);
    for k = 1:n-1
        yc(1,1) = ycAGO(1,1);
        yc(1,k+1) = c * exp(-a*k);                              %yc即预测值
    end
    s = sum((ysh - yc).^2);                                     %s即离差平方和
    M(1,i) = m;
    A(1,i) = a;
    U(1,i) = u;
    C(1,i) = c;
    S(1,i) = s;
    F(1,i) = m;
    F(2,i) = s;
end
Y = sortrows(F',2)                                           %依据残差大小进行排序
yc
plot(M,S)                                         %绘制离方平均和与权重系数的关系图

结果
由结果可知，当权重系数为0.500的时候，离差平方和最小，所以选择权重系数为0.500（下面的离差平方和相等的权重系数也可以）

## 若不满足级比检验的要求，则应该对原始数据进行适当的平移：Y（0k） = Y（0k） + C，直到满足要求为止。`clear;clc;

ysh = [0.5997,0.597,0.559,0.542,0.509,0.437,0.409,0.395,0.391,0.393333333,0.371666667]; %输入的数据
n = length(ysh);
AGO = zeros(1,n);
AGO(1,1) = ysh(1,1);
for k = 2:n
AGO(1,k) = AGO(1,k-1) + ysh(1,k);
end

for j = -10:0.01:10
aysh(n) = ysh(n) + j;
aS(n) = aysh(n)/aysh(n-1);
for i = 2:n
if exp(-2/(n+1))<=aS(n)<=exp(-2/(n+1))
i
j
disp(‘满足’)
else
j;
end
end
end

从上述代码的结果中找到适当的平移数C##

选择权重系数为0.500，进行灰色模型的构建，具体代码如下所示：

clear;clc;
ysh = [0.851,0.66733333,0.603,0.779,0.857666667,0.374666667,0.669,0.813333333,0.252,0.536333333,0.190666667];
n = length(ysh);
AGO = zeros(1,n);
AGO(1,1) = ysh(1,1); 
for k = 2:n
    AGO(1,k) = AGO(1,k-1) + ysh(1,k);
end                                 
z = zeros(1,n-1);
for k = 1:n - 1
    z(1,k) = 0.5*AGO(1,k+1) + 0.5*AGO(1,k);            %权值选择为0.5
end
B = zeros(n-1,2);
Y = zeros(n-1,1);
for k = 1:n-1
    B(k,1) = -z(1,k);
    B(k,2) = 1;
    Y(k,1) = ysh(1,k+1);
end
N = zeros(2,1);
N = inv(B'*B)*B'*Y;
a = N(1,1);
u = N(2,1);
c1 = 0;c2 = 0;
for k = 2:n
    c1 = c1 + ysh(1,k)*exp(-a*(k-1));
    c2 = c2 + exp(-2*a*(k-1));
end
c = ((ysh(1,1)-u/a)/(1-exp(a)) + c1)/((1-exp(a))^(-2)+c2);
ycAGO = zeros(1,n);
for k = 0:n-1
    ycAGO(1,k+1) = c*(1-exp(a))^(-1)*exp(-a*k)+u/a;
end
yc = zeros(1,n);
for k = 1:n-1
    yc(1,1) = ycAGO(1,1);
    yc(1,k+1) = c * exp(-a*k);
end
s = sum((ysh - yc).^2);
yc                           %预测值
H = yc(1:11);         
%计算残差序列
epsilon = ysh - yc;

%法一：相对残差Q检验
%计算相对误差序列
delta = abs(epsilon./ysh);
%计算相对误差Q
disp('相对残差Q检验：')
Q = mean(delta)

%法二：方差比C检验
disp('方差比C检验：')
C = std(epsilon, 1)/std(ysh, 1)

%法三：小误差概率P检验
S1 = std(ysh, 1);
tmp = find(abs(epsilon - mean(epsilon))< 0.6745 * S1);
disp('小误差概率P检验：')
P = length(tmp)/n

%绘制曲线图
t1 = 1:11;
t2 = 1:11;

plot(t1, ysh,'ro'); hold on;
plot(t2, yc, ':b*');
xlabel('xx'); ylabel('yy');
legend('实际值','预测值');

grid on;
    

结果如下图所示

可以根据三个检验的值来判断模型的精度，具体如下图所示：

参考文献> 1.杨旭. 改进的灰色预测GM(1,1)模型的MATLAB实现[J]. 江苏科技信息, 2014(7):69-70.
2.组合灰色模型在黑土区玉米产量预测中的应用