plot_fit2.m

用于模拟时变非平稳的ARMA过程

clear;tfpm;

re_im= 'i';
mo_no= 'n';
alpha= 1/2;
beta = 1/2;
rho= 1-log(12);

N= 256;
M= 3;
L= 2;
MSEm= [];
meanMm= [];
meanLm= [];
resultname= sprintf('data/%04d/%s%d%d%s%s.mat', N, 'ord2C', M, L, re_im, mo_no);
load(resultname)
sprintf('%s%4d', resultname, MM)

   
%Ns= (2:(1+Mmax))'*(2*(0:Lmax)+1)-1;%FOR TFAR!!!!!!!!!!!!!
Ns= (2*(1:Mmax))'*(2*(0:Lmax)+1)-1;%FOR TFARMA!!!!!!!!!!!!!

eeprime= log(EEprime(:, :, 1, :));

MDLhist= zeros(Mmax, Lmax+1);
BIChist= zeros(Mmax, Lmax+1);
AIChist= zeros(Mmax, Lmax+1);
GIChist= zeros(Mmax, Lmax+1);
MIChist= zeros(Mmax, Lmax+1);

for mm= 1:MM
   eprime= eeprime(:, :, 1, mm);
   MDL= eprime + (log(N+1)+rho)  *(Ns+1/2)/N;
   BIC= eprime + log(N)          *(Ns    )/N;
   AIC= eprime + 2               *(Ns    )/N;
   GIC= eprime + (2+rho)         *(Ns    )/N;
   MIC= eprime + (1+rho+log(N)/2)*(Ns    )/N;

   [MMDL, LMDL] = find(MDL==min(min(MDL)));LMDL= LMDL-1;
   [MBIC, LBIC] = find(BIC==min(min(BIC)));LBIC= LBIC-1;
   [MAIC, LAIC] = find(AIC==min(min(AIC)));LAIC= LAIC-1;
   [MGIC, LGIC] = find(GIC==min(min(GIC)));LGIC= LGIC-1;
   [MMIC, LMIC] = find(MIC==min(min(MIC)));LMIC= LMIC-1;

%   [N M L mm MMDL, LMDL, MBIC, LBIC, MAIC, LAIC, MGIC, LGIC, MMIC, LMIC]
%   [N M L mm Mmdl, Lmdl, Mbic, Lbic, Maic, Laic, Mgic, Lgic, Mmic, Lmic]

   MDLhist(MMDL, LMDL+1)= MDLhist(MMDL, LMDL+1) + 1;
   BIChist(MBIC, LBIC+1)= BIChist(MBIC, LBIC+1) + 1;
   AIChist(MAIC, LAIC+1)= AIChist(MAIC, LAIC+1) + 1;
   GIChist(MGIC, LGIC+1)= GIChist(MGIC, LGIC+1) + 1;
   MIChist(MMIC, LMIC+1)= MIChist(MMIC, LMIC+1) + 1;
end;

