国产无码综合区,色欲AV无码国产永久播放,无码天堂亚洲国产AV,国产日韩欧美女同一区二区

【多目標進化優(yōu)化】MOPSO 原理與代碼實現(xiàn)

2年前作者：青年有志分類：Toy博客閱讀(26)違法舉報

這篇具有很好參考價值的文章主要介紹了【多目標進化優(yōu)化】MOPSO 原理與代碼實現(xiàn)。希望對大家有所幫助。如果存在錯誤或未考慮完全的地方，請大家不吝賜教，您也可以點擊"舉報違法"按鈕提交疑問。

?? 博主相信： 有足夠的積累，并且一直在路上，就有無限的可能！??！

?????個人主頁： 青年有志的博客

?? Gitee 源碼地址： https://gitee.com/futurelqh/Multi-objective-evolutionary-optimization

前言

前驅(qū)知識

粒子群優(yōu)化算法 PSO：https://blog.csdn.net/qq_46450354/article/details/127464089
Pareto 最優(yōu)解集： https://blog.csdn.net/qq_46450354/article/details/127917026

粒子群優(yōu)化算法 PSO

pbest： 粒子本身經(jīng)歷過的最優(yōu)位置

gbest： 粒子群整體經(jīng)歷過的最優(yōu)位置

算法思路： 在單目標優(yōu)化中，通過自我認知 pbest 以及社會認知 gbest 來計算出每個個體下一次移動的速度，從而更新個體所處的位置 x，并更新個體 pbest，以及群體的 gbest，更新的依據(jù)直接通過比較目標函數(shù)值 function value 的大小即可。

mopso,多目標進化優(yōu)化,python,開發(fā)語言

mopso,多目標進化優(yōu)化,python,開發(fā)語言

PSO 算法流程圖

mopso,多目標進化優(yōu)化,python,開發(fā)語言

PSO 引出 MOPSO

PSO 在用于單目標優(yōu)化過程中，由于只有一個函數(shù)，兩個個體之間可以通過直接比較大小來判斷好壞，從而判斷是否更新 pbest，gbest。
在多目標問題中存在多個目標值，比如求兩個目標函數(shù)的最小值，有兩個個體 ${P_1} = (8, 15), P_2 = (10, 12)$ ， ${P_1}$ 的第一個目標函數(shù)值小于 ${P_2}$ ，而第二個目標函數(shù)值大于 $P_2$ ，如何判斷誰更好呢？

問題：

如何更新多目標中的 pbest ？
如何更新多目標中的 gbest ？
對于多目標優(yōu)化問題，輸出的結果為一個集合，即 rep (非支配集/歸檔集) ，此集合的大小如何維護？
并且 PSO 很容易陷入局部最優(yōu)，如何避免呢？

注： rep 集，在部分文獻中會稱為 Archive 集，均為同一種，即當前迭代中的的非支配集，也稱為歸檔集，為避免混淆，本文同一使用 rep 。

MOPSO 所解決的正是這些個問題，從而將單目標 PSO，變?yōu)榱四軌蚋咝蠼舛嗄繕藘?yōu)化的算法。下面我們來看看這些問題是如何解決的

1. 如何更新多目標中的 pbest ？

分為三種情況進行更新：

① 當新產(chǎn)生的個體 i 在每個目標函數(shù)值上都比 pbest 更優(yōu)越，則更新 pbest （i 支配 pbest）
② 當新產(chǎn)生的個體 i 在每個目標函數(shù)值上都比 pbest 更劣，則不更新（pbest 支配 i）
③ 當新產(chǎn)生的個體 i 在部分目標函數(shù)值上比 pbest 優(yōu)越，則隨機按照一定的概率進行更新（i 與 pbest 互不支配）

2. 如何更新多目標中的 gbest ？

如何在所有互不支配的個體當中選出這個 leader 個體呢？這就需要引入自適應網(wǎng)格算法了

自適應網(wǎng)格算法

算法思路： 將所有個體劃分在不同的網(wǎng)格里面，再計算網(wǎng)格的密度，密度越小的網(wǎng)格，粒子越稀疏，個體被選擇的概率越大（從而保證粒子的分布性）。具體步驟如下：

mopso,多目標進化優(yōu)化,python,開發(fā)語言

Step1： 遍歷集合種群中的所有個體，分別在每一個目標函數(shù)上找出所有個體中的最小值和最大值，如圖所示，獲得目標函數(shù) $f_1$ 上的最小值 $min_1$ 、最大值 $max_1$ ，以及目標函數(shù) $f_2$ 上的最小值 $min_2$ 、最大值 $max_2$ 。為了保證邊界粒子 $x_1、x_5$ 也能在網(wǎng)格內(nèi)，需要進行延長操作，把最大最小值根據(jù)一定比例稍微延長一點，最后獲得網(wǎng)格的邊界 $LB_1、UB_1、LB_2、UB_2$

mopso,多目標進化優(yōu)化,python,開發(fā)語言

