- A robot uses a range sensor that can measure ranges from 0m and 3m. For simplicity, assume that the actual ranges are distributed uniformly in this interval. Unfortunately, the sensor can be faulty. When the sensor is faulty, it constantly outputs a range below 1m, regardless of the actual range in the sensor's measurement cone. We know that the prior probability for a sensor to be faulty is $p = 0.01$.
Suppose the robot queried it's sensor N times, and every single time the measurement value is below 1m. What is the posterior probability of a sensor fault, for $N = 1, 2, ..., 10$. Formulate the corresponding probabilistic model.
$ P_N = pr(\mathtt{faulty} \mid d^{<1m}_N))$ is the probability that the sensor is faulty given that there have been N instances of a distance less than 1m. $p$ is the probability that the sensor is faulty.
Using Bayes Rule: $$ P_N = \frac{pr(d^{<1m}_N \mid \mathtt{faulty})p}{pr(d^{<1m}_N)}$$
We already know a couple values:
$ pr(d^{<1m}_N \mid \mathtt{faulty}) = 1$
$ p = 0.01 $
One might be tempted to say probability of getting range < 1m is 1/3 (area under uniform probability distribution between 0 and 1), but that assumes the sensor is not faulty. we should use rule of total probability:
$$ pr(d^{<1m}_1) = pr(d^{<1m}_1 \mid \mathtt{not \, faulty}) * pr(\mathtt{not \, faulty}) + pr(d^{<1m}_1 \mid \mathtt{faulty}) * pr(\mathtt{faulty})$$
$$ pr(d^{<1m}_1) = 1/3 * 0.99 + 1.0 * 0.01 = 0.34 $$
So, for $N=1$ $$P_1 = \frac{1.0 \times 0.01}{0.34}$$
P_1 = 0.01 / .34
P_1
There is a 2.9% probability that the sensor is faulty if we receive one distance measurement < 1.
For $N = 2$ we need to calculate the total probability of receiving two measurements less than 1m:
$$ pr(d^{<1m}_2) = pr(d^{<1m}_2 \mid \mathtt{not \, faulty}) * pr(\mathtt{not \, faulty}) + pr(d^{<1m}_2 \mid \mathtt{faulty}) * pr(\mathtt{faulty})$$
Let's assume the robot is moving around and all distance measurements are independent:
$$ pr(d^{<1m}_2 \mid \mathtt{not \, faulty}) = (1/3)^2 $$
$$ pr(d^{<1m}_2) = (1/3)^2 * 0.99 + 1 * 0.01 $$
$$ P_2 = \frac{1.0 \times 0.01}{(1/3)^2 * 0.99 + 1 * 0.01} $$
P_2 = 0.01 / ((1/3)**2 + 0.01)
P_2
The probability of a faulty sensor is 8.2% if we receive two measurements.
Now we can see the pattern:
$$ P_n = \frac{1.0 \times 0.01}{(1/3)^n \times 0.99 + 1 \times 0.01} $$
for n in range(1, 11):
print("{}: {:.3f}".format(n, 0.01/(0.33**n + 0.01)))
The probability that the sensor is faulty increases rapidly with an increasing number of measurements.
On the other hand, if the robot and item to which distance is being measured are both stationary, then the subsequent distance measurements are not independent - we would expect to obtain the same distance measurement over and over. In this scenario:
$$ pr(d^{<1m}_n \mid \mathtt{not \, faulty}) = 1/3 $$
The probability is independent of the number of measurements n. the probability that the sensor is faulty will be 2.9% regardless of the number of measurements taken.
Suppose we live in a place where days are either sunny, cloudy, or rainy. The weather transition function is a Markov chain with the following transition table:
Table 2-1
| Today's Weather | Tomorrow Sunny | Tomorrow Cloudy | Tomorrow Rainy |
|---|---|---|---|
| Sunny | 0.8 | 0.2 | 0 |
| Cloudy | 0.4 | 0.4 | 0.2 |
| Rainy | 0.2 | 0.6 | 0.2 |
Suppose Day1 is sunny.
What is the probability of this sequence of days: Day2 = cloudy, Day3 = cloudy, Day4 = rainy
$p = 0.2 \times 0.4 \times 0.4 \times 0.6$
p = 0.2*0.4*0.4*0.6
"{:.4}".format(p)
Write a simulator that can randomly generate sequences of weather from state transition function.
import random
def weather(day1="sunny", num_days=1000000):
today = day1
tomorrow = None
sunny = 0
cloudy = 0
rainy = 0
total_days = 0
prob = {"sunny": [0.8, 1.0], "cloudy": [0.4, 0.8], "rainy": [0.2, 0.8]}
for day in range(num_days):
rnd = random.uniform(0, 1)
if rnd <= prob[today][0]:
tomorrow = "sunny"
sunny = sunny + 1
elif rnd <= prob[today][1]:
tomorrow = "cloudy"
cloudy = cloudy + 1
else:
tomorrow = "rainy"
rainy = rainy + 1
total_days = total_days + 1
today = tomorrow
print("Probability Sunny: {}".format(sunny/total_days))
print("Probability Cloudy: {}".format(cloudy/total_days))
print("Probability Rainy: {}".format(rainy/total_days))
Use the simulator to determine the stationary distribution of the Markov chain -- the probability that any given day will be sunny, cloudy, or rainy.
weather()
weather(day1="cloudy")
weather(day1="rainy")
Calculate the stationary probability distribution using a closed form solution.
Using the total probability theorem:
$$\tag{2-1} p(x) = \sum_{y}p(x \mid y)p(y)$$
$p(s) = $ probability sunny
$p(c) = $ probability cloudy
$p(r) = $ probability rainy
$$ p(s) = p(s \mid s)p(s) + p(s \mid c)p(c) + p(s \mid r)p(r) $$ $$ p(c) = p(c \mid s)p(s) + p(c \mid c)p(c) + p(c \mid r)p(r) $$ $$ p(r) = p(r \mid s)p(s) + p(r \mid c)p(c) + p(r \mid r)p(r) $$
And from the definition of probability: $ p(s) + p(c) + p(r) = 1 $
Substituting values from the weather transition function:
$$ p(s) = 0.8p(s) + 0.4p(c) + 0.2p(r) $$ $$ p(c) = 0.2p(s) + 0.4p(c) + 0.6p(r) $$ $$ p(r) = 0p(s) + 0.2p(c) + 0.2p(r) $$
Collecting like terms:
$$ -0.2p(s) + 0.4p(c) + 0.2p(r) = 0 $$ $$ 0.2p(s) - 0.6p(c) + 0.6p(r) = 0 $$ $$ 0p(s) + 0.2p(c) - 0.8p(r) = 0 $$
Probabilities must sum to 1:
$$ p(s) + p(c) + p(r) = 1 $$
Using only two of the total probability equations (the third equation is a linear combination of the first two) and the fact that the probabilities sum to 1, we get a set of 3 equations with three unknowns that are represented by the following matrix:
$$ \begin{bmatrix} -0.2 & 0.4 & 0.2 & 0\\ 0.2 & -0.6 & 0.6 & 0\\ 1 & 1 & 1 & 1\\ \end{bmatrix} $$
We'll use the Python sympy package to solve for $p(s)$, $p(c)$, and $p(r)$.
import sympy
aug = sympy.Matrix([
[0, 0.2, -0.8, 0],
[0.2, -0.6, 0.6, 0],
[1, 1, 1, 1]
])
s, c, r = sympy.symbols('s, c, r')
stationary = sympy.linsolve(aug, [s, c, r])
stationary
Calculate the entropy of the stationary distribution.
$$ H_p(x) = -\sum_x p(x) \, log_2 \, p(x) $$
$$ H_p(x) = 0.6429 \times log_2(0.6429) + 0.2857 \times log_2(0.2857 + 0.0714 \times log_2(0.0714) $$
from math import log
probs = stationary.args[0]
H = (probs[0]*log(probs[0], 2) + probs[1]*log(probs[1], 2) + probs[2]*log(probs[2], 2)) * -1
H
Calculate the probability of yesterday's weather given today's weather.
Bayes rule: $$ p(x \mid y) = \frac{p(y \mid x) \, p(x)}{p(y)} $$
$ p(s_t) = $ probability that today is sunny
$ p(s_y) = $ probability that yesterday was sunny
$ p(c_t) = $ probability that today is cloudy
$ p(c_y) = $ probability that yesterday was cloudy
$ p(r_t) = $ probability that today is rainy
$ p(r_y) = $ probability that yesterday was rainy
Note that $ p(s_t) = p(s_y) = p(s) = $ stationary probability that any particular day is sunny.
Note that $ p(c_t) = p(c_y) = p(c) = $ stationary probability that any particular day is cloudy.
Note that $ p(r_t) = p(r_y) = p(r) = $ stationary probability that any particular day is rainy.
# Stationary probabilities
p_sun = probs[0]
p_cld = probs[1]
p_rny = probs[2]
For today and yesterday being sunny: $$ p(s_y \vert s_t) = \frac{p(s_t \vert s_y) p(s_y)}{p(s_t)} $$ $$ p(s_y) = p(s_t) $$ $$ p(s_y \vert s_t) = p(s_t \vert s_y) = 0.8 $$
Similarly for today and yesterday both being cloudy or rainy $$ p(c_y \vert c_t) = p(c_t \vert c_y) = 0.4 $$ $$ p(r_y \vert r_t) = p(r_t \vert r_y) = 0.2 $$
Yesterday sunny, today cloudy: $$ p(s_y \vert c_t) = \frac{p(c_t \vert s_y)p(s_y)}{p(c_t)} $$
sun_yest_cld_today = 0.2 * p_sun / p_cld
sun_yest_cld_today
Yesterday sunny, today rainy: $$ p(s_y \vert r_t) = \frac{p(r_t \vert s_y)p(s_y)}{p(r_t)} $$
sun_yest_rny_today = 0 * p_sun / p_rny
sun_yest_rny_today
Yesterday cloudy, today sunny: $$ p(c_y \vert s_t) = \frac{p(s_t \vert c_y)p(c_y)}{p(s_t)} $$
cld_yest_sun_today = 0.4 * p_cld / p_sun
cld_yest_sun_today
Yesterday cloudy, today rainy: $$ p(c_y \vert r_t) = \frac{p(r_t \vert c_y)p(c_y)}{p(r_t)} $$
cld_yest_rny_today = 0.2 * p_cld / p_rny
cld_yest_rny_today
Yesterday rainy, today sunny: $$ p(r_y \vert s_t) = \frac{p(s_t \vert r_y)p(r_y)}{p(s_t)} $$
rny_yest_sun_today = 0.2 * p_rny / p_sun
rny_yest_sun_today
Yesterday rainy, today cloudy: $$ p(r_y \vert c_t) = \frac{p(c_t \vert r_y)p(r_y)}{p(c_t)} $$
rny_yest_cld_today = 0.6 * p_rny / p_cld
rny_yest_cld_today
Final Results
| Today's Weather | Yesterday Sunny | Yesterday Cloudy | Yesterday Rainy |
|---|---|---|---|
| Sunny | 0.8 | 0.1778 (8/45) | 0.0222 (2/9) |
| Cloudy | 0.45 | 0.4 | 0.15 |
| Rainy | 0 | 0.8 | 0.2 |
Add seasons to model. State transition function above would apply only to Summer. Other transition functions would apply to Winter, Spring, and Fall. Would this violate the Markov property of this process?
As long as we added the current date as a state variable and formulated the trnasition function such that the season and correct transitition probabilities could be selected using the date and today's weather, the Markov property would not be violated. But if the set of state variables does not incldue the date, then the Markov property is violated.
Stacy Irwin, 8 July 2018