# MethodError: no method matching Float64 when using ForwardDiff.jl despite me pre-allocating arrays to take floating point

**URL:** <https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859>\
**Category:** Optimization (Mathematical)\
**Created:** [November 7, 2022, 7:59am UTC](https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859 "2022-11-07T07:59:00Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![sumant28](https://avatars.discourse-cdn.com/v4/letter/s/a698b9/32.png) [@sumant28](https://discourse.julialang.org/u/sumant28)\
**Post date:** [November 7, 2022, 7:59am UTC](https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859/1 "2022-11-07T07:59:00Z")

</div>

```julia
using Optimization, OptimizationOptimJL, OptimizationBBO, OptimizationMOI, ForwardDiff, ModelingToolkit, Random, OptimizationOptimisers, Statistics, Plots

function MSEg(u, p)
  X = p.X
  y = p.y
  h = u[1]
  beta = u[2:end]
  N = length(y)
  Xbeta = X*beta

  krnl = u -> exp.(-(u .* u) ./ 2) ./ sqrt(2*π)

  yhat = Array{Float64, 2}(undef, N, 1)
  weight_i = Array{Float64, 2}(undef, N, 1)
  weight_sum = Array{Float64, 2}(undef, N, 1)

  for i = 1:N
    u = (Xbeta .- Xbeta[i]) ./ h
    Ku = krnl(u)
    weight_i[i] = Ku[i]
    weight_sum[i] = sum(Ku)
    yhat[i] = sum(Ku .* y) ./ sum(Ku)
  end

  residual = (y .- yhat) ./ (1 .- (weight_i ./ weight_sum))
  SBeta = sum(residual.^2)
end

Random.seed!(1)
N = 1000
X = -1 .+ 2 .* rand(N, 2)
beta = [1;3]
Xbeta = X*beta
y = Xbeta.^2 .* (1 .+ sin.(Xbeta)) .+ randn(N, 1)

y_ = y .- X[:,1]
X_ = X[:,2]
beta0 = X_'*X_\X_'*y_ 
beta0 = [1; beta0]

x_ = X*beta0
hx = median(abs.(x_ .- median(x_)))/0.6745*(4/3/N)^0.2
hy = median(abs.(y .- median(y)))/0.6745*(4/3/N)^0.2
h = sqrt(hy*hx)

param0 = [h;beta0]
u = param0

struct smol
  X
  y
end

p = smol(X, y)

cons(res, u, p) = (res .= u[2])

lb = [0;-Inf;-Inf]
ub = [Inf;Inf;Inf]

optf = OptimizationFunction(MSEg, Optimization.AutoForwardDiff(), cons = cons)
prob = OptimizationProblem(optf, u, p, lb = lb, ub = ub, lcons = 1, ucons = 1)
sol = solve(prob, NelderMead())

```

I had received this error in the past when I tried to use a for loop to populate elements of arrays and matrices. I had not had trouble in the past with this error in a while but now it is coming back again and I cannot understand why

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [November 7, 2022, 8:11am UTC](https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859/2 "2022-11-07T08:11:45Z")

</div>

> despite me pre-allocating arrays to take floating point

This is the problem. You need:

```julia
function MSEg(u::AbstractVector{T}, p) where {T}
    # ...
     yhat = Array{T, 2}(undef, N, 1)
    #...

```

Your function must not assume that `u` is a `Vector{Float64}`.

---

<div class="post-metadata">

**Author:** ![sumant28](https://avatars.discourse-cdn.com/v4/letter/s/a698b9/32.png) [@sumant28](https://discourse.julialang.org/u/sumant28)\
**Post date:** [November 7, 2022, 9:31am UTC](https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859/3 "2022-11-07T09:31:44Z")

</div>

> [@odow](#):
>
> `yhat = Array{T, 2}(undef, N, 1)`

I am still confused. I changed the function below and get the same error when running the script

```julia
function MSEg(u::AbstractVector{T}, p) where {T}
  typeof(u)
  X = p.X
  y = p.y
  h = u[1]
  beta = u[2:end]
  N = length(y)
  Xbeta = X*beta

  krnl = u -> exp.(-(u .* u) ./ 2) ./ sqrt(2*π)

  # yhat = Array{Float64, 2}(undef, N, 1)
  # weight_i = Array{Float64, 2}(undef, N, 1)
  # weight_sum = Array{Float64, 2}(undef, N, 1)

  yhat = Array{T, 2}(undef, N, 1)
  weight_i = Array{T, 2}(undef, N, 1)
  weight_sum = Array{T, 2}(undef, N, 1)

  for i = 1:N
    u = (Xbeta .- Xbeta[i]) ./ h
    Ku = krnl(u)
    weight_i[i] = Ku[i]
    weight_sum[i] = sum(Ku)
    yhat[i] = sum(Ku .* y) ./ sum(Ku)
  end
  residual = (y .- yhat) ./ (1 .- (weight_i ./ weight_sum))
  SBeta = sum(residual.^2)
end

```

---

<div class="post-metadata">

**Author:** ![ChrisRackauckas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/chrisrackauckas/32/77_2.png) [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)\
**Post date:** [November 7, 2022, 12:46pm UTC](https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859/4 "2022-11-07T12:46:09Z")

</div>

> [@sumant28](#):
>
> ```julia
> yhat = Array{Float64, 2}(undef, N, 1)
> weight_i = Array{Float64, 2}(undef, N, 1)
> weight_sum = Array{Float64, 2}(undef, N, 1)
> 
> ```

These won’t be Float64, so the type is incorrect.

---

<div class="post-metadata">

**Author:** ![odow](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/odow/32/28685_2.png) [@odow](https://discourse.julialang.org/u/odow)\
**Post date:** [November 7, 2022, 7:07pm UTC](https://discourse.julialang.org/t/methoderror-no-method-matching-float64-when-using-forwarddiff-jl-despite-me-pre-allocating-arrays-to-take-floating-point/89859/5 "2022-11-07T19:07:31Z")

</div>

You can simplify your `MSEg` code a lot (and there were a bunch of other problems in your code, like the wrong optimizer). How about this:

```julia
using Optimization
import OptimizationMOI
import ForwardDiff
import Ipopt
import Random
using Statistics

kernel(u) = exp(-(u * u) / 2) / sqrt(2π)

function MSEg(u::AbstractArray{T}, p) where {T}
    X, y = p.X, p.y
    h, beta = u[1], u[2:end]
    Xbeta = X * beta
    residuals = zero(T)
    for i in 1:length(y)
        Ku = kernel.((Xbeta .- Xbeta[i]) ./ h)
        yhat = sum(Ku .* y) ./ sum(Ku)
        residuals += ((y[i] - yhat) / (1 - (Ku[i] / sum(Ku))))^2
    end
    return residuals
end

Random.seed!(1)
N = 1000
X = -1 .+ 2 .* rand(N, 2)
beta = [1;3]
Xbeta = X*beta
y = Xbeta.^2 .* (1 .+ sin.(Xbeta)) .+ randn(N, 1)

y_ = y .- X[:,1]
X_ = X[:,2]
beta0 = X_' * X_ \ X_' * y_ 
beta0 = [1; beta0]

x_ = X * beta0
hx = median(abs.(x_ .- median(x_))) / 0.6745 * (4 / 3 / N)^0.2
hy = median(abs.(y .- median(y))) / 0.6745 * (4 / 3 / N)^0.2
h = sqrt(hy*hx)

struct smol
  X
  y
end

p = smol(X, y)

cons(res, u, p) = (res .= u[2])

optf = OptimizationFunction(MSEg, Optimization.AutoForwardDiff(), cons = cons)
prob = OptimizationProblem(
  optf, 
  [h; beta0], 
  p; 
  lb = [0, -Inf, -Inf],
  ub = [Inf, Inf, Inf],
  lcons = [1.0],
  ucons = [1.0],
)
sol = solve(prob, Ipopt.Optimizer())

```