Adding measurement components back to a measurement after iteratively solving for a value

[Boxylmer]

Hey! I've been using Measurements since I started using Julia, and I've finally run into a project where I need to be able to get components back out via uncertainty_components. I've considered a number of ways to recover the actual partial derivatives, but none are too attractive in terms of effort and how they would work applied generically to any problem.

I'm wondering if, upon calculating the iterative solution, I could go ahead and calculate the partial derivatives of the solution either with finite differences or some other method, then add them back in to the der field of the new measurement I construct. Is this doable?

To reiterative my problem a bit more concretely: here's what I'm doing now and also why it isn't working

1. I've defined some function of measurements x, y, z where f(x, y, z) = a
2. Later I have some other function of measurements u, v, w where g(u, v, w) = b
3. Finally, I have some downstream processing where (among other things) h(a, b) = c, and c is the value I care about.
4. I want to know how much x, y, z, u, v, and w contribute to the error of c, but the problem is, that in the function f, I solve a iteratively, then later go back and analytically solve for the uncertainty of a by solving the variance formula for whatever direct expression was originally available. The result of this is that I had to define a "new" measurement for a, since we originally found it by guessing (through some optimizer), which means it has no partial derivative / history of the pathway.
5. If I can use finite difference partial derivatives for a by varying x, y, and z slightly, then plug their partials back into the der field with appropriate (val, err, tag) keys, then I can send a on its merry way and it's business as usual...

...Or so I think.

Does this problem make sense, and is there an obvious way around it in the measurments library? Or am I not seeing another better solution in finding error contributions?

For clarity, I already am able to get the overall uncertainty and have experimentally confirmed that these uncertainties are correct. The issue is in finding what contributes to uncertainty the most.

[giordano]

Thanks for the interesting question. If I understand it correctly, the Measurement object a which comes from f(x, y, z) doesn't have any link with x, y, and z? If so, you can try using @uncertain. If we're lucky, a = @uncertain f(x, y, z) should do for you what you suggested in point 5.

[Boxylmer]

Thanks for your advice! And yep! Basically, the act of solving for a iteratively breaks the chain of partial derivatives that tie the resulting calculated measurement object back to the fundamental measurements that were used to create it. My question effectively asks if we can somehow manually rebuild that chain.

So I have to ask, what kind of deep magic is the @uncertain macro doing that lets it not only get the overall error, but also get the uncertainty components? If that works, it would save us a great deal of effort, as one of my colleagues is looking into figuring out how to do automatic implicit differentiation on our formula, haha. Also, apologies if this was in the documentation regarding @uncertain also tracking the uncertainty_components! I figured it wouldn't be able to.

For anyone else with this kind of problem where you can't get the components but can get the overall error: My alternative "brute force" solution is to simply solve the calculation again with every variable's uncertainty except for one set to zero. This way, the variance formula will be simplified to only contain that component, and the resulting error will simply be the square root of the variables contribution. For example: if you're starting with the variable x for some function f(x, y, x), the result is sigma_f = sqrt(sigma_x^2 * (df/dx)^2), take this value squared over the true uncertainty using all variables.

[giordano]

So I have to ask, what kind of deep magic is the @uncertain macro doing that lets it not only get the overall error, but also get the uncertainty components?

In @uncertain f(x, y, z) it numerically computes the partial derivatives of f with respect to x, y, and z and builds the resulting Measurement object accordingly. See also the examples in the documentation

[Boxylmer]

I cannot believe I missed this. Thank you Mosè! Heres a MWE demonstrating this working just in case anyone else sees this in the future:

using Optim
using Measurements

f(x, y, z) = - x * y + z*x + z*x*y +sin(y)
# say we know f, x, y but not z
f_known = 1 ± 0.2
x_known = 23 ± 2.4
y_known = 2 ± 0.01

function solve_for_z(x_meas, y_meas, f_meas)
    f_error(z) = (f(x_meas, y_meas, z[1]) - f_meas)^2
    z_solved = Optim.optimize(f_error, [1.], GradientDescent(); autodiff = :forward)
    return Optim.minimizer(z_solved)[1]
end

# Note that running this throws an error as Optim currently cannot handle measurement types
# solve_for_z(x_known, y_known, f_known)  
z = @uncertain solve_for_z(x_known, y_known, f_known)

x_key = (x_known.val, x_known.err, x_known.tag)
y_key = (y_known.val, y_known.err, y_known.tag)
f_key = (f_known.val, f_known.err, f_known.tag)
components = uncertainty_components(z)
contributions = (components[x_key], components[y_key], components[f_key]).^2 ./ z.err^2

println(contributions)
println(sum(contributions))  # if we accounted for all error sources, this will sum to 1

[giordano]

The only caveat is that I'm using Calculus.jl to compute the numerical derivates: so the result makes sense as long as Calculus.jl is able to compute a meaningful numerical derivative. Whether that's possible or not, it's out of my control. Pull request #82 will switch to FiniteDifferences.jl instead of Calculus.jl, but in principle this shouldn't make much difference for the end-user.

[Boxylmer]

Thanks for the tip! I forgot to ask one quick question: will the finite differencing used in calculating the partials through @Uncertain affect the base calculations at all? I know I can forwarddiff through the @Uncertain macro without the typical issues seen before, but I was unsure if there was any additional performance degradation in the actual values themselves.

[giordano]

will the finite differencing used in calculating the partials through @Uncertain affect the base calculations at all?

The value is computed normally, Calculus.jl is used to compute the derivative/gradient of the function and build the error. You can just read the code, it's rather simple:

Measurements.jl/src/math.jl

Lines 150 to 170 in ba56443

	macro uncertain(expr::Expr)
	f = esc(expr.args[1]) # Function name
	n = length(expr.args) - 1
	if n == 1
	a = esc(expr.args[2]) # Argument, of Measurement type
	return quote
	let x = measurement($a)
	result($f(x.val), Calculus.derivative($f, x.val), x)
	end
	end
	else
	a = expr.args[2:end] # Arguments, as an array of expressions
	args = :([]) # Build up array of arguments
	[push!(args.args, :(measurement($(esc(a[i]))))) for i=1:n] # Fill the array
	argsval =:([]) # Build up the array of values of arguments
	[push!(argsval.args, :($(args.args[i]).val)) for i=1:n] # Fill the array
	return :( result($f($argsval...),
	Calculus.gradient(x -> $f(x...), $argsval),
	$args) )
	end
	end

I was unsure if there was any additional performance degradation in the actual values themselves.

Are you talking about speed performance (it'll be slower) or "goodness" performance (the value shouldn't be any different, only the error part would be slightly different compared to not using @uncertain, if that's possible, but usually within sqrt(eps) relative error)