Левенберг-Марквардт в пакетах Нелинейных наименьших квадратов Юлии - PullRequest
Новая
       17

Левенберг-Марквардт в пакетах Нелинейных наименьших квадратов Юлии

1 голос
/ 16 мая 2019

У меня есть проблема параметрической оценки, которую можно решить с помощью нелинейной оптимизации наименьших квадратов.У меня есть аналитическая модель с двумя неизвестными параметрами x[0] и x[1].У меня есть вектор выборок измеренных данных, с помощью которого я пишу cost function. 

function cost_function(x)
    measured = data_read()   # This is where I read my measured samples
    model = analytical_model(x)  # I calculate the values from my 
                                   analytical model method
    residual = abs.(measured - model)
    return residual
end

In LsqFit package в Julia, который имеет реализацию Levenberg-Marquardt (LM),есть только метод curve_fit, который принимает model (analytical_model() здесь), xdata, p, ydata и передает эту функцию оптимизатору LM.(то есть) model(xdata,p) - ydata становится остатком.Но для меня мои measured и model являются комплексными числами, и именно поэтому я должен вернуть abs.().

Я пробовал 

import LsqFit
result = LsqFit.levenberg_marquardt(cost_function, rand(2))

, но для этого требуется якобиан моего cost_function(x) в качестве другого аргумента.Я не знаю якобиана и хочу, чтобы мой оптимизатор рассчитал для меня якобиан, используя приближение прямой разницы.Есть ли способ сделать это в Джулии?

1 Ответ

2 голосов
/ 16 мая 2019

Должен работать следующий вызов: 

LsqFit.lmfit(cost_function, [0.5, 0.5], Float64[])

(пожалуйста, дайте мне знать, если это сработало, [0.5, 0.5] является примером отправной точки)

Вот пример: 

julia> function cost_function(x)
           measured = [2:2:20;] .+ 0.5im
           model = [1:10;] .* x[1] .+ x[2] .* im
           residual = abs.(measured .- model)
           return residual
       end
cost_function (generic function with 1 method)

julia> LsqFit.lmfit(cost_function, [0.5, 0.5], Float64[])
LsqFit.LsqFitResult{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,1}}([2.0, 0.5], [5.18665e-8, 5.36734e-8, 5.65567e-8, 6.03625e-8, 6.49286e-8, 7.01067e-8, 7.57714e-8, 8.18218e-8, 8.81784e-8, 9.47796e-8], [0.000658313 -0.00846345; 0.00131663 -0.00846342; … ; 0.00592487 -0.00846286; 0.00658318 -0.00846272], true, Float64[])

julia> LsqFit.lmfit(cost_function, [0.5, 0.5], Float64[], autodiff=:forward)
LsqFit.LsqFitResult{Array{Float64,1},Array{Float64,1},Array{Float64,2},Array{Float64,1}}([2.0, 0.5], [4.44089e-16, 8.88178e-16, 1.77636e-15, 1.77636e-15, 1.77636e-15, 3.55271e-15, 3.55271e-15, 3.55271e-15, 3.55271e-15, 3.55271e-15], [-1.0 -0.0; -2.0 -0.0; … ; -9.0 -0.0; -10.0 -0.0], true, Float64[])

и мы видим, что верное решение [2.0, 0.5] возвращается.

...