import numpy as np
from itertools import chain
%matplotlib inline
import matplotlib.pyplot as plt
res = []
for _ in xrange(50000):
kids = np.random.randint(1,10,100)
boys = np.array([np.random.randint(0, 2, k).sum() for k in kids])
girls = kids-boys
boy_brothers = [[b-1]*b for b in boys]
girl_brothers = [[b]*g for b, g in zip(boys, girls)]
res.append(np.array(list(chain(*girl_brothers))).mean() - np.array(list(chain(*boy_brothers))).mean())
print np.array(res).mean()
0.0117380724123
plt.hist(res, bins=100)
plt.show()