Standard error of stationary point of a 2nd order polynomial whose coefficient std errors are known I have a 3-variable full quadratic polynomial
$$
f(x,y,z) = c_0 + c_1 x + c_2 y + c_3 z + c_4 x y + c_5 x z + c_6 y z + c_7 x^2 + c_8 y^2 + c_9 z^2
$$
Each of the coefficients $c_i$ are known with a standard error $s_i$.
Ignoring the standard errors $s_i$, I have one single stationary point, with values $(x_0, y_0, z_0)$ (thus satisfying $f'(x_0,y_0,z_0) = 0$), but I wonder: how to deduce, from the standard errors of the coefficients, the standard errors for $x_0$, $y_0$ and $z_0$ ? If it helps, the errors in the coefficients are (1) very small, and (2) uncorrelated.
I can see that this problem can be reformulated as "transferring" the errors from the coefficients $c_i$ toward the stationary point $(x_0,y_0,z_0)$ in the system of first order equations resulting from $f'(x_0,y_0,z_0) = 0$:
\begin{cases} 
c_1 + c_4 y_0 + c_5 z_0 + 2 c_7 x_0 = 0 \\
c_2 + c_4 x_0 + c_6 z_0 + 2 c_8 y_0 = 0 \\
c_3 + c_5 x_0 + c_6 y_0 + 2 c_9 z_0 = 0
\end{cases} 
Edit:
I have tried to simplify Aleco's expression but started to get the jitters and then flew away after seeing the expression of $C^{-1}$:
$$
\frac{1}{2 (c_9 c_4^2-c_5 c_6 c_4+c_6^2 c_7+c_8 (c_5^2-4 c_7 c_9))} %\left( 
\begin{array}{ccc}
c_6^2-4 c_8 c_9 & 2 c_4 c_9-c_5 c_6 & 2 c_5 c_8-c_4 c_6 \\
2 c_4 c_9-c_5 c_6 & c_5^2-4 c_7 c_9 & 2 c_6 c_7-c_4 c_5 \\
2 c_5 c_8-c_4 c_6 & 2 c_6 c_7-c_4 c_5 & c_4^2-4 c_7 c_8 \end{array} %\right)
$$
What I ended up doing is hypothesizing further that the coefficients were normally distributed. Then I drew a sample of 1000 coefficients for the $c_i$ from $N(c_i,s_i)$, and computed the resulting stationary point with each of the 1000 set of coefficients. These stationary points were also normally distributed, and I was finally able to deduce from their distribution their standard error.
 A: Matrix algebra may help.  
Manipulate the steady-state equations as follows:
Define the $3\times 3$ coefficient matrix
$$C = \left [\begin{matrix}
2 c_7 & c_4 & c_5 \\
c_4 & 2 c_8 & c_6 \\
c_5 & c_6 & 2 c_9 \\
\end{matrix} \right]
$$
and note that it is symmetric, $C'=C$, and that the diagonal elements are uncorrelated with the off-diagonal elements.
Define the $3\times 1$ coefficient column vector
$$a = \left [\begin{matrix}
c_1 \\
c_2  \\
c_3 \\
\end{matrix} \right]
$$ 
and note that it contains elements that are uncorrelated with the elements of the $C$ matrix.
Finally, define the $3\times 1$ variables column vector
$$w_0 = \left [\begin{matrix}
x_0 \\
y_0  \\
z_0 \\
\end{matrix} \right] $$
Then the steady state can be expressed as
$$ Cw_0 + a = \mathbf 0$$
Assuming that $C$ is invertible we have
$$ w_0 =-C^{-1}a $$
Since you are talking about standard deviations, all these are random variables - random vectors.
The co-variance matrix of $w_0$ is defined as
$$\text {Cov}(w_0) = E(w_0w_0') - E(w_0)(E(w_0))' $$
$$ = E\Big(-C^{-1}a\Big)\Big(-C^{-1}a\Big)'-E\Big(-C^{-1}a\Big)\left[E\Big(-C^{-1}a\Big)\right]'$$
$$ = E\Big(C^{-1}aa'C^{-1}\Big)-E\Big(C^{-1}a\Big)\left[E\Big(C^{-1}a\Big)\right]'$$
where we have used the fact that the inverse of a symmetric matrix is also symmetric.
If uncorrelatedness between $C$ and $a$ can be extended to independence we will further have
$$\text {Cov}(w_0) = E\Big(C^{-1}aa'C^{-1}\Big)-E\left(C^{-1}\right)E(a)E(a)'E\left(C^{-1}\right)$$
$$ = E\Big(C^{-1}aa'C^{-1}\Big)+E\left(C^{-1}\right)(-E(a)E(a)')E\left(C^{-1}\right)$$
and using the analogous expression for the covariance matrix for $a$ as for $w_0$
$$ = E\Big(C^{-1}aa'C^{-1}\Big)+E\left(C^{-1}\right)\left[\text {Cov}(a) -E(aa')\right]E\left(C^{-1}\right)$$
$$\Rightarrow \text {Cov}(w_0)= E\left(C^{-1}\right)\text {Cov}(a)E\left(C^{-1}\right) + E\Big(C^{-1}aa'C^{-1}\Big) -E\left(C^{-1}\right)E(aa')E\left(C^{-1}\right)$$
Now you should work on the $C^{-1}$ matrix, using various methods, properties, decompositions and shortcuts that exist for the inverse of a symmetric matrix.  Remember that you are only interested in the diagonal elements of $ \text {Cov}(w_0)$. The term that may give you trouble in terms of calculating its expected value is  $E\Big(C^{-1}aa'C^{-1}\Big)$, because even under independence, it may contain quotients of random variables, in which case even if they are independent, the expected value of the quotient does not equal the quotient of the expected values. But you must arrive there to find out. Then, approximations and accepting some bias.
