Julia：使用CurveFit包非線性

Question

爲什麼下面的代碼不起作用？

xa = [0 0.200000000000000 0.400000000000000 1.00000000000000 1.60000000000000 1.80000000000000 2.00000000000000 2.60000000000000 2.80000000000000 3.00000000000000 3.80000000000000 4.80000000000000 5.00000000000000 5.20000000000000 6.00000000000000 6.20000000000000 7.40000000000000 7.60000000000000 7.80000000000000 8.60000000000000 8.80000000000000 9.00000000000000 9.20000000000000 9.40000000000000 10.0000000000000 10.6000000000000 10.8000000000000 11.2000000000000 11.6000000000000 11.8000000000000 12.2000000000000 12.4000000000000]; 

ya = [-0.183440428023042 -0.131101157495126 0.0268875670852843 0.300000000120000 0.579048247883555 0.852605831272159 0.935180993484717 1.13328608090532 1.26893326843583 1.10202945535186 1.09201137189664 1.14279083803453 0.72 0.909735376251797 0.417067545528244 0.460107770989798 -0.516307074859654 -0.333994077331822 -0.504124744955962 -0.945794320817293 -0.915934553082780 -0.975458595671737 -1.09943707404275 -1.11254211607374 -1.29739980589100 -1.23440439602665 -0.953807504156356 -1.12240274852172 -0.609284630192522 -0.592560286759450 -0.402521296049042 -0.510090363150962]; 

x0 = vec(xa) 
y0 = vec(ya) 
fun(x,a) = a[1].*sin(a[2].*x - a[3]) 
a0 = [1,2,3] 
eps = 0.000001 
maxiter=200 
coefs, converged, iter = CurveFit.nonlinear_fit(x0 , fun , a0 , eps, maxiter) 
y0b = fit(x0) 
Winston.plot(x0, y0, "ob", x0, y0b, "r-", linewidth=3) 

Error: LoadError: MethodError: convert has no method matching convert(::Type{Float64}, ::Array{Float64,1}) This may have arisen from a call to the constructor Float64(...), since type constructors fall back to convert methods. Closest candidates are: call{T}(::Type{T}, ::Any) convert(::Type{Float64}, !Matched::Int8) convert(::Type{Float64}, !Matched::Int16)

while loading In[269], in expression starting on line 8 in nonlinear_fit at /home/jmarcellopereira/.julia/v0.4/CurveFit/src/nonlinfit.jl:75

Answer 1

fun函數必須返回一個Float64類型的剩餘價值r，在數據的每個迭代中計算，如下：

r = y - fun(x, coefs) 

所以你的功能y=a1*sin(x*a2-a3)將被定義爲：

fun(x,a) = x[2]-a[1]*sin(a[2]*x[1] - a[3]) 

其中：

x[2] is a value of 'y' vector 
    x[1] is a value of 'x' vector 
    a[...] is the set of parameters 

函數fun必須返回一個單獨的Float64，所以操作符不能是'dot ve rsion'（.*）。

通過調用nonlinear_fit功能，第一個參數必須是一個數組NX2，用含有x N個值和所述第二的第一列，包含y N個值，所以必須將兩者連接起來的載體x和y在兩列數組：

xy = [x y] 

最後，調用函數：

coefs, converged, iter = CurveFit.nonlinear_fit(xy , fun , a0 , eps, maxiter) 

的swering您關於返回係數的評論是不正確的：

的y = 1 * sin (x * a2-a3)是諧波功能，所以係數從函數調用返回，取決於嚴重的參數a0（「初始猜測爲每個配件參數」），您將作爲第三個參數發送（與maxiter=200_000）：

a0=[1.5, 1.5, 1.0] 
coefficients: [0.2616335317043578, 1.1471991302529982,0.7048665905560775] 

a0=[100.,100.,100.] 
coefficients: [-0.4077952060368059, 90.52328921205392, 96.75331155303707] 

