Examples are taken from http://ipl.physics.harvard.edu/wp-uploads/2013/03/PS3_Error_Propagation_sp13.pdf and used on MCMC to show how the answers work
Example: suppose you measure the height H of a door and get 2.00 ± 0.03 m. This means that H = 2.00 m and δH = 0.03 m. The door has a knob which is a height h = 0.88 ± 0.04 m from the bottom of the door. Then the distance from the doorknob to the top of the door is Q = H − h = 1.12 m. What is the uncertainty in Q?
Q = 1.12 ± 0.05 m
In [2]:
import numpy as np
import pymc as mc
In [2]:
H = mc.Normal('H', 2.00, (0.03)**-2)
h = mc.Normal('h', 0.88, (0.04)**-2)
@mc.deterministic()
def Q(H=H, h=h):
return H-h
model = mc.MCMC((H,h,Q))
In [3]:
model.sample(1e4, burn=100, burn_till_tuned=True)
# mc.Matplot.plot(model)
# mc.Matplot.plot(Q)
print(Q.summary())
print("MCMC gives {0:.2f} +/- {1:.2f}, analytic gives {2} +/- {3}".format(np.mean(Q.trace()), np.std(Q.trace()), 1.12, 0.05))
Example: a bird flies a distance d = 120 ± 3 m during a time t = 20.0 ± 1.2 s. The average speed of the bird is v = d/t = 6 m/s. What is the uncertainty of v?
0.39 m/s.
In [4]:
d = mc.Normal('d', 123, (3)**-2)
t = mc.Normal('t', 20.0, (1.2)**-2)
@mc.deterministic()
def v(d=d, t=t):
return d/t
model = mc.MCMC((d, t, v))
In [5]:
model.sample(1e5, burn=1000, burn_till_tuned=True)
print(v.summary())
print("MCMC gives {0:.2f}, analytic gives {1}".format(np.std(v.trace()), 0.39))
Example: the period of an oscillation is measured to be T = 0.20 ± 0.01 s. Thus the frequency is f = 1/T = 5 Hz. What is the uncertainty in f? Answer: the percent uncertainty in T was 0.01/0.20 = 5%. Thus the percent uncertainty in f is also 5%, which means that δf = 0.25 Hz. So f = 5.0 ± 0.3 Hz (after rounding).
f = 5.0 ± 0.3 Hz
In [6]:
T = mc.Normal('T', 0.20, (0.01)**-2)
@mc.deterministic()
def f(T=T):
return 1/T
model = mc.MCMC((T, f))
In [7]:
model.sample(1e4, burn=100, burn_till_tuned=True)
In [8]:
print(f.summary())
print("MCMC gives {0:.1f} +/- {1:.1f}, analytic gives {2} +/- {3}".format(np.mean(f.trace()), np.std(f.trace()),
5.0, 0.3))
Example: a ball is tossed straight up into the air with initial speed v0 = 4.0 ± 0.2 m/s. After a time t = 0.60±0.06 s, the height of the ball is y = v0t− 1 2 gt2 = 0.636 m. What is the uncertainty of y? Assume g = 9.80 m/s2 (no uncertainty in g).
Thus y would be properly reported as 0.6 ± 0.4 m.
In [9]:
g = 9.80
t = mc.Normal('t', 0.60, (0.06)**-2)
v0 = mc.Normal('v0', 4.0, (0.2)**-2)
@mc.deterministic()
def h(t=t, v0=v0):
return v0*t - 0.5*g*t**2
model = mc.MCMC((t, v0, h))
In [10]:
model.sample(1e5, burn=100, burn_till_tuned=True)
In [11]:
print(h.summary())
print("MCMC gives {0:.1f} +/- {1:.1f}, analytic gives {2} +/- {3}".format(np.mean(h.trace()), np.std(h.trace()),
0.6, 0.4))
In [12]:
A = mc.Normal('A', 3.6, (0.2)**-2)
B = mc.Normal('B', 3.3, (0.3)**-2)
@mc.deterministic()
def D(A=A, B=B):
return A-B
model = mc.MCMC((A, B, D))
In [13]:
model.sample(1e5, burn=100, burn_till_tuned=True)
In [14]:
print(D.summary())
print("MCMC gives {0:.1f} +/- {1:.1f}, analytic gives {2} +/- {3}".format(np.mean(D.trace()), np.std(D.trace()),
0.3, 0.4))
# mc.Matplot.plot(model)
In [15]:
# https://www.lhup.edu/~dsimanek/scenario/errorman/propagat.htm
Example: An angle is measured to be 30° ±0.5°. What is the error in the sine of this angle?
Solution: Use your electronic calculator. The sine of 30° is 0.5; the sine of 30.5° is 0.508; the sine of 29.5° is 0.492. So if the angle is one half degree too large the sine becomes 0.008 larger, and if it were half a degree too small the sine becomes 0.008 smaller. (The change happens to be nearly the same size in both cases.) So the error in the sine would be written ±0.008.
The size of the error in trigonometric functions depends not only on the size of the error in the angle, but also on the size of the angle. A one half degree error in an angle of 90° would give an error of only 0.00004 in the sine.
In [3]:
a = mc.Normal('a', 30.00, (0.5)**-2)
@mc.deterministic()
def sina(a=a):
return np.sin(np.deg2rad(a))
model = mc.MCMC((a, sina))
In [14]:
model.sample(1e5, burn=5000, burn_till_tuned=False)
# mc.Matplot.plot(model)
# mc.Matplot.plot(Q)
print(sina.summary())
print("MCMC gives {0:.3f} +/- {1:.3f}, analytic gives {2} +/- {3}".format(np.mean(sina.trace()), np.std(sina.trace()), 0.5, 0.008))
for i in range(4):
model.sample(1e5, burn=5000, burn_till_tuned=False)
mc.Matplot.summary_plot(sina)
Two data quantities, X and Y, are used to calculate a result, R = XY. X = 38.2 ± 0.3 and Y = 12.1 ± 0.2. What is the error in R?
The product rule requires fractional error measure. The fractional error in X is 0.3/38.2 = 0.008 approximately, and the fractional error in Y is 0.017 approximately. Adding these gives the fractional error in R: 0.025. Multiplying this result by R gives 11.56 as the absolute error in R, so we write the result as R = 462 ± 12. Note that once we know the error, its size tells us how far to round off the result (retaining the first uncertain digit.) Note also that we round off the error itself to one, or at most two, digits. This is why we could safely make approximations during the calculations of the errors.
In [18]:
X = mc.Normal('X', 38.2, (0.3)**-2)
Y = mc.Normal('Y', 12.1, (0.2)**-2)
@mc.deterministic()
def R(X=X, Y=Y):
return X*Y
model = mc.MCMC((X, Y, R))
model.sample(1e5, burn=5000, burn_till_tuned=False)
# mc.Matplot.plot(model)
# mc.Matplot.plot(Q)
print(R.summary())
print("MCMC gives {0:.3f} +/- {1:.3f}, analytic gives {2} +/- {3}".format(np.mean(R.trace()), np.std(R.trace()), 462, 12))
In [ ]: