Computer programs For S-Distributions Workshop Main Page a. Calculating the fixing parameter "τ(m,ρ)" ( in the sense that it fixes the mixture to be density). Search The body function of the root and suitable r , as it explained in the paper, the program are in the file named as “Tao.txt” and “FoundR”. The main program is in the file” Positive Root” which contain these parts: Program Environment: MATLAB 7.8.0(R2009a) Inputs: “N, P”, which are “m, ρ”, the mixture order parameter and the right side skewness parameter, respectively in the paper. Outputs: Rmin: which is the fixing parameter "τ(m,ρ)" .
a1. For the positive Root:
function [ Rmin ]= Positive_Root(N,p,rstep)
L=[]; for i= 1:N+1 l=(0:1/i:1); l=l(1:end); L=[L,[l;(0:i);zeros(size(l))+i-1]]; end
index_p=1; index_k=2; index_N=3;
[b m n]=unique(L(index_p,:)); [a in]=sort(L(index_p,:));
QL=L(:,m); NQ=L(:,in); PNQ=NQ(index_p,:)>0; NQL=NQ(:,PNQ); L_in=find(QL(index_p,:)<p);
is=0; [X R tao]=FoundR(NQL,QL,p,N,rstep,is);
Xmin=X; Rmin=R; MXmin=[X]; MRmin=[R]; for is=1:size(L_in,2)-2 [x r]=FoundR(NQL,QL,p,N,rstep,is); MXmin=[x MXmin]; MRmin=[r MRmin]; if r<Rmin Rmin=r; end end % end
a2. Tao: function tao=Tao(n,k,p) Rn=inline('sin((k.*pi)/(n+1)) ./ sin(p.*pi - (k.*pi)/(n+1))'); tao=Rn(k,n,p); end
a3. FoundR:
function [X R tao]=FoundR(NQL,QL,p,N,rstep,is)
index_p=1; index_k=2; index_N=3;
L_in=find(QL(index_p,:)<p); lp=QL(index_p,L_in(end)); U_in=find(QL(index_p,:)>p); if is==0 U_in=find(NQL(index_p,:)==lp);%QL(index_p,U_in(1))); else U_in=find(NQL(index_p,:)==QL(index_p,L_in(end-is))); end up=NQL(index_p,U_in(end)+1); Up=NQL(:,U_in(end)); Q=NQL(:,find(NQL(index_p,:)<up)); T=[]; for n=1:N if mod(sum(Q(index_N,:)==n),2)==1 T=[T n]; end end
n=Up(index_N); k=Up(index_k); r=rstep:rstep:10; tao=repmat(Tao(n,k,p),[size(r)]); rtao=r./tao; SRT=zeros(size(r)); for i=1:size(T,2) SRT=SRT+rtao.^T(i); end SRT=SRT.*(SRT<1); [srt indexr]=max(SRT); R=r(indexr); X=rtao(indexr); tao=tao(1); end
b. MLE program for parameters estimation:
Programs Environment: MATLAB 7.8.0(R2009a)
For calculating MLE's parameters, firstly we have to define the body function of the Log likelihood for the S distribution, which have three parts in the files named as (f1.txt, avg_f1.txt, f2.txt) respectively. Then in the program named as, main7.txt , we estimate the parameter by using MLE method, this program has these parts: Program Environment: MATLAB 7.8.0(R2009a) Inputs: finp1='tr.xls'; which is the matrix of "m " and " ρ" and the fixing parameter "τ(m,ρ)" ( in the sense that it fixes the mixture to be density). finp2='AMEX.xls'; the file of data that we can used dat: for reading data from excel file m: which is the mixture order parameter Outputs: The MLE's estimated parameters, "rh, ah, bh, rhoh " which are " r, α, λ, ρ" respectively, in the paper.
b1. f1.txt:
function y=f1(t,n,rho,a,b) y=0; for j=1:(n+1) y=y+nchoosek(n+1,j).*b.^(j-1).*t.^(a*j-a).*sin(j.*pi.*rho); end; end
b2. avg_f1.txt:
function y=avg_f1(t,m,rho,r,a,b)
y=0; x=r*sin(pi*rho);
b3. f2.txt:
function y=f2(t,m,rho,r,a,b,tr,nrho) if nargin==7 nrho=size(tr,1)+1; end; if(rho<=0||rho>=1||b<=0||r<=0||a<=0||a>=2||r>=tr(rho*nrho,m)) y=1e7; else y=-sum(log(avg_f1(t,m,rho,r,a,b))); end; end
b4. main7.txt:
clear all;clc; finp1='tr.xls'; finp2='AMEX.xls'; tr=xlsread(finp1,1,'B2:K10'); tr(tr>=1000)=10; dat= xlsread(finp2,1,'D1:D1200'); [nr_tr nc_tr]=size(tr); m=2; tic; for i=1:nr_tr rho=i/(nr_tr+1); lf=@(th) f2(dat,m,rho,th(1),th(2),th(3),tr); [x fx]= patternsearch(lf,[min(tr(i,m)/2,.6),.4,.4],[],[],[],[],[0,0,0],[tr (rho*10,m),2,10000000000]) A(i,:)=[x,rho,fx]; fprintf('i = %3d\t,\ttime : %4d seconds\n',i,round(toc)); end; [lmin,imin]=min(A(:,5)); th_hat=A(imin,:); rh=th_hat(1); ah=th_hat(2); bh=th_hat(3); rhoh=th_hat(4); th_hat
|
|