I randomly decided to write a program that tests the Birthday Problem/Paradox.

The gist of the problem/paradox is that if you have a group of people in a room, how many people are required for there to be a 50% chance that 2 people will have the same birthday? Near 100% chance? It turns out only 23 people need to be in a room for there to be a 50% chance, and with 70 people, there's essentially a 100% chance (99.9%). This is assuming a even distribution of birthdays. This was hard for me to believe, so I decided to test it!

The main things I want critique of:

• generate-random-birthdays was originally much simpler. I had it defined simply as (repeat n-birthdays (new-random-birthday rand-gen)), but the laziness seemed to interfere with the random generator. It was returning a list of the correct size, but where every element was the same (which obviously destroyed the point of the test). A fold seemed like the obvious next tool, and it works, but it seems rather ugly for such a simple task. Is there an easier was to achieve a strict repeat?

• Is my use of a transient in generate-random-birthdays justified? I timed it with and without transients, and using them I got around a 10% speed increase. This seems somewhat significant. This is really the first time I've ever used them, so is there anything I need to watch out for that would make their use less favorable (besides a bit of code-bloat)?

• count-matching-birthdays has a very simple definition, but I suspect the (count (distinct birthdays)) part is expensive. I can't profile it though as VisualVM is broken on my Surface. The simplicity of it though has me suspicious (since when is the best way ever the simplest?). Is there a more efficient way of counting the number of duplicates in a list? Really, I only need to know whether or not there are any matches, not a specific number, but in the future I may way to modify this test to actually count the number of matches.

• (if (> (count-matching-birthdays rand-birthdays) 0) (inc n-dupes) n-dupes)) stinks to me. This is connected to the point above. Is there a more succinct way of saying this?

• This is more of a side question, but why when I tried using transients in test-n-people-range did it slow it down by a fair amount (almost as much as it sped up the other function)? Are there specific circumstances transients should be used for?

• Of course, anything else you notice.

The seemingly unrelated functions at the top were just inlined from my personal library. If they seem useless here it's because using new-rand-gen prevents me from needing to import the random generator, and random-int allows for deterministic behavior, which the standard functions are disappointingly lacking.

(ns bits.birthday-problem
(:import [java.util Random]))

(defn new-rand-gen
([seed] (Random. seed)))

(defn random-int
"Returns a random integer (actually a long) between min (inclusive) and max (exclusive)."
^long [^long min ^long max ^Random rand-gen]
(+ min (.nextInt rand-gen (- max min))))

(defn iterate-many [x repetitions iterate-f]
(nth
(iterate iterate-f x)
repetitions))

(defn new-random-birthday [rand-gen]
(random-int 0 365 rand-gen))

(defn generate-random-birthdays
"Generates a list of random birthdays."
[n-birthdays rand-gen]
(persistent!
(reduce (fn [acc _]
(conj! acc (new-random-birthday rand-gen)))
(transient [])
(range n-birthdays))))

(defn count-matching-birthdays
"Returns the number of matching birthdays in a list."
[birthdays]
(- (count birthdays) (count (distinct birthdays))))

(defn percentage-matching-birthdays
"Returns the percentage of tests that contained matching birthdays for the given number of people."
[n-tests n-people rand-gen]
(let [n-matches (iterate-many 0 n-tests
(fn [n-dupes]
(let [rand-birthdays (generate-random-birthdays n-people rand-gen)]
(if (> (count-matching-birthdays rand-birthdays) 0)
(inc n-dupes)
n-dupes))))]
(double (/ n-matches n-tests))))

(defn test-n-people-range
"Applies percentage-matching-birthdays to a range of different numbers of people.
Returns a 2D list in the form [[n-of-people perc-matching-birthdays]]"
[n-tests-per min-people max-people rand-gen]
(reduce
(fn [acc n-people]
(conj acc [n-people (percentage-matching-birthdays n-tests-per n-people rand-gen)]))
[]
(range min-people (inc max-people))))


Example:

(clojure.pprint/pprint
(test-n-people-range 100000 2 365 (new-rand-gen)))


Produces (after 29.9 minutes!):

