非線性最小二乘問題的方法


1.簡介和定義............................... 1
2.設計方法.................................................. 5
  2.1.最陡下降法. ..................... 7
  2.2.牛頓法. ........................................................ 8
  2.3.線搜索...................................................................................... 9
  2.4.信賴域和阻尼方法........................................ 11
3.非線性最小二乘問題........................................ 17
  3.1.高斯-牛頓法........................................ 20
  3.2. Levenberg–Marquardt方法........................................ .24
  3.3.鮑威爾的狗腿法........................................ 29
  3.4.混合方法:LM和擬牛頓.........................................34
  3.5. L–M方法的割線形式........................................ 40
  3.6.狗腿法的一個正割版本........................................ 45
  3.7.最后的評論........................................ 47
附錄........................................ 50
參考資料........................................ 55
索引........................................57

1.引言和定義

在本手冊中,我們考慮以下問題

   定義1.1.  最小二乘問題
  查找x,一個局部最小化器,用於1)

 范例1.1. 最小二乘問題的重要來源是數據擬合。 例如,請考慮以下所示的數據點(t1,y1),...,(tm,ym

圖1.1  數據點{(ti,yi)}(用+標記)和模型M(x,t)(用實線標記)

此外,我們給出了擬合模型,

 模型取決於參數x = [x1,x2,x3,x4]T。 我們假設存在一個x,因此

 {εi}是數據坐標上的(測量)誤差,假定像“白噪聲”一樣。

對於x的任何選擇,我們都可以計算殘差

 對於最小二乘擬合,將參數確定為殘差平方和的最小值x。 可以看出這是定義1.1中n = 4形式的問題。在圖1.1中用實線顯示了M(x,t)的圖。

最小二乘問題是更常見問題的一個特殊變體:給定函數F:IRn→IR,找到參數F,該參數給出該所謂的目標函數成本函數的最小值。

定義1.2  全局最小化器

 一般而言,這個問題很難解決,我們僅介紹解決以下簡單問題的方法:找到F的局部極小值,這是一個自變量矢量,在某個區域內給出了F的最小值,其大小由δ給出,其中 δ是一個小的正數。

定義1.3 本地最小化器

 在本介紹的其余部分中,我們將討論優化中的一些基本概念,第2章簡要介紹了為一般成本函數找到局部最小化器的方法。 有關更多詳細信息,請參閱Frandsen等人(2004年)。 在第3章中,我們提供了針對最小平方問題專門調整的方法。

我們假設成本函數F是可微的且非常平滑,以至於隨后的泰勒展開式是有效的,2

 其中g是梯度,

 H是Hessian(在數學中,Hessian矩陣或Hessian是標量值函數或標量場的二階偏導數的平方矩陣。 它描述了許多變量的函數的局部曲率。)

如果x是局部極小值並且||h||足夠小,則我們找不到具有較小F值的點x + h。 將此觀察值與(1.4a)結合,我們得到

定理1.5 局部最小化器的必要條件。

如果x是局部極小值,則

 對於滿足必要條件的參數,我們使用特殊名稱:

定義1.6 固定點。 如果

 那么xs被認為是F的固定點。

因此,局部最小化器也是固定點,但是局部最大化器也是如此。 既不是局部最大化器也不是局部最小化器的固定點稱為鞍點。 為了確定給定的固定點是否是局部極小值,我們需要在泰勒級數(1.4a)中包括二階項。 插入xs我們看到

 根據Hessian的定義(1.4c),可以得出任何H是對稱矩陣。 如果我們要求Hs是正定的,則其特征值大於某個值δ> 0(請參閱附錄A),並且

 這表明,對於足夠小的(1.7)右側的第三項將由第二項占主導。 這個詞是肯定的,所以我們得到

定理1.8 局部最小化器的充分條件

假設xs是一個固定點,並且F''(xs)是正定的。 那么xs是一個局部最小化器。

如果Hs是負定數,則xs是局部最大化器。 如果Hs是不確定的(即,它同時具有正和負特征值),則xs是一個鞍點。

2.設計方法

  所有用於非線性優化的方法都是迭代的:從起點x0開始,該方法會產生一系列向量x1,x2,...,(希望地)收斂為給定函數的局部極小值x,請參見定義1.3大多數方法都有強制降低條件的措施

 這會阻止收斂到最大化器,也使收斂到鞍點的可能性降低。 如果給定的函數具有多個最小化器,則結果將取決於起點x0。 我們不知道會找到哪個最小化器。 它不一定是最接近x0的最小化器。

在許多情況下,該方法產生的向量在兩個明顯不同的階段收斂到最小化器。 當x0遠離解時,我們希望該方法產生迭代,並朝x穩定地移動。 在迭代的這個“全局階段”,如果除了第一步之外的錯誤沒有增加,我們就感到滿意,即

 其中ek表示當前錯誤,

 在迭代的最后階段,xk接近x ∗,我們想要更快的收斂速度。 我們區分

線性收斂(Linear convergence):

 二次收斂(Quadratic convergence):

 超線性收斂(Superlinear convergence):

 本講義中介紹的方法是下降方法,在迭代的每個步驟中都滿足下降條件(2.1)當前迭代的第一步在於

  1.找到下降方向hd(如下所述),然后
  2.找到一個步長,使F值大大降低。

因此,下降法的概述是

 考慮F值沿着從x開始並沿h方向的一半線的變化。 從泰勒展開式(1.4a)中,我們看到

 我們說,如果F(xh)是α在α= 0時的遞減函數,則h是下降方向。這導致以下定義。

 定義2.6. 下降方向

 如果不存在這樣的h,則F'(x)= 0,這表明x在這種情況下是固定的。 否則,我們必須選擇α,即在hd給出的方向上距x應該走多遠,以便使目標函數的值減小。 一種方法是找到(近似)

 該過程稱為行搜索,將在第2.3節中討論。 但是,首先,我們將介紹兩種計算下降方向的方法。

2.1.The Steepest Descent method最陡下降法

從(2.5)可以看出,當我們以正α執行步長αh時,函數值的相對增益就可以滿足

 其中θ是向量h和F'(x)之間的夾角。 這表明如果θ=π,即使用以下公式給出的最陡下降方向hsd,我們將獲得最大增益率:

 基於(2.8)的方法(即算法2.4中的hd = hsd)被稱為最速下降法或梯度法。 下降方向的選擇(局部)是“最佳”,我們可以將其與精確的線搜索(2.7)結合起來。 這樣的方法會收斂,但是最終收斂是線性的,並且通常很慢。 Frandsen等人(2004年)中的示例表明,具有精確線搜索和有限計算機精度的最速下降方法如何無法找到二次多項式的極小值。 但是,對於許多問題,該方法在迭代過程的初始階段具有相當好的性能。

諸如此類的考慮導致了所謂的混合方法,顧名思義,該方法基於兩種不同的方法。 一種在初始階段很好,例如梯度法,另一種在最后階段很好的方法,例如牛頓法。 請參閱下一節。 混合方法的主要問題是在適當時在兩種方法之間切換的機制。

2.2. 牛頓法 Newton’s Method

我們可以從x*是一個固定點的條件中得出這種方法。 根據定義1.6,它滿足F'(x*)=0。這是一個非線性方程組,根據泰勒展開式

 我們推導出牛頓的方法:找到hn作為解

 並計算下一個迭代

 假設H是正定的,則它是非奇異的(暗示(2.9a)具有唯一解),並且對於所有非零u,。 因此,通過在(2.9a)的兩邊乘以h> n,我們得到

 表明hn是下降方向:它滿足定義2.6中的條件。

牛頓的方法在迭代的最后階段非常好,其中x接近x ∗。 可以證明(參見Frandsen等人(2004)),如果解中的Hessian是正定的(滿足定理1.8中的充分條件),並且我們在x周圍的區域內的位置是F''( x)是正定的,則我們得到二次收斂(在(2.3)中定義)。 另一方面,如果x處在F''(x)到處都是負定且存在固定點的區域,則基本牛頓法(2.9)會(正交)收斂於該固定點,即 最大化器。 我們可以通過要求所有步驟都沿下降方向來避免這種情況。

我們可以基於牛頓法和最速下降法建立一種混合方法。 根據(2.10),如果F00(x)是正定的,則牛頓階躍保證為下坡,因此該混合算法的中心部分的草圖可以為

 在此,hsd是最陡的下降方向,並且通過線搜索找到α; 請參閱第2.3節。 Cholesky方法(請參閱附錄A)是檢查矩陣的正定性的一個很好的工具,該方法成功后也可用於求解所討論的線性系統。 因此,確定性檢查幾乎是免費的。

在第2.4節中,我們介紹了一些方法,其中搜索方向hd和步長α的計算是同時進行的,並給出了(2.11)的版本,而沒有行搜索。 這種混合方法可能非常有效,但是幾乎從未使用過。 原因是它們需要F''(x)的實現,對於復雜的應用程序問題,這是不可用的。 取而代之的是,我們可以根據逐漸接近H* = F''(x*)的一系列矩陣使用所謂的擬牛頓法。 在第3.4節中,我們介紹了這種方法。 另請參閱Frandsen等人(2004年)的第5章。

2.3. 線搜索Line Search

給定點x和下降方向h。 下一個迭代步驟是從x向h方向移動。 為了找出移動的距離,我們研究了給定函數沿x沿h方向從x的一半線的變化,

 2.1(α)的行為示例如圖2.1所示。

 圖2.1. 成本函數沿搜索線的變化。

我們h是下降的方向,可以確保

 表示如果α足夠小,我們滿足下降條件(2.1),它等於

 通常,我們會對α進行初步猜測,例如,使用牛頓法,α= 1。 圖2.1說明了可能出現的三種不同情況

  1 ◦ α很小,以至於目標函數的值增益很小。 α應增加。
  2 ◦ α太大:ϕ(α)≥ϕ(0), 減小α以滿足下降條件(2.1)。
  3 α接近ϕ(α)的極小值 1)。 接受此α值。

    1)更准確地說:local的最小局部最小值。 如果我們將α增加到超過圖2.1中所示的間隔,則很可能我們接近F的另一個局部最小值。