figure(1);clf
%subplot(1, 2, 1);
   bar([sum(MDLhist')' sum(AIChist')'], 1, 'grouped')
   axis([1 Mmax 0 100])
   colormap(gray)
   grid
%   ylabel('percentage')
   xlabel('Mhat')
   set(gca, 'XTick', 1:8)
   set(gca, 'XTicklabel', ['M1'; 'M2'; 'M3'; 'M4'; 'M5'; 'M6'; 'M7'; 'M8'])
   title('a')
%print -deps2 -tiff -r300 /users/mjachan/tex/prj/TF/figures/ordhist1.eps

figure(2);clf
%subplot(1, 2, 2);
   bar([sum(MDLhist); sum(AIChist)]', 1, 'grouped')
   axis([1 Lmax+2 0 100])
   colormap(gray)
   grid
%   ylabel('percentage')
   xlabel('Lhat')
   set(gca, 'XTick', 1:6)
   set(gca, 'XTicklabel', ['L0'; 'L1'; 'L2'; 'L3'; 'L4'; 'L5'])
   title('b')
%print -deps2 -tiff -r300 /users/mjachan/tex/prj/TF/figures/ordhist2.eps












N= 256;
L= 2;
MSEm= [];
corrMm= [];
corrLm= [];
for M= 2:5
   resultname= sprintf('data/%04d/%s%d%d%s%s.mat', N, 'ord2C', M, L, re_im, mo_no);
   load(resultname)
   sprintf('%s%4d', resultname, MM)

%   Ns= (2:(1+Mmax))'*(2*(0:Lmax)+1)-1;%FOR TFAR!!!!!!!!!!!!!
   Ns= (2*(1:Mmax))'*(2*(0:Lmax)+1)-1;%FOR TFARMA!!!!!!!!!!!!!
   
   eeprime= log(EEprime(:, :, 1, :));

   MDLhist= zeros(Mmax, Lmax+1);
   BIChist= zeros(Mmax, Lmax+1);
   AIChist= zeros(Mmax, Lmax+1);
   GIChist= zeros(Mmax, Lmax+1);
   MIChist= zeros(Mmax, Lmax+1);

   for mm= 1:MM
      eprime= eeprime(:, :, 1, mm);
      MDL= eprime + (log(N+1)+rho)  *(Ns+1/2)/N;
      BIC= eprime + log(N)          *(Ns    )/N;
      AIC= eprime + 2               *(Ns    )/N;
      GIC= eprime + (2+rho)         *(Ns    )/N;
      MIC= eprime + (1+rho+log(N)/2)*(Ns    )/N;

      [MMDL, LMDL] = find(MDL==min(min(MDL)));LMDL= LMDL-1;
      [MBIC, LBIC] = find(BIC==min(min(BIC)));LBIC= LBIC-1;
      [MAIC, LAIC] = find(AIC==min(min(AIC)));LAIC= LAIC-1;
      [MGIC, LGIC] = find(GIC==min(min(GIC)));LGIC= LGIC-1;
      [MMIC, LMIC] = find(MIC==min(min(MIC)));LMIC= LMIC-1;

%      [N M L mm MMDL, LMDL, MBIC, LBIC, MAIC, LAIC, MGIC, LGIC, MMIC, LMIC]
%      [N M L mm Mmdl, Lmdl, Mbic, Lbic, Maic, Laic, Mgic, Lgic, Mmic, Lmic]

      MDLhist(MMDL, LMDL+1)= MDLhist(MMDL, LMDL+1) + 1;
      BIChist(MBIC, LBIC+1)= BIChist(MBIC, LBIC+1) + 1;
      AIChist(MAIC, LAIC+1)= AIChist(MAIC, LAIC+1) + 1;
      GIChist(MGIC, LGIC+1)= GIChist(MGIC, LGIC+1) + 1;
      MIChist(MMIC, LMIC+1)= MIChist(MMIC, LMIC+1) + 1;
   end;
   
   if(MM==50)
      AIChist= 2*AIChist;
      MDLhist= 2*MDLhist;
   end;
   
   MDLMmarg= sum(MDLhist');
   AICMmarg= sum(AIChist');
   corrMm= [corrMm [MDLMmarg(M); AICMmarg(M)]];
end;

figure(3);clf;hold on
%subplot(1, 2, 1);hold on;
   plot(corrMm(1, :), 'Linewidth', 2, 'Linestyle', '-', 'Color', 'k')
   plot(corrMm(2, :), 'Linewidth', 2, 'Linestyle', '--', 'Color', 'k')
   axis([1 4 0 100])
   box;grid
%   ylabel('percentage')
   xlabel('M')
   set(gca, 'XTick', 1:8)
   set(gca, 'XTicklabel', ['M2'; 'M3'; 'M4'; 'M5'])
   title('c')
%print -deps2 -tiff -r300 /users/mjachan/tex/prj/TF/figures/ordcorr1.eps


   
   
N= 256;
M= 2;
MSEm= [];
corrLl= [];
for L= 0:3
   resultname= sprintf('data/%04d/%s%d%d%s%s.mat', N, 'ord2C', M, L, re_im, mo_no);
   load(resultname)
   sprintf('%s%4d', resultname, MM)

   Ns= (2*(1:Mmax))'*(2*(0:Lmax)+1)-1;%FOR TFARMA!!!!!!!!!!!!!
   
   eeprime= log(EEprime(:, :, 1, :));

   MDLhist= zeros(Mmax, Lmax+1);
   BIChist= zeros(Mmax, Lmax+1);
   AIChist= zeros(Mmax, Lmax+1);
   GIChist= zeros(Mmax, Lmax+1);
   MIChist= zeros(Mmax, Lmax+1);

   for mm= 1:MM
      eprime= eeprime(:, :, 1, mm);
      MDL= eprime + (log(N+1)+rho)  *(Ns+1/2)/N;
      BIC= eprime + log(N)          *(Ns    )/N;
      AIC= eprime + 2               *(Ns    )/N;
      GIC= eprime + (2+rho)         *(Ns    )/N;
      MIC= eprime + (1+rho+log(N)/2)*(Ns    )/N;

      [MMDL, LMDL] = find(MDL==min(min(MDL)));LMDL= LMDL-1;
      [MBIC, LBIC] = find(BIC==min(min(BIC)));LBIC= LBIC-1;
      [MAIC, LAIC] = find(AIC==min(min(AIC)));LAIC= LAIC-1;
      [MGIC, LGIC] = find(GIC==min(min(GIC)));LGIC= LGIC-1;
      [MMIC, LMIC] = find(MIC==min(min(MIC)));LMIC= LMIC-1;

%      [N M L mm MMDL, LMDL, MBIC, LBIC, MAIC, LAIC, MGIC, LGIC, MMIC, LMIC]
%      [N M L mm Mmdl, Lmdl, Mbic, Lbic, Maic, Laic, Mgic, Lgic, Mmic, Lmic]

      MDLhist(MMDL, LMDL+1)= MDLhist(MMDL, LMDL+1) + 1;
      BIChist(MBIC, LBIC+1)= BIChist(MBIC, LBIC+1) + 1;
      AIChist(MAIC, LAIC+1)= AIChist(MAIC, LAIC+1) + 1;
      GIChist(MGIC, LGIC+1)= GIChist(MGIC, LGIC+1) + 1;
      MIChist(MMIC, LMIC+1)= MIChist(MMIC, LMIC+1) + 1;
   end;

   if(MM==50)
      AIChist= 2*AIChist;
      MDLhist= 2*MDLhist;
   end;

   MDLLmarg= sum(MDLhist);
   AICLmarg= sum(AIChist);
   corrLl= [corrLl [MDLLmarg(L+1); AICLmarg(L+1)]];
end;

figure(4);clf;hold on
%subplot(1, 2, 2);hold on;
   plot(corrLl(1, :), 'Linewidth', 2, 'Linestyle', '-', 'Color', 'k')
   plot(corrLl(2, :), 'Linewidth', 2, 'Linestyle', '--', 'Color', 'k')
   axis([1 4 0 100])
   box;grid
%   ylabel('percentage')
   xlabel('L')
   set(gca, 'XTick', 1:4)
   set(gca, 'XTicklabel', ['L0'; 'L1'; 'L2'; 'L3'])
   title('d')
%print -deps2 -tiff -r300 /users/mjachan/tex/prj/TF/figures/ordcorr2.eps