[[2 0.00285]
[3 0.00812]
[4 0.01773]
[5 0.02732]
[6 0.04091]
[7 0.05785]
[8 0.07359]
[9 0.09435]
[10 0.1169]
[11 0.14066]
[12 0.16749]
[13 0.19378]
[14 0.22449]
[15 0.25166]
[16 0.28343]
[17 0.31648]
[18 0.34677]
[19 0.37685]
[20 0.41249]
[21 0.44335]
[22 0.47167]
[23 0.50717]
[24 0.53855]
[25 0.56909]
[26 0.6001]
[27 0.62695]
[28 0.65704]
[29 0.67904]
[30 0.70723]
[31 0.72846]
[32 0.75329]
[33 0.77335]
[34 0.79474]
[35 0.81418]
[36 0.83126]
[37 0.84894]
[38 0.86519]
[39 0.87792]
[40 0.89141]
[41 0.90191]
[42 0.91346]
[43 0.92318]
[44 0.93354]
[45 0.94125]
[46 0.94852]
[47 0.95511]
[48 0.95993]
[49 0.96638]
[50 0.96969]
[51 0.97444]
[52 0.97813]
[53 0.98111]
[54 0.98373]
[55 0.98646]
[56 0.98826]
[57 0.99004]
[58 0.99134]
[59 0.99313]
[60 0.99443]
[61 0.99487]
[62 0.99618]
[63 0.99674]
[64 0.99709]
[65 0.99772]
[66 0.99793]
[67 0.99835]
[68 0.9987]
[69 0.99905]
[70 0.99925]
[71 0.99931]
[72 0.99941]
[73 0.99957]
[74 0.99956]
[75 0.99982]
[76 0.99983]
[77 0.99982]
[78 0.99986]
[79 0.99993]
[80 0.99992]
[81 0.99996]
[82 0.99996]
[83 0.99996]
[84 0.99995]
[85 0.99998]
[86 0.99998]
[87 0.99998]
[88 0.99999]
[89 0.99999]
[90 1.0]
[91 0.99999]
[92 0.99999]
[93 0.99999]
[94 0.99998]
[95 1.0]
[96 1.0]
[97 1.0]
[98 1.0]
[99 1.0]
[100 1.0]
[101 1.0]
[102 1.0]
[103 1.0]
[104 1.0]
[105 1.0]
[106 1.0]
[107 1.0]
[108 1.0]
[109 1.0]
[110 1.0]
[111 1.0]
[112 1.0]
[113 1.0]
[114 1.0]
[115 1.0]
[116 1.0]
[117 1.0]
[118 1.0]
[119 1.0]
[120 1.0]
[121 1.0]
[122 1.0]
[123 1.0]
[124 1.0]
[125 1.0]
[126 1.0]
[127 1.0]
[128 1.0]
[129 1.0]
[130 1.0]
[131 1.0]
[132 1.0]
[133 1.0]
[134 1.0]
[135 1.0]
[136 1.0]
[137 1.0]
[138 1.0]
[139 1.0]
[140 1.0]
[141 1.0]
[142 1.0]
[143 1.0]
[144 1.0]
[145 1.0]
[146 1.0]
[147 1.0]
[148 1.0]
[149 1.0]
[150 1.0]
[151 1.0]
[152 1.0]
[153 1.0]
[154 1.0]
[155 1.0]
[156 1.0]
[157 1.0]
[158 1.0]
[159 1.0]
[160 1.0]
[161 1.0]
[162 1.0]
[163 1.0]
[164 1.0]
[165 1.0]
[166 1.0]
[167 1.0]
[168 1.0]
[169 1.0]
[170 1.0]
[171 1.0]
[172 1.0]
[173 1.0]
[174 1.0]
[175 1.0]
[176 1.0]
[177 1.0]
[178 1.0]
[179 1.0]
[180 1.0]
[181 1.0]
[182 1.0]
[183 1.0]
[184 1.0]
[185 1.0]
[186 1.0]
[187 1.0]
[188 1.0]
[189 1.0]
[190 1.0]
[191 1.0]
[192 1.0]
[193 1.0]
[194 1.0]
[195 1.0]
[196 1.0]
[197 1.0]
[198 1.0]
[199 1.0]
[200 1.0]
[201 1.0]
[202 1.0]
[203 1.0]
[204 1.0]
[205 1.0]
[206 1.0]
[207 1.0]
[208 1.0]
[209 1.0]
[210 1.0]
[211 1.0]
[212 1.0]
[213 1.0]
[214 1.0]
[215 1.0]
[216 1.0]
[217 1.0]
[218 1.0]
[219 1.0]
[220 1.0]
[221 1.0]
[222 1.0]
[223 1.0]
[224 1.0]
[225 1.0]
[226 1.0]
[227 1.0]
[228 1.0]
[229 1.0]
[230 1.0]
[231 1.0]
[232 1.0]
[233 1.0]
[234 1.0]
[235 1.0]
[236 1.0]
[237 1.0]
[238 1.0]
[239 1.0]
[240 1.0]
[241 1.0]
[242 1.0]
[243 1.0]
[244 1.0]
[245 1.0]
[246 1.0]
[247 1.0]
[248 1.0]
[249 1.0]
[250 1.0]
[251 1.0]
[252 1.0]
[253 1.0]
[254 1.0]
[255 1.0]
[256 1.0]
[257 1.0]
[258 1.0]
[259 1.0]
[260 1.0]
[261 1.0]
[262 1.0]
[263 1.0]
[264 1.0]
[265 1.0]
[266 1.0]
[267 1.0]
[268 1.0]
[269 1.0]
[270 1.0]
[271 1.0]
[272 1.0]
[273 1.0]
[274 1.0]
[275 1.0]
[276 1.0]
[277 1.0]
[278 1.0]
[279 1.0]
[280 1.0]
[281 1.0]
[282 1.0]
[283 1.0]
[284 1.0]
[285 1.0]
[286 1.0]
[287 1.0]
[288 1.0]
[289 1.0]
[290 1.0]
[291 1.0]
[292 1.0]
[293 1.0]
[294 1.0]
[295 1.0]
[296 1.0]
[297 1.0]
[298 1.0]
[299 1.0]
[300 1.0]
[301 1.0]
[302 1.0]
[303 1.0]
[304 1.0]
[305 1.0]
[306 1.0]
[307 1.0]
[308 1.0]
[309 1.0]
[310 1.0]
[311 1.0]
[312 1.0]
[313 1.0]
[314 1.0]
[315 1.0]
[316 1.0]
[317 1.0]
[318 1.0]
[319 1.0]
[320 1.0]
[321 1.0]
[322 1.0]
[323 1.0]
[324 1.0]
[325 1.0]
[326 1.0]
[327 1.0]
[328 1.0]
[329 1.0]
[330 1.0]
[331 1.0]
[332 1.0]
[333 1.0]
[334 1.0]
[335 1.0]
[336 1.0]
[337 1.0]
[338 1.0]
[339 1.0]
[340 1.0]
[341 1.0]
[342 1.0]
[343 1.0]
[344 1.0]
[345 1.0]
[346 1.0]
[347 1.0]
[348 1.0]
[349 1.0]
[350 1.0]
[351 1.0]
[352 1.0]
[353 1.0]
[354 1.0]
[355 1.0]
[356 1.0]
[357 1.0]
[358 1.0]
[359 1.0]
[360 1.0]
[361 1.0]
[362 1.0]
[363 1.0]
[364 1.0]
[365 1.0]]


