声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 2526|回复: 4

[共享资源] 一阶加权局域法,在线等

[复制链接]
发表于 2015-7-28 16:55 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
求吕金虎版的第五章一阶加权局域法的Matlab编程代码
回复
分享到:

使用道具 举报

发表于 2015-9-23 10:02 | 显示全部楼层
个人收集到两个一阶加权局域发的代码,一致没机会使用,你可以试试看是否好用

该程序用加权一阶局域法对数据进行进行一步预测
  1. %skyhawk
  2. clear all;

  3. m=6;     %嵌入维数
  4. N=80;    %预测后N个点

  5. A=load('kj.txt');
  6. P=26; % 北空的平均循环周期=26

  7. whl=A(:,4);
  8. [whsl,lll]=size(whl);  

  9. % lmd_1=lyapunov(m,m,whl,whsl);%求lyapunov指数
  10. % lmd_mm=lmd_1(m);
  11. for j=1:whsl            
  12.     whlsj(j)=whl(j);
  13. end   

  14. fch=0;
  15. for i=whsl-N+1:whsl         %预测后N个点
  16.     [lmd_m,idx,min_d,idx1,min_d1]=lyapunov(m,whlsj,i-1,P);  
  17.     [y(i),z(i)]=pre_by_lya(m,lmd_m,whlsj,i-1,idx,min_d);%预测第i+1个点  
  18.    
  19.     fch=fch+(y(i)-whl(i))*(y(i)-whl(i));
  20. %     fch=fch+(z(i)-whl(i))*(z(i)-whl(i));
  21. %     clear whlsj;

  22.     iii=whsl-i     %显示进度
  23. end

  24. fch=sqrt(fch)/N

  25. % for i=whsl-N+1:whsl
  26. %     p(i-(whsl-N+1)+1)=y(i);
  27. %     q(i-(whsl-N+1)+1)=z(i);
  28. %     w(i-(whsl-N+1)+1)=whl(i);
  29. % end

  30. % kk=1:N;
  31. % plot(kk,p,'r',kk,w)

  32. yyy=[whl,y'];
  33. save('kjyc.txt','yyy','-ASCII');

  34. kk=1:whsl;
  35. plot(kk,whl,'b',kk,y,'r')
复制代码


发表于 2015-9-23 10:03 | 显示全部楼层
该程序利用一阶局域加权法进行混沌时间序列的预测

  1. %AOLMM多步预测函数
  2. function [FChaosPredict] = FunctionChaosPredict(Data,N,mtbp,deltaT,tao,d,MaxStep)
  3. %Data是一维信号时间序列,N是信号数据长度,mtbp,deltaT,tao,d分别是重构相空间的平均时间序列、采样周期、时延及嵌入维
  4. roll=Data;%取横摇数据
  5. M = N - (d - 1)*tao;
  6. for i = 1 : M
  7.     for j = 1 : d
  8.         MatrixX(i,j) = roll(i + (j - 1)*tao);
  9.     end
  10. end
  11. %计算相空间中第M点与各点的距离
  12. for j = 1 : (M - 1)
  13.          Dis(j) = norm(MatrixX(M,:) - MatrixX(j,:),2);
  14. end
  15. %排序计算相空间中第M点的(m+1)个参考邻近点
  16. for i = 1 : (d + 1)
  17.     NearDis(i) = Dis(i);
  18.     NearPos(i) = i;
  19. end
  20. for i = (d + 2) : (M - mtbp)
  21.     for j = 1 : (d + 1)
  22.         if (abs(i-j)>mtbp) %& abs(i-j)<10*mtbp  
  23.         if(Dis(i) < NearDis(j))
  24.             NearDis(j) = Dis(i);
  25.             NearPos(j) = i;
  26.             break;
  27.         end
  28.         end
  29.     end
  30. end
  31. SortedDis = sort(NearDis);
  32. MinDis = SortedDis(1);
  33. %计算第M点的(m+1)个参考邻近点的权P[i]
  34. SumP = 0;
  35. for i = 1 : (d + 1)
  36.     P(i) = exp(-NearDis(i)/MinDis);
  37.     SumP = SumP + P(i);
  38. end
  39. P = P/SumP;
  40. %用最小二乘法计算a[],b[]
  41. for step=1:1:MaxStep
  42.     aCoe1 = 0;
  43.     aCoe2 = d;
  44.     bCoe1 = 0;
  45.     bCoe2 = 0;
  46.     e = 0;
  47.     f = 0;
  48.     for i = 1 : (d + 1)
  49.         aCoe1 = aCoe1 + P(i)*sum(MatrixX(NearPos(i),:));
  50.         bCoe1 = bCoe1 + P(i)*(MatrixX(NearPos(i),:)*MatrixX(NearPos(i),:)');
  51.         e = e + P(i)*(MatrixX(NearPos(i) + step,:)*MatrixX(NearPos(i),:)');
  52.         f = f + P(i)*sum(MatrixX(NearPos(i) + step,:));
  53.     end
  54.     bCoe2 = aCoe1;
  55.     CoeMatrix = [aCoe1,bCoe1;aCoe2,bCoe2];
  56.     ResultMatrix = [e;f];
  57.     abResult = pinv(CoeMatrix)*ResultMatrix;
  58.     a = abResult(1);
  59.     b = abResult(2);
  60. for j = 1 : d
  61. %     MatrixX(M + step,j) = a + b*MatrixX(M,j); %以历史上相近点的演化规律作为中心点的演化规律以中心点为基准进行预报
  62.      
  63.     MatrixX(M + step,j) = 0;
  64.     for i = 1 : (d + 1)
  65.         MatrixX(M + step,j) = MatrixX(M + step,j) + P(i)*(a + b*MatrixX(NearPos(i),j)); %以历史上相近点的演化加权和直接作为中心点的演化点进行预报
  66.     end
  67. end
  68. %误差修正
  69.     if M-tao+step+(d-1)*tao < N+1
  70.         for j=1:d-1
  71.             err(j)=MatrixX(M + step,j)-roll(M+step+(j-1)*tao);
  72.         end
  73.         ppp=1:d-1;ttt=err;neterr=newrbe(ppp,ttt);xxx=2:d;errp=sim(neterr,xxx);
  74.         PredictedData(step) = MatrixX(M + step,d) - errp(d-1);
  75.         roll(N+step)=PredictedData(step);
  76. else
  77.     PredictedData(step) = MatrixX(M + step,d);
  78.    end%        roll(N+k)=PredictedData(k);
  79. FChaosPredict(step) = PredictedData(step);
  80. % FChaosPredict(step) = MatrixX(M + step,d);
  81. end
复制代码
发表于 2015-9-23 10:03 | 显示全部楼层
在附送一个加权零阶局域法的代码

  1. function Forcast=WZOLL1(data)
  2. % Weighted zero-order local law,加权零阶局域预测法
  3. % data为重构后的时间序列
  4. % Xi中为临域的值
  5. [row,col]=size(data);
  6. % 确定临域中的各个点到基准点的欧氏距离
  7. Xk=data(row,:);%Xk为中心点
  8. Distance=zeros(1,row-1);
  9. for i=1:row-1
  10. Mid=data(i,:);
  11. Minus=Xk-Mid;
  12. D=0;
  13. for j=1:col
  14. D=D+Minus(j)^2;
  15. end
  16. Distance(i)=sqrt(D);
  17. end
  18. % 找出欧式距离最小的M个点作为中心点的最近临域其值放在Xi中
  19. M=2*col+1;
  20. Xi=zeros(M,col);
  21. % 确定Distance中的最小值
  22. dm=min(Distance);
  23. % D中为距离值
  24. D=zeros(1,M);
  25. for i=1:M
  26. Midmin=min(Distance);
  27. D(i)=Midmin;
  28. [~,d]=find(Distance==Midmin);
  29. d=d(1);
  30. Xi(i,:)=data(d,:);
  31. Distance(d)=nan;
  32. end
  33. %计算权重P
  34. p=zeros(1,M);
  35. for i=1:M
  36. fenzi=exp(dm-D(i));
  37. fenmu=0;
  38. for j=1:M
  39. fenmu=fenmu+exp(dm-D(j));
  40. end
  41. p(i)=fenzi/fenmu;
  42. end
  43. Forcast=zeros(1,col);
  44. for i=1:M-1
  45. Forcast=Forcast+p(i)*Xi(i+1,:);
  46. end
  47. % 处理奇异值 和零阶局域法计算的结果对比
  48. % 计算加权零阶局域法预测出来的值
  49. Test=mean(Xi,1);
  50. if Test(col)/2>Forcast(col)
  51. Forcast=mean(data(row-2:row,:));
  52. end
  53. end
复制代码
 楼主| 发表于 2015-9-23 12:50 | 显示全部楼层
happy 发表于 2015-9-23 10:03
在附送一个加权零阶局域法的代码

太感谢了,你会万福的!
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-10 18:30 , Processed in 0.063606 second(s), 17 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表