Matlab语言分析

筱风听雨

最近试了一下灰色模型，这个代码是别人基础上改的，不知道为啥C总是1？
%建立符号变量a(发展系数)和b(灰作用量)
syms a b;
c = [a b]';

%原始数列 A
A = [3182.8203
3228.9778
3270.196
3303.9235
3348.3485
3394.5033
3457.8879
3509.7971
3575.7637
3632.3111
3683.8811
3746.2196
3810.6418
3864.6374
3913.0927
3966.0405
4015.4459
4062.5406
4105.4635
4150.3338
4191.2074
4231.1742
4148.5447
4185.7676
4222.4273
4254.2255
4283.5667
4311.2439
4339.1287
4368.4125
4400.1038
4432.1581
4462.2489
4488.4367
4503.9321
4522.1468
4542.1607
4565.6316
4592.2644
4622.0636
4647.5728
4666.1318];
n = length(A);

%对原始数列 A 做累减得到数列 B
for i = 2:n
    H(i) = A(i) - A(i - 1);
end
H(1) = [];


%对原始数列 A 做紧邻均值生成
for i = 2:n
    C(i) = (A(i) + A(i-1))/2;
end
C(1) = [];

%构造数据矩阵
D = [-C; C.^2];
Y = H; Y = Y';

%使用最小二乘法计算参数 a(发展系数)和b(灰作用量)
c = inv(D*D')*D*Y;
c = c';
a = c(1); b = c(2);

%得到预测出的数据
F = []; F(1) = A(1);
for i = 2:(n+n)
    F(i) = (a*A(1))/(b*A(1)+(a - b*A(1))*exp(a*(i-1)));
end

disp('预测数据为：');
F

%绘制曲线图
t1 = 0:n-1;
t2 = 0:2*n-1;
plot(t1, A, 'ro'); hold on;
plot(t2, F);
xlabel('年份'); ylabel('人口数');
legend('实际数量','预测数量');
title('人口预测');
grid on;
disp('预测数据为：');
F

%模型检验

H = F(1:20);
%计算残差序列
epsilon = A - H;

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

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

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