1 | initial version |

**First:**
The probability for $X>40$, where $X$ is binomially distributed as in the problem, is given by

```
p = sum( [ binomial(64,k) for k in [41..64] ] ) / 2^64
```

which is exactly:

```
sage: p
302210786238405829/18446744073709551616
sage: p.n()
0.0163828795494116
```

which is also the number from the posted code:

```
sage: pSupKmax
0.0163828795494116
```

It is also excluded, that the "more than $40$ heads" in the text translates as $X\ge 40$, since the "other number" is

```
sage: ( sum( [ binomial(64,k) for k in [40..64] ] ) / 2^64 ).n()
0.0299705947834996
```

**Second:**
I do not know which is the "rule of thumb", so i try to guess. The random variable $X$ is binomially distributed, its mean is $64\cdot \frac 12 = 32$, its variance is $64\cdot \frac 12\cdot\left(1- \frac 12\right) = 16=4^2$, so the variable that can be approximated with a *standard* normally distributed variable $Y$ is $(X-32)/4$,
$$ Y\approx \frac 14(X-32)\ .$$
Let $\Phi$ be the distribution of $Y$.
Then in the statistics one "knows" the following values for $\Phi$, corresponding to "confidence intervals" of often usage:
68-95-99 rule
\begin{align}
\mathbb P(\ |X-32| \le 1\cdot 4\ )& \approx \mathbb P( |Y| <1) \approx 0.6827\ ,\\\\
\mathbb P(\ |X-32| \le 2\cdot 4\ )& \approx \mathbb P( |Y| <2 ) \approx 0.9545\ ,\\\\
\mathbb P(\ |X-32| \le 3\cdot 4\ )& \approx \mathbb P( |Y| <3) \approx 0.9974\ .\\\\
\end{align}
In our case, the problem wants "by coincidence" the number that we can easily extract from the second line.
We have quickly
$$\mathbb P(Y>2) = \frac 12(\mathbb P(Y<-2)+\mathbb P(Y>2)) =\frac 12(1-\mathbb P(|Y|<2))\approx \frac 12 \cdot 0.0455
=0.02275\ .$$
This last value correpsonds well to the value `1-Phi(2) = 0.0227501319482`

in the following sage dialog:

```
sage: T = RealDistribution( 'gaussian', 1 )
sage: for k in [0..3]:
....: print "Phi(%s) = %s" % ( k, T.cum_distribution_function(k) )
....:
Phi(0) = 0.5
Phi(1) = 0.841344746069
Phi(2) = 0.977249868052
Phi(3) = 0.998650101968
sage: for k in [0..3]:
....: print "1-Phi(%s) = %s" % ( k, 1-T.cum_distribution_function(k) )
....:
1-Phi(0) = 0.5
1-Phi(1) = 0.158655253931
1-Phi(2) = 0.0227501319482
1-Phi(3) = 0.00134989803163
```

Note: In statistics, there is one more trick used when approximating some "binomial $X$" by some "standard normal $Y$" as in our case. The random variable $X$ is discrete, so "more than $40$" is both $X>40$ and $X\ge 41$. When we translate this in terms of $Y$ it is better to translate using the middle point $40.5$. And we get:

```
sage: 1 - T.cum_distribution_function( (40.5-32)/4 )
0.016793306448448786
```

which fits nicely in the "reality", which is `0.0163828795494116`

. This important experience can be done only by using something like sage, so they should introduce it in the school, since is the only free "thing" doing the job (without pain)... !

**Third**
Sorry, i could not reproduce the error, at any rate, the following works for me...

```
sage: show( 'Using normal law integration 1-normalCDFnum :'
....: , 1 - T.cum_distribution_function( (40.5-32)/4 ) )
\newcommand{\Bold}[1]{\mathbf{#1}}\verb|Using|\phantom{\verb!xx!}\verb|normal|\phantom{\verb!x!}\verb|law|\phantom{\verb!x!}\verb|integration|\phantom{\verb!x!}\verb|1-normalCDFnum|\phantom{\verb!x!}\verb|:| 0.0167933064484
sage: show( 'Using normal law integration 1-normalCDFnum :'
....: , ( 1 - T.cum_distribution_function( (40.5-32)/4 ) ).n() )
\newcommand{\Bold}[1]{\mathbf{#1}}\verb|Using|\phantom{\verb!xx!}\verb|normal|\phantom{\verb!x!}\verb|law|\phantom{\verb!x!}\verb|integration|\phantom{\verb!x!}\verb|1-normalCDFnum|\phantom{\verb!x!}\verb|:| 0.0167933064484488
```

Copyright Sage, 2010. Some rights reserved under creative commons license. Content on this site is licensed under a Creative Commons Attribution Share Alike 3.0 license.