Extremely new to Julia, so please pardon any obvious oversights on my end
I am trying to estimate a piecewise likelihood function through optimization. I have the code functional in R, but have begun translating it to Julia in the hopes of faster estimation, for eventual bootstrapping
Here is the current block of code that I am trying (v and x are already as 1000x1 vectors elsewhere defined elsewhere):
function est(a,b)
function pwll(v,x)
if v>4
ILL=pdf(Poisson(exp(a+b*x)), v)
elseif v==4
ILL=pdf(Poisson(exp(a+b*x)), 4)+pdf(Poisson(exp(a+b*x)),3)+pdf(Poisson(exp(a+b*x)),2)
else v==0
ILL=pdf(Poisson(exp(a+b*x)), 1)+pdf(Poisson(exp(a+b*x)), 0)
end
return(ILL)
end
ILL=pwll.(v, x)
function fixILL(x)
if x==0
x=0.00000000000000001
else
x=x
end
end
ILL=fixILL.(ILL)
LILL=log10.(ILL)
LL=-1*LILL
return(sum(LL))
end
using Optim
params0=[1,1]
optimize(est, params0)
And the error message(s) I am getting are:
ERROR: InexactError: Int64(NaN)
Stacktrace:
[1] Int64(x::Float64)
@ Base ./float.jl:788
[2] x_of_nans(x::Vector{Int64}, Tf::Type{Int64}) (repeats 2 times)
@ NLSolversBase ~/.julia/packages/NLSolversBase/kavn7/src/NLSolversBase.jl:60
[3] NonDifferentiable(f::Function, x::Vector{Int64}, F::Int64; inplace::Bool)
@ NLSolversBase ~/.julia/packages/NLSolversBase/kavn7/src/objective_types/nondifferentiable.jl:11
[4] NonDifferentiable(f::Function, x::Vector{Int64}, F::Int64)
@ NLSolversBase ~/.julia/packages/NLSolversBase/kavn7/src/objective_types/nondifferentiable.jl:10
[5] promote_objtype(method::NelderMead{Optim.AffineSimplexer, Optim.AdaptiveParameters}, x::Vector{Int64}, autodiff::Symbol, inplace::Bool, args::Function)
@ Optim ~/.julia/packages/Optim/tP8PJ/src/multivariate/optimize/interface.jl:63
[6] optimize(f::Function, initial_x::Vector{Int64}; inplace::Bool, autodiff::Symbol, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
@ Optim ~/.julia/packages/Optim/tP8PJ/src/multivariate/optimize/interface.jl:86
[7] optimize(f::Function, initial_x::Vector{Int64})
@ Optim ~/.julia/packages/Optim/tP8PJ/src/multivariate/optimize/interface.jl:83
[8] top-level scope
@ ~/Documents/Projects/ki_new/peicewise_ll.jl:120
I understand that it seems the error is coming from the function to be optimized being non-differentiable. A fairly direct translation works well in R, using the built in optim() function.
Can anyone provide any insight?
I have tried the above code displayed above, with multiple variations. The function to be optimized is functional, I am struggling with the optimization (the issues of which may stem from function being inefficiently written)
Here's an adapted version of your code which produces a solution:
using Distributions, Optim
function pwll(v, x, a, b)
d = Poisson(exp(a+b*x))
if v > 4
return pdf(d, v)
elseif v == 4
return pdf(d, 4) + pdf(d, 3) + pdf(d, 2)
else
return pdf(d, 1) + pdf(d, 0)
end
end
fixILL(x) = iszero(x) ? 1e-17 : x
est(a, b, v, x) = sum(-1 .* log10.(fixILL.(pwll.(v, x, a, b))))
v = 4; x = 0.5 # Defining these here as they are not given in your post
obj(input; v = v, x = x) = est(input[1], input[2], v, x)
optimize(obj, [1.0, 1.0])
I have no idea whether this is correct of course, check this against some sort of known result if you can.