Step2： 假設設置的參數(shù) n = 3，將網(wǎng)格平均切割成 3 X 3 = 9 個格子

mopso,多目標進化優(yōu)化,python,開發(fā)語言

Step3： 摘選出其中有個體的網(wǎng)格，并計算這些網(wǎng)格中存儲的個體總數(shù)，如：① = 1，② = 1，③ = 1，④ = 2

mopso,多目標進化優(yōu)化,python,開發(fā)語言

Step4： 由前面獲得的個體數(shù)根據(jù)公式計算這些網(wǎng)格的被選概率，個體越少被選概率越大，歸一化所有網(wǎng)格的被選概率和為 1。

Step5： 根據(jù)前面計算的概率，采取輪盤賭方法選擇一個網(wǎng)格，該網(wǎng)格中的個體作為 gBest 備選集，在備選集中采用隨機選擇的方法選擇其中一個個體作為全局最優(yōu) gBest

3. 如何維護 rep （非支配集/歸檔集) ？

同樣可以采取自適應網(wǎng)格算法。即假設在加入新的非支配個體后，rep 庫現(xiàn)在的長度是 $m (m > N)$ ，此時我們需要刪除 $m ? N$ 個個體來保證 rep 不溢出。此時只需要通過自適應網(wǎng)格法從 rep 庫中選擇出一個最差勁個體然后從rep庫中刪除，如此循環(huán) $m ? N$ 次，問題迎刃而解

MOPSO 進行了三輪篩選

首先，根據(jù)支配關系進行第一輪篩選，將個體自身與上一次相比的劣解去除，即更新 pbest
將所有個體根據(jù)支配關系進行第二輪篩選，將劣解去除，并加入到 rep 歸檔集中，并計算個體所處的網(wǎng)格位置
最后，若歸檔集數(shù)量超過了存檔閥值，則根據(jù)自適應網(wǎng)格進行篩選刪除，直到閥值限額為止，重新進行網(wǎng)格劃分。

4. 如何避免算法陷入局部最優(yōu)？

Step1： 對當前的粒子 x，根據(jù)當前迭代次數(shù) It 和最大迭代次數(shù) MaxIt 以及突變率 m 計算擾亂算子 p。

$-\frac{ ( It - 1 ) }{ ( MaxIt - 1 )} )^ {( 1 / m )}$

Step2： 獲取隨機數(shù) rand，如果 rand < p，隨機選擇粒子的某個決策變量進行變異，其余決策變量不變；若 rand 不小于 p，則保持不變。

Step3： 如果粒子需要變異，首先隨機選擇到粒子 x 的第 j 個決策變量 x(j)，然后根據(jù)擾亂算子計算變異范圍 d，

$d = p ? (Va r M a x ? Va r M in);$

式中， $Va r M a x$ 和 $Va r min$ 分別是規(guī)定的決策變量的最大值和最小值。

Step4： 由 d 計算變異的上下界，上界 $U B = x (j) + d$ ，下界 $L B = x (j) ? d$ 。同時注意越界處理：
$UB=min(UB,Var_{max}),LB=max(LB,Var_{min})$

Step5： 最后根據(jù)變異的范圍 $[L B, U B]$ 隨機獲取新的 $x (j) = u ni f r n d (l b, u b)$ ，即獲得區(qū)間 $(l b, u b)$ 中的一個隨機數(shù)。

Step6： 在粒子變異后，要比較變異后的粒子是否更優(yōu)秀，如果變異后的粒子更優(yōu)秀則更新粒子，否則不更新；同樣還要比較新的粒子和粒子的 pBest，若新的粒子比 pBest 優(yōu)秀則更新 pBest，否則不更新。

算法流程

mopso,多目標進化優(yōu)化,python,開發(fā)語言

綜上可知 MOPSO 算法主要是針對以下一些問題進行研究：

如何選擇 pbest。對于多目標來說兩個粒子的對比，并不能對比出哪個好一些。如果粒子的每個目標都要好的話，則該粒子更優(yōu)。若有些更好，有些更差的話，就無法嚴格的說哪個好些，哪個差一些。
如何選擇 gbest。對于多目標來說，最優(yōu)的個體有很多個。該如何選擇，涉及到最優(yōu)個體的存檔、存檔的管理等問題（多樣性/分布性）
速度更新公式優(yōu)化；
位置更新處理；
增加擾亂算子即變異算子的選擇問題，主要是為了解決 PSO 快速收斂到局部最優(yōu)的問題，增加擾亂可以使得收斂到全局最優(yōu)
增加一些其他的操作來改善收斂性和多樣性

補充

上述解決方案為較簡單的方式，還可以有其他很多方式來解決這三類問題，從而實現(xiàn)不同的算法。

MOPSO 做出的改進方向如下：

速度更新公式的優(yōu)化，如：引入了一個收縮因子等
解決算法快速收斂陷入局部最優(yōu)的問題，如：變異操作
算法收斂性和多樣性的實現(xiàn)方式，如本文的自適應網(wǎng)格算法；
MOPSO 的存檔方法

算法執(zhí)行圖

算法在 ZDT1 測試函數(shù)下的執(zhí)行結果

mopso,多目標進化優(yōu)化,python,開發(fā)語言

利用綜合指標 （HV、IGD） 測試算法在測試函數(shù)上的收斂性與分布性，在迭代完成后加上下面兩段代碼即可：

zdt1.mat：為測試函數(shù) zdt1 的真是 Pareto 解集，在文章開頭給出的 Gitee 地址中有

load('zdt1.mat'); % 真實 Pareto 解集
HVResult = HV([rep.Cost]', zdt1)
IGDResult = IGD([rep.Cost]', zdt1)

mopso,多目標進化優(yōu)化,python,開發(fā)語言

代碼實現(xiàn)

mopso.mlx

Problem Definition

CostFunction = @(x) ZDT(x);      % Cost Function

nVar = 5;             % Number of Decision Variables

VarSize = [1 nVar];   % Size of Decision Variables Matrix

VarMin = 0;          % Lower Bound of Variables
VarMax = 1;          % Upper Bound of Variables

MOPSO Parameters

MaxIt = 200;           % Maximum Number of Iterations

nPop = 200;            % Population Size

nRep = 100;            % Repository Size

w = 0.5;              % Inertia Weight
wdamp = 0.99;         % Intertia Weight Damping Rate
c1 = 1;               % Personal Learning Coefficient
c2 = 2;               % Global Learning Coefficient

nGrid = 7;            % Number of Grids per Dimension
alpha = 0.1;          % Inflation Rate

beta = 2;             % Leader Selection Pressure
gamma = 2;            % Deletion Selection Pressure

mu = 0.1;             % Mutation Rate

Initialization

empty_particle.Position = [];
empty_particle.Velocity = [];
empty_particle.Cost = [];
empty_particle.Best.Position = [];
empty_particle.Best.Cost = [];
empty_particle.IsDominated = [];
empty_particle.GridIndex = [];
empty_particle.GridSubIndex = [];

pop = repmat(empty_particle, nPop, 1);

for i = 1:nPop
    
    pop(i).Position = unifrnd(VarMin, VarMax, VarSize); % 產(chǎn)生范圍內(nèi)的連續(xù)隨機數(shù)
    
    pop(i).Velocity = zeros(VarSize);
    
    pop(i).Cost = CostFunction(pop(i).Position);
    
    
    % Update Personal Best
    pop(i).Best.Position = pop(i).Position;
    pop(i).Best.Cost = pop(i).Cost;
    
end

% Determine Domination
pop = DetermineDomination(pop);

rep = pop(~[pop.IsDominated]); % 選出當前非支配個體

Grid = CreateGrid(rep, nGrid, alpha); % 定義網(wǎng)格

for i = 1:numel(rep)
    rep(i) = FindGridIndex(rep(i), Grid); % 給每個個體分配網(wǎng)格位置
end

MOPSO Main Loop

for it = 1:MaxIt
    
    for i = 1:nPop
        
        leader = SelectLeader(rep, beta); % 輪盤賭選擇 gbest
        
        pop(i).Velocity = w*pop(i).Velocity ...
            +c1*rand(VarSize).*(pop(i).Best.Position-pop(i).Position) ...
            +c2*rand(VarSize).*(leader.Position-pop(i).Position);
        
        % 更新位置
        pop(i).Position = pop(i).Position + pop(i).Velocity;
        
        % 越界判斷
        pop(i).Position = max(pop(i).Position, VarMin);
        pop(i).Position = min(pop(i).Position, VarMax);
        
        % 計算新的目標值
        pop(i).Cost = CostFunction(pop(i).Position);
        
        if Dominates(pop(i), pop(i).Best) % 如果 i 支配 Best，則更新 Best
            pop(i).Best.Position = pop(i).Position;
            pop(i).Best.Cost = pop(i).Cost;
            
        elseif Dominates(pop(i).Best, pop(i)) % 若 Best 支配 i ，則不更新
            % Do Nothing
            
        else % 若互不支配，則產(chǎn)生隨機數(shù)判斷是否更新
            if rand<0.5
                pop(i).Best.Position = pop(i).Position;
                pop(i).Best.Cost = pop(i).Cost;
            end
        end
        
    end
    
    % Add Non-Dominated Particles to REPOSITORY
    rep = [rep
         pop(~[pop.IsDominated])]; %#ok
    
    % Determine Domination of New Resository Members
    rep = DetermineDomination(rep);
    
    % Keep only Non-Dminated Memebrs in the Repository
    rep = rep(~[rep.IsDominated]);
    
    % Update Grid
    Grid = CreateGrid(rep, nGrid, alpha);

    % Update Grid Indices
    for i = 1:numel(rep)
        rep(i) = FindGridIndex(rep(i), Grid);
    end
    
    % Check if Repository is Full
    if numel(rep)>nRep
        
        Extra = numel(rep)-nRep;
        for e = 1:Extra
            rep = DeleteOneRepMemebr(rep, gamma);
        end
        
    end
    
    % Plot Costs
    figure(1);
    PlotCosts(pop, rep);
    pause(0.01);
    
    % Show Iteration Information
    disp(['Iteration ' num2str(it) ': Number of Rep Members = ' num2str(numel(rep))]);
    
    % Damping Inertia Weight
    w = w*wdamp;
    
end

ZDT.m：測試函數(shù) ZDT1

function z = ZDT(x)
% Function: z = ZDT(x)
%
% Description:  此系列的測試函數(shù)總共有 6 個，ZDT1 ~ ZDT6, 屬于數(shù)值型測試函數(shù)，
%               當前實現(xiàn)的測試函數(shù)為 ZDT1
%
%
% Syntax:
%   
%
% Parameters:
%   x：決策變量
%
% Return:
%   z：函數(shù)結果 f1, f2
%
%                  Young99
%         Revision:    Data: 
%*************************************************************************

    n = numel(x);

    f1 = x(1);
    
    g = 1+9/(n-1)*sum(x(2:end));
    
    h = 1-sqrt(f1/g);
    
    f2 = g*h;
    
    z = [f1
       f2];

end

DetermineDomination.m：判斷是否為非支配個體

function pop = DetermineDomination(pop)
% Function: pop = DetermineDomination(pop)
%
% Description: 計算個體之間的支配關系，若被其他個體支配，將 IsDominated 設置為 true
%
%
% Syntax:
%   
%
% Parameters:
%   pop：種群
%
% Return:
%   pop：進行非支配標記后的種群
%
%                  Young99
%         Revision:1.0     Data: 2022-12-07
%*************************************************************************

    nPop = numel(pop);
    
    % 初始化每個個體均為非支配個體
    for i = 1:nPop
        pop(i).IsDominated = false;
    end
    
    % 兩兩個體比較
    for i = 1:nPop-1
        for j = i+1:nPop
            
            if Dominates(pop(i), pop(j))
               pop(j).IsDominated = true;
            end
            
            if Dominates(pop(j), pop(i))
               pop(i).IsDominated = true;
            end
            
        end
    end

end

Dominates.m：兩兩個體之間比較大小

function b = Dominates(x, y)
% Function: b = Dominates(x, y)
%
% Description: 比較兩個個體的函數(shù)值
%
%
% Syntax:
%   
%
% Parameters:
%   pop：種群
%
% Return:
%   b：比較結果，為布爾型
%
%                  Young99
%         Revision:1.0     Data: 2022-12-07
%*************************************************************************
    if isstruct(x)
        x = x.Cost;
    end
    
    if isstruct(y)
        y = y.Cost;
    end

    b = all(x <= y) && any(x<y);

end

CreateGrid.m：創(chuàng)建自適應網(wǎng)格

function Grid = CreateGrid(pop, nGrid, alpha)
% Function: Grid = CreateGrid(pop, nGrid, alpha)
%
% Description: 為每一個目標函數(shù)，選出 pop 中個體目標函數(shù)值的最小及最大值
%              將最小最大值利用 alpha 參數(shù)擴展一定的范圍，作為當前維度網(wǎng)格的上下限
%
%
% Syntax:
%   
%
% Parameters:
%   pop：rep 即為當前輪的非支配個體
%   nGrid：網(wǎng)格的維度大小
%   alpha：自定義參數(shù)，用于調(diào)整網(wǎng)格的最大最小上下線
%
% Return:
%   Grid：定義后的網(wǎng)格
%
%                  Young99
%         Revision:1.0     Data: 2022-12-07
%*************************************************************************

    c = [pop.Cost];
    
    cmin = min(c, [], 2);
    cmax = max(c, [], 2);
    
    dc = cmax-cmin;
    cmin = cmin-alpha*dc;
    cmax = cmax+alpha*dc;
    
    nObj = size(c, 1);
    
    empty_grid.LB = [];
    empty_grid.UB = [];
    Grid = repmat(empty_grid, nObj, 1);
    
    for j = 1:nObj
        
        cj = linspace(cmin(j), cmax(j), nGrid+1);
        
        Grid(j).LB = [-inf cj];
        Grid(j).UB = [cj +inf];
        
    end

end

FindGridIndex.m：為每個個體確定在網(wǎng)格中的位置

function particle = FindGridIndex(particle, Grid)
% Function: particle = FindGridIndex(particle, Grid)
%
% Description: 為 rep 中的每個個體劃分所在的網(wǎng)格，其中 
%              GridSubIndex 表示每維目標所對應的網(wǎng)格位置
%              GridIndex 表示將 GridSubIndex 歸并為一個數(shù)來表示
%
%
% Syntax:
%   
%
% Parameters:
%   particle：當前非支配集中的一個個體
%   Grid：已經(jīng)劃分好的網(wǎng)格，用于自適應網(wǎng)格算法
%
% Return:
%   particle：劃分到對應網(wǎng)格后的個體
%
%                  Young99
%         Revision:1.0     Data: 2022-12-07
%*************************************************************************
    nObj = numel(particle.Cost);
    
    nGrid = numel(Grid(1).LB);
    
    particle.GridSubIndex = zeros(1, nObj);
    
    for j = 1:nObj
        
        particle.GridSubIndex(j) = ...
            find(particle.Cost(j)<Grid(j).UB, 1, 'first');
        
    end
    
    % 將高維轉化為一維的空間索引
    particle.GridIndex = particle.GridSubIndex(1);
    for j = 2:nObj
        particle.GridIndex = particle.GridIndex-1;
        particle.GridIndex = nGrid*particle.GridIndex;
        particle.GridIndex = particle.GridIndex+particle.GridSubIndex(j);
    end
    
end

SelectLeader.m：利用網(wǎng)格聚集密度選擇 gBest

function leader = SelectLeader(rep, beta)
% Function: leader = SelectLeader(rep, beta)
%
% Description: 利用輪盤賭選出其中一個網(wǎng)格，網(wǎng)格聚集密度越小，越容易被選中，
%              再在選中的網(wǎng)格中隨機選出一個個體作為 leader 即 gbest
%
%
% Syntax:
%   
%
% Parameters:
%   rep：當前迭代的非支配個體
%   beta：自定義參數(shù)
%
% Return:
%   leader：被選中的 gbest 
%
%                  Young99
%         Revision:1.0     Data: 2022-12-07
%*************************************************************************


    % Grid Index of All Repository Members
    GI = [rep.GridIndex];
    
    % Occupied Cells
    OC = unique(GI);
    
    % Number of Particles in Occupied Cells
    N = zeros(size(OC));
    for k = 1:numel(OC)
        N(k) = numel(find(GI == OC(k)));
    end
    
    % Selection Probabilities
    P = exp(-beta*N);  % 注意這是有負數(shù)，表示個體數(shù)越多的網(wǎng)格被選擇的概率越小
    P = P/sum(P); % 轉化為和為 1 的概率
    
    % Selected Cell Index
    sci = RouletteWheelSelection(P);
    
    % Selected Cell
    sc = OC(sci);
    
    % Selected Cell Members
    SCM = find(GI == sc);
    
    % Selected Member Index，random selection
    smi = randi([1 numel(SCM)]);
    
    % Selected Member
    sm = SCM(smi);
    
    % Leader
    leader = rep(sm);

end

RouletteWheelSelection.m：輪盤賭

function i = RouletteWheelSelection(P)
% 輪盤賭用于選擇某個網(wǎng)格
    r = rand;
    
    C = cumsum(P);
    
    i = find(r <= C, 1, 'first');

end

Mutation.m：變異操作

function particle = Mutation(particle, pm, VarMin, VarMax)
% Function: particle = Mutation(particle, pm, VarMin, VarMax)
%
% Description: 隨機選擇一個決策變量進行變異，避免陷入局部最優(yōu)
%
%
% Syntax:
%   
%
% Parameters:
%   particle：個體決策變量
%   pm：擾亂算子
%   VarMin：變量最小值
%   VarMax：變量最大值
%
% Return:
%   particle：變異后的個體
%
%                  99Young99
%         Revision:1.0     Data: 2022-12-09
%*************************************************************************


    nVar = numel(particle);
    p = randi([1 nVar]); % 隨機選取一個變異位置 p
    
    x = particle(p);

    d = pm * (VarMax - VarMin);
    ub = x + d; 
    lb = x - d;
    ub = min(ub, VarMax);
    lb = max(lb, VarMin);

    x = unifrnd(lb, ub, 1);

    particle(p) = x;

end

DeleteOneRepMemebr.m：當 rep 超出閾值，進行刪除操作

function rep = DeleteOneRepMemebr(rep, gamma)
% Function: rep = DeleteOneRepMemebr(rep, gamma)
%
% Description: 利用輪盤賭選出其中一個網(wǎng)格，網(wǎng)格聚集密度越大，越容易被選中，
%              再在選中的網(wǎng)格中隨機選出一個個體刪除
%
%
% Syntax:
%   
%
% Parameters:
%   rep：當前迭代的非支配個體
%   beta：自定義參數(shù)
%
% Return:
%   leader：被選中的 gbest 
%
%                  Young99
%         Revision:1.0     Data: 2022-12-07
%*************************************************************************
    % Grid Index of All Repository Members
    GI = [rep.GridIndex];
    
    % Occupied Cells
    OC = unique(GI);
    
    % Number of Particles in Occupied Cells
    N = zeros(size(OC));
    for k = 1:numel(OC)
        N(k) = numel(find(GI == OC(k)));
    end
    
    % Selection Probabilities
    P = exp(gamma*N);
    P = P/sum(P);
    
    % Selected Cell Index
    sci = RouletteWheelSelection(P);
    
    % Selected Cell
    sc = OC(sci);
    
    % Selected Cell Members
    SCM = find(GI == sc);
    
    % Selected Member Index
    smi = randi([1 numel(SCM)]);
    
    % Selected Member
    sm = SCM(smi);
    
    % Delete Selected Member
    rep(sm) = [];

end

PlotCosts.m：畫圖可視化

function PlotCosts(pop, rep)

    pop_costs = [pop.Cost];
    plot(pop_costs(1, :), pop_costs(2, :), 'ko');
    hold on;
    
    rep_costs = [rep.Cost];
    plot(rep_costs(1, :), rep_costs(2, :), 'r*');
    
    xlabel('1^{st} Objective');
    ylabel('2^{nd} Objective');
    
    grid on;
    
    hold off;

end

HV.m

function [Score,PopObj] = HV(PopObj,PF)
% <metric> <max>
% Hypervolume

%------------------------------- Reference --------------------------------
% E. Zitzler and L. Thiele, Multiobjective evolutionary algorithms: A
% comparative case study and the strength Pareto approach, IEEE
% Transactions on Evolutionary Computation, 1999, 3(4): 257-271.
%------------------------------- Copyright --------------------------------
% Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for
% research purposes. All publications which use this platform or any code
% in the platform should acknowledge the use of "PlatEMO" and reference "Ye
% Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform
% for evolutionary multi-objective optimization [educational forum], IEEE
% Computational Intelligence Magazine, 2017, 12(4): 73-87".
%--------------------------------------------------------------------------

% 
% 這段代碼是用來計算超體積指標的。它包含了一個 HV 函數(shù)，用于計算超體積指標，并輸出其值和輸入種群的結果。
% 
% HV 函數(shù)的輸入?yún)?shù)包括種群結果 PopObj 和理論最優(yōu)結果 PF。它首先會將種群結果進行歸一化，然后根據(jù)超體積指標的定義，計算種群結果與理論最優(yōu)結果的相似度。
% 
% 最后，HV 函數(shù)會輸出計算得到的超體積指標值和歸一化后的種群結果。


    % Normalize the population according to the reference point set
    [N,M]  = size(PopObj);
    fmin   = min(min(PopObj,[],1),zeros(1,M));
    fmax   = max(PF,[],1);
    PopObj = (PopObj-repmat(fmin,N,1))./repmat((fmax-fmin)*1.1,N,1);
    PopObj(any(PopObj>1,2),:) = [];
    % The reference point is set to (1,1,...)
    RefPoint = ones(1,M);
    if isempty(PopObj)
        Score = 0;
    elseif M < 4
        % Calculate the exact HV value
        pl = sortrows(PopObj);
        S  = {1,pl};
        for k = 1 : M-1
            S_ = {};
            for i = 1 : size(S,1)
                Stemp = Slice(cell2mat(S(i,2)),k,RefPoint);
                for j = 1 : size(Stemp,1)
                    temp(1) = {cell2mat(Stemp(j,1))*cell2mat(S(i,1))};
                    temp(2) = Stemp(j,2);
                    S_      = Add(temp,S_);
                end
            end
            S = S_;
        end
        Score = 0;
        for i = 1 : size(S,1)
            p     = Head(cell2mat(S(i,2)));
            Score = Score + cell2mat(S(i,1))*abs(p(M)-RefPoint(M));
        end
    else
        % Estimate the HV value by Monte Carlo estimation
        SampleNum = 1000000;
        MaxValue  = RefPoint;
        MinValue  = min(PopObj,[],1);
        Samples   = unifrnd(repmat(MinValue,SampleNum,1),repmat(MaxValue,SampleNum,1));
        if gpuDeviceCount > 0
            % GPU acceleration
            Samples = gpuArray(single(Samples));
            PopObj  = gpuArray(single(PopObj));
        end
        for i = 1 : size(PopObj,1)
            drawnow();
            domi = true(size(Samples,1),1);
            m    = 1;
            while m <= M && any(domi)
                domi = domi & PopObj(i,m) <= Samples(:,m);
                m    = m + 1;
            end
            Samples(domi,:) = [];
        end
        Score = prod(MaxValue-MinValue)*(1-size(Samples,1)/SampleNum);
    end
