%------------------------------主程序
%粒子群優化算法在參數最小二乘估計中的應用
%陳毓靈
%Post2Ling@hotmail.com
%function [Pgd,Fgd] = PSO(Pop,nIter,Wmax)
%function PSO
%points=50; %樣本點數
clear;
vect=3; %待估參數向量維數

Pop=50; %種群大小
Vmax=20; %最大速度
Vmin=-20; %最小速度
nIter=150; %迭代次數
Wmax=1.2; %最大慣性因子
Wmin=0.0; %最小慣性因子
Wn=1.0; %慣性因子衰減指數
FaiI=1.2; %自身學習率
FaiG=1.2; %社會學習率
global Pgd;
global Fgd;
%輸入樣本 X Y
%X=51+rand(50,1)*50;
%Y = 100+(200./(1+exp(0.3.*(X-75))))+(rand(50,1)-0.5)*20;


%粒子群
%當前位置vec維 個體歷史最優位置vec維 全局最優位置vec維 當前值 個體歷史最優值 全局最優值
X = (rand(Pop,vect) - 0.5).*20;             %當前位置
%Pid = (rand(Pop,vect) - 0.5).*20;           %個體歷史最優
Pid = X;
Pgd = Pid(1,:);                            %全局最優


% 更新適應值
F = Fitness(X); %當前值
Fid = F; %個體歷史最優值
Fgd = F(1); %全局最優值
% 更新全局最優位置和適應值
for j=1:Pop
    if Fid(j)<Fgd
        Fgd=F(j);
        Pgd = Pid(j,:);
    end
end

%粒子群速度
V=(rand(Pop,vect)-0.5)*10;

for i=1:nIter

    % 更新慣性因子
    W=((nIter-i)^Wn)/(nIter^Wn)*(Wmax-Wmin)+Wmin;

    % 更新速度
    V=W.*V + FaiI.*rand.*ones(Pop,vect)*(Pid-X)+FaiG.*rand.*ones(Pop,vect)*(ones(Pop,1)*(Pgd-X));
    % 限制速度
    V(find(V>Vmax))=Vmax;
    V(find(V<Vmin))=Vmin;

    % 更新位置
    X=X + V;

    % 更新適應值
    F = Fitness(X);
    %P(:,vec*3+1)=Fitness(P(:,1:vec),X,Y);

    % 更新個體最優位置和適應值
    for j=1:Pop
        if F(j) < Fid(j);
            Fid(j)=F(j);
            Pid(j,:) = X(j,:);
        end
    end

    % 更新全局最優位置和適應值
    for j=1:Pop
        if Fid(j)<Fgd
            Fgd=Fid(j);
            Pgd = Pid(j,:);
        end
    end
Rec(i) = Fgd
end
plot(Rec);
%輸出結果
% [P(1,vec*2+1:vec*3),P(1,vec*3+3)]
%[P(1,vec*2+1),P(1,vec*2+2),P(1,vec*2+3)/300,P(1,vec*2+4)]

%BP=51:100;
%HR = P(1,vec*2+1)+(P(1,vec*2+2)./(1+exp((P(1,vec*2+3)./300).*(BP-P(1,vec*2+4)))));%+(rand(1,50)-0.5)*10;

%h=plot(X,Y,'r*',BP,HR,'b');
%xlabel('x')
%ylabel('Logistic')
%title('4參數Logistic函數圖像與樣本比較 (p1=100;p2=200;p3=0.3;p4=75)')
%legend('樣本','估計函數')