% G-P 算法求關聯維(輸入時間序列數據)% 使用平臺 - Matlab6.5 / Matlab7.0% 作者：陸振波，海軍工程大學% 歡迎同行來信交流與合作，更多文章與程序下載請訪問我的個人主頁% 電子郵件：luzhenbo@yahoo.com.cn% 個人主頁：http://luzhenbo.88uu.com.cnclcclear allclose all%--------------------------------------------------------------------------% 產生 Lorenz 時間序列% dx/dt = sigma*(y-x)% dy/dt = r*x - y - x*z% dz/dt = -b*z + x*ysigma = 16;             % Lorenz方程參數r = 45.92;               b = 4;          y = [-1;0;1];           % 起始點 (3x1 的列向量)h = 0.01;               % 積分時間步長k1 = 30000;             % 前面的迭代點數k2 = 5000;              % 后面的迭代點數X = LorenzData(y,h,k1+k2,sigma,r,b);X = X(k1+1:end,1);      % 時間序列(列向量)%--------------------------------------------------------------------------% G-P算法計算關聯維rr = 0.5;Log2R = -6:rr:0;        % log2(r)R = 2.^(Log2R);t = 10;                 % 時延dd = 1;                 % 嵌入維間隔D = 2:dd:10;            % 嵌入維    p = 50;                 % 限制短暫分離，大于序列平均周期(不考慮該因素時 p = 1)   tic Log2Cr = log2(CorrelationIntegral(X,t,D,R,p));       % 輸出每一行對應一個嵌入維toc%--------------------------------------------------------------------------% 結果作圖figureplot(Log2R,Log2Cr','k.-'); axis tight; grid on; hold on;xlabel('log2(r)'); ylabel('log2(C(r))');title(['Lorenz, length = ',num2str(k2)]);%--------------------------------------------------------------------------% 最小二乘擬合Linear = [3:9];                       % 線性似合區域[A,B] = LM1(Log2R,Log2Cr,Linear);     % 最小二乘求斜率for i = 1:length(D)    Y = polyval([A(i),B(i)],Log2R(Linear),1);    plot(Log2R(Linear),Y,'r');    endhold off;%--------------------------------------------------------------------------% 求梯度Slope = diff(Log2Cr,1,2)/rr;           % 求梯度xSlope = Log2R(1:end-1);               % 梯度所對應的log2(r)figure;plot(xSlope,Slope','k.-'); axis tight; grid on;xlabel('log2(r)'); ylabel('slope');title(['Lorenz, length = ',num2str(k2)]);%--------------------------------------------------------------------------% 關聯維曲線figure;plot(D,A,'k.-'); grid on; axis tight;xlabel('m'); ylabel('Correlation Dimension');title(['Lorenz, length = ',num2str(k2)]);