end

function S = Slice(pl,k,RefPoint)
    p  = Head(pl);
    pl = Tail(pl);
    ql = [];
    S  = {};
    while ~isempty(pl)
        ql  = Insert(p,k+1,ql);
        p_  = Head(pl);
        cell_(1,1) = {abs(p(k)-p_(k))};
        cell_(1,2) = {ql};
        S   = Add(cell_,S);
        p   = p_;
        pl  = Tail(pl);
    end
    ql = Insert(p,k+1,ql);
    cell_(1,1) = {abs(p(k)-RefPoint(k))};
    cell_(1,2) = {ql};
    S  = Add(cell_,S);
end

function ql = Insert(p,k,pl)
    flag1 = 0;
    flag2 = 0;
    ql    = [];
    hp    = Head(pl);
    while ~isempty(pl) && hp(k) < p(k)
        ql = [ql;hp];
        pl = Tail(pl);
        hp = Head(pl);
    end
    ql = [ql;p];
    m  = length(p);
    while ~isempty(pl)
        q = Head(pl);
        for i = k : m
            if p(i) < q(i)
                flag1 = 1;
            else
                if p(i) > q(i)
                    flag2 = 1;
                end
            end
        end
        if ~(flag1 == 1 && flag2 == 0)
            ql = [ql;Head(pl)];
        end
        pl = Tail(pl);
    end  
end

function p = Head(pl)
    if isempty(pl)
        p = [];
    else
        p = pl(1,:);
    end
end

function ql = Tail(pl)
    if size(pl,1) < 2
        ql = [];
    else
        ql = pl(2:end,:);
    end
end

function S_ = Add(cell_,S)
    n = size(S,1);
    m = 0;
    for k = 1 : n
        if isequal(cell_(1,2),S(k,2))
            S(k,1) = {cell2mat(S(k,1))+cell2mat(cell_(1,1))};
            m = 1;
            break;
        end
    end
    if m == 0
        S(n+1,:) = cell_(1,:);
    end
    S_ = S;     
end

IGD.m

function Score = IGD(PopObj,PF)

% PopObej：算法求得的非支配解集
% PF：真實 Pareto 解集

% <metric> <min>
% Inverted generational distance

%------------------------------- Reference --------------------------------
% C. A. Coello Coello and N. C. Cortes, Solving multiobjective optimization
% problems using an artificial immune system, Genetic Programming and
% Evolvable Machines, 2005, 6(2): 163-190.
%------------------------------- Copyright --------------------------------
% Copyright (c) 2018-2019 BIMK Group. You are free to use the PlatEMO for
% research purposes. All publications which use this platform or any code
% in the platform should acknowledge the use of "PlatEMO" and reference "Ye
% Tian, Ran Cheng, Xingyi Zhang, and Yaochu Jin, PlatEMO: A MATLAB platform
% for evolutionary multi-objective optimization [educational forum], IEEE
% Computational Intelligence Magazine, 2017, 12(4): 73-87".
%--------------------------------------------------------------------------

    % 通過 pdist 計算 PF 與 PopObj 兩兩歐式距離，通過 min 取每一行的最小值作為列向量
    % 最后再取所有距離的平均值，越小說明分布性和收斂性越

    Distance = min(pdist2(PF,PopObj),[],2); 
    Score    = mean(Distance);
