由於工作內容接觸到點雲標定，需要用到最小二乘法，所以特意花了點時間研究lm演算法，但是由於大學的高等數學忘得差不多了，所以本文從最基本的一些數學概念開始；
在最優化演算法中，都是要求乙個函式的極小值，每一步迭代中，都要求目標函式值是下降的，而信賴域法，顧名思義，就是從初始點開始，先假設乙個可以信賴的最大位移，然後在以當前點為中心，以為半徑的區域內，通過尋找目標函式的乙個近似函式（二次的）的最優點，來求解得到真正的位移。在得到了位移之後，再計算目標函式值，如果其使目標函式值的下降滿足了一定條件，那麼就說明這個位移是可靠的，則繼續按此規則迭代計算下去；如果其不能使目標函式值的下降滿足一定的條件，則應減小信賴域的範圍，再重新求解。
雅可比矩陣幾乎在所有的最優化演算法中都有提及，因此我們很有必要了解一下其具到底是什麼，關於這一點，下方截圖說的很清楚；
從上面可以了解，雅可比矩陣實際上就是一階偏導數所組成的矩陣，其列數由未知引數個數決定，其行數由我們提供的輸入引數組決定；
需要注意的是,對於lm演算法，可以具體到下種形式：
其中，r是殘差；
lm演算法的關鍵是用模型函式 f 對待估引數向量p在其領域內做線性近似，忽略掉二階以上的導數項，從而轉化為線性最小二乘問題，它具有收斂速度快等優點。
lm演算法需要對每乙個待估引數求偏導，所以，如果你的擬合函式 f 非常複雜，或者待估引數相當地多，那麼就不適合使用lm演算法了，可以使用powell演算法，powell演算法不需要求導。
需要說明的是，這是非線性無約束的問題，如果待估引數是有約束的，暫時還沒有涉及到這個領域；
就是從初始點開始，先假設乙個可以信賴的最大位移，然後在以當前點為中心，以為半徑的區域內，通過尋找目標函式的乙個近似函式（二次的）的最優點，來求解得到真正的位移。在得到了位移之後，再計算目標函式值，如果其使目標函式值的下降滿足了一定條件，那麼就說明這個位移是可靠的，則繼續按此規則迭代計算下去；如果其不能使目標函式值的下降滿足一定的條件，則應減小信賴域的範圍，再重新求解。
在使用levenberg-marquart時，先設定乙個比較小的μ值，當發現目標函式反而增大時，將μ增大使用梯度下降法快速尋找，然後再將μ減小使用牛頓法進行尋找。
6.阻尼係數的調整
當阻尼係數足夠大時，使演算法更接近最速下降法，所以在殘差沒有明顯變化時可以使用；當阻尼係數足夠小時，演算法更接近高斯牛頓演算法，此時迭代速度更快；
有演算法精度ep和上一次殘差e，當eep時，lamda = lamda*5，當lamda < ep時，lamda = lamda；
**如下：
% 計算函式f的雅克比矩陣
syms a b y x real;
f=a*cos(b*x) + b*sin(a*x)
jsym=jacobian(f,[a b])
data_1=[ 0, 0.2, 0.4, 0.6, 0.8, 1.0, 1.2, 1.4, 1.6, 1.8, 2.0, 2.2, 2.4, 2.6, 2.8, 3.0,3.2, 3.4, 3.6, 3.8, 4.0, 4.2, 4.4, 4.6, 4.8, 5.0, 5.2, 5.4, 5.6, 5.8, 6.0, 6.2 ];
obs_1=[102.225 ,99.815,-21.585,-35.099, 2.523,-38.865,-39.020, 89.147, 125.249,-63.405, -183.606, -11.287,197.627, 98.355, -131.977, -129.887, 52.596, 101.193,5.412, -20.805, 6.549, -40.176, -71.425, 57.366, 153.032,5.301, -183.830, -84.612, 159.602, 155.021, -73.318, -146.955];
% 2. lm演算法
% 初始猜測初始點
a0=100; b0=100;
y_init = a0*cos(b0*data_1) + b0*sin(a0*data_1);
% 資料個數
ndata=length(obs_1);
% 引數維數
nparams=2;
% 迭代最大次數
n_iters=60;
% lm演算法的阻尼係數初值
lamda=0.1;
%lm演算法的精度
ep=100
% step1: 變數賦值
updatej=1;
a_est=a0;
b_est=b0;
% step2: 迭代
for it=1:n_iters
if updatej==1
% 根據當前估計值，計算雅克比矩陣，雅可比矩陣只需要在第一次迴圈時計算一次就好
j=zeros(ndata,nparams); % 雅可比矩陣的行數由原始輸入資料個數決定，列數由待估引數個數決定
for i=1:length(data_1)
j(i,:)=[cos(b_est*data_1(i))+data_1(i)*b_est*cos(a_est*data_1(i)) -sin(b_est*data_1(i))*a_est*data_1(i)+sin(a_est*data_1(i)) ]; % 雅可比矩陣由偏導組成
end% 根據當前引數，得到函式值
y_est = a_est*cos(b_est*data_1) + b_est*sin(a_est*data_1);
% 計算誤差
d=obs_1-y_est;
% 計算（擬）海塞矩陣
h=j'*j;
% 若是第一次迭代，計算誤差
if it==1
e=dot(d,d); % 可以認為e是初始值計算所估誤差
endend
% 根據阻尼係數lamda混合得到h矩陣
h_lm=h+(lamda*eye(nparams,nparams));
% 計算步長dp，並根據步長計算新的可能的\引數估計值
dp=inv(h_lm)*(j'*d(:))
%求誤差大小
g = j'*d(:);
a_lm=a_est+dp(1); % 在初始值上加上所求步長，作為新的評估引數
b_lm=b_est+dp(2);
% 計算新的可能估計值對應的y和計算殘差e
y_est_lm = a_lm*cos(b_lm*data_1) + b_lm*sin(a_lm*data_1);
d_lm=obs_1-y_est_lm
e_lm=dot(d_lm,d_lm) % 這個值後面主要用於和上一次誤差進行比對，從而調整阻尼係數
% 根據誤差，決定如何更新引數和阻尼係數
if e_lm
LM演算法詳解
殘差函式f x 為非線性函式，對其一階泰勒近似有 這裡的j是殘差函式f的雅可比矩陣，帶入損失函式的 令其一階導等於0，得 這就是 裡常看到的normal equation。lm是對高斯牛頓法進行了改進，在求解過程中引入了阻尼因子 2.1 阻尼因子的作用 2.2 阻尼因子的初始值選取 乙個簡單的策略就...
LM螞蟻聚類演算法
實驗表明 工蟻能在幾小時內將分散在蟻穴各處的大小不同的螞蟻屍體聚成幾類，小的蟻堆通過吸引螞蟻積攢更多的屍體來逐漸變大，這種正反饋會導致蟻堆逐漸越積越大，以達到聚類資料的目的。deneubourg 等人提出了一種基本模型 basic model 簡稱bm 用來解釋螞蟻屍體堆積成螞蟻墓的行為，並模擬實現...
LM演算法的C 實現
這是乙個資料擬合的例子，並沒有採用物件導向的設計方法是使能更好的理解lm演算法的流程，簡約而不簡單。演算法詳細過程不多介紹。程式中用到opencv庫中的矩陣類mat。例 pragma once include include opencv2 core core.hpp pragma comment ...