hurst_expo.m

hurst parameter: 給入一序列的資類

function hurst_expo(x)  
%
% hurst_expo(pnl) will plot a pnl loglog plot.
%
%
hurst_base=x;
H_estimate=zeros(size(x));
stdx=zeros(size(x));
avgx=zeros(size(x));
Rx=zeros(size(x));


for i=21:length(hurst_base)%basic number of data 
    disp(i);
   
    if i>100 && mod(i,2)~=0
        continue
    elseif i>1000 && mod(i,20)~=0
        continue
    elseif i>10000 && mod(i,100)~=0
        continue
    end
     periodN=fix(length(hurst_base)/i);
    hurst_elem=nan;
    for k=1:periodN %cal very segement
        hurst_elem(k)=hurst_expo_segement(hurst_base(1+i*(k-1):i*k));
    end
    hurst_para(i)=sum(hurst_elem)/length(hurst_elem);
    fithurst(i,1)=hurst_para(i);
    xaxis(i)=i;
    fitx(i,1)=i;
end

%%fitting program
hurst_para=f_zero2nan(hurst_para);
%fithurst(:,1)=f_zero2nan(fithurst(:,1));
Hi=1;

Hperiod=fix(length(hurst_para)/100);
for i=fix(length(hurst_para)/10): Hperiod:length(hurst_para)
    %'Lower',[0,0],'Upper',[Inf,max(hurst_para(1:i))], 'Startpoint',hurst_para(1)
    s=fitoptions('Method','NonlinearLeastSquares','StartPoint',[0 0],'Lower',[0 0]);
    f = fittype('a*x^h','options',s);%,'problem',h,'options',s);
    
    %s = fitoptions('Method','SmoothingSpline');
    %options.lower=[-inf 0];
    %s.upper=[inf 1];
 
    [fitting_func,gof2] = fit(fitx(1:i,1),fithurst(1:i,1),f);
    if Hi==1 
        H_estimate(1:i)=fitting_func.h;
        H_great=fitting_func.h;
        A_great=fitting_func.a;
        fitting_func_great=fitting_func;
    else
        H_estimate(i-Hperiod:i)=fitting_func.h;
        if fitting_func.h > H_great
            H_great=fitting_func.h;
             A_great=fitting_func.a;
            fitting_func_great=fitting_func;
        end
    end
    Hi=Hi+1;
end
Haxis=1:length(H_estimate);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%plot all figure

ax(1)=subplot(2,2,2);
%plot(fitting_func,xaxis,hurst_para);
plot(xaxis,A_great.*xaxis.^H_great,xaxis,hurst_para);
legend('boxoff');
title('hurst fitting');
grid on;
ax(2)=subplot(2,2,4);
plot(Haxis,f_zero2nan(H_estimate));
title('H estimate');
grid on;
ax(3)=subplot(2,2,[1 3]);
loglog(xaxis,hurst_para,'r',xaxis,A_great.*xaxis.^H_great);    

grid on;
title('hurst');        
            
function one_segement=hurst_expo_segement(x)  
%
% hurst_expo(pnl) will plot a pnl loglog plot.
%still going... how to estimate H
%
hurst_base=x;
hurst_elem=zeros(size(x));
stdx=zeros(size(x));
avgx=zeros(size(x));
Rx=zeros(size(x));

   stdx=std(hurst_base,1);
   avgx=sum(hurst_base)/length(hurst_base);

    Rx_base=hurst_base-avgx;
   
        for m=1:length(hurst_base)
            if m==1 
                Rx_elem(m)=Rx_base(m);
            else
                Rx_elem(m)=Rx_base(m)+Rx_elem(m-1);
            end
        end
        Rx=max(Rx_elem)-min(Rx_elem);
  
    one_segement=Rx/stdx;