end

Reference

部分理論來源于網(wǎng)絡，如有侵權請聯(lián)系刪除文章來源地址http://www.zghlxwxcb.cn/news/detail-743406.html

??小手手動起來，點個贊再走唄 ~??

到了這里，關于【多目標進化優(yōu)化】MOPSO 原理與代碼實現(xiàn)的文章就介紹完了。如果您還想了解更多內(nèi)容，請在右上角搜索TOY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關文章，希望大家以后多多支持TOY模板網(wǎng)！

本文來自互聯(lián)網(wǎng)用戶投稿，該文觀點僅代表作者本人，不代表本站立場。本站僅提供信息存儲空間服務，不擁有所有權，不承擔相關法律責任。如若轉載，請注明出處：如若內(nèi)容造成侵權/違法違規(guī)/事實不符，請點擊違法舉報進行投訴反饋，一經(jīng)查實，立即刪除！

分享到：

領支付寶紅包贊助服務器費用

智能優(yōu)化算法：多目標粒子群優(yōu)化算法（MOPSO）
目錄一、粒子群優(yōu)化算法（Particle Swarm Optimization，PSO）二、多目標優(yōu)化Pareto支配三、多目標粒子群優(yōu)化算法（Multiple Objective Particle Swarm Optimization，MOPSO）四、參考文獻：粒子群算法(PSO)[1]是一種群智能優(yōu)化算法，于1995年Kennedy和Eberhart提出，其靈感來源于鳥群的協(xié)作覓食行
2023年04月19日
瀏覽(34)
五種多目標優(yōu)化算法（MOAHA、MOGWO、NSWOA、MOPSO、NSGA2）性能對比，包含6種評價指標，9個測試函數(shù)（提供MATLAB代碼）
1.1MOAHA 1.2MOGWO 1.3NSWOA 1.4MOPSO 1.5NSGA2 為了測試5種算法的性能將其求解9個多目標測試函數(shù)（zdt1、zdt2 、zdt3、 zdt4、 zdt6 、Schaffer、?Kursawe 、Viennet2、 Viennet3），其中Viennet2 與Viennet3的目標數(shù)為3，其余測試函數(shù)的目標數(shù)為2，并采用6種評價指標（IGD、GD、HV、Coverage、Spread、Spacing）進
2024年02月21日
瀏覽(95)
五種多目標優(yōu)化算法（MOJS、MOGWO、NSWOA、MOPSO、NSGA2）性能對比，包含6種評價指標，9個測試函數(shù)（提供MATLAB代碼）
1.1MOJS 1.2MOGWO 1.3NSWOA 1.4MOPSO 1.5NSGA2 為了測試5種算法的性能將其求解9個多目標測試函數(shù)（zdt1、zdt2 、zdt3、 zdt4、 zdt6 、Schaffer、?Kursawe 、Viennet2、 Viennet3），其中Viennet2 與Viennet3的目標為3，其余測試函數(shù)的目標為2，并采用6種評價指標（IGD、GD、HV、Coverage、Spread、Spacing）進行評
2024年02月21日
瀏覽(101)
【無功優(yōu)化】基于多目標差分進化算法的含DG配電網(wǎng)無功優(yōu)化模型【IEEE33節(jié)點】（Matlab代碼實現(xiàn)）
???????? 歡迎來到本博客 ???????? ??博主優(yōu)勢： ?????? 博客內(nèi)容盡量做到思維縝密，邏輯清晰，為了方便讀者。 ?? 座右銘：行百里者，半于九十。 ?????? 本文目錄如下： ?????? 目錄 ??1 概述 ??2 運行結果 ? ??3?參考文獻 ??4 Matlab代碼、數(shù)據(jù)、文
2023年04月16日
瀏覽(92)
進化算法——多目標優(yōu)化
所有的實際優(yōu)化問題都是多目標的，如果不是顯式的至少也是隱式的。接下來討論多目標優(yōu)化問題（MOP）如何修改進化算法。實際的優(yōu)化問題包含多個目標，那些目標常常互相沖突。例如：在購買汽車時，我們可能想要車最舒適并且花錢最少。最舒適的汽車太貴，最便宜的汽
2023年04月13日
瀏覽(21)
【單目標優(yōu)化算法】海鷗優(yōu)化算法（Matlab代碼實現(xiàn)）
?? ?? ?? ?? 歡迎來到本博客 ?? ?? ?? ?? ?? 博主優(yōu)勢： ?? ?? ??博客內(nèi)容盡量做到思維縝密，邏輯清晰，為了方便讀者。 ? 座右銘：行百里者，半于九十。 ?? ?? ?? 本文目錄如下： ?? ?? ?? 目錄 ??1 概述 ??2 運行結果 ??3 參考文獻 ??4 Matlab代碼實現(xiàn) 編
2024年02月03日
瀏覽(28)
【單目標優(yōu)化算法】雜草優(yōu)化算法（Matlab代碼實現(xiàn)）
?? ?? ?? ?? 歡迎來到本博客 ?? ?? ?? ?? ?? 博主優(yōu)勢： ?? ?? ??博客內(nèi)容盡量做到思維縝密，邏輯清晰，為了方便讀者。 ? 座右銘：行百里者，半于九十。 ?? ?? ?? 本文目錄如下： ?? ?? ?? 目錄 ??1 概述 ??2 運行結果 ??3 參考文獻 ??4 Matlab代碼實現(xiàn) 雜
2024年02月11日
瀏覽(29)
多目標優(yōu)化NSGA-II的實現(xiàn)（MATLAB完整代碼）
由于歷史原因，沒有整理好完整的代碼，所以在【多目標優(yōu)化NSGA-II的實現(xiàn)和測試（MATLAB實現(xiàn)）】中只放了部分代碼。現(xiàn)在已經(jīng)整理好了代碼，此部分的代碼測試內(nèi)容為：ZDT1、ZDT2、ZDT3、ZDT4、ZDT6。目錄主要內(nèi)容代碼模塊其他內(nèi)容運行注意事項 ?代碼 nsga2_test nsga2_main get_
2024年02月05日
瀏覽(27)
【經(jīng)濟調(diào)度】基于多目標宇宙優(yōu)化算法優(yōu)化人工神經(jīng)網(wǎng)絡環(huán)境經(jīng)濟調(diào)度研究（Matlab代碼實現(xiàn)）
????????? 歡迎來到本博客 ???????? ??博主優(yōu)勢： ?????? 博客內(nèi)容盡量做到思維縝密，邏輯清晰，為了方便讀者。 ?? 座右銘：行百里者，半于九十。 ?????? 本文目錄如下： ?????? 目錄 ??1 概述 ??2 運行結果 ??3?參考文獻 ??4 Matlab代碼及數(shù)據(jù) 基于
2024年02月16日
瀏覽(24)
目標檢測--邊框回歸損失函數(shù)SIoU原理詳解及代碼實現(xiàn)
1.1 原理有關IoU損失函數(shù)，如(GIoU, DIoU, CIoU)沒有考慮到真實框與預測框框之間的方向，導致收斂速度較慢，對此SIoU引入真實框和預測框之間的向量角度，重新定義相關損失函數(shù)，具體包含四個部分：（1）角度損失(Angle cost)，定義如下 Λ = 1 ? 2 ? sin ? 2 ( arcsin ? ( c h σ ) ?
2024年01月18日
瀏覽(18)

<sub id="lvslj"><kbd id="lvslj"></kbd></sub>