精確線搜索是產生一系列α1,α2...的迭代過程。 目的是找到在(2.7)中定義的真正的最小化子αe,並且算法在迭代的αs滿足時停止

 其中τ是一個小的正數。 在迭代中,我們可以基於計算的值對ϕ(α)的變化使用近似值。

 有關詳細信息,請參見Frandsen等人(2004年)的第2.5 – 2.6節。

精確的行搜索會浪費大量的計算時間:當x遠離x ∗時,搜索方向h可能會偏離x ∗ −x方向,因此無需非常精確地找到the的真實最小值。 這就是所謂的“軟線搜索”的背景,如果它不屬於上面列出的1或2類別,則我們接受該值。 我們使用降級條件(2.1)的嚴格版本,即

 這確保了我們不會處於情況2中。 情況1對應於點(α,ϕ(α))太接近起始切線,我們補充了條件

 如果對α的初始猜測滿足這兩個標准,那么我們將其視為αs。 否則,我們必須按照草圖進行迭代以進行精確的線搜索。 有關詳細信息,請參見Frandsen等人(2004年)的2.5節。

2.4. 信任區域和阻尼方法Trust Region and Damped Methods

假設我們在當前迭代x的附近有一個F行為的模型L,

 其中和矩陣是對稱的。 本節的基本思想可以推廣到模型的其他形式,但是在本手冊中,我們只需要(2.14)中給出的L形式。 通常,模型是F圍繞x的二階泰勒展開,就像(1.4a)右側的前三個項或L(h)可能是該展開的近似值。 通常,只有當h足夠小時,這種模型才是好的。 我們將在確定步長h的過程中引入包括這方面的兩種方法,步長h是下降方向,可以在算法2.4中與α= 1一起使用。

在信任區域方法中,我們假設我們知道一個正數∆,以使模型在半徑為∆的球(以x為中心)內足夠精確,並將步長確定為

 在阻尼方法中,該步驟確定為

 其中阻尼參數µ≥0。認為會懲罰較大的步長

基於這些方法之一的算法2.4的中心部分具有以下形式

如果步長h滿足下降條件(2.1),則這對應於α= 1。 否則,α= 0,即我們不動2)但是,我們不會陷入x的困境(除非x = x ∗):通過對∆或µ進行適當的修改,我們的目標是在下一個迭代步驟中獲得更好的運氣。

 2)這些方法有多種版本,其中包括適當的線搜索以查找具有較小F值的點x +αh,並且在線搜索過程中收集的信息用於∆或µ的更新。 對於許多問題,此類版本使用較少的迭代步驟,但使用更多的函數值累積量。

由於對於h足夠小,L(h)被認為是F(x + h)的良好近似值,因此步失敗的原因是h太大,應將其減小。 此外,如果接受該步驟,則有可能使用新迭代中的較大步驟,從而減少我們達到x*之前所需的步驟數。

具有計算步驟的模型的質量可以通過所謂的增益比來評估

 即功能值實際減少與預測減少之間的比率。 通過構造,分母為正,如果台階不是下坡,則分子為負–太大,應減小。 使用信任區域方法,我們通過半徑∆的大小監視步長。 以下更新策略被廣泛使用,

 因此,如果,我們決定使用較小的步長,而%> 34表示可能使用較大的步長。 信任區域算法對閾值0.25和0.75,除數p1 = 2或因數p2 = 3的微小變化不敏感,但是選擇數字p1和p2以使Δ值不發生振盪很重要。

在阻尼方法中,%的較小值表示我們應該增加阻尼系數,從而增加大步幅的代價。 %的大值表示對於計算出的h,L(h)非常接近F(x + h),並且可以減小阻尼。 以下是一種廣泛使用的策略,它與(2.19)類似,最初是由Marquardt(1963)提出的,

 同樣,該方法對閾值0.25和0.75或數字p1 = 2和p2 = 3的微小變化不敏感,但重要的是選擇數字p1和p2,以使µ值不會振盪。 經驗表明,跨越閾值0.25和0.75的不連續變化會引起“波動”(如第27頁的示例3.7所示),這會減慢收斂速度,我們在Nielsen(1999)中證明,以下策略在總體上勝過 (2.20),

 將因數ν初始化為ν=2。請注意,一系列連續的故障導致µ值迅速增加。 這兩個更新公式如下所示。

 

 圖2.2. Mar的策略(2.20)(虛線)用ν= 2(實線)將µ更新為(2.21)。

2.4.1. 步驟計算Computation of the step.

 在阻尼方法中,該步被計算為函數的固定點

 

 

 這意味着hdm是解決

 

 

 從(2.14)中L(h)的定義中我們看到,這等效於

 

 

 我是單位矩陣 如果µ足夠大,則對稱矩陣B + µI是正定的(如附錄A所示),然后從定理1.8可以看出hdmL的極小值

