​

一、 BP神經網路預測演算法簡介

說明：1.1節主要是概括和幫助理解考慮影響因素的BP神經網路演算法原理，即常規的BP模型訓練原理講解（可根據自身掌握的知識是否跳過）。1.2節開始講基於歷史值影響的BP神經網路預測模型。

使用BP神經網路進行預測時，從考慮的輸入指標角度，主要有兩類模型：

1.1 受相關指標影響的BP神經網路演算法原理

如圖一所示，使用MATLAB的newff函式訓練BP時，可以看到大部分情況是三層的神經網路（即輸入層，隱含層，輸出層）。這裡幫助理解下神經網路原理：
1）輸入層：相當於人的五官，五官獲取外部資訊，對應神經網路模型input埠接收輸入資料的過程。
2）隱含層：對應人的大腦，大腦對五官傳遞來的資料進行分析和思考，神經網路的隱含層hidden Layer對輸入層傳來的資料x進行對映，簡單理解為一個公式hiddenLayer_output=F(w*x+b)。其中，w、b叫做權重、閾值引數，F()為對映規則，也叫啟用函式，hiddenLayer_output是隱含層對於傳來的資料對映的輸出值。換句話說，隱含層對於輸入的影響因素資料x進行了對映，產生了對映值。
3）輸出層：可以對應為人的四肢，大腦對五官傳來的資訊經過思考（隱含層對映）之後，再控制四肢執行動作（向外部作出響應）。類似地，BP神經網路的輸出層對hiddenLayer_output再次進行對映，outputLayer_output=w *hiddenLayer_output+b。其中，w、b為權重、閾值引數，outputLayer_output是神經網路輸出層的輸出值（也叫模擬值、預測值）（理解為，人腦對外的執行動作，比如嬰兒拍打桌子）。
4）梯度下降演算法：通過計算outputLayer_output和神經網路模型傳入的y值之間的偏差，使用演算法來相應調整權重和閾值等引數。這個過程，可以理解為嬰兒拍打桌子，打偏了，根據偏離的距離遠近，來調整身體使得再次揮動的胳膊不斷靠近桌子，最終打中。

再舉個例子來加深理解：

圖一所示BP神經網路，具備輸入層、隱含層和輸出層。BP是如何通過這三層結構來實現輸出層的輸出值outputLayer_output，不斷逼近給定的y值，從而訓練得到一個精準的模型的呢？

從圖中串起來的埠，可以想到一個過程：坐地鐵，將圖一想象為一條地鐵線路。王某某坐地鐵回家的一天：在input起點站上車，中途經過了很多站（hiddenLayer），然後發現坐過頭了（outputLayer對應現在的位置），那麼王某某將會根據現在的位置離家（目標Target）的距離（誤差Error），返回到中途的地鐵站（hiddenLayer）重新坐地鐵（誤差反向傳遞，使用梯度下降演算法更新w和b），如果王某某又一次發生失誤，那麼將再次進行這個調整的過程。

從在嬰兒拍打桌子和王某某坐地鐵的例子中，思考問題：BP的完整訓練，需要先傳入資料給input，再經過隱含層的對映，輸出層得到BP模擬值，根據模擬值與目標值的誤差，來調整引數，使得模擬值不斷逼近目標值。比如（1）嬰兒受到了外界的干擾因素（x），從而作出反應拍桌（predict），大腦不斷的調整胳膊位置，控制四肢拍準（y、Target）。（2）王某某上車點（x），過站點（predict），不斷返回中途站來調整位置，到家（y、Target）。

在這些環節中，涉及了影響因素資料x，目標值資料y（Target）。根據x，y，使用BP演算法來尋求x與y之間存在的規律，實現由x來對映逼近y，這就是BP神經網路演算法的作用。再多說一句，上述講的過程，都是BP模型訓練，那麼最終得到的模型雖然訓練準確，但是找到的規律（bp network）是否準確與可靠呢。於是，我們再給x1到訓練好的bp network中，得到相應的BP輸出值（預測值）predict1，通過作圖，計算Mse，Mape，R方等指標，來對比predict1和y1的接近程度，就可以知道模型是否預測準確。這是BP模型的測試過程，即實現對資料的預測，並且對比實際值檢驗預測是否準確。
​
圖一 3層BP神經網路結構圖

1.2 基於歷史值影響的BP神經網路

以電力負荷預測問題為例，進行兩種模型的區分。在預測某個時間段內的電力負荷時：

一種做法，是考慮t時刻的氣候因素指標，比如該時刻的空氣溼度x1，溫度x2，以及節假日x3等的影響，對t時刻的負荷值進行預測。這是前面1.1所說的模型。

另一種做法，是認為電力負荷值的變化，與時間相關，比如認為t-1，t-2，t-3時刻的電力負荷值與t時刻的負荷值有關係，即滿足公式y(t)=F(y(t-1),y(t-2),y(t-3))。採用BP神經網路進行訓練模型時，則輸入到神經網路的影響因素值為歷史負荷值y(t-1),y(t-2),y(t-3)，特別地，3叫做自迴歸階數或者延遲。給到神經網路中的目標輸出值為y(t)。

二、麻雀演算法

優化問題是科學研究和工程實踐領域中的熱門問題。智慧優化演算法大多是受到人類智慧、生物群體社會性或自然現象規律的啟發，在解空間內進行全域性優化。麻雀演算法於2020年由薛建凱[1]首次提出，是基於麻雀種群的覓食和反捕食行為的一種新型智慧優化演算法。

