```
In [1]:
```from __future__ import print_function, division
import thinkbayes2
import thinkplot
import numpy as np
from scipy import stats
%matplotlib inline

`prior`

=$x$, and the likelihoods of the data, `prob_given_sick`

=$P(fever|sick)$ and `prob_given_not`

=$P(fever|not sick)$

```
In [2]:
```def update(prior, prob_given_sick, prob_given_not):
suite = thinkbayes2.Suite()
suite['sick'] = prior * prob_given_sick
suite['not sick'] = (1-prior) * prob_given_not
suite.Normalize()
return suite['sick']

```
In [3]:
```prior = 0.1
prob_given_sick = 0.9
prob_given_not = 0.3
post = update(prior, prob_given_sick, prob_given_not)
post

```
Out[3]:
```

`prob_given_sick`

=$P(fever|sick)$ and $t = $ `prob_given_not`

=$P(fever|not sick)$, but we think they are uniformly distributed and independent.

```
In [4]:
```dist_s = thinkbayes2.Beta(1, 1)
dist_t = thinkbayes2.Beta(1, 1)
dist_s.Mean(), dist_t.Mean()

```
Out[4]:
```

```
In [5]:
```n = 1000
ss = dist_s.Sample(n)
ts = dist_t.Sample(n)

Just checking that the samples have the right distributions:

```
In [6]:
```thinkplot.Cdf(thinkbayes2.Cdf(ss))
thinkplot.Cdf(thinkbayes2.Cdf(ts))

```
Out[6]:
```

Now computing the posteriors:

```
In [7]:
```posts = [update(prior, s, t) for s, t in zip(ss, ts)]

Here's what the distribution of values for $x\prime$ looks like:

```
In [8]:
```cdf = thinkbayes2.Cdf(posts)
thinkplot.Cdf(cdf)

```
Out[8]:
```

And here's the mean:

```
In [9]:
```cdf.Mean()

```
Out[9]:
```

This result implies that if our prior probability for $x$ is 0.1, and then we learn that the patient has a fever, we should be uncertain about $x\prime$, and this distribution describes that uncertainty. It says that the fever probably has little predictive power, but might have quite a lot.

The mean of this distribution is a little higher than the prior, which suggests that our priors for $s$ and $t$ are not neutral with respect to updating $x$. It's surprising that the effect is not symmetric, because our beliefs about $s$ and $t$ are symmetric. But then again, we just computed an arithmetic mean on a set of probabilities, which is a bogus kind of thing to do. So maybe we deserve what we got.

```
In [10]:
```dist_s = thinkbayes2.Beta(1, 1)
dist_t = thinkbayes2.Beta(1.75, 1)
n = 10000
ss = dist_s.Sample(n)
ts = dist_t.Sample(n)
thinkplot.Cdf(thinkbayes2.Cdf(ss))
thinkplot.Cdf(thinkbayes2.Cdf(ts))

```
Out[10]:
```

```
In [11]:
```posts = [update(prior, s, t) for s, t in zip(ss, ts)]
np.array(posts).mean()

```
Out[11]:
```

```
In [12]:
```def prob_sick(x, s, t):
return x * s / (x * s + (1-x) * t)

```
In [13]:
```dist_s = thinkbayes2.Beta(1, 1)
dist_t = thinkbayes2.Beta(1, 1)
n = 10000
ss = dist_s.Sample(n)
ts = dist_t.Sample(n)

```
In [14]:
```x = 0.1
probs = [prob_sick(x, s, t) for s, t in zip(ss, ts)]

```
In [15]:
```np.array(probs).mean()

```
Out[15]:
```

```
In [16]:
```cohort = np.random.random(len(probs)) < probs

```
In [17]:
```cohort.mean()

```
Out[17]:
```

April 6, 2016

Suppose

- t is known to be 0.5
- x is known to be 0.1
- s is equally likely to be 0.2 or 0.8, but we don't know which

If we take the average value of s and compute p = p(sick|fever), we would consider p to be the known quantity 0.1

```
In [26]:
```s = 0.5
t = 0.5
x = 0.1
p = x * s / (x * s + (1-x) * t)
p

```
Out[26]:
```

```
In [36]:
```t = 0.8
p1 = x * s / (x * s + (1-x) * t)
p1, 5/77

```
Out[36]:
```

```
In [35]:
```t = 0.2
p2 = x * s / (x * s + (1-x) * t)
p2, 5/23

```
Out[35]:
```

```
In [29]:
```(p1 + p2) / 2

```
Out[29]:
```

```
In [34]:
```def logodds(p):
return np.log(p / (1-p))
logodds(p1) - logodds(0.1), logodds(p2) - logodds(0.1)

```
Out[34]:
```

```
In [37]:
```p_mix = (p1 + p2) / 2
logodds(p_mix) - logodds(0.1)

```
Out[37]:
```

```
In [30]:
```import random
def simulate_patient():
s = 0.5
t = random.choice([0.2, 0.8])
x = 0.1
p = x * s / (x * s + (1-x) * t)
return random.random() < p
simulate_patient()

```
Out[30]:
```

In this simulation, the fraction of patients with fever who turn out to be sick is close to 0.141

```
In [31]:
```patients = [simulate_patient() for _ in range(10000)]
sum(patients) / len(patients)

```
Out[31]:
```

```
In [153]:
```x = 0.1
s = 0.5
t1 = 0.2
t2 = 0.8

```
In [154]:
```import pandas as pd
d1 = dict(feverbad=0.5, fevergood=0.5)
d2 = dict(sick=x, notsick=1-x)
d3 = dict(fever='t', notfever='1-t')
iterables = [d1.keys(), d2.keys(), d3.keys()]
index = pd.MultiIndex.from_product(iterables, names=['first', 'second', 'third'])
df = pd.DataFrame(np.zeros(8), index=index, columns=['prob'])

```
In [155]:
```t_map = dict(fevergood=t2, feverbad=t1)
for label1, p1 in d1.items():
t = t_map[label1]
for label2, p2 in d2.items():
for label3, p3 in d3.items():
if label2 == 'sick':
p = p1 * p2 * 0.5
else:
p = p1 * p2 * eval(p3)
df.prob[label1, label2, label3] = p
df

```
Out[155]:
```

```
In [157]:
```df.prob[:, 'sick', 'fever'].sum() / df.prob[:, :, 'fever'].sum()

```
Out[157]:
```

If we know fevergood, P(sick|fever) is

```
In [158]:
```p_sick_fevergood = df.prob['fevergood', 'sick', 'fever'].sum() / df.prob['fevergood', :, 'fever'].sum()
p_sick_fevergood

```
Out[158]:
```

If we know feverbad, P(sick|fever) is

```
In [149]:
```p_sick_feverbad = df.prob['feverbad', 'sick', 'fever'].sum() / df.prob['feverbad', :, 'fever'].sum()
p_sick_feverbad

```
Out[149]:
```

If we think there's a 50-50 chance of feverbad

```
In [159]:
```(p_sick_fevergood + p_sick_feverbad) / 2

```
Out[159]:
```

If fevergood, here's the fraction of all patients with fever:

```
In [143]:
```p_fever_fevergood = df.prob['fevergood', :, 'fever'].sum()
p_fever_fevergood

```
Out[143]:
```

If fever bad, here's the fraction of all patients with fever

```
In [144]:
```p_fever_feverbad = df.prob['feverbad', :, 'fever'].sum()
p_fever_feverbad

```
Out[144]:
```

```
In [146]:
```p_feverbad_fever = p_fever_feverbad / (p_fever_feverbad + p_fever_fevergood)
p_feverbad_fever

```
Out[146]:
```

And if we compute the weighted sum of the two possible worlds

```
In [150]:
```p_feverbad_fever * p_sick_feverbad + (1-p_feverbad_fever) * p_sick_fevergood

```
Out[150]:
```

```
In [ ]:
```