Going back to the OP, CoupledFields is based on kernel regression.
If yᵢ = g(xᵢ), then it follows from kernel ridge regression that:
∇g(xᵢ) = ∇K(xᵢ, ⋅) (Gₓ + 10ᵇnI)⁻¹ Y
where Gₓ is the gram matrix of X, and ∇K() is the gradient function of a kernel function.
In CoupledFields, for gradvecfield([a b], ...), a is a smoothness parameter (that scales an auto-estimate of kernel width), and b is a ridge parameter (as in the equation above).
Using a=1.0 should typically work well.
Here’s an example. The black line is the true value, the gradient function is estimated using gradvecfield([a -7], X, Y, GaussianKP(X)) for 3 values of a.
For high-dimensional problems, the average product of jacobians can be used to find an active subspace.