How to find largest lyapunov for Mackey Glass data?

조회 수: 10 (최근 30일)
Zahra kamkar
Zahra kamkar 2015년 8월 16일
답변: alireza ghavami 2023년 8월 7일
Hello.
I want to find largest Lyapunov for Mackey Glass data. How can I do that?

답변 (2개)

Zahra kamkar
Zahra kamkar 2015년 8월 22일
편집: Star Strider 2015년 8월 22일
Hello
I have another question. I used a Lyapunov Exponent Matlab code to find the largest Lyapunov result. The code seems to be correct but I don't know what should be the value of Maximum Lag in this code!
here is the code:
% close all;
clear all;
y= xlsread('MackeyGlass.xlsx');
warning off
tic
%___________________________Defining lags for y____________________________
%%Here is the my problem
maxlag=3;
y=y(:);
[~,nyc]=size(y);
yL=lagmatrix(y,1:maxlag);
y=y(maxlag+1:end,1);
yL=yL(maxlag+1:end,:);
[ryL cyL]=size(yL);
%_________________ Regressors Up to degree 3_______________________________
X1=yL;
num1=0;
X2ij=[];
for i=1:cyL
for j=i:cyL
X2ij=[X2ij yL(:,i).*yL(:,j)];
Indexij(num1+1,1)=i;
Indexij(num1+1,2)=j;
num1=num1+1;
end
end
num2=0;
X3ijk=[];
for i=1:cyL
for j=i:cyL
for k=j:cyL
X3ijk=[ X3ijk yL(:,i).*yL(:,j).*yL(:,k)];
Indexijk(num2+1,1)=i;
Indexijk(num2+1,2)=j;
Indexijk(num2+1,3)=k;
num2=num2+1;
end
end
end
X=[ones(ryL,1) X1 X2ij X3ijk];
beta =inv (X'*X)*X'*y;
betaX1=beta(2:cyL+1);
betaX2ij=beta(cyL+2:cyL+1+num1);
betaX3ijk=beta(cyL+2+num1:cyL+1+num1+num2);
for i=1:cyL
Inij1=find(Indexij(:,1)==i);
Inij2=find(Indexij(:,2)==i);
Inijk1=find(Indexijk(:,1)==i);
Inijk2=find(Indexijk(:,2)==i);
Inijk3=find(Indexijk(:,3)==i);
Df1(:,i)=betaX1(i,1)*ones(ryL,1);
Df210=zeros(ryL,1);
for h21=1:length(Inij1)
Df210=betaX2ij(Inij1(h21),1)*yL(:,Indexij(Inij1(h21),2))+Df210;
end
Df21(:,i)=Df210;
Df220=zeros(ryL,1);
for h22=1:length(Inij2)
Df220=betaX2ij(Inij2(h22),1)*yL(:,Indexij(Inij2(h22),1))+Df220;
end
Df22(:,i)=Df220;
Df310=zeros(ryL,1);
for h31=1:length(Inijk1)
Df310=betaX3ijk(Inijk1(h31),1)*yL(:,Indexijk(Inijk1(h31),2))...
.*yL(:,Indexijk(Inijk1(h31),3))+Df310;
end
Df31(:,i)=Df310;
Df320=zeros(ryL,1);
for h32=1:length(Inijk2)
Df320=betaX3ijk(Inijk2(h32),1)*yL(:,Indexijk(Inijk2(h32),1))...
.*yL(:,Indexijk(Inijk2(h32),3))+Df320;
end
Df32(:,i)=Df320;
Df330=zeros(ryL,1);
for h33=1:length(Inijk3)
Df330=betaX3ijk(Inijk3(h33),1)*yL(:,Indexijk(Inijk3(h33),1))...
.*yL(:,Indexijk(Inijk3(h33),2))+Df330;
end
Df33(:,i)=Df330;
end
Df=Df1+Df21+Df22+Df31+Df32+Df33;
%QR decomposition
M=floor(ryL/50);
Lyap=[Df(1,:);eye(maxlag-1,maxlag)]-[Df(1,:);eye(maxlag-1,maxlag)];
for bl=1:50
LAMBDA=[Df(1,:);eye(maxlag-1,maxlag)]-[Df(1,:);eye(maxlag-1,maxlag)];
Q0=eye(maxlag);
for i=(bl-1)*M+1:bl*M;
[Q,R]=qr([Df(i,:);eye(maxlag-1,maxlag)]*Q0);
LAMBDA1=logm(R);
LAMBDA=LAMBDA+real(LAMBDA1);
Q0=Q;
end
Lyap=LAMBDA+Lyap;
end
Lyap=Lyap/(M*50);
t=toc
%_________________________________END______________________________________
  댓글 수: 5
Star Strider
Star Strider 2015년 8월 22일
It is still a heuristic exercise. Use the delay that gives you the result you want. I haven’t analysed chaotic time series in a while, but I’ve not discovered any particular method for determining the best delay without simply trying different ones to see which worked best according to whatever criteria you set. This Scholarpedia article on the Lyapounov exponent may provide you with some guidance.
alireza ghavami
alireza ghavami 2023년 8월 7일
Hi dear ,
i just followed your conversation with ms.Zahra,
since Im preparing for my M.Sc thesis, is it possible to help me too?
i need to find the largest Lyapunov exponent for my EEG signal as well as other non-linear features.
could you please help me with the code?

댓글을 달려면 로그인하십시오.


alireza ghavami
alireza ghavami 2023년 8월 7일
Hi dear Zahra ,
i just found your problem by searching in google since im trying to find the way to calculate or find the largest Lyapunov exponent in my EEG signal,
could you please help me where can i find the solution?
i tried in youtube, aparat and also anywhere i thought i could find but no result!...
i consider it has long time ago you made your question but if its possible for you to help me, i will be grateful.
here's my gmail address: alirezagh2222@gmail.com
thanks alot

카테고리

Help CenterFile Exchange에서 Matrix Computations에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by