scipy.optimize.leastsq無法適合簡單模型

Question

我一直試圖使用python的scipy.optimize庫來估計模型的參數，但目前爲止沒有成功。我試圖使用使用levenberg-marquardt算法的scipy.optimize.leastsq。不幸的是，即使我將初始參數的猜測設置得非常接近最佳擬合，它總是無法找到我的模型函數的最小值。實際上，它總是返回最初猜測的參數。所以，我認爲我做錯了什麼。我的模型是一個簡單的圓，但爲了使事情更簡單，只有半徑是實際參數，數據中圓的中心是已知的並且是硬編碼的。數據是一個10x10像素的浮點圖像，中心爲5,5和半徑爲4的圓。實際上，數據是使用我試圖擬合的模型生成的。所以，存在一個完美契合。這裏是我的代碼：爲了去除任何數據的依賴性，並允許代碼安裝scipy任何機器上運行開箱

import math 
import numpy 
import scipy.optimize 

# ======================================================================== 

g_data_width = 10 
g_data_height = 10 
g_xo = 5.0 
g_yo = 5.0 

def evaluate_model01(x,y,r): 
    x2 = x*x 
    y2 = y*y 
    r2 = r*r 
    v = 0.0 
    if(x2 + y2 <= r2): 
     v = 20.0 
    return v 

def model01(params,data_o): 
    data_r = numpy.zeros(g_data_height*g_data_width) 
    r = params[0] 
    for y in range(g_data_height): 
     for x in range(g_data_width): 
      xnew = x - g_xo 
      ynew = y - g_yo 
      data_r[y*g_data_width+x] = math.fabs(data_o[y,x]-evaluate_model01(xnew,ynew,r)) 
    return data_r 

# ======================================================================== 

g_data_o = numpy.array([[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], 
         [ 0, 0, 0, 0, 0, 20, 0, 0, 0, 0], 
         [ 0, 0, 0, 20, 20, 20, 20, 20, 0, 0], 
         [ 0, 0, 20, 20, 20, 20, 20, 20, 20, 0], 
         [ 0, 0, 20, 20, 20, 20, 20, 20, 20, 0], 
         [ 0, 20, 20, 20, 20, 20, 20, 20, 20, 20], 
         [ 0, 0, 20, 20, 20, 20, 20, 20, 20, 0], 
         [ 0, 0, 20, 20, 20, 20, 20, 20, 20, 0], 
         [ 0, 0, 0, 20, 20, 20, 20, 20, 0, 0], 
         [ 0, 0, 0, 0, 0, 20, 0, 0, 0, 0]],dtype=numpy.float32) 
g_params = numpy.array([8.0]) 
print(scipy.optimize.leastsq(model01,g_params,args=(g_data_o),full_output=1)) 

# ======================================================================== 

我硬編碼的數據。我不太明白的是model01函數應該返回。根據文件，它應該返回一個數組。一個什麼樣的數組？在我的代碼中，我假設我必須爲每個數據點返回一個殘差數組。那是對的嗎？我的數據是二維數組，因爲它是一個圖像，但我的殘差是一個平坦的二維殘差數組。這可以嗎？有人能告訴我我究竟做錯了什麼嗎？有人可以修改和修復代碼嗎？正如我上面提到的，代碼應該在安裝有scipy和numpy的任何機器上運行。如果scipy.optimize.leastsq不能實現我想要實現的功能，那麼能否使用levenberg-marquardt算法推薦一些適合模型的其他庫？

Answer 1

恐怕你不能用leastsq解決你的特殊問題。 leastsq是FORTRAN庫minipack的環繞，它調用MINPACK的lmdif和lmder算法。重要的是，它基於最小二乘目標函數的雅可比矩陣和Hessian矩陣。你的目標函數不具有光滑的衍生物，由於：

if(x2 + y2 <= r2): 
    v = 20.0 

和

math.fabs(......) 

，所以leastsq總是返回初始啓動參數。

您應該嘗試使用一些不需要梯度/衍生的方法，例如Powell的方法fmin_powell或Nelder-Mead fmin。

關於您的model101()正在發生什麼。它首先製作一個width*height的大小的一維數組，稱爲data_r。然後它遍歷g_data，計算每個元素的math.fabs(data_o[y,x]-evaluate_model01(xnew,ynew,r))並將該值放入一維數組data_r。最後，它返回date_r。

你的解釋是正確的，model101()爲你的二維數據返回平坦的1D殘差陣列。