17-PSO优化-逼近最小二乘法结果（23个决策变量）

Halcom

17-PSO优化-逼近最小二乘法结果
百度网盘链接：https://pan.baidu.com/s/1KqoDNIlIXH8_0Jw1j-ghcg
具体链接在halcom.cn论坛，联系人QQ：3283892722
该论坛是一个学习交流平台，我会逐一的和大家分享学习。
欢迎大家录制视频，你可在论坛进行打赏分享。
视频专用播放器：http://halcom.cn/forum.php?mod=viewthread&tid=258&extra=page%3D1\

例如一个庞大的数组，Jj = 120x23, 目标ee = 120x1，则需要求解最小二乘拟合结果x，理论上
x = pinv(Jj)*ee
% x=inv(Jj'*Jj)*Jj'*ee

此时，x是一个23x1的未知数矩阵，如果采用生物智能算法进行计算，那么如何准确的计算这个x呢？
我们知道生物智能算法，例如pso算法，对于未知数个数较少的时候，寻优是比较靠谱的，然后未知数是23个，基本没法准确的寻优了。
（1）不假思索，直接写出适应度函数，进行PSO寻优

1. function y = fun(x)
   
2. Jj=[-300.300……0];
   
3. ee=[6.316300……000];
   
4. % Jj*x = ee
   
5. % 误差最小
   
6. y = sum( abs( Jj*x' - ee ) );

复制代码

然而我们采用2改进的粒子群算法进行计算时，误差y特别大，根本不可接受。那么此时我们该如何做呢？
（2）改进适应度函数求解：

1. function [x, s] = fun2(x, popmin, popmax)
   
2. % 上下限
   
3. lb = popmin*ones(length(x), 1);
   
4. ub = popmax*ones(length(x), 1);
   
5. 
   
6. options = optimoptions( 'fmincon', 'Algorithm', 'sqp', 'Display', 'off' );
   
7. [x, s] = fmincon( @(x)fun(x), x, [],[],[],[], lb, ub, [], options );

复制代码

PS：此方法类比于PSO优化BP神经网络，神经网络通常被优化的为权值和阈值，而且是初始的权值和阈值，那么由神经网络得到的反馈最优权值和阈值可作为第二次迭代的输入；则PSO优化-逼近最小二乘法结果通过fun2函数即可解开谜团。
附上主程序：

1. % PSO 参数
   
2. c1 = 1.4995;  
   
3. c2 = 1.4995;
   
4. Vmin = -1;
   
5. Vmax = 1;
   
6. maxiter = 10;  % 迭代次数
   
7. sizepop = 10;  % 种群数量
   
8. popmin = -1.5;  popmax = 2; % x
   
9. nvar = 23;   % 23个未知量
   
10. %% 初始化种群
    
11. for i=1:sizepop
    
12.     pop(i,:) = popmin + (popmax-popmin).*rand(1,nvar);
    
13.     [pop(i,:), fitness(i)] = fun2( pop(i,:), popmin, popmax );
    
14.     V(i,1:nvar) = Vmin + (Vmax-Vmin).*rand(1,nvar);
    
15. end
    
16. % 记录一组最优值
    
17. [bestfitness,bestindex]=min(fitness);
    
18. zbest=pop(bestindex,:);   % 全局最佳
    
19. gbest=pop;                % 个体最佳
    
20. fitnessgbest=fitness;     % 个体最佳适应度值
    
21. fitnesszbest=bestfitness; % 全局最佳适应度值
    
22. wmax = 0.9;  wmin = 0.4;  % 权重
    
23. %% 迭代寻优
    
24. for i=1:maxiter
    
25.     disp(['当前迭代次数：', num2str(i)])
    
26.     for j=1:sizepop
    
27.         % 自适应权重1
    
28. %         w = wmin + (wmax-wmin)*(fitnessgbest(j)-min(fitness))/( max(fitness)-min(fitness) );
    
29.         % 自适应权重2
    
30. %         w = wmin - (wmax-wmin)*(fitnessgbest(j)-min(fitness))/( mean(fitness)-min(fitness) );
    
31.         % 自适应权重3
    
32.         if fitnessgbest(j)<=mean(fitness)
    
33.             w = wmin - (wmax-wmin)*(fitnessgbest(j)-min(fitness))/( mean(fitness)-min(fitness) );
    
34.         else
    
35.             w = wmax;
    
36.         end
    
37.         % 速度更新
    
38.         V(j,:) = w*V(j,:) + c1*rand*(gbest(j,:) - pop(j,:)) + c2*rand*(zbest - pop(j,:));
    
39.         V(j,:) = lb_ub(V(j,:), Vmin, Vmax);
    
40.         
    
41.         % 个体更新
    
42.         pop(j,:) = pop(j,:) + 0.5 * V(j,:);
    
43.         pop(j,:) = lb_ub(pop(j,:), popmin, popmax);
    
44.         
    
45.         % 适应度更新
    
46.         [pop(j,:), fitness(j)] = fun2(pop(j,:), popmin, popmax);
    
47.         
    
48.         % 比较  个体间比较
    
49.         if fitness(j)<fitnessgbest(j)
    
50.             fitnessgbest(j) = fitness(j);
    
51.             gbest(j,:) = pop(j,:);
    
52.         end
    
53.         if fitness(j)<bestfitness
    
54.             bestfitness = fitness(j);
    
55.             zbest =  pop(j,:);
    
56.         end
    
57.         
    
58.     end
    
59.    
    
60.     fitness_iter(i) = bestfitness;
    
61. end

复制代码