$$
\left\{\begin{array}{l}
s_{k+1}=\boldsymbol{A}_{k} s_{k}+\boldsymbol{B}_{k} \boldsymbol{u}_{k} \\
y_{k+1}=f\left(s_{k+1}\right)
\end{array}\right.
$$
其中, $s_{k} = (s_{1}, s_{2}, · · · , s_{n})^{T}$為當前狀態, 代表優化問題的一個候選解; $\boldsymbol{A}_{k}$ 或$\boldsymbol{B}_{k}$ 為狀態轉移矩陣, 可以看成是算法中算子所產生的的變換效果; $\boldsymbol{u}_{k}$ 為當前狀態及歷史狀態的函數, 可以看成是控制變量; $f(·)$為目標函數或者評價函數。
2.2 基本狀態轉移算法的理論實現
狀態轉移算法中專門設計了局部、全局和啟發式搜索算子，其中全局搜索算子是為了保證有一定概率使產生的候選解集形成的鄰域中包含全局最優解； 局部搜索算子具有在較小鄰域進行精細搜索的功能;；啟發式搜索算子用來產生具有一定潛在價值的候選解，避免搜索的盲目性。 另一方面，狀態轉移算法中設計了特定的選擇與更新策略， 用來適應不同優化問題的需要。 更進一步地， 狀態轉移算法中實施了智能化策略， 比如全局和局部算子的交替調用，可以很大程度上縮短尋找全局最優解的時間和避免陷入局部最優； 采樣方式的使用,可以選擇有代表性的候選解， 避免窮舉鄰域內的所有解，將極大地縮短搜索時間。一般地，狀態轉移算法的基本流程如下:
i ) 從當前最好解出發，通過某種狀態變換算子生成具有某種特定性質的鄰域；
ii ) 以某種采樣機制從該鄰域內采集一定數目的樣本；
iii ) 根據選擇與更新策略，比較當前最好解與采集樣本中最好解的優劣， 並以某種更新機制更新當前最好解； iv ) 返回步驟i )， 並交替使用各種狀態變換算子， 直到滿足某種終止條件。
$$
s_{k+1}=s_{k}+\alpha \frac{1}{n\left\|s_{k}\right\|_{2}} \boldsymbol{R}_{r} s_{k}{(2)}
$$
其中$\alpha$為旋轉因子；$\boldsymbol{R}_{r} \in \mathbb{R}^{n \times n}$為隨機矩陣，其元素取值於[-1,1]之間均勻分布；$\|\cdot\|_{2}$表示取向量的模。旋轉變換能夠以當前解$s_{k}$為原點，在給定半徑$\alpha$內進行局部搜索。
$$
s_{k+1}=x_{k}+\beta \boldsymbol{R}_{t} \frac{s_{k}-s_{k-1}}{\left\|s_{k}-s_{k-1}\right\|_{2}}
$$
其中$\beta$為平移因子；$\boldsymbol{R}_{t} \in \mathbb{R}^{n \times n}$是一個隨機向量，其元素取值於[0,1]之間均勻分布。平移變換沿着當前解$s_{k}$和前一個歷史解$s_{k-1}$所形成的的直線上進行啟發式搜索。
$$
s_{k+1}=s_{k}+\gamma \boldsymbol{R}_{e} s_{k}
$$
其中$\gamma$為伸縮因子；$\boldsymbol{R}_{e} \in \mathbb{R}^{n \times n}$是一個非零元素取值服從正態分布的隨機對角矩陣。從式子中不難理解，伸縮變換能夠以正態分布的概率，將當前解$s_{k}$的每個維度上的值伸縮到實軸上的任意位置，具有全局搜索的能力。
$$
s_{k+1}=s_{k}+\delta \boldsymbol{R}_{a} s_{k}
$$
其中$\delta$為坐標因子；$\boldsymbol{R}_{a} \in \mathbb{R}^{n \times n}$是一個非零元素取值服從正態分布的稀疏對角矩陣，在基本狀態轉移算法中，$\boldsymbol{R}_{a}$僅有一個位置為非零元素。坐標變換能夠從當前解$s_{k}$出發，隨機挑選一個維度進行搜索，能夠增強單一維度的搜索能力。
2.2.2 鄰域、采樣、選擇與更新
鄰域：
由當前狀態$s_{k}$出發，通過狀態變換算子能夠產生候選狀態$s_{k+1}$。不難看出，由於隨機性，某種特定的算子能夠反復產生不相等但是具有同質性的一系列候選狀態。這些候選解將構成一個集合，形成一個鄰域。理論上，鄰域中具有無數個候選解，但在實際中，必須采用一定的采樣機制選出有限個候選解。
采樣：
對於鄰域中的具有同質性的候選解，以隨機策略采樣$SE$個候選解，產生一個大小有限且可控的集合$State$作為實際當中的狀態候選集合。$SE$可理解為“搜索強度”、“采樣力度”或者“樣本個數”。
選擇與更新：
在算法運行的階段中，考慮$State$中的$SE$個候選解和當前最優解，采取“貪婪策略”來進行選擇與更新最優解。每次選取具有最優適應度函數值的解來更新當前最優解。
2.2.3 交替輪換
為了實現可控、最優、全局、穩定的搜索，具有不同搜索功能的狀態變換算子被交替輪換地調用。這種設計思路，從一開始就清晰地規划了優化中各種可能需要的搜索形式，並將其抽象封裝為功能各異的狀態變換算子。於是，在解決實際問題的時候，不同搜索功能可以有選擇地被強化或者削弱，以適應不同場景的具體需要。進一步改進優化算法的思路也就從漫無目的地構思新的動物算法，變成了如何為已知的狀態變換算子設計合適的調用策略。
2.2.4 參數設置
在基本的狀態轉移算子中，參數設置為：$\alpha_{max}=1$，$\alpha_{min}=1e-4$，$\beta=1$，$\gamma=1$，$\delta=1$，$SE=30$，$fc=2$。其中平移、伸縮和坐標因子的大小都保持恆定，旋轉因子以$1/fc$為衰減率，在最大值$\alpha_{max}$和最小值$\alpha_{min}$之間周期變化。
如圖為基本狀態轉移算法的MATLAB源碼的全部文件。本部分將首先對各個M文件給出概括性的功能描述，然后針對具體代碼給出給出解釋與分析。
3.1 概括性的功能描述
Test_sta.m為算法的調用示例；STA.m為封裝好的基本狀態轉移算法；initialization.m為解的初始化過程；fitness.m為適應度計算與比較；Get_Test_Functions.m為各類常見優化問題目標函數的匯總封裝；rotate.m，expand.m，axesion.m為對應狀態變換算子的外層邏輯；op_rotate.m，op_expand.m，op_axesion.m，op_translate.m分別是四種狀態變換算子的矩陣運算實現。
3.2 源碼解析
Test_sta.m ：
clear all
clc
currentFolder = pwd;
addpath(genpath(currentFolder))
% parameter setting
warning('off')
SE = 30; %搜索強度
Dim = 30;%問題維度定義
Range = repmat([-5.12;5.12],1,Dim);%搜索范圍定義
Iterations = 1e3;% 最大迭代次數定義
tic
[Best,fBest,history] = STA(@Rastrigin,SE,Dim,Range,Iterations);
%調用sta.m，以SE的搜索強度，迭代最多Iterations次，在搜索范圍Range內，解決變量數為Dim個的典型優化問題Rastrigin函數。並返回最優解Best，最優值fBest以及搜索歷史history。
toc
history
semilogy(history)
STA.m ：
%參數初始化
alpha_max = 1;
alpha_min = 1e-4; alpha = alpha_max; beta = 1; gamma = 1; delta = 1; fc = 2;
% 初始化狀態集合State State = initialization(SE,Dim,Range); [Best,fBest] = fitness(funfcn,State);%計算適應度函數，返回當前最優解Best和最優值fBest
counter = 0; %迭代過程 for iter = 1:Iterations if alpha < alpha_min alpha = alpha_max; end oldfBest = fBest; [Best,fBest] = expand(funfcn,Best,SE,Range,beta,gamma); [Best,fBest] = rotate(funfcn,Best,SE,Range,alpha,beta); [Best,fBest] = axesion(funfcn,Best,SE,Range,beta,delta); %終止條件判斷 if norm(oldfBest-fBest) < 1e-8 % 計算矩陣范數表示兩代最優值之間的誤差 counter = counter + 1; if counter > 10 % can be changed break; end else counter = 0; end fprintf('iter=%d ObjVal=%g\n',iter,fBest); history(iter) = fBest; alpha = alpha/fc;%更新參數alpha的值 end
initialization.m ：
function State=initialization(SE,Dim,Range) Pop_Lb = Range(1,:);%Range的第一行，全是-5.12 Pop_Ub = Range(2,:);%Range的第二行，全是5.12 State = rand(SE,Dim).*repmat(Pop_Ub-Pop_Lb,SE,1)+repmat(Pop_Lb,SE,1);%生成初始狀態集合State，大小為Dim*SE
fitness.m :
function [Best,fBest] = fitness(funfcn,State)
% calculate fitness
SE = size(State,1);
fState = zeros(SE,1);
for i = 1:SE
fState(i) = feval(funfcn,State(i,:));
end
[fGBest, g] = min(fState);%返回最小值fGBest和其索引g
fBest = fGBest;
Best = State(g,:);%最值對應那一行的解
op_rotate.m :
function y=op_rotate(Best,SE,alpha) n = length(Best); y = repmat(Best',1,SE) + alpha*(1/n/(norm(Best)+eps))*reshape(unifrnd(-1,1,SE*n,n)*Best',n,SE);%旋轉算子的矩陣運算實現 y = y';%轉置