麻雀搜尋演算法的具體步驟描述以及公式介紹：

構建麻雀種群：
​

其中，d表示待優化問題的維數，n表示麻雀種群的數量。所有麻雀種群的適應度函式可以表示成如下形式：

​

其中，Fx表示適應度函式值。

麻雀演算法中的麻雀具有兩大類分別是發現者和加入者，發現者負責為整個種群尋找食物併為加入者提供覓食的方向，因此，發現者的覓食搜尋範圍要比加入者的覓食搜尋範圍大。在每次迭代過程中，發現者按照公式（3）進行迭代。

​

其中，t表示當前迭代次數，Xij表示第i個麻雀種群在第j維中的位置資訊,阿爾法表示的0到1的隨機數，itermax表示最大迭代次數，Q表示一個服從正態分佈的隨機數，L是一個1*d並且元素全為1的矩陣，R2屬於0-1表示麻雀種群位置的預警值，ST屬於0.5-1表示麻雀種群位置的安全值。

當R2<ST時表示 預警值小於安全值，此時覓食環境中沒有捕食者，發現者可以進行廣泛搜尋操作；當R2>ST時意味著種群中有部分麻雀已經發現捕食者，並向種群中的其他麻雀發出預警，所有麻雀都需要飛往安全區域進行覓食。

在覓食過程中，部分加入者會時刻監視發現者，當發現者發現更好的食物，加入者會與其進行爭奪，若成功，會立即獲得該發現者的食物，否則加入者按照公式（4）進行位置更新。

​

其中，XP表示目前發現者所發現的最優位置，Xworst表示當前全域性最差的位置，A表示其元素隨機賦值為1或-1的1*d的矩陣並且滿足一下關係：

L仍然是一個1*d並且元素全為1的矩陣。當i>n/2時這表明第i個加入者沒有獲得食物，處於飢餓狀態，此時需要飛往其他地方進行覓食，以獲得更多的能量。

在麻雀種群中，意識到危險的麻雀數量佔總數的10%到20%，這些麻雀的位置是隨機產生的，按照公式（5）對意識到危險的麻雀的位置進行不斷更新。

​

其中，Xbest表示當前全域性最優位置，是服從標準正態分佈的隨機數用來作為步長控制引數，貝塔是一個屬於-1到1的隨機數，fi表示當前麻雀個體的適應度值，fg表示全域性最佳適應度值，fw表示全域性最差適應度值，像左耳朵一樣的這個是讀"一不洗諾"嗎？"一不洗諾"表示一個避免分母為0的常數。當fi>fg時表示此時麻雀處於種群邊緣，極易受到捕食者的攻擊，當fi=fg時表示處於種群中間的麻雀也受到了危險，此時需要靠近其他麻雀以減少被捕食的風險。

三、部分程式碼

function [FoodFitness,FoodPosition,Convergence_curve]=SSA(N,Max_iter,lb,ub,dim,fobj)
if size(ub,1)==1
    ub=ones(dim,1)*ub;
    lb=ones(dim,1)*lb;
end
Convergence_curve = zeros(1,Max_iter);
%Initialize the positions of salps
SalpPositions=initialization(N,dim,ub,lb);
FoodPosition=zeros(1,dim);
FoodFitness=inf;
%calculate the fitness of initial salps
for i=1:size(SalpPositions,1)
    SalpFitness(1,i)=fobj(SalpPositions(i,:));
end
[sorted_salps_fitness,sorted_indexes]=sort(SalpFitness);
for newindex=1:N
    Sorted_salps(newindex,:)=SalpPositions(sorted_indexes(newindex),:);
end
FoodPosition=Sorted_salps(1,:);
FoodFitness=sorted_salps_fitness(1);
%Main loop
l=2; % start from the second iteration since the first iteration was dedicated to calculating the fitness of salps
while l<Max_iter+1
    c1 = 2*exp(-(4*l/Max_iter)^2); % Eq. (3.2) in the paper
    for i=1:size(SalpPositions,1)
        SalpPositions= SalpPositions';
        if i<=N/2
            for j=1:1:dim
                c2=rand();
                c3=rand();
                %%%%%%%%%%%%% % Eq. (3.1) in the paper %%%%%%%%%%%%%%
                if c3<0.5 
                    SalpPositions(j,i)=FoodPosition(j)+c1*((ub(j)-lb(j))*c2+lb(j));
                else
                    SalpPositions(j,i)=FoodPosition(j)-c1*((ub(j)-lb(j))*c2+lb(j));
                end
                %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
            end
        elseif i>N/2 && i<N+1
            point1=SalpPositions(:,i-1);
            point2=SalpPositions(:,i);
            SalpPositions(:,i)=(point2+point1)/2; % % Eq. (3.4) in the paper
        end
        SalpPositions= SalpPositions';
    end
    for i=1:size(SalpPositions,1)
        Tp=SalpPositions(i,:)>ub';Tm=SalpPositions(i,:)<lb';SalpPositions(i,:)=(SalpPositions(i,:).*(~(Tp+Tm)))+ub'.*Tp+lb'.*Tm;
        SalpFitness(1,i)=fobj(SalpPositions(i,:));
        if SalpFitness(1,i)<FoodFitness
            FoodPosition=SalpPositions(i,:);
            FoodFitness=SalpFitness(1,i); 
        end
    end
    Convergence_curve(l)=FoodFitness;
    l = l + 1;
end

四、模擬結果

​

​

​

五、參考文獻及程式碼私信博主

《基於BP神經網路的寧夏水資源需求量預測》

​