Which matches the expected output almost exactly (although it reached 100% with less people than expected in this test).

I've since changed count-matching birthdays to:

(defn count-matching-birthdays
"Returns the number of matching birthdays in a list."
[birthdays]
(- (count birthdays) (count (into #{} birthdays))))


It's 70% faster! I guess just putting the birthdays into a set actually makes the most sense, in retrospect. Using transients in test-n-people-range now seems to in fact cause a speed increase. Not sure why, unless my timing test yesterday was botched.

Replacing the use of iterate-many with a fold gave a non-negligible speed increase. I think I'll redefine iterate-many in terms of a fold instead of interate/nth.

You can simplify most of your functions as follows:

(defn random-int [min max]
(+ min (rand-int (- max min))))

(defn new-random-birthday []
(random-int 0 365))

(defn generate-random-birthdays [n]
(repeatedly n new-random-birthday))

(defn percentage-matching-birthdays [n-tests n-people]
(let [n-matches (->> #(generate-random-birthdays n-people)
(repeatedly n-tests)
(remove (partial apply distinct?))
count)]
(double (/ n-matches n-tests))))

(defn test-n-people-range [n-tests-per min-people max-people]
(mapv
(fn [n-people]
[n-people (percentage-matching-birthdays n-tests-per n-people)])
(range min-people (inc max-people))))


I've tested this with

(percentage-matching-birthdays 1000 23)
;0.505


..., much as expected.

The above generates the random integers in native Clojure, hence can dispense with explicit reference to a Java randomiser. I've also got rid of the optimisations, which I think are premature. But the main change is to use apposite core functions wherever possible.

I think your iterate-many ...

(defn iterate-many [x n f]
(nth (iterate f) n))


... is wrong. I think it ought to be ...

(defn iterate-many [x n f]
(take n (iterate f x)))


But I haven't used or tested it.

• I want my random functions to take a random generator as an argument. I like having deterministic behavior. My iterate-many also behaves as I'd like. I only want it to return a final value; not all the values leading up to it. I don't want it to act as reductions. – Carcigenicate Nov 9 '16 at 0:43
• I'll upvote this since this is faster, but I'm not sure I really want to accept it since it disregards certain aspects of the question, and appears to give skewed results for some reason. – Carcigenicate Nov 9 '16 at 0:47
• Also, you have take taking 3 parameters. I think you intended x to be an argument to iterate. – Carcigenicate Nov 9 '16 at 0:50
• @Carcigenicate I've tried a few of my results against the theoretical probabilities, and they look at first sight at least as close as yours. – Thumbnail Nov 9 '16 at 2:03
• OK, I'll check again. – Carcigenicate Nov 9 '16 at 2:03

On second thoughts ...

The problem splits into two pieces:

• generating a seeded random sequence of birthdays and
• chopping up a sequence and classifying the results.

It follows that

• The sequence of birthdays doesn't have to know what is to be done to it.
• The slicing and dicing can apply to any endless sequence of whatever-you-like.

This way, we avoid having to pass mutable seeded random number generators through the layers of function calls. We just pass the immutable sequences they generate (off which we take as little or as much as we like). We can test our classification reporting stuff on any sequences, as deterministically as we like.

To generate a seeded sequence of birthdays:

(defn birthday-seq [seed]
(let [gen (Random. seed)]
(repeatedly #(.nextInt gen 365))))


For example,

(take 10 (seeded-random-birthday-seq 1))
;(75 313 142 343 309 189 144 166 58 228)


This is deterministic: do it again, you get the same result.

Now for slicing a sequence and classifying the slices:

(defn slice-and-classify [slice-count slice-size classifier coll]
(->> coll
(partition slice-size)
(take slice-count)
(map classifier)
frequencies))


For example,

(slice-and-classify 20 7 (comp odd? #(apply + %)) (range))
;{true 10, false 10}


What we want is proportions:

(defn normalise [m]
(let [vs (vals m)
sum (apply + vs)]
(zipmap (keys m) (map #(double (/ % sum)) vs))))


For example,

(normalise {:one 1, :two 2, :three 3, :four 4})
;{:four 0.4, :three 0.3, :two 0.2, :one 0.1}


We could roll up normalise into slice-and-classify, but they seem clearer apart.

We've now got some useful functions that we can plug together to do the necessary. Let's just see if sample size 23 does tend to produce about half the samples with duplicates:

(slice-and-classify 100000 23 (partial apply distinct?) (birthday-seq 42))
;{false 50592, true 49408}


Yes, it does.

The above was put together with no thought for performance. As you've noticed, the hot spot is the classifier function argument to slice-and-classify, for which (partial apply distinct?) is pretty slow.

If you look at the source code for distinct?, you'll find it doesn't use transients, presumably because contains? isn't available. But count is available, and still works provided you remember that the transient is mutable:

(defn unique? [[x & xs]]
(loop [[y & ys :as coll] xs, seen (transient #{x})]
(or (not coll)
(let [n (count seen) ; remember this before we mutate seen
seen (conj! seen y)]
(and (> (count seen) n)
(recur ys seen))))))


This is modified to take a single collection argument, thereby dispensing with the need to apply it.

(slice-and-classify 10000 23 unique? (birthday-seq 42))
;{false 5092, true 4908}


It works, but - though faster than (partial apply distinct?) - it's not much faster than your #(= (count (into #{} %)) (count %)) or the #(= (count (set %)) (count %)) variant.

To seriously improve performance, I'd try the following:

• Generate the birthdays as an endless sequence of Java vectors-of-long, cut to size.
• To test for the distinctness of its elements, loop through the vector, accumulating a Java BitSet.

With appropriate type hints/casts, this should avoid a lot of boxing and unboxing.

• I tried to improve the performance, as suggested. Doing so made little perceptible difference :(. – Thumbnail Nov 14 '16 at 14:30