Assume that I have a resistor with an unknown resistance R. I measure the current trough it, I, and the voltage across it, U. I can then estimate R: $$ \hat{R} = \frac{U}{I} $$ Let us assume I get the following result: $$ \hat{R} = \frac{5V}{100mA} = 50 \Omega $$ But how certain is this estimate? Let us assume that the uncertainty in U is: $U \pm \delta U$ where $\delta U = 1mV$ and that the uncertainty in I is: $I \pm \delta I$ where $\delta I = 1mA$.
How do I calculate the the uncertainty in R ($R \pm \delta R$)?
Taylor approximation: $$ R \approx R_0 + \frac{\partial R}{\partial U}dU + \frac{\partial R}{\partial I}dI $$ Definition of variance: $$ Var[R] = E[(R - \mu_R)^2] = E[(R - R_0)^2] = E[(\frac{\partial R}{\partial U}dU + \frac{\partial R}{\partial I}dI)^2] \approx E[(\frac{\partial R}{\partial U}dU)^2 + (\frac{\partial R}{\partial I}dI)^2] $$ where second order terms ($dUdI$) have been omitted. Now $$ E[dU]=E[U - \mu_U]= \sigma_U $$ and $$ E[dI]=E[I - \mu_I]= \sigma_I $$ Therefore: $$ \sigma_R^2 \approx (\frac{\partial R}{\partial U}\sigma_U)^2 + (\frac{\partial R}{\partial I}\sigma_I)^2 = (\frac{\sigma_U}{I})^2 + (\frac{U\sigma_I}{I^2})^2 = (\frac{5e-4}{0.1})^2 + (\frac{5*5e-4}{0.1^2})^2 $$ $$ \sigma_R \approx 5e-4\sqrt{10^2+500^2} = 0.25 $$ Is my derivation above correct?
However one practical problem remains; R is NOT normally distributed! In fact no simple analytical expression exists for the distribution of R. The same will be true for many other problems.
How can I find a 95 % interval around $\mu$ for an unknown distribution where I only know $\mu$ and $\sigma$?
I suspect the answer is that I can not. If so the error propagation approach above seem fundamentally flawed for real world problems?
However I ran a Monte Carlo simulation and found:
R_min=4.951e+01,R_max=5.050e+01
which is almost the same as predicted by $R=\mu_R \pm 2*\sigma_R$! Just dumb luck? Or does the Central limit theorem come into action here?
Here is my code:
Here is my code:
import numpy as np
import matplotlib.pyplot as plt
mu_U = 5
sigma_U = 5e-4
mu_I = 0.1
sigma_I = 5e-4
N = 100*1000
U = np.random.normal(mu_U, sigma_U, N)
I = np.random.normal(mu_I, sigma_I, N)
R = U/I
R_min = np.percentile(R, 2.5)
R_max = np.percentile(R, 97.5)
print(f"R_min={R_min:.3e},R_max={R_max:.3e}")
plt.hist(R, bins=1000)
plt.show()
(1)