求救：哪位大侠有改进的CC算法，感激不尽
最近用吕金虎那本书里介绍的C-C算法求嵌入维和时间延迟，然后代入G-P算法中求关联维，并用Chaos Toolbox Ver.2.0工具箱里的C-C、G-P程序仿真运行了下，结果始终不尽如意。我用logistics混沌序列进行验证，具体过程如下：首先当u=4时，我产生40000个混沌序列点，然后代入C-C算法，画出了三条曲线图如下：file:///C:/DOCUME~1/ME/LOCALS~1/Temp/msohtml1/01/clip_image002.gif
蓝线----smean(t)
绿点---delt_s
红*-----s_cor 由图求的delt_s第一个极小值点：4
s_cor最小值点在9
即：Tw=9,tau=4.
然后我代入G_P算法算出m=2:10下的ln_C~ln_r，如图：
file:///C:/DOCUME~1/ME/LOCALS~1/Temp/msohtml1/01/clip_image002.gif
但是画出的不同嵌入维下ln_C~ln_r的斜率始终不趋于稳定。
并且下面几条曲线用最小二乘法拟合的时候，输出均为Inf，不知道为什么。
起初以为数据格式不够，我也试着改了下，还是不好使。后来看到一篇论文里用改进的C-C算法，就是在关联积分计算过程中引入了权衡计算精度与速度的可调参数。效果还不错，不知哪位大侠有没有相关程序，发给小妹一份，不胜感激。
还有在论坛上看到用fft求最佳延迟时间，如果哪位大侠有程序也发给小妹看一下。谢谢，谢谢啊~~~
我的邮zhl20062159@126.com
logistic 序列产生程序：
clear all;
N=40000;
% initial conditions:
x(1)=0.2;
% parameter:
% r=3.8;% chaos
% r=2.0;% steady state
% r=3.2; % period 2
r=4; % period 4?
for n=1:N
x(n+1)=r*x(n)*(1-x(n));
end
plot(x,'k')
title('logistic map')
xlabel('n')
ylabel('x(n)')
figure
plot(x(2:N),x(1:N-1),'k.','markersize',5)
title('logistic map')
xlabel('x(n+1)')
ylabel('x(n)')
fid=fopen('logistic.txt','w');
fprintf(fid,'%f\n',x);
fclose(fid);
C-C算法：
function =C_CMethod(data)
%this function calculate time delay and embedding demension with C-C
% 该函数用C-C算法来计算延迟时间和嵌入维数
%Method,which proved by H.S.Kim
%skyhawk&flyinghawk
% %****************调试程序段****************************
clear all;
data1=load('E:\chaos\4.15\Chaos Toolbox Ver.2.0_zhu\Chaos Toolbox Ver.2.0\Main\logistic.txt');
data=data1(1:20:length(data1)); %%数据采样
% max(data)
% min(data)
% for i=1:length(data) %%数据的归一化
% data(i)=(data(i)-min(data))/(max(data)-min(data));
% end
% %************************************************
N=length(data) % 数据向量长度
max_d=20; %the maximum value of the time delay 时间延迟的最大值
sigma=std(data) %calcute standard deviation s_d 计算标准偏差s_d
for t=1:max_d
t
s_t=0;
delt_s_s=0;
for m=2:5
m
s_t1=0;
for j=1:4%zhu:对r的循环
j
r=sigma*j/2;
data_d=disjoint(data,N,t); %将时间序列分解成t个不相交的时间序列
=size(data_d)
s_t3=0;
for i=1:t %zhu:对t个序列分别计算关联积分C
i
Y=data_d(i,:); %zhu:取data_d中的第i行，即第i个子序列
C_1(i)=correlation_integral(Y,N_d,r); %计算C(1,N_d,r,t)zhu:即长度为N_d的一维序列
X=reconstitution(Y,N_d,m,t); %相空间重构，zhu:新的序列为总长度为N_d产生的m维序列
N_r=N_d-(m-1)*t; %zhu:t=tao 新的序列个数为N_r
C_I(i)=correlation_integral(X,N_r,r); %计算C(m,N_r,r,t)
s_t3=s_t3+(C_I(i)-C_1(i)^m); %对t个不相关的时间序列求和
end
s_t2(j)=s_t3/t; %计算S(m,N,r,t)
s_t1=s_t1+s_t2(j); %对rj求和
end
delt_s_m(m)=max(s_t2)-min(s_t2); %求delt S(m,t)
delt_s_s=delt_s_s+delt_s_m(m); %delt S(m,t)对m求和
s_t0(m)=s_t1;
s_t=s_t+s_t0(m); %S对m求和
end
s(t)=s_t/16;
delt_s(t)=delt_s_s/4;
s_cor(t)=delt_s(t)+abs(s(t));
end
t=1:max_d;
plot(t,s,t,delt_s,'.',t,s_cor,'*')
fid=fopen('s(t).txt','w');
fprintf(fid,'%f\n',s(t));
fclose(fid);
fid=fopen('delt_s(t).txt','w');
fprintf(fid,'%f\n',delt_s(t));
fclose(fid);
fid=fopen('s_cor(t).txt','w');
fprintf(fid,'%f\n',s_cor(t));
fclose(fid);
% 寻找时间延迟tau：即Sdeltmean第一个极小值点对应的t
for i=2:length(delt_s)-1
if delt_s(i)<delt_s(i-1)&delt_s(i)<delt_s(i+1)
tau=i;
break;
end
end
% 寻找时间延迟tau：即Sdeltmean第一个极小值点对应的t
% 寻找时间窗口tw：即Scor最小值对应的t
for i=1:length(s_cor)
if s_cor(i)==min(s_cor)
tw=i;
break;
end
end
% 寻找时间窗口tw：即Scor最小值对应的t
G-P算法：
function =G_P(data,N,tau,min_m,max_m,ss)
% the function is used to calculate correlation dimention with G-P algorithm
% 计算关联维数的G－P算法
% data:the time series 时间序列
% N: the length of the time series 时间序列长度
% tau: the time delay 时间延迟
% min_m:the least embedded dimention m 最小的嵌入维数
% max_m:the largest embedded dimention m 最大的嵌入维数
% ss:the stepsize of r r的步长 !!!我的理解—ss—r的循环数目!!!
%skyhawk
%%%%%%%%%%%zhu test %%%%%%%%%%%%%%%%%%%%%
clear all;
data1=load('E:\chaos\4.15\Chaos Toolbox Ver.2.0_zhu\Chaos Toolbox Ver.2.0\Main\logistic.txt');
data=data1(1:20:length(data1));
% max(data)
% min(data)
% for i=1:length(data) %%数据的归一化
% data(i)=(data(i)-min(data))/(max(data)-min(data));
% end
min_m=2;
max_m=10;
N=length(data)
tau=4;
ss=20;
%%%%%%%%%zhu test %%%%%%%%%%%%%%%%%%%%%
for m=min_m:max_m
m
Y=reconstitution(data,N,m,tau); %reconstitute state space
M=N-(m-1)*tau; %the number of points in state space
for i=1:M-1
% i
for j=i+1:M
d(i,j)=max(abs(Y(:,i)-Y(:,j)));%calculate the distance of each two
end %points in state space计算状态空间中每两点之间的距离
end
max_d=max(max(d)) %the max distance of all points 得到所有点之间的最大距离
d(1,1)=max_d;
min_d=min(min(d)) %the min distance of all points 得到所有点间的最短距离
delt=(1.2*max_d-min_d)/ss; %the stepsize of r 得到r的步长
for k=1:ss
k
r=min_d+k*delt;
C(k)=correlation_integral(Y,M,r); %calculate the correlation integral
ln_C(m,k)=log(C(k)); %lnC(r)
ln_r(m,k)=log(r); %lnr
fprintf('%d/%d/%d/%d\n',k,ss,m,max_m);
fprintf('%d/%d/%d\n',k,ss,m);
end
plot(ln_r(m,:),ln_C(m,:));
xlabel('ln_r');
ylabel('ln_C(r)');
hold on;
end
%%%%%%%%%%%%%%%%%%% 画 Correlation dimensional versus embedding dimension.
figure
S=[];ii=1;
for mm=min_m:max_m
x=1:floor(length(ln_C(mm,:)));
pp=polyfit(ln_r(mm,x),ln_C(mm,x),1) %%曲线拟合，最小二乘多项式拟合，n=1，为直线拟合。
S(ii)=pp(1)%%% 返回斜率。
ii=ii+1;
end
MM=min_m:max_m
plot(MM,S,'*');
xlabel('嵌入维m');
ylabel('关联维D');
%%%%%%%%%%%%%%%%%%% 画 Correlation dimensional versus embedding dimension.
fid=fopen('lnr.txt','w');
fprintf(fid,'%f %f\n',ln_r);
fclose(fid);
fid = fopen('lnC.txt','w');
fprintf(fid,'%f %f\n',ln_C);
fclose(fid);
页:
[1]