Turing type error: expected Float64, got ForwardDiff.Dual

Specific Domains Probabilistic programming
, ,
1

Hi All,

I’m trying to do a multivariable linear regression y = \Phi(A, x) using Turing, where y \in \mathbb R^N, x \in \mathbb R^N, and the function \Phi is linear in the unknown parameters A \in \mathbb R^{N \times N}. For my particular problem, the function \Phi = (\phi_1, \cdots, \phi_n)^T needs to be constructed row-by-row.

It’s straightforward to build a Turing model to infer the unknown matrix A from a dataset containing the inputs x and outputs y. However, when I run the sampler, I find the following error:

ERROR: TypeError: in typeassert, expected Float64, got ForwardDiff.Dual{Nothing,Float64,9}

Below is a minimal example to reproduce this error. I also found that the error does not occur if all the rows of \Phi are constructed in a single step.

I checked this related previous post, but the solutions there do not seem to solve the problem.

Does anybody have a suggestion about how to fix this issue?

Thanks!

using Distributions
using LinearAlgebra
using StatsBase

using Turing

function Φ(A, x)
    N = size(x,1)
    temp = zeros(N)

    for i = 1:N
        temp[i] =  dot(A[i, :], x)
    end

    return temp
end

N = 5
number_of_samples = 100

# Generate data
A = rand(Normal(), N,N) + I(N)
xdata = rand(LogNormal(0,1), (N, number_of_samples) )
ydata = hcat([Φ(A, x) for x in eachcol(xdata)]...)

# Build model
@model function regression(x, y)
    N = size(x, 1)
    number_of_samples = size(x, 2)

    # Priors
    σ ~ truncated(Normal(0, 100), 0, Inf)
    A ~ filldist(Normal(0,1), N, N)

    # Equations
    for k = 1: number_of_samples
        # mu  =  A*x[:, k]
         mu =  Φ(A, x[:, k])
        y[:, k] ~ MvNormal(mu, sqrt(σ))
    end
end

# Inference
model = regression(xdata, ydata)
chain = sample(model, NUTS(0.65), 300)
2

This line:

will create an array that is of type Float64, and I suspect that x contains dual number because NUTS needs derivatives with respect to it. Doing temp = zeros(eltype(x), N) might fix the problem.

An even easier option might be to replace the loop with a comprehension, which will take care of types for you:

temp = [dot(A[i, :], x) for i in 1:N]
1 Like
3

Why not just write

temp = A * x

?
That’s both simpler, clearer, and potentially much faster.

3 Likes
4

This, too, will be slow. Why not just

ydata = Φ(A, xdata)

?

1 Like
5

Unfortunately, your first suggestion temp = zeros(eltype(x), N) does not work (do you understand why?)
The second suggestion temp = [dot(A[i, :], x) for i in 1:N] solves the problem.
Thanks!

6

For my particular problem, I cannot write the model as a multiplication of a matrix and a vector.
Note the code I pasted is only a minimal case that reproduces the error I was finding.
Thanks for answering.

7

This does not work because the inner product in Φ will not correctly broadcast.

8

Yes, I just tried it. It’s because NUTS also wants the gradient with respect to A, and so A can contain dual numbers. This fixes it:

temp = zeros(promote_type(eltype(x), eltype(A)), N)

Of course, the list comprehension is an easier option when possible to use.

1 Like
9

Or you can use the super secret special function Base.promote_eltype :shushing_face:

julia> Base.promote_eltype(rand(Float64, 3), rand(Int, 3))
Float64
3 Likes
10

That’s the thing with minimal reproducible examples, they can be difficult to make.

But, I think that in order for something to be a minimal example, it shouldn’t be possible to look at it and then dramatically simplify it further.