示例2.1 在阻尼牛頓法中,模型L(h)由c = F'(x)和B = F''(x)給出,而(2.22)的形式為

 即在接近最陡下降方向的方向上的短步。 另一方面,如果μ很小,則hdn接近牛頓步距hn。 因此,我們可以認為阻尼牛頓法是最速下降法和牛頓法的混合體。

我們將在第3.2節中返回阻尼方法。

在信任區域方法中,步驟htr是約束優化問題的解決方案,

 

 

 對此問題進行任何詳細討論都超出了本手冊的范圍(請參閱Madsen等人(2004年)或Nocedal and Wright(1999年)中的4.1節。我們只想提及一些特性。

如果(2.14)中的矩陣B是正定的,則L的無約束最小化子為

 

 

 如果該值足夠小(如果滿足hT h≤∆2),則這是所需的步驟htr。 否則,約束將處於活動狀態,並且問題更加復雜。 使用與第11頁類似的論點,我們可以看到我們不必計算(2.23)的真實解,並且在第3.3和3.4節中,我們提供了兩種計算htr近似值的方法。

最后,在B為正定的情況下,我們在阻尼方法和信任區域方法之間存在兩個相似之處:如果無約束極小子在信任區域之外,則可以證明(Madsen等人(2004年,定理2.11 ))存在一個λ> 0使得

 

 

 通過對該方程重新排序並將其與(2.22)進行比較,我們可以看到htr與使用阻尼參數µ =λ計算的阻尼階躍hdm相同。 另一方面,還可以證明(Frandsen等人(2004年)中的定理5.11),如果我們針對給定的µ≥0計算hdm,則

 

 

 即hdm等於對應於信任區域半徑∆ = ||hdm||的htr。 因此,這兩種方法密切相關,但是對於給出相同步驟的∆值和µ值之間的聯系,沒有簡單的公式。

3. 非線性最小二乘問題NON-LINEAR LEAST SQUARES PROBLEMS

在本講義的其余部分,我們將討論非線性最小二乘問題的方法。 給定一個向量函數f:,m≥n。 我們要最小化||f(x)||或等效地找到

 

 最小二乘問題可以通過一般的優化方法來解決,但我們將介紹更有效的特殊方法。 在許多情況下,即使不需要二階導數,它們也比線性收斂,有時甚至是二次收斂要好。 在本章方法的描述中,我們將需要F的導數的公式:假設f具有連續的二階偏導數,我們可以將其泰勒展開式寫為

 

 其中是雅可比行列式。 這是一個矩陣,其中包含功能組件的一階偏導數,

 

 關於F:,它來自(3.1b)中的第一個公式,

 那1)

 

 1)如果在定義(3.1b)中未使用因子1/2,則在許多表達式中我們將得到令人討厭的因子2。

 

 因此,梯度(1.4b)為

 

 我們還需要F的Hessian。從(3.3)可以看出,位置(j,k)的元素是

 

 顯示

 

 示例3.1 (3.1)的最簡單情況是f(x)的形式為

 

給出了向量和矩陣。 我們說這是一個線性最小二乘問題。 在這種情況下,對於所有x,J(x)= −A,從(3.4a)中我們看到

 

 

 對於確定為所謂正態方程的解的x,該值為零。

 

 

 問題可以以以下形式寫

 

 

 或者我們可以通過正交變換來解決它:找到一個正交矩陣Q使得

 

 

 其中是上三角。 通過系統中的反向替換找到解決方案2)

 

 

 2)像up:q這樣的表達式用於表示元素ui為i = p,的子向量... ,q. 矩陣A的第i行和第j列分別表示為Ai ,:和A,J,

 

 

 該方法比通過正態方程的解決方案更准確

在MATLAB中,假定數組A和b分別包含矩陣A和向量b。 然后,命令A \ b返回通過正交變換計算的最小二乘解。

正如手冊標題所暗示的,我們假定f是非線性的,因此將不詳細討論線性問題。 我們參考Madsen和Nielsen(2002)中的第2章,或Golub和Van Loan(1996)中的5.2節。

示例3.2 在示例1.1中,我們看到了由數據擬合引起的非線性最小二乘問題。 另一個應用是解決非線性方程組,

 

 

 我們可以使用牛頓-拉夫森(Newton-Raphson)的方法:根據初始猜測x0,通過以下算法計算x1,x2,...,該算法基於尋找h使得f(x + h)= 0並忽略術語O(||h||2 )在(3.2a)中,

 

 

 這里,雅可比定律J由(3.2b)給出。 如果J(x ∗)是非奇異的,則該方法具有二次最終收斂性,即,如果小,則 但是,如果xk遠離x ∗,那么我們就有可能走得更遠。

我們可以重新構造問題,使其能夠使用本章將要介紹的所有“工具”:(3.6)的解決方案是由(3.1)定義的函數F的全局最小值,

 

 

 因為如果f(x)= 0,則F(x*)= 0且F(x)>0。例如,我們可以將(3.6)中的近似解的更新替換為

 

 

 其中,通過對函數進行線搜索來找到αk。 作為一個具體示例,我們將考慮以下問題,該問題取自Powell(1970),

 

 

 x = 0是唯一的解決方案。 雅可比式是

 

 

 在解決方案上是唯一的。

如果我們取並使用上述算法進行精確行搜索,則迭代收斂到,這不是解決方案。 另一方面,很容易看出,算法(3.6)給出的迭代為,其中,即我們對解具有線性收斂性。 在許多示例中,我們將返回到該問題,以查看不同的方法如何處理它。

3.1. 高斯-牛頓法The Gauss–Newton Method

此方法是我們將在下一部分中介紹的非常有效的方法的基礎。 它基於向量函數各成分的已實現一階導數。 在特殊情況下,它可以給出二次收斂,就像牛頓法進行一般優化一樣,請參見Frandsen等人(2004年)。

高斯-牛頓法基於x附近的f的線性近似(f的線性模型):對於小||h||,我們從泰勒展開式(3.2)中得出:

 

 

 將其插入F的定義(3.1)中,我們看到

 

 

 (其中f = f(x)且J = J(x))。 高斯-牛頓步長hgn將L(h)最小化,

 

 

 可以很容易地看出L的梯度和Hessian為

 

 

 與(3.4a)的比較表明L'(0)= F'(x)。 此外,我們看到矩陣L''(h)獨立於h。 它是對稱的,並且如果J具有最高秩,即,如果列是線性獨立的,則L''(h)也是正定的,請參閱附錄A。這意味着L(h)具有唯一的極小值,可以通過求解來找到

 

 

 這是F的下降方向,因為

 

 

 

 因此,我們可以在算法2.4中將hgn用於hd。 典型的步驟是

 

 

 通過行搜索找到α。 經典的高斯-牛頓法在所有步驟中都使用α= 1。 可以證明,采用線搜索的方法具有瓜爾特事先收斂性,條件是

  a){x | F(x)≤F(x0)}有界,並且
  b)雅可比J(x)在所有步驟中都具有最高等級

在第二章中,我們看到牛頓的優化方法具有二次收斂性。 高斯-牛頓法通常不是這種情況。