a0=[1.2, 0.5, 0.5] 
coefficients: [1.192007321713507, 0.49426296880933257, 0.19863645732313934] 

我認爲你得到的結果是諧波，作爲圖：

其中：用Gadfly生成

blue line: 
f1(xx)=0.2616335317043578*sin(xx*1.1471991302529982-0.7048665905560775) 

yellow line: 
f2(xx)=1.192007321713507*sin(xx*0.49426296880933257-0.19863645732313934) 

pink line: 
f3(xx)=-0.4077952060368059*sin(xx*90.52328921205392-96.75331155303707) 

blue dots are your initial data. 

圖爲：

plot(layer(x=x,y=y,Geom.point),layer([f1,f2,f3],0.0, 15.0,Geom.line)) 

_{與朱莉婭版0.4測試。3}

Answer 2

從文件：

we are trying to fit the relationship fun(x, a) = 0

所以，如果你想找到的a元素的方式：每個xi,yi在[x0 y0] =>a[1].*sin(a[2].*xi - a[3])==yi，那麼正確的做法是：

fun(xy,a) = a[1].*sin(a[2].*xy[1] - a[3])-xy[2]; 
xy=hcat(x0,y0); 
coefs,converged,iter = CurveFit.nonlinear_fit(xy,fun,a0,eps,maxiter); 

Answer 3

rapz所以你救我一命）。 居然沒有注意到正弦函數的諧波特性等的初始值是至關重要的最後係數。解釋很優秀的，我推薦後如曲線擬合教程GitHub上爲「juleiros」第一代碼和我一樣的頁面，那些教程是值得很多推廣的語言 - 去了其他同事誰也開始在語言中。

再次感謝成倍的言論。我正在慢慢馴服朱莉婭;）。

Answer 4

我發現LsqFit package有點簡單易用，只需要首先定義模型，並與你的約會「適應它」：使用LsqFit沒有解決的問題

using DataFrames, Plots, LsqFit 

xa   = [0 0.200000000000000 0.400000000000000 1.00000000000000 1.60000000000000 1.80000000000000 2.00000000000000 2.60000000000000 2.80000000000000 3.00000000000000 3.80000000000000 4.80000000000000 5.00000000000000 5.20000000000000 6.00000000000000 6.20000000000000 7.40000000000000 7.60000000000000 7.80000000000000 8.60000000000000 8.80000000000000 9.00000000000000 9.20000000000000 9.40000000000000 10.0000000000000 10.6000000000000 10.8000000000000 11.2000000000000 11.6000000000000 11.8000000000000 12.2000000000000 12.4000000000000]; 
ya   = [-0.183440428023042 -0.131101157495126 0.0268875670852843 0.300000000120000 0.579048247883555 0.852605831272159 0.935180993484717 1.13328608090532 1.26893326843583 1.10202945535186 1.09201137189664 1.14279083803453 0.72 0.909735376251797 0.417067545528244 0.460107770989798 -0.516307074859654 -0.333994077331822 -0.504124744955962 -0.945794320817293 -0.915934553082780 -0.975458595671737 -1.09943707404275 -1.11254211607374 -1.29739980589100 -1.23440439602665 -0.953807504156356 -1.12240274852172 -0.609284630192522 -0.592560286759450 -0.402521296049042 -0.510090363150962]; 
x0   = vec(xa) 
y0   = vec(ya) 
xbase  = collect(linspace(minimum(x0),maximum(x0),100)) 
p0   = [1.2, 0.5, 0.5]              # initial value of parameters 
fun(x0, p) = p[1] .* sin.(p[2] .* x0 .- p[3])          # definition of the model 
fit  = curve_fit(fun,x0,y0,p0)            # actual fitting job 
yFit  = [fit.param[1] * sin(fit.param[2] * x - fit.param[3]) for x in xbase] # building the fitted values 
# Plotting.. 
scatter(x0, y0, label="obs") 
plot!(xbase, yFit, label="fitted") 

注闕從Gomiero強調初始條件的依賴..