Jacobian returns zero matrix

[rakshith95]

Hello,
I’m working with images, and I’m seeing some strange behaviour on Jacobian computation, so I suspect my function formulation is incorrect.

I have two patches in an image, and my function F which returns a vector of differences in intensities between the patches.
The difference vector which is returned is computed as image[ y ] - image[ g(y, x) ] where x are the variables, and y are some coordinates.

Now, if I try computing the jacobian of F wrt the variables of x with ForwardDiff, my Jacobian is a matrix of all zeros, implying that there is no relation between the function vector, and the input variables.

Is the way I am formulating this incorrect?

The function code:

function F_photometric(x̂::AbstractVector, A_data::AbstractVector, u_data::AbstractVector, imagebgr::AbstractVecOrMat, region1_coordinates::AbstractMatrix, image_dims::AbstractVector)
    # A_diff = x̂[1:4] - A_data
    # u_diff = x̂[5:8] - u_data
    â = reshape(x̂[1:4],2,2)
    
    region1_coordinates_int = Int.(round.(region1_coordinates[1:2,:]))
    region2_coordinates = x̂[7:8] .+  â*(region1_coordinates_int .- x̂[5:6])  
    region2_coordinates_int = Int.(round.(region2_coordinates))
    
    region1_indices = Array{Int64}(undef, size(region1_coordinates)[2])
    region2_indices = Array{Int64}(undef, size(region2_coordinates)[2])
    for j=1:size(region1_coordinates)[2]
        region1_indices[j] = convert_2Dto1D_index(image_dims,region1_coordinates_int[:,j])
        region2_indices[j] = convert_2Dto1D_index(image_dims,region2_coordinates_int[:,j])
    end

    region1_bgr = imagebgr[region1_indices]
    region2_bgr = imagebgr[region2_indices]
    photo_diff = region1_bgr - region2_bgr
    # photo_diff_norm = reduce(vcat, photo_diff)

    # print(norm(photo_diff),"\t")
    return photo_diff
end
  
Jₒ = ForwardDiff.jacobian(x̂ -> F(x̂, A_data, u_data, image_bgr, region1, img_dims), x̂)

[mcabbott]

Can't run this, but the output seems to depend smoothly on imagebgr alone, not the first argument, in which case zero is correct.

[rakshith95]

Yeah, sorry, this wasn't meant to be a MWE, just wanted to show my function.

Can't run this, but the output seems to depend smoothly on imagebgr alone, not the first argument, in which case zero is correct.

Yeah, I suspected it was this. The coordinates of imagebgr are functions of my variables, but is that the output itself doesn't 'see' that?

[mcabbott]

If you perturb x slightly, the output is usually unchanged. There are contexts where you can think of the gradient of x -> round(x) as being a comb of delta functions, but almost always zero. AD isn't that smart, it isn't going to turn indexing into interpolation at special points.

You may need to replace the function with one which already interpolates smoothly as x is changed.

[rakshith95]

Ah, that makes perfect sense, thank you. I think I have what I need