一阶加权局域法，在线等

混蛋理论

求吕金虎版的第五章一阶加权局域法的Matlab编程代码

happy

个人收集到两个一阶加权局域发的代码，一致没机会使用，你可以试试看是否好用

该程序用加权一阶局域法对数据进行进行一步预测

1. %skyhawk
   
2. clear all;
   
3. 
   
4. m=6;     %嵌入维数
   
5. N=80;    %预测后N个点
   
6. 
   
7. A=load('kj.txt');
   
8. P=26; % 北空的平均循环周期＝26
   
9. 
   
10. whl=A(:,4);
    
11. [whsl,lll]=size(whl);  
    
12. 
    
13. % lmd_1=lyapunov(m,m,whl,whsl);%求lyapunov指数
    
14. % lmd_mm=lmd_1(m);
    
15. for j=1:whsl            
    
16.     whlsj(j)=whl(j);
    
17. end   
    
18. 
    
19. fch=0;
    
20. for i=whsl-N+1:whsl         %预测后N个点
    
21.     [lmd_m,idx,min_d,idx1,min_d1]=lyapunov(m,whlsj,i-1,P);  
    
22.     [y(i),z(i)]=pre_by_lya(m,lmd_m,whlsj,i-1,idx,min_d);%预测第i+1个点  
    
23.    
    
24.     fch=fch+(y(i)-whl(i))*(y(i)-whl(i));
    
25. %     fch=fch+(z(i)-whl(i))*(z(i)-whl(i));
    
26. %     clear whlsj;
    
27. 
    
28.     iii=whsl-i     %显示进度
    
29. end
    
30. 
    
31. fch=sqrt(fch)/N
    
32. 
    
33. % for i=whsl-N+1:whsl
    
34. %     p(i-(whsl-N+1)+1)=y(i);
    
35. %     q(i-(whsl-N+1)+1)=z(i);
    
36. %     w(i-(whsl-N+1)+1)=whl(i);
    
37. % end
    
38. 
    
39. % kk=1:N;
    
40. % plot(kk,p,'r',kk,w)
    
41. 
    
42. yyy=[whl,y'];
    
43. save('kjyc.txt','yyy','-ASCII');
    
44. 
    
45. kk=1:whsl;
    
46. plot(kk,whl,'b',kk,y,'r')

复制代码

happy

该程序利用一阶局域加权法进行混沌时间序列的预测

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

复制代码

happy

在附送一个加权零阶局域法的代码

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

复制代码

混蛋理论

happy 发表于 2015-9-23 10:03
在附送一个加权零阶局域法的代码

太感谢了，你会万福的！