马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
- %**************************************************************************
- % 子函数功能:利用C-C法求延迟时间tau和时间窗口tw
- %程序原理说明: 1.吕金虎—混沌时间序列分析及其应用—P67-69
- % 创作者:何&*(桂林电子科技大学) QQ:458689834
- % 时间: 2015-5-26-00:00:00
- %I/O口说明:
- % data: 输入时间序列
- % max_lag: 最大时间延迟
- % tau: 计算得到的延迟时间
- % tw: 时间窗口
- %**************************************************************************
- function [tau,tw]=C_CMethod_he(data,max_lag)
- M=5; %嵌入最大维数,可认为修改!
- R=4; %最大半径搜索点!
- sigma=std(data); %时间序列标准差!
- S_mean1=zeros(M-1,max_lag); %P69页公式3.31,式子右侧j从1到4(R)的和!
- S_mean2=zeros(1,max_lag); %P69页公式3.31的左边(j和m取完值的和)!
- S_delta=zeros(M-1,max_lag); %P69页公式3.32右侧的三角形S(m,t)!
- S_delta_mean=zeros(1,max_lag); %对S_delta求均值,公式3.32的左侧标识符!
- S_cor=zeros(1,max_lag); %P69(3.33)!
- for t=1:max_lag
- t
- Y=Disjoint(data,t); %将原始序列分成t个不相交的序列(P68公式3.27的代码)!
- [m0,n0]=size(Y); %m0其实在这里没有用,但是n0其实等于t的值!
- L1=floor(length(data)/t); %t个不相交序列的长度(向下取整)!
- for m=2:M %嵌入维数循环!
- B1=[]; %B1用于存储P69页公式3.30右侧大括号内的S(m,rj,t),当m和t固定时,B1是一行数据,所以可以求出最大值和最小值!
- for r=1:R %半径循环开始!
- ri=sigma*r/2; %P69页公式3.31上面那行文字和数字的代码!
- A1=0; %A1用于存储P68页公式3.28中右侧中括号内的差值(该差值有t个)!
- for t1=1:n0 %分别对不同不相交序列进行求解关联积分值!
- Y1=Y(:,t1); %取出第t1个不相交的时间序列!
- X1=Reconstitution(Y1,1,t); %对其进行嵌入维数为1,时间延迟为t的相空间重构!
- C_sm1rt=Correlation_integral(X1,L1,ri); %调用关联积分子函数求公式3.28中括号-号后面的值!
-
- X2=Reconstitution(Y1,m,t); %对其进行嵌入维数为m,时间延迟为t的相空间重构!
- L2=L1-(m-1)*t; %重构后相空间中相点的个数!
- C_smrt=Correlation_integral(X2,L2,ri); %调用关联积分子函数求公式3.28中括号-号前面式子的值!
- A1=A1+ C_smrt-C_sm1rt^m; %两个关联积分值做差!
- end
- B1=[B1 A1/t]; %r循环一圈得到3.29左边的值,这里S(m,r,t)中r变量可以消除,得到二维(m,t)矩阵,因为r被遍历了!
- end
- S_delta(m-1,t)=max(B1)-min(B1); %因为m从2开始,所以m-1;然后该式子得到3.30左边的值!一一对应!
- S_mean1(m-1,t)=sum(B1); %P69页公式3.31,式子右侧j从1到4(R)的和!
- end
- S_mean2(t)=sum(S_mean1(:,t))/(M-1)/R; %公式3.31的值!M=5,R=4!
- S_delta_mean(t)=sum(S_delta(:,t))/(M-1) %公式3.32的值!
- S_cor(t)=S_delta_mean(t)+abs(S_mean2(t)); %公式3.33的值!
- end
- figure
- plot(S_mean2,'r');hold on;
- plot(S_delta_mean,'b');hold on;
- plot(S_cor,'k-');legend('S_mean2','S_delta_mean','S_cor');grid on; %做出三条曲线!
- %--------------------------------------------------------------------------
- % 从三个统计指标中取出时间延迟—tau和时间窗宽度—tw
- %--------------------------------------------------------------------------
- tw=min(S_cor); %公式3.33下面的话!
- for t=2:max_lag
- if S_mean(t-1)*S_mean(t)<0
- tau=t;
- end
- end
- %--------------------------------------------------------------------------
- %若您发现有问题或不懂地方,请及时通知本人,发生邮件或加QQ都行!
复制代码
|