So, we had Secret Santa at work.

We’re 8 people.

We each took turns and pulled a small piece of paper from a bowl with a name on it.

The only rule : If you pull your name, you have to put the piece of paper back in the bowl and try again.Let’s call the people A, B, C, D, E, F, G, H, which is also the order in which they picked their piece of paper.

We did the gift exchange last night.

A was F’s secret santa.

B was E’s secret santa.

C was D’s secret santa.

D was C’s secret santa.

E was B’s secret santa.

F was A’s secret santa.

G was H’s secret santa.

H was G’s secret santa.See what happened? We made couples.

A and F were each other’s secret santa.

B and E were each other’s secret santa.

C and D were each other’s secret santa.

G and H were each other’s secret santa.What is the probability of this happening and how do you calculate it?

**Answer**

The total number of assignments among 2n people, where nobody is assigned to themselves, is d(2n) = (2n)!(1/2 – 1/6 + \cdots + (-1)^k/k! + \cdots + 1/(2n)!). (These are called derangements.) The value is very close to (2n)! / e.

If they correspond to perfect pairings, then they are a product of disjoint transpositions. This implies their cycle structure is of the form

(a_{11}a_{12})(a_{21}a_{22})\cdots(a_{n1}a_{n2}).

The number of distinct such patterns is the order of the group of all permutations of the 2n names divided by the order of the stabilizer of the pattern. A stabilizing element may swap any number of the pairs and it may also permute the n! pairs, whence there are 2^n n! stabilizing elements. Therefore there are

p(2n) = \frac{(2n)!}{2^n n!}

such pairings.

Since all such perfect pairings are derangements, and all derangements are equally likely, the chance equals

\frac{p(2n)}{d(2n)} = \frac{1}{2^n n!(1 – 1/2 + 1/6 – \cdots + (-1)^k/k! + \cdots + 1/(2n)!)} \approx \frac{e}{2^n n!}.

For 2n=8 people the exact answer therefore is 15/2119 \approx 0.00707881 while the approximation is e/(2^4\, 4!) \approx 0.00707886: they agree to five significant figures.

To check, this `R`

simulation draws a million random permutations of eight objects, retains only those that are derangements, and counts those that are perfect pairings. It outputs its estimate, the standard error of the estimate, and a Z-score to compare it to the theoretical value. Its output is

```
p.hat se Z
0.006981031 0.000137385 -0.711721705
```

The small Z-score is consistent with the theoretical value. (These results would be consistent with any theoretical value between 0.0066 and 0.0073.)

```
paired <- function(x) crossprod(x[x] - 1:length(x))==0
good <- function(x) sum(x==1:length(x)) == 0
n <- 8
set.seed(17)
x <- replicate(1e6, sample(1:n, n))
i.good <- apply(x, 2, good)
i.paired <- apply(x, 2, paired)
n.deranged <- sum(i.good)
k.paired <- sum(i.good & i.paired)
p.hat <- k.paired / n.deranged
se <- sqrt(p.hat * (1-p.hat) / n.deranged)
(c(p.hat=p.hat, se=se, Z=(p.hat - 15/2119)/se))
```

**Attribution***Source : Link , Question Author : hermann , Answer Author : whuber*