為此,我們比較了兩種方法中使用的搜索方向,

 

 

 我們已經在(3.8)處指出兩個右側相同,但是從(3.4b)和(3.8)中我們看到系數矩陣不同:

 

 

 因此,如果f(x)= 0,則x的接近x*,並且我們也用高斯-牛頓法得到二次收斂。 如果函數{fi}的曲率小或小,我們可以期望超線性收斂,但是總的來說我們必須期望線性收斂。 值得注意的是,F(x)的值控制着收斂速度。

示例3.3 考慮n = 1的簡單問題,m = 2由

 

 

 它遵循

 

 

 所以x = 0是F的固定點。現在,

這表明如果λ<1,則F''(0)> 0,因此x = 0是局部最小化器-實際上,它是全局最小化器。

雅可比式是

 

 

 xk的經典高斯-牛頓法給出

 

 

 現在,如果且xk接近零,則

 

 

 因此,如果|λ| <1,我們有線性收斂。 如果λ<-1,則經典高斯牛頓法無法找到最小化器。 例如,當λ= − 2且x0 = 0.1時,我們得到了迭代器看似混亂的行為,

 

 

 最后,如果λ= 0,則

 

 

 即,我們一步找到解決方案。 原因是在這種情況下,f是線性函數。

示例3.4 對於示例1.1中的數據擬合問題,雅可比矩陣的第i行為

 

 

 如果問題是一致的,那么只要x1 1與x2 顯着不同,則帶線搜索的高斯-牛頓法將具有二次最終收斂性。 如果,則,並且高斯牛頓法失敗。

如果一個或多個測量誤差較大,則f(x ∗)的分量較大,這可能會降低收斂速度。

在MATLAB中,我們可以提供一個非常緊湊的函數來計算f和J:假設x保持當前迭代,並且m×2數組ty保持數據點的坐標。 以下函數返回分別包含f(x)和J(x)的f和J。

 

 

 示例3.5 考慮示例3.2中的問題,其中f(x ∗)= 0,其中f:。 如果我們使用Newton-Raphson的方法來解決此問題,則典型的迭代步驟是

 

 

 用於最小化的高斯-牛頓法具有典型步驟

注意,J(x)是一個方矩陣,並且我們假定它不是奇異的。 然后存在,並由此得出hgn = hnr。 因此,將高斯-牛頓法應用於例3.2中的鮑威爾問題時,將遇到與該例中的牛頓-拉夫森法相同的麻煩

這些示例表明,無論是否進行線搜索,高斯-牛頓法都可能會失敗。 盡管如此,它在許多應用中仍具有相當好的性能,盡管通常它僅具有線性收斂性,而與牛頓方法中實現二次導數的二次收斂性相反。

在3.2和3.3節中,我們給出了兩種具有較高全局性能的方法,而在3.4節中,我們對第一種方法進行了修改,以實現最終的超線性收斂。

3.2. The Levenberg–Marquardt Method

Levenberg(1944)和后來的Marquardt(1963)建議使用阻尼高斯-牛頓法,請參見第2.4節。 通過對(3.9)的以下修改來定義步驟hlm

 

 

 在此,J = J(x),f = f(x)。 阻尼參數µ有幾個作用:

  a)對於所有μ> 0,系數矩陣為正定的,這可確保hlm為下降方向,比照(3.10)。

  b)對於較大的µ,我們得到

    

     即沿最陡的下降方向的一小步。 如果當前的迭代距離解還很遠,那么這很好。

 

 

   c)如果µ很小,則當x接近x ∗時,hlm'hgn是迭代最后階段的一個好步驟。 如果F(x)= 0(或很小),那么我們可以(幾乎)獲得二次最終收斂。

因此,阻尼參數會影響步長的方向和大小,從而導致我們做出一種無需特定線搜索的方法。 初始µ值的選擇應與中元素的大小有關,例如

 

 

其中τ由用戶選擇3)在迭代過程中,μ的大小可以按照2.4節所述進行更新。 更新由增益比控制

 

 

     3)該算法對τ的選擇不是很敏感,但是根據經驗,如果x0被認為是對x ∗的良好近似,則應使用較小的值,例如τ= 10-6。 否則,使用τ= 10-3甚至τ= 1。

 

 

 分母是線性模型(3.7b)預測的增益,

 

 

 注意,均為正,因此保證為正。

%的大值表示L(hlm)是F(x + hlm)的良好近似值,我們可以減小µ,以便下一個Levenberg-Marquardt步驟更接近於Gauss-Newton步驟。 如果%小(甚至可能為負),則L(hlm)的近似值很差,我們應該增加µ的雙重目的是接近最陡的下降方向並減小步長。 這些目標可以通過不同的方式實現,請參見第14頁和下面的示例3.7。

該算法的停止標准應反映出,在全局最小化器中,我們有,因此我們可以使用

 

 

 其中ε1是用戶選擇的較小的正數。 另一個相關標准是,如果x的變化很小,則停止;

 

 

 該表達式給出了從漸進變化的相對步長ε2(當|x|大時)到絕對步長(如果x接近0)的變化。最后,像在所有迭代過程中一樣,我們需要防范無限循環的保護措施,

 

 

 ε2和kmax也由用戶選擇。

最后兩個標准生效,例如,如果ε1選擇得太小以至於舍入誤差的影響會很大。 通常,這將使F的實際增益與線性模型預測的增益(3.7b)之間的符合性差,從而導致µ在每個步驟中都增大。 增大μ的策略(2.21)表示在這種情況下,μ增長很快,導致khlmk小,並且該過程將被(3.15b)停止。

該算法總結如下。

示例3.6 通過比較(3.9)和正規方程(3.5),我們可以看到hgn只是線性問題的最小二乘解

 

 

 同樣,L-M方程(3.13)是線性問題的正規方程

 

 

 

 

 

 如示例3.1所述,通過正交變換可以找到最准確的解決方案。 但是,解hlm只是迭代過程中的一個步驟,不需要非常精確地計算,並且由於通過正態方程的解是“便宜”的,因此通常使用此方法。

示例3.7 我們已對示例1.1和3.4中的數據擬合問題使用了算法3.16。 圖1.1表示x1和x2均為負,並且

滿足這些條件。 此外,我們在表達式(3.14)中將τ= 10−3用於µ0,將停止條件由(3.15)給出,其中ε12= 10−8,kmax =200。該算法在使用x進行62次迭代之后停止 。 性能如下所示; 注意對數縱坐標軸。

這個問題不一致,因此我們可以期望線性最終收斂。 最后的7個迭代步驟表明收斂性更好(超線性)。 推論是,是ti的緩慢變化函數,具有“隨機”符號,因此對(3.12)中的“遺忘項”的貢獻幾乎抵消了 出來。 這種情況發生在許多數據擬合應用程序中。

 

 圖3.2a。 L-M方法適用於示例1.1中的擬合問題。

為了進行比較,圖3.2b顯示了使用更新策略(2.20)時的性能。 從步驟5到步驟68,我們看到µ的每個減小都緊隨其后增大,並且梯度的范數具有粗糙的行為。 這減慢了收斂速度,但是最后階段如圖3.2a所示。

 

 圖3.2b。 使用更新策略的效果(2.20)

示例3.8 圖3.3從示例3.2和3.5展示了應用於鮑威爾問題的算法3.16的性能。 起點是x0 = [3,1]>,由(3.14)中的τ= 1給出的µ0,在停止標准(3.15)中,我們使用ε1=ε2= 10-15,kmax = 100。

圖3.3。 L-M方法適用於鮑威爾的問題。

迭代似乎在步驟22和30之間停止。這是(幾乎)奇異雅可比矩陣的影響。 在那之后,似乎出現了線性收斂。 迭代由“保護”在點x = [-3.82e-08,-1.38e-03]T處停止。 與示例3.2相比,這是對x = 0的更好近似,但是我們希望能夠做得更好。 參見示例3.10和3.17。

3.3.Powell’s Dog Leg Method鮑威爾的狗腿法

作為Levenberg–Marquardt方法,此方法適用於高斯–牛頓和最陡的下降方向的組合。 現在,可通過信任區域的半徑進行顯式控制,請參見第2.4節。 鮑威爾的名字與算法有關,因為他提出了如何找到由(2.23)定義的htr的近似值。

給定。 在當前迭代x高斯-牛頓步長hgn是線性系統的最小二乘解

 

 可以通過求解正規方程來計算

 

 最陡的下降方向為

 

 這是一個方向,而不是一個步驟,要了解我們應該走多遠,我們看一下線性模型

 

 α的這個函數對於

 

 現在,對於從當前點x采取的步驟,我們有兩個候選項:a =αhsd和b = hgn。 當信任區域的半徑為∆時,鮑威爾建議使用以下策略選擇步驟。 該策略的最后一種情況如圖3.4所示。

 

 

 圖3.4  信任區域和狗腿步4)

 

狗腿這個名字取自高爾夫:“狗腿洞”處的球道的形狀是從x(發球點)到a的終點到hdl(洞)的終點的線。 鮑威爾是一位敏銳的高爾夫球手!

使用上面定義的a和b,以及,我們可以寫

 

 

 我們為該二階多項式求根,請注意對於β→−∞,ψ→+∞;。因此,ψ在] 0,1 [中具有一個負根和一個根。 我們尋找后者,最精確的計算由

 

 

 與L-M方法一樣,我們可以使用增益比

 

 

 監視迭代。 同樣,L是線性模型

 

 

 在L-M方法中,我們使用%來控制阻尼參數的大小。 在這里,我們使用它來控制信任區域的半徑∆。 %的大值表示線性模型良好。 我們可以增加∆從而采取更長的步驟,並且它們將更接近高斯-牛頓方向。 如果%小(甚至可能為負),則我們減小∆,這意味着更小的步長,更接近最陡的下降方向。 下面我們總結一下算法。

我們有以下幾點評論。

  1◦初始化。 x0和∆0應由用戶提供。
  2◦我們使用停止條件(3.15)加上kf(x)k∞≤ε3,反映出在m = n的情況下f(x ∗)= 0,即非線性方程組。
  3◦如果m = n,則用“ =”代替“'”,參見(3.6),並且我們不使用繞等式(3.18a)繞行的方法; 參見例3.9。

  4◦與(3.20a)中的三種情況相對應,我們可以證明

    

 

 

   5◦策略(2.19)用於更新信任區半徑。

  6◦額外的停止條件。 如果,則在下一步驟中肯定會滿足(3.15b)。

 

 

 示例3.9 在示例3.6中,我們簡要討論了步驟hlm的計算,並認為我們也可以通過法線方程(3.13)對其進行計算。 如果μ不是很小,則矩陣條件良好,並且舍入誤差不會有太大影響。

狗腿法的目的是在非線性方程組上也表現出色,即(3.17)是線性方程組的平方系統

 

 

 溶液h = hnr,則采用牛頓-拉夫森步驟,請參見示例3.2。 Ja cobian J可能是病態的(甚至是奇異的),在這種情況下,舍入誤差往往是解決方案的主導。 如果我們使用(3.18a)來計算hgn,則會使此問題惡化。

在immoptibox的實現方法中,針對這些問題計算了(3.17)的解決方案。 如果J(x)的列不是線性獨立的,則最小二乘解h不是唯一的,並且hgn被計算為具有最小范數的h。 此計算的一些詳細信息在附錄B中給出。

示例3.10 圖3.5說明了從示例3.2和3.8中以起點x0 = [3,1]>,∆0 = 1且由ε1=ε2= 10-15,ε3給出的停止准則,適用於鮑威爾問題的Dog Leg方法的性能 = 10-20,kmax = 100。

 

 

 圖3.5  狗腿法適用於鮑威爾的問題。

由於梯度小,迭代在37個步驟后停止,並返回x = [-2.41·10-35,1.26·10-9]>,這非常近似於x ∗ =0。如圖3.3所示, 最終收斂是線性的(由奇異的J(x ∗)引起),但是比Marquardt方法要快得多。

示例3.11 對於示例1.1、3.4和3.7中的數據擬合問題,我們使用了算法3.21。 如例3.7所示,我們使用起點x0 = [−1,-2,1,-1]>,取∆0 = 1且由ε123= 10-8,kmax = 200給出的停止標准 該算法在經過30個迭代步驟后以x'[−4,-5,4,-4]>停止。 性能說明如下。 如圖3.3所示,我們注意到非常快的最終收斂速度。

 

 

 圖3.6。 狗腿法適用於示例1.1中的擬合問題。

最后兩個例子似乎表明,“狗腿法”比“ Levenberg-Marquardt法”要好得多。 當最小二乘問題來自非線性方程組時,這是正確的。 狗腿法目前被認為是求解非線性方程組的最佳方法。

對於一般的最小二乘問題,狗腿法與LM方法具有相同的劣勢:如果F(x ∗)= 0,則最終收斂可以期望是線性的(且緩慢)。對於給定的問題和給定的開始猜測 x0不可能事先說出兩種方法中哪一種更快。

3.4. 混合方法:LM和擬牛頓A Hybrid Method: L–M and Quasi–Newton

1988年,Madsen提出了一種混合方法,將L–M方法(如果F(x ∗)= 0則為二次收斂,否則為線性收斂)與Quasi5)-Newton方法相結合,即使F(x ∗)也具有超線性收斂。 =0。迭代從LM方法的一系列步驟開始。 如果性能表明F(x ∗)明顯非零,則我們改用准牛頓法以獲得更好的性能。 我們可能會發現有跡象表明,最好改回到LM方法,因此也有一種機制。

 

 

   5)源自拉丁語:“准” =“幾乎”。 有關准牛頓法的一般介紹,請參見Frandsen等人(2004年)的第5章。

如果條件滿足,則切換到擬牛頓法

 

 

 在三個連續成功的迭代步驟中得到滿足。 這被解釋為表明我們正在接近F *(x *)= 0且F(x *)明顯非零的x *。 如結合(3.12)所述,這可能導致緩慢的線性收斂。

擬牛頓法基於在當前迭代x處具有近似於Hessian F''(x)的B,並且通過求解可以找到步驟hqn

 

 這是牛頓方程(2.9a)的近似值。

近似值B通過BFGS策略進行更新,請參見Frandsen等人(2004年)的5.10節:近似值矩陣序列中的每個B都是對稱的(與任何F''(x)一樣)並且是正定的。 這樣可確保hqn為“下坡”,參見(2.10)。 我們從對稱的正定矩陣B0 = I開始,BFGS更新由要添加到當前B的2級矩陣組成。Madsen(1988)使用以下版本,由Al-Baali和Fletcher(1985)提出 ,

 

 其中J = J(x),Jnew = J(xnew)。 如前所述,電流B是正定的,並且僅在h> y> 0時才改變。在這種情況下,可以證明新的B也是正定的。

擬牛頓法在迭代的全局階段並不穩健。 它不能保證下降。 在解x ∗處,我們有F'(x ∗)= 0,並且迅速減小的值指示出良好的最終收斂性。 如果這些標准值的下降速度不夠快,那么我們將切換回L–M方法。

該算法總結如下。 它調用輔助函數LMstep和QNstep,實現這兩種方法。

 

 我們有以下幾點評論:

  1◦初始化。 可以通過(3.14)找到µ0。 停止標准由(3.15)給出。

  2◦點表示我們還傳輸了f和J等的當前值,因此我們不必為相同的x重新計算它們。

  3◦請注意,L-M和擬牛頓步驟都為Hessian矩陣的近似提供了信息。

下面給出了兩個輔助功能:

 

 

 我們對LMstep和QNstep函數有以下說明:

  4◦像算法3.16一樣,增益比%也用於更新µ。

  5◦指示可能是時候切換方法了。 在算法3.25開始時,參數計數被初始化為零。

  6 ◦  在三個連續的迭代步驟中都滿足了(3.22),所有這些步驟的%> 0,即x在這些步驟的每一個中都發生了變化。

  7◦我們將擬牛頓法與信任區域法結合在一起,並對邊界有效的情況進行了簡單處理,請參見第15f頁。 在從L–M方法切換時,將∆初始化為

  8◦未顯示:∆通過(2.19)更新。

   9◦在算法的這一部分中,我們着重於使F0接近零,因此我們接受F值的輕微增加,例如δ=√εM,其中εM是計算機的單位舍入。

  10◦梯度下降得不夠快。

示例3.12 注意,在更新公式(3.24)中,y的計算涉及乘積。 這意味着我們必須存儲先前的雅可比矩陣。 相反,我們可以使用

 

 在更新公式中,但Madsen(1988)發現(3.24)表現更好。 Madsen(1988)中沒有包括Q-N步驟中的信任區域方法,但是在immoptibox函數nlshybrid的開發過程中,發現這種想法可以提高性能。 它減少了徒勞嘗試Q-N步驟的次數,即10°的條件立即返回到L-M方法。

示例3.13 對於示例3.7和3.8中討論的問題,此混合方法的性能不會優於算法3.16。 在后一種情況下(參見圖3.3),F(x)→0,並且永遠不會滿足注釋5◦處的切換條件。 在前一種情況下,F(x ∗)明顯為非零,但是-如示例3.7中所述-簡單的LM方法具有所需的超線性最終收斂。

為了證明算法3.25的效率,我們考慮改進的Rosen brock問題,請參見Frandsen等人(1999年)的示例5.5,由f給出:

 

 

 

 其中可以選擇參數λ。的極小值是,其中

下面,我們給出一些值λ的算法3.16和3.25的結果。 在所有情況下,我們使用x0 = [-1.2,1]>,初始阻尼參數µ0由(3.14)中的τ= 10-3定義,並且(ε1,ε2,kmax)=(10-10,10-14, 200)中的停止條件(3.15)。

 在前兩種情況下,λ太小而不能真正影響迭代,但是對於較大的λ值,我們看到混合方法比簡單的Levenberg-Marquardt算法要好得多-特別是在獲得的精度方面。

在圖3.7中,我們說明了兩種算法在λ= 104的情況下的性能。

 

 圖3.7.  Levenberg–Marquardt方法(左)和混合方法(右)

使用L–M方法后,所有步驟均不執行。 15未能改善目標功能; µ迅速增加,並且在步驟no滿足停止標准(3.15b)。 23.對於混合方法,從步驟Nos開始,嘗試了幾種使用擬牛頓法的方法。 5、11和17。最后一次嘗試成功,並且在22個步驟之后,迭代(3.15a)停止。

 3.5. L–M方法的割線版本A Secant Version of the L–M Method

 本手冊中討論的方法假定向量函數f是可微的,即雅可比行列式

 

 

 

 存在。 在許多實際的優化問題中,碰巧我們無法給出J中元素的公式,例如,因為f是由“黑匣子”給出的。 LM方法的割線版本適用於此類問題。

 最簡單的解決方法是將J(x)替換為通過數值微分獲得的矩陣B:第(i,j)個元素通過有限差分近似法近似

 

 

 

 其中ej是第j坐標方向上的單位矢量,而δ是適當小的實數。 使用這種策略,每個迭代x都需要對f進行n + 1次評估,並且由於δ可能比距離kx-x ∗ k小得多,因此我們獲得的關於f的全局行為的信息不會比僅評估時要多 f(x)。 我們想要更好的效率。

 示例3.14 令m = n = 1並考慮一個非線性方程

 

 

 

 對於這個問題,我們可以將牛頓-拉夫森算法(3.6)編寫為

 

 

 

 如果我們無法實現f'(x),則可以將其近似為

 

 

 

 δ選擇得適當小。 一般而言,我們可以用(3.29)替換為

 

 

 

 假設我們已經知道xprev和f(xprev)。 然后我們可以通過要求固定系數b(近似於f'(x))

 

 

 

 這給出,通過選擇b,我們將(3.30)視為割線方法,例如參見Elden等人(2004)的第70f頁。

 割線法相對於牛頓-拉夫森法的另一種有限差分近似法的主要優點是,每個迭代步驟只需要一個函數求值,而不是兩個。

 現在,考慮f的線性模型(3.7a):

 

 

 

 

 

 我們將替換為

 

 

 

 其中B是J(x)的當前近似值。 在下一個迭代步驟中,我們需要Bnew,以便

 

 

 

 特別地,我們希望該模型對h = x-xnew相等,即

 

 

 

 這給了我們Bnew的m·n個未知元素中的m個方程,因此我們需要更多條件。 Broyden(1965)建議用(3.31a)補充

 

 

 

 很容易驗證是否滿足條件(3.31a–b)

 

 

 

 注意,在n = 1的情況下,條件(3.31a)對應於割線條件(3.30b)。我們說這種方法是廣義割線方法。 具有此修改的算法3.16的中心部分的簡要草圖具有以下形式

 

 

 

 鮑威爾表明,如果向量x0,x1,x2,...的集合收斂到x ∗,並且步驟{hk≡xk-xk-1}的集合滿足{hk-n + 1,.....的條件。 對於每個k≥n,……,hk}是線性獨立的(它們跨越整個IRn),則不管B0的選擇如何,近似集{Bk}都收斂到J(x ∗)。

 然而,實際上,經常會發生前n個步驟未覆蓋整個IRn的情況,並且存在以下風險:在某些迭代步驟之后,當前B與真實的Jacobian矩陣的近似值如此之差,即-B> f (x)甚至不是下坡方向。 在這種情況下,x將保持不變,而µ會增加。 近似值B發生了變化,但是仍然可能是一個較差的近似值,從而導致µ的進一步增加等。最終,盡管x可能與x相距甚遠,但由於hslm太小而無法滿足(3.15b),因此該過程停止了。 

 已經提出了許多策略來克服這個問題,例如通過有限的差異偶爾對B進行重新計算。 在下面的Algo rithm 3.34中,我們用循環的,按坐標方式進行的一系列更新來補充由迭代過程確定的更新:令h表示當前步驟,令j表示當前坐標號。 如果h和ej之間的角度θ為“大”,則我們將計算J的第j列的有限差分近似值。更具體地說,如果

 

 

 

 實驗表明,(相當悲觀的)選擇γ= 0.8提供了良好的性能。 通過這種選擇,我們可以預期每個迭代步驟(幾乎)需要向量函數f的兩次評估。

 現在我們准備介紹該算法。 阻尼參數µ的監控與算法3.16相同,為清楚起見,我們在演示中將其省略。

 

 

 

 我們有以下幾點評論:

  1◦初始化。 輸入x0,輸入B0或由(3.28)計算。 停止標准(3.15)中的參數和(3.28)中使用的步距δ也是輸入值。
  2◦Cf(3.33)。 mod(j,n)是除以n后的余數。
  3◦步驟η由xj = 0給出,

 

 

             

  4◦僅在滿足下降條件(2.1)的情況下更新迭代x,而近似值B在每一步中都更新。 因此,當f(x)不變時,近似梯度g也可能改變。

示例3.15 我們對例3.13中的修正Rosenbrock問題使用了算法3.34,其中λ=0。如果我們使用與該例相同的起點和終點准則,並且差近似(3.28)中的δ= 10−7,我們發現 29個迭代步驟后的解,涉及f(x)的53個評估。 為了進行比較,“真正的” LM算法僅需17個步驟,這意味着總共對f(x)和J(x)進行了18次評估。

 我們還對測試1.1、3.7和3.11中的數據擬合問題使用了割線算法。 在δ= 10-7且起點和停止條件與示例3.7相同的情況下,迭代在94個步驟后被(3.15a)停止,涉及f(x)的總共192次評估。 為了進行比較,算法3.16需要62個迭代步驟。

 這兩個問題表明算法3.34是健壯的,但它們也說明了一條一般的經驗法則:如果有梯度信息可用,通常使用它是值得的。

 在許多應用中,數字m和n很大,但是每個函數fi(x)僅取決於x中的一些元素。 在那種情況下,大多數為零,我們說J(x)是一個稀疏矩陣。 在Levenberg-Marquardt方程(3.13)的解中,有一些利用稀疏性的有效方法,例如,Nielsen(1997)。 但是,在更新公式(3.32)中,向量h和u中的所有元素通常都不為零,因此Bnew將是一個密集矩陣。 討論如何解決這個問題超出了本手冊的范圍。 我們指的是Gill等人(1984)和Toint(1987)。

3.6. 狗腿法的割線版本 A Secant Version of the Dog Leg Method

 使用割線近似於雅可比行列的想法當然也可以與3.3節中的“狗腿法”結合使用。 在本節中,我們將考慮m = n的特殊情況,即在非線性方程組的解中。 Broyden(1965)不僅給出了定義3.32中的公式,

 

 

 

 用於更新近似雅可比行列式。 他還給出了一個公式,用於更新雅可比行列的近似逆D'J(x)-1。 公式是

 

 

 

 其中h和y在(3.35a)中定義。

 對於這些矩陣,最陡的下降方向hsd和高斯–牛頓步長(在這種情況下與牛頓步長相同,參見示例3.5),hgn(3.18)近似為

 

 

 可以輕松修改算法3.21以使用這些近似值。 初始B = B0可以通過差分近似(3.28)求出,D0計算為B-1 0。 很容易表明,那么當前的B和D滿足BD =I。步長參數α由(3.19)找到,而J(x)由B代替。

像L–M方法的割線版本一樣,此方法需要額外的更新,以使B和D保持與當前Jacobian及其逆的良好近似。 我們發現圍繞(3.33)討論的策略在這種情況下也很好用。 還應該提到的是,(3.35b)中的分母可以為零或非常小。 如果

 

 

 那么D不更新,而是按D = B-1計算。

用(3.35)進行的每次更新都會“耗費” 10n2觸發器6),通過(3.36)計算兩個階躍向量,再通過(3.19)來計算α,則需要6n2觸發器。 因此,使用Dog Leg方法的無梯度版本的每個迭代步驟的成本約為16n2觸發器加上f(xnew)的評估。 對於比較而言,采用算法3.21的每一步的成本約為23n3 + 6n2觸發器,加上f(xnew)和J(xnew)的評估。 因此,對於較大的n值,無梯度版本每步便宜。 但是,迭代步驟的數量通常要大得多,如果可以使用雅可比矩陣,則梯度版本通常會更快。

 

 

  6)一個“翻牌”是兩個浮點數之間的簡單算術運算。

示例3.16 我們在Rosenbrock函數f:上使用了算法3.21和無梯度Dog Leg方法,由

 

 

 參見示例3.13。 該函數具有一個根,x ∗ = [1,1]>,在兩種方法中,我們都使用起始點x0 = [-1.2,1]>和ε1=ε2= 10-12,kmax = 100作為停止准則 (3.15),而(3.28)中的δ= 10-7。 Algo rithm 3.21在17個迭代步驟后,即在對f及其雅可比行列式進行18次評估后,在求解處停止。 割線版本也停止了解決方案。 這需要28個迭代步驟和總共49次f評估。

3.7。 結束語Final Remarks

我們討論了許多用於解決非線性最小二乘問題的算法。 它們全部出現在任何好的程序庫中,並且可以通過Internet地址上的GAMS(可用數學軟件指南)找到實現。

本手冊中的示例是在MATLAB中計算的。 程序可在工具箱immoptibox中獲得,該工具箱可從以下網址獲得:

最后,應該提到的是,有時問題的重新制定可以使其更容易解決。 我們將通過示例來說明此主張,其中涉及可能也適用於您的問題的想法。示例3.17 在示例3.2、3.8和3.10的鮑威爾問題中,變量x2僅作為x22出現。 我們可以引入新的變量 ,問題的形式為:找使得f(z ∗)= 0,其中

 

 對於所有z,該雅可比行列都是非奇異的。 在停止准則中起點為z0 = [3,1]>,τ= 10-16和ε1=ε2= 10-15的LM算法3.16

(3.15)使用z'[-1.40e-25,9.77e-25]> 3個步驟后停止。 這是z ∗ = 0的良好近似值。
示例3.18 可以將示例1.1、3.7和3.11中的數據擬合問題重新編寫為僅具有兩個參數x1和x2:我們可以用以下形式編寫模型:

 

 其中,對於給定的x,找到向量c = c(x)∈IR2作為線性問題的最小二乘解

 其中。 與示例1.1中一樣,函數f由由(E)i行給出

定義,導致

 

 其中

行給出。 與示例1.1中一樣,函數f由f

定義,導致

 

 可以證明雅可比矩陣是

 其中,對於任何向量u,我們定義對角矩陣[u] = diag(u),並且

 

 

 

 算法3.16具有與示例3.7相同的較差的開始猜測,x0 = [-1,-2]>,τ= 10-3且ε1=ε2= 10-8得出解x'[-4,-5]> 經過13個迭代步驟; 大約是4參數模型所需步數的15。這種方法可以推廣到其中某些參數線性出現的任何模型。 它具有可分離的最小二乘法的名稱,例如在Nielsen(2000)和Golub and Pereyra(2003)中進行了討論。

示例3.19 最后一個例子說明了一個最小二乘問題的常見困難:通常,當問題按比例縮放以便所有(非零)
最小時,算法最有效。 具有相同的數量級。
考慮所謂的邁耶問題

 

 

 

 

 替代公式是

 

 然后,對於第一個公式,在175個迭代步驟之后,迭代將由(3.15b)停止,而在按比例縮放的情況下,在88個步驟之后,迭代將由(3.15a)停止。

A.對稱的正定矩陣

 

 如果A = AT,即aij = aji對於所有i,j,矩陣A∈IRn×n是對稱的

 

 此類定理的一些有用屬性在下面的定理A.2中列出。 可以通過將幾乎所有關於線性代數和數值線性代數的教科書中的定理結合起來找到證明。 在本附錄的最后,我們給出了該定理的一些實際含義。

 

 顯示A為正半定數。 如果m≥n並且J中的列線性獨立,則x ≠ 0⇒y ≠ 0且y> y>0。因此,在這種情況下,A為正定。

從下面的(A.3)中可以立即得出

 

對於任何μ∈IR。 將其與定理A.2中的2◦相結合,我們可以看到,如果A是對稱的且為半正定的且µ> 0,則矩陣A + µI也是對稱的,並且可以保證為正定的。

 

   定理A.2。 令A∈IRn×n是對稱的,令A = LU,其中L是單位下三角矩陣,U是上三角矩陣。 此外,令{(λj,vj)} nj = 1表示A的本征解,即

 

然后

1◦特征值是實數,λj∈IR,特征向量{vj}構成IRn的正交基。 2◦以下語句等效:a)A是正定(正半定)

b)所有λj> 0(λj≥0)

c)所有uii> 0(uii≥0)。

如果A為正定,則

3◦LU分解在數值上是穩定的。

4◦U = DL>,D = diag(uii)。 5◦A = C> C,即Cholesky分解。 C∈IRn×n是上三角。

對稱矩陣A的條件數為

 

如果A為正(半)定且µ> 0,則

 

這是µ的遞減函數。

最后,對定理A.2和實踐細節作了一些說明:一個單位較低的三角矩陣L的特征是'ii = 1,而j> i的ij = 0。 請注意,LU分解A = LU無需旋轉。 還要注意,點4◦–5◦在LU-和Cholesky分解之間具有以下關系
B.最小范數最小二乘解

 

顯示

 

可以通過以下算法直接計算Cholesky因式分解(即,沒有中間結果L和U),該算法包括正定性檢驗。

 

該算法的“成本”約為觸發器。 一旦計算出C,就可以通過正向和反向替換來求解系統Ax = b

 

分別。 每個步驟的費用約為n2觸發器

B.最小范數最小二乘解
考慮最小二乘問題:給定 且m≥n且b∈IR m,求h∈IR n使得

 

被最小化。 為了對此進行分析,我們將使用A的奇異值分解(SVD),

 

矩陣U∈IRm×m和V∈IRn×n正交

 

σj是奇異值。 p是A的秩,等於子空間R(A)⊆IRm的尺寸(所謂的A范圍),該子空間包含A列的每種可能的線性組合。

令{uj} mj = 1和{vj} nj = 1分別表示U和V中的列。 由於矩陣是正交的,因此向量形成兩個正交集,即

 

從(B.1)和(B.2)中可以得出

 

{uj}和{vj}可以分別用作IRm和IRn中的正交基,因此我們可以編寫

 

 通過使用(B.3),我們看到

 

 並借助(B.2)

 

 最小化

 

 因此,最小二乘解可以表示為

 

 當p <n時,最小二乘解具有n-p個自由度:ηp+ 1,...,ηn是任意的。 與關於(B.5)的討論類似,我們看到當ηp+ 1 =··=ηn= 0時kh ∗ k最小。對應於自由參數選擇的解決方案是所謂的最小范數解,

 

 重新制定源自(B.4)和(B.2)。

R EFERENCES
1. M. Al-Baali and R. Fletcher (1985): Variational Methods for Non-Linear Least  Squares. J. Opl. Res. Soc. 36 , No. 5, pp 405–421.
2. C.G. Broyden (1965): A Class of Methods for Solving Nonlinear Simultaneous  Equations. Maths. Comp. 19 , pp 577–593.
3. J.E. Dennis, Jr. and R.B. Schnabel (1983): Numerical Methods for Uncon  strained Optimization and Nonlinear Equations, Prentice Hall.
4. L. Elden, L. Wittmeyer-Koch and H.B. Nielsen (2004):  ´
Introduction to Numer ical Computation – analysis and M ATLAB illustrations , to appear.
5. P.E. Frandsen, K. Jonasson, H.B. Nielsen and O. Tingleff (2004):  Unconstrained Optimization, 3rd Edition, IMM, DTU. Available as  http://www.imm.dtu.dk/courses/02611/uncon.pdf
6. P.E. Gill, W. Murray, M.A. Saunders and M.H. Wright (1984): Sparse Matrix  Methods in Optimization, SIAM J.S.S.C. 5 , pp 562–589.
7. G. Golub and V. Pereyra, (2003): Separable Nonlinear Least Squares: The  Variable Projection Method and its Applications , Inverse Problems, 19 , pp R1–  R26.
8. G. Golub and C.F. Van Loan (1996): Matrix Computations , 3rd edition. John  Hopkins Press.
9. K. Levenberg (1944): A Method for the Solution of Certain Problems in Least  Squares. Quart. Appl. Math. 2 , pp 164–168.
10. K. Madsen (1988): A Combined Gauss-Newton and Quasi-Newton Method  for Non-Linear Least Squares. Institute for Numerical Analysis (now part of  IMM), DTU. Report NI-88-10.
11. K. Madsen and H.B. Nielsen (2002): Supplementary Notes for 02611 Opti mization and Data Fitting , IMM, DTU. Available as  http://www.imm.dtu.dk/courses/02611/SN.pdf
12. K. Madsen, H.B. Nielsen and J. Søndergaard (2002): Robust Subroutines for  Non-Linear Optimization. IMM, DTU. Report IMM-REP-2002-02. Available  as http://www.imm.dtu.dk/ km/F  package/robust.pdf
13. K. Madsen, H.B. Nielsen and O. Tingleff (2004): Optimization with Con straints, 2nd Edition. IMM, DTU. Available as  http://www.imm.dtu.dk/courses/02611/ConOpt.pdf
14. D. Marquardt (1963): An Algorithm for Least Squares Estimation on Nonlinear  Parameters. SIAM J. A PPL . M ATH . 11 , pp 431–441.
15. H.B. Nielsen (1997): Direct Methods for Sparse Matrices. IMM, DTU. Avail able as http://www.imm.dtu.dk/ hbn/publ/Nsparse.ps
16. H.B. Nielsen (1999): Damping Parameter in Marquardt’s Method. IMM, DTU.  Report IMM-REP-1999-05. Available as  http://www.imm.dtu.dk/ hbn/publ/TR9905.ps
17. H.B. Nielsen (2000): Separable Nonlinear Least Squares. IMM, DTU. Report  IMM-REP-2000-01. Available as  http://www.imm.dtu.dk/ hbn/publ/TR0001.ps
18. J. Nocedal and S.J. Wright (1999): Numerical Optimization , Springer.
19. M.J.D. Powell (1970): A Hybrid Method for Non-Linear Equations. In P. Rabi nowitz (ed): Numerical Methods for Non-Linear Algebraic Equations , Gordon  and Breach. pp 87ff.
20. P.L. Toint (1987): On Large Scale Nonlinea

 

 
顯示A為正半定數。 如果m≥n並且J中的列線性獨立,則x = 0⇒y = 0且y> y>0。因此,在這種情況下,A為正定。


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM