gasanfis.m

使用自适应神经模糊推理系统的方法预测时间序列

%xuebao07:Forecasting Coalmine Gas Concentration Based on ANFIS
epoch_n=10;%训练次数为10（LJ）
mf_n=2;
ddm=5;%嵌入维数

load gas_data.dat % Sample data of Gas concentration
num =(1:2000)';
ts = gas_data(1:2200, 1);
sample=zeros(2000,5);
trn_data = zeros(1500, 5);
chk_data = zeros(500, 5);

delay=4;%延迟时间

start = 1;
sample(:, 1) = ts(start:start+2000-1); %1：2000
start = start + delay;%=1+4=5
sample(:, 2) = ts(start:start+2000-1); %5：2004
start = start + delay;%5+4=9
sample(:, 3) = ts(start:start+2000-1); %9：2008
start = start + delay;%9+4=13
sample(:, 4) = ts(start:start+2000-1); %13：2012
start = start + delay;%13+4=17
sample(:, 5) = ts(start:start+2000-1); %17：2016

trn_data=sample(1:1500, 1:5);
chk_data=sample(1501:2000, 1:5);

mf_type = 'gbellmf';
in_fismat = genfis1(trn_data, mf_n, mf_type);

ss = 0.1;
mf_type = 'gbellmf';
tic;
[trn_out_fismat trn_error step_size chk_out_fismat chk_error] =	anfis(trn_data, in_fismat, [epoch_n nan ss nan nan], [1,1,1,1], chk_data);
toc;

figure('name', 	['ANFIS: Gas Concentration Prediction (file: ' mfilename ')'],	'number', 'off');

subplot(221);
tmp = [trn_error chk_error];
plot(tmp);
title('Error Curves');
axis([0 epoch_n min(tmp(:)) max(tmp(:))]);
legend('Training Error', 'Checking Error');

subplot(222);
plot(step_size);
title('Step Sizes');

subplot(223);
input = [trn_data(:, 1:4); chk_data(:, 1:4)];
anfis_output = evalfis(input, trn_out_fismat);
index = 1:2000;
outputt=[trn_data(:, 5); chk_data(:, 5)];

plot(num(index), [outputt anfis_output]);
axis([min(num(index)) max(num(index)) min(ts(index)) max(ts(index))]);
title('Desired and ANFIS Outputs')

subplot(224);

diff = outputt-anfis_output;
plot(num(index), abs(diff));
axis([min(num(index)) max(num(index)) min(abs(diff)) max(abs(diff))]);
title('Prediction Errors')

figure('name',	['ANFIS: Gas Concentration Prediction (file: ' mfilename ')'],	'number', 'off');

in_range = getfis(trn_out_fismat, 'inrange');
mf_type = getfis(trn_out_fismat, 'inmftypes');
mf_para = getfis(trn_out_fismat, 'inmfparams');
in_n = getfis(trn_out_fismat, 'numinputs');
in_mf_n = getfis(trn_out_fismat, 'NumInputMfs');
start = [0 cumsum(in_mf_n)] + 1;

subplot(221);
input_index = 1;
curr_in_range = in_range(input_index, :);
curr_mf_type = mf_type(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
curr_mf_para = mf_para(start(input_index): 	start(input_index)+in_mf_n(input_index)-1, :);
x = linspace(curr_in_range(1), curr_in_range(2), 101);
mf = evalmmf(x, curr_mf_para, curr_mf_type);
plot(x', mf');
axis([curr_in_range 0 1.2]);
title('Final MFs on Input 1, x(t - 18)');

subplot(222);
input_index = 2;
curr_in_range = in_range(input_index, :);
curr_mf_type = mf_type(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
curr_mf_para = mf_para(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
x = linspace(curr_in_range(1), curr_in_range(2), 101);
mf = evalmmf(x, curr_mf_para, curr_mf_type);
plot(x', mf');
axis([curr_in_range 0 1.2]);
title('Final MFs on Input 2, x(t - 12)');

subplot(223);
input_index = 3;
curr_in_range = in_range(input_index, :);
curr_mf_type = mf_type(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
curr_mf_para = mf_para(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
x = linspace(curr_in_range(1), curr_in_range(2), 101);
mf = evalmmf(x, curr_mf_para, curr_mf_type);
plot(x', mf');
axis([curr_in_range 0 1.2]);
title('Final MFs on Input 3, x(t - 6)');

subplot(224);
input_index = 4;
curr_in_range = in_range(input_index, :);
curr_mf_type = mf_type(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
curr_mf_para = mf_para(start(input_index): start(input_index)+in_mf_n(input_index)-1, :);
x = linspace(curr_in_range(1), curr_in_range(2), 101);
mf = evalmmf(x, curr_mf_para, curr_mf_type);
plot(x', mf');
axis([curr_in_range 0 1.2]);
title('Final MFs on Input 4, x(t)');

figure;

subplot(2,2,1); plot(tmp);
title('Error Curves');
axis([0 epoch_n min(tmp(:)) max(tmp(:))]);
xlabel('Epochs');
ylabel('Root-mean-squared Error')

legend('Training Error', 'Checking Error');

trn_m=anfis_output(1:1500)';
tst_m=anfis_output(1501:2000)';
trn_d=trn_m-trn_data(:, ddm)';
tst_d=tst_m-chk_data(:, ddm)';
trn_error=0;
tst_error=0;
for i=1:1500
    trn_error=trn_error+trn_d(i)*trn_d(i);
end
for i=1:500
    tst_error=tst_error+tst_d(i)*tst_d(i);
end
trn_error=sqrt(trn_error/1500);
t_e=trn_error/std(trn_data(:, ddm)');
tst_error=sqrt(tst_error/500);
c_e=tst_error/std(chk_data(:, ddm)');
errorr=[trn_error; tst_error; t_e; c_e]
figure;
subplot(312)
plot([1:1:1500], [trn_data(:, ddm) trn_m'])
figure;
subplot(312)
plot([1:1:500], [chk_data(:, ddm) tst_m'])
figure;
subplot(312);

plot(num(index), abs(diff));
axis([min(num(index)) max(num(index)) min(abs(diff)) max(abs(diff))]);