home *** CD-ROM | disk | FTP | other *** search
/ ftp.pasteur.org/FAQ/ / ftp-pasteur-org-FAQ.zip / FAQ / puzzles / archive / probability < prev    next >
Encoding:
Internet Message Format  |  1996-04-26  |  49.0 KB

  1. Received: from MIT.EDU (SOUTH-STATION-ANNEX.MIT.EDU [18.72.1.2]) by bloom-picayune.MIT.EDU (8.6.13/2.3JIK) with SMTP id OAA03476; Sat, 20 Apr 1996 14:53:15 -0400
  2. Received: from [199.164.164.1] by MIT.EDU with SMTP
  3.     id AA07977; Sat, 20 Apr 96 14:13:31 EDT
  4. Received: by questrel.questrel.com (940816.SGI.8.6.9/940406.SGI)
  5.     for news-answers-request@mit.edu id LAA25285; Sat, 20 Apr 1996 11:14:26 -0700
  6. Newsgroups: rec.puzzles,news.answers,rec.answers
  7. Path: senator-bedfellow.mit.edu!bloom-beacon.mit.edu!spool.mu.edu!howland.reston.ans.net!europa.eng.gtefsd.com!uunet!questrel!chris
  8. From: chris@questrel.questrel.com (Chris Cole)
  9. Subject: rec.puzzles Archive (probability), part 31 of 35
  10. Message-Id: <puzzles/archive/probability_745653851@questrel.com>
  11. Followup-To: rec.puzzles
  12. Summary: This is part of an archive of questions
  13.  and answers that may be of interest to
  14.  puzzle enthusiasts.
  15.  Part 1 contains the index to the archive.
  16.  Read the rec.puzzles FAQ for more information.
  17. Sender: chris@questrel.questrel.com (Chris Cole)
  18. Reply-To: archive-comment@questrel.questrel.com
  19. Organization: Questrel, Inc.
  20. References: <puzzles/archive/Instructions_745653851@questrel.com>
  21. Date: Wed, 18 Aug 1993 06:06:50 GMT
  22. Approved: news-answers-request@MIT.Edu
  23. Expires: Thu, 1 Sep 1994 06:04:11 GMT
  24. Lines: 1312
  25. Xref: senator-bedfellow.mit.edu rec.puzzles:25020 news.answers:11540 rec.answers:1940
  26. Apparently-To: news-answers-request@mit.edu
  27.  
  28. Archive-name: puzzles/archive/probability
  29. Last-modified: 17 Aug 1993
  30. Version: 4
  31.  
  32.  
  33. ==> probability/amoeba.p <==
  34. A jar begins with one amoeba.  Every minute, every amoeba
  35. turns into 0, 1, 2, or 3 amoebae with probability 25%
  36. for each case ( dies, does nothing, splits into 2, or splits 
  37. into 3).  What is the probability that the amoeba population
  38. eventually dies out?
  39.  
  40. ==> probability/amoeba.s <==
  41. If p is the probability that a single amoeba's descendants will die
  42. out eventually, the probability that N amoebas' descendents will all
  43. die out eventually must be p^N, since each amoeba is independent of
  44. every other amoeba.  Also, the probability that a single amoeba's
  45. descendants will die out must be independent of time when averaged
  46. over all the possibilities.  At t=0, the probability is p, at t=1 the
  47. probability is 0.25(p^0+p^1+p^2+p^3), and these probabilities must be
  48. equal.  Extinction probability p is a root of f(p)=p.  In this case,
  49.  p = sqrt(2)-1.
  50.  
  51. The generating function for the sequence P(n,i), which gives the
  52. probability of i amoebas after n minutes, is f^n(x), where f^n(x) ==
  53. f^(n-1) ( f(x) ), f^0(x) == x .  That is, f^n is the nth composition
  54. of f with itself.
  55.  
  56. Then f^n(0) gives the probability of 0 amoebas after n minutes, since
  57. f^n(0) = P(n,0). We then note that:
  58.  
  59.     f^(n+1)(x) = ( 1 + f^n(x) + (f^n(x))^2 + (f^n(x))^3 )/4
  60.  
  61. so that if f^(n+1)(0) -> f^n(0) we can solve the equation.
  62.  
  63. The generating function also gives an expression for the expectation
  64. value of the number of amoebas after n minutes. This is d/dx(f^n(x))
  65. evaluated at x=1. Using the chain rule we get f'(f^(n-1)(x))*d/dx(f^(n-1)(x))
  66. and since f'(1) = 1.5  and f(1) = 1, we see that the result is just
  67. 1.5^n, as might be expected.
  68.  
  69. ==> probability/apriori.p <==
  70. An urn contains one hundred white and black balls.  You sample one hundred
  71. balls with replacement and they are all white.  What is the probability
  72. that all the balls are white?
  73.  
  74. ==> probability/apriori.s <==
  75. This question cannot be answered with the information given.
  76.  
  77. In general, the following formula gives the conditional probability
  78. that all the balls are white given you have sampled one hundred balls
  79. and they are all white:
  80.  
  81.     P(100 white | 100 white samples) =
  82.  
  83.             P(100 white samples | 100 white) * P(100 white) 
  84.         -----------------------------------------------------------
  85.         sum(i=0 to 100) P(100 white samples | i white) * P(i white)
  86.  
  87. The probabilities P(i white) are needed to compute this formula.  This
  88. does not seem helpful, since one of these (P(100 white)) is just what we
  89. are trying to compute.  However, the following argument can be made:
  90. Before the experiment, all possible numbers of white balls from zero to
  91. one hundred are equally likely, so P(i white) = 1/101.  Therefore, the
  92. odds that all 100 balls are white given 100 white samples is:
  93.  
  94.     P(100 white | 100 white samples) =
  95.  
  96.         1 / ( sum(i=0 to 100) (i/100)^100 ) =
  97.  
  98.         63.6%
  99.  
  100. This argument is fallacious, however, since we cannot assume that the urn
  101. was prepared so that all possible numbers of white balls from zero to one
  102. hundred are equally likely.  In general, we need to know the P(i white)
  103. in order to calculate the P(100 white | 100 white samples).  Without this
  104. information, we cannot determine the answer.
  105.  
  106. This leads to a general "problem": our judgments about the relative
  107. likelihood of things is based on past experience.  Each experience allows
  108. us to adjust our likelihood judgment, based on prior probabilities.  This
  109. is called Bayesian inference.  However, if the prior probabilities are not
  110. known, then neither are the derived probabilities.  But how are the prior
  111. probabilities determined?  For example, if we are brains in the vat of a
  112. diabolical scientist, all of our prior experiences are illusions, and
  113. therefore all of our prior probabilities are wrong.
  114.  
  115. All of our probability judgments indeed depend upon the assumption that
  116. we are not brains in a vat.  If this assumption is wrong, all bets are
  117. off.
  118.  
  119. ==> probability/bayes.p <==
  120. One urn contains black marbles, and the other contains white or black
  121. marbles with even odds.  You pick a marble from an urn; it is black;
  122. you put it back; what are the odds that you will draw a black marble on
  123. the next draw?  What are the odds after n black draws?
  124.  
  125. ==> probability/bayes.s <==
  126. Every time you draw a black marble, you throw out (from your
  127. probability space) half of those possible urns that contain both
  128. colors.  So you have 1/2^n times as many ways to have a white marble in
  129. the urn after n draws, all black, as at the start.  But you have
  130. exactly the same number of ways to have both marbles black.  The
  131. numbers (mixed cases vs. all-black cases) go as 1:1, 1:2, 1:4, 1:8,...
  132. and the chance of having a white marble in the urn goes as 1/2, 1/3,
  133. 1/5, 1/9, ..., 1/(1+2^(n-1)), hence the odds of drawing a white marble
  134. on the nth try after n-1 consecutive drawings of black are
  135.  
  136. 1/4    the first time
  137. 1/6    the second time
  138. 1/10   the third time
  139. ...
  140. 1/(2+2^n)  the nth time
  141.  
  142. ==> probability/birthday/line.p <==
  143. At a movie theater, the manager announces that they will give a free ticket
  144. to the first person in line whose birthday is the same as someone who has
  145. already bought a ticket.  You have the option of getting in line at any
  146. time.  Assuming that you don't know anyone else's birthday, that birthdays
  147. are distributed randomly throughtout the year, etc., what position in line
  148. gives you the greatest chance of being the first duplicate birthday?
  149.  
  150. ==> probability/birthday/line.s <==
  151. Suppose you are the Kth person in line.  Then you win if and only if the
  152. K-1 people ahead all have distinct birtdays AND your birthday matches
  153. one of theirs.  Let 
  154.  
  155. A = event that your birthday matches one of the K-1 people ahead
  156. B = event that those K-1 people all have different birthdays
  157.  
  158. Then 
  159.  
  160. Prob(you win) = Prob(B) * Prob(A | B)
  161.  
  162. (Prob(A | B) is the conditional probability of A given that B occurred.)
  163.  
  164. Now let P(K) be the probability that the K-th person in line wins,
  165. Q(K) the probability that the first K people all have distinct
  166. birthdays (which occurs exactly when none of them wins).  Then
  167.  
  168. P(1) + P(2) + ... + P(K-1) + P(K) = 1 - Q(K)
  169. P(1) + P(2) + ... + P(K-1)        = 1 - Q(K-1)
  170.  
  171. P(K) = Q(K-1) - Q(K)   <--- this is what we want to maximize.
  172.  
  173. Now if the first K-1 all have distinct birthdays, then assuming
  174. uniform distribution of birthdays among D days of the year,
  175. the K-th person has K-1 chances out of D to match, and D-K+1 chances
  176. not to match (which would produce K distinct birthdays).  So 
  177.  
  178. Q(K) = Q(K-1)*(D-K+1)/D = Q(K-1) - Q(K-1)*(K-1)/D
  179. Q(K-1) - Q(K) = Q(K-1)*(K-1)/D = Q(K)*(K-1)/(D-K+1)
  180.  
  181. Now we want to maximize P(K), which means we need the greatest K such
  182. that P(K) - P(K-1) > 0.  (Actually, as just given, this only
  183. guarantees a local maximum, but in fact if we investigate a bit
  184. farther we'll find that P(K) has only one maximum.)
  185. For convenience in calculation let's set K = I + 1.  Then
  186.  
  187. Q(I-1) - Q(I) = Q(I)*(I-1)/(D-I+1)
  188. Q(I) - Q(I+1) = Q(I)*I/D
  189.  
  190. P(K) - P(K-1) = P(I+1) - P(I)
  191.               = (Q(I) - Q(I+1)) - (Q(K-2) - Q(K-1))
  192.               = Q(I)*(I/D - (I-1)/(D-I+1))
  193.  
  194. To find out where this is last positive (and next goes negative), solve
  195.  
  196. x/D - (x-1)/(D-x+1) = 0
  197.  
  198. Multiply by D*(D+1-x) both sides:
  199.  
  200. (D+1-x)*x - D*(x-1) = 0
  201. Dx + x - x^2 - Dx + D = 0
  202. x^2 - x - D = 0
  203.  
  204. x = (1 +/- sqrt(1 - 4*(-D)))/2    ... take the positive square root
  205.   = 0.5 + sqrt(D + 0.25)
  206.  
  207. Setting D=365 (finally deciding how many days in a year!),
  208.  
  209. desired I = x = 0.5 + sqrt(365.25) = 19.612 (approx).
  210.  
  211. The last integer I for which the new probability is greater then the old
  212. is therefore I=19, and so K = I+1 = 20.  You should try to be the 20th
  213. person in line.
  214.  
  215. Computing your chances of actually winning is slightly harder, unless
  216. you do it numerically by computer.  The recursions you need have already
  217. been given.
  218.  
  219. -- David Karr (karr@cs.cornell.edu)
  220.  
  221.  
  222.  
  223.  
  224. ==> probability/birthday/same.day.p <==
  225. How many people must be at a party before you have even odds or better
  226. of two having the same bithday (not necessarily the same year, of course)?
  227.  
  228. ==> probability/birthday/same.day.s <==
  229. 23.
  230.  
  231. See also:
  232.     archive entry "coupon"
  233.  
  234. ==> probability/cab.p <==
  235. A cab was involved in a hit and run accident at night.  Two cab companies,
  236. the Green and the Blue, operate in the city.  Here is some data:
  237.  
  238.     a)  Although the two companies are equal in size, 85% of cab 
  239. accidents in the city involve Green cabs and 15% involve Blue cabs.
  240.  
  241.     b)  A witness identified the cab in this particular accident as Blue.
  242. The court tested the reliability of the witness under the same circumstances
  243. that existed on the night of the accident and concluded that the witness 
  244. correctly identified each one of the two colors 80% of the time and failed
  245. 20% of the time.
  246.  
  247. What is the probability that the cab involved in the accident was 
  248. Blue rather than Green?
  249.  
  250. If it looks like an obvious problem in statistics, then consider the
  251. following argument:
  252.  
  253. The probability that the color of the cab was Blue is 80%!  After all,
  254. the witness is correct 80% of the time, and this time he said it was Blue!
  255.  
  256. What else need be considered?  Nothing, right?
  257.  
  258. If we look at Bayes theorem (pretty basic statistical theorem) we 
  259. should get a much lower probability.  But why should we consider statistical
  260. theorems when the problem appears so clear cut?  Should we just accept the 
  261. 80% figure as correct?
  262.  
  263. ==> probability/cab.s <==
  264. The police tests don't apply directly, because according to the
  265. wording, the witness, given any mix of cabs, would get the right
  266. answer 80% of the time.  Thus given a mix of 85% green and 15% blue
  267. cabs, he will say 20% of the green cabs and 80% of the blue cabs are
  268. blue.  That's 20% of 85% plus 80% of 15%, or 17%+12% = 29% of all the
  269. cabs that the witness will say are blue.  Of those, only 12/29 are
  270. actually blue.  Thus P(cab is blue | witness claims blue) = 12/29.
  271. That's just a little over 40%.
  272.  
  273. Think of it this way... suppose you had a robot watching parts on a
  274. conveyor belt to spot defective parts, and suppose the robot made a
  275. correct determination only 50% of the time (I know, you should
  276. probably get rid of the robot...).  If one out of a billion parts are
  277. defective, then to a very good approximation you'd expect half your
  278. parts to be rejected by the robot.  That's 500 million per billion.
  279. But you wouldn't expect more than one of those to be genuinely
  280. defective.  So given the mix of parts, a lot more than 50% of the
  281. REJECTED parts will be rejected by mistake (even though 50% of ALL the
  282. parts are correctly identified, and in particular, 50% of the
  283. defective parts are rejected).
  284.  
  285. When the biases get so enormous, things starts getting quite a bit
  286. more in line with intuition.
  287.  
  288. For a related real-life example of probability in the courtroom see
  289. People v. Collins, 68 Cal 2d319 (1968).
  290.  
  291. ==> probability/coupon.p <==
  292. There is a free gift in my breakfast cereal. The manufacturers say that
  293. the gift comes in four different colors, and encourage one to collect
  294. all four (& so eat lots of their cereal). Assuming there is an equal
  295. chance of getting any one of the colors,  what is the expected number
  296. of boxes I must consume to get all four?  Can you generalise to n
  297. colors and/or unequal probabilities?
  298.  
  299. ==> probability/coupon.s <==
  300. The problem is well known under the name of "COUPON COLLECTOR PROBLEM".
  301. The answer for n equally likely coupons is
  302. (1)        C(n)=n*H(n)    with H(n)=1+1/2+1/3+...+1/n.
  303. In the unequal probabilities case, with p_i the probability of coupon i,
  304. it becomes
  305. (2)        C(n)=int_0^infty [1-prod_{i=1}^n (1-exp(-p_i*t))] dt
  306. which reduces to (1) when p_i=1/n (An easy exercise).
  307. In the equal probabilities case  C(n)~n log(n). For a Zipf law,
  308. from (2), we have C(n)~n log^2(n).
  309.  
  310. A related problem is the "BIRTHDAY PARADOX" familiar to people
  311. interested in hashing algorithms: With a party of 23 persons,
  312. you are likely (i.e. with probability >50%) to find two with
  313. the same birthday. The non equiprobable case was solved by:
  314.     M. Klamkin and D. Newman, Extensions of the birthday
  315.     surprise, J. Comb. Th. 3 (1967), 279-282.
  316.  
  317. ==> probability/darts.p <==
  318. Peter throws two darts at a dartboard, aiming for the center.  The
  319. second dart lands farther from the center than the first.  If Peter now
  320. throws another dart at the board, aiming for the center, what is the
  321. probability that this third throw is also worse (i.e., farther from 
  322. the center) than his first?  Assume Peter's skilfulness is constant.
  323.  
  324. ==> probability/darts.s <==
  325. Since the three darts are thrown independently,
  326. they each have a 1/3 chance of being the best throw.  As long as the
  327. third dart is not the best throw, it will be worse than the first dart.
  328. Therefore the answer is 2/3.
  329.  
  330. Ranking the three darts' results from A (best) to C
  331. (worst), there are, a priori, six equiprobable outcomes.
  332.  
  333. possibility #    1    2    3    4    5    6
  334. 1st throw    A    A    B    B    C    C
  335. 2nd throw    B    C    A    C    A    B
  336. 3rd throw    C    B    C    A    B    A
  337.  
  338. The information from the first two throws shows us that the first
  339. throw will not be the worst, nor the second throw the best.  Thus
  340. possibilities 3, 5 and 6 are eliminated, leaving three equiprobable
  341. cases, 1, 2 and 4.  Of these, 1 and 2 have the third throw worse than
  342. the first; 4 does not.  Again the answer is 2/3.
  343.  
  344. ==> probability/derangement.p <==
  345. 12 men leave their hats with the hat check.  If the hats are randomly
  346. returned, what is the probability that nobody gets the correct hat?
  347.  
  348. ==> probability/derangement.s <==
  349. Suppose we are handing out hats to n people.  First we start with all
  350. the possible outcomes.  Then we subtract off those that assign the
  351. right hat to a given person, for each of the n people.  But this
  352. double-counts each outcome that assigned 2 hats correctly, so we have
  353. to add those outcomes back in.  But now we've counted one net copy of
  354. each outcome with 3 correct hats in our set, so we have to subtract
  355. those off again.  But now we've taken away each 4-correct-hat outcome
  356. once too often, and have to put it back in, and so forth ... until we
  357. add or subtract the outcome that involves all n people getting the
  358. correct hats.
  359.  
  360. Putting it all in probabilities, the measure of the original set is 1.
  361. For a given subset of k people, the probability that they all get
  362. their correct hats is (n-k)!/n!, while there are (n choose k) such
  363. subsets of k people altogether.  But then
  364.  
  365.    (n choose k)*(n-k)!/n! = (n!/((n-k)!*k!))*(n-k)!/n! = 1/k!
  366.  
  367. is the total probability measure we get by counting each subset of k
  368. people once each.  So we end up generating the finite series
  369.  
  370.    1 - 1/1! + 1/2! - 1/3! +- ... +/- 1/n!
  371.  
  372. which is of course just the first n+1 terms of the Taylor series
  373. expansion for f(x) = e^x centered at 0 and evaluated at -1, which
  374. converges to 1/e quite fast.  You can compute the exact probability for
  375. any n >= 1 simply by rounding n!/e to the nearest whole number, then
  376. dividing again by n!.  Moreover I think you will find you are always
  377. rounding down for odd n and rounding up for even n.
  378.  
  379. For example,
  380.  
  381.    12! / e = 176214840.95798...
  382.  
  383. which is within 0.05 (absolute error, not relative) of the correct
  384. intermediate result, 176214841.
  385.  
  386. Another fact is that if you set the probability of no matching hats
  387. equal to m/n!, then m is an integer divisible by (n-1).  That's
  388. because the number of ways you can give hat #2 to person #1 is exactly
  389. the same as the number of ways you can give that person hat #3,
  390. likewise for hat #4, and so forth.
  391.  
  392. -- David Karr (karr@cs.cornell.edu)
  393.  
  394. ==> probability/family.p <==
  395. Suppose that it is equally likely for a pregnancy to deliver
  396. a baby boy as it is to deliver a baby girl.  Suppose that for a
  397. large society of people, every family continues to have children
  398. until they have a boy, then they stop having children.
  399. After 1,000 generations of families, what is the ratio of males
  400. to females?
  401.  
  402. ==> probability/family.s <==
  403. The ratio will be 50-50 in both cases.  We are not killing off any
  404. fetuses or babies, and half of all conceptions will be male, half
  405. female.  When a family decides to stop does not affect this fact.
  406.  
  407. ==> probability/flips/once.in.run.p <==
  408. What are the odds that a run of one H or T (i.e., THT or HTH) will occur
  409. in n flips of a fair coin?
  410.  
  411. ==> probability/flips/once.in.run.s <==
  412. References:
  413.     John P. Robinson, Transition Count and Syndrome are Uncorrelated, IEEE
  414.     Transactions on Information Theory, Jan 1988.
  415.  
  416. First we define a function or enumerator P(n,k) as the number of length
  417. "n" sequences that generate "k" successes.  For example,
  418.  
  419.      P(4,1)= 4  (HHTH, HTHH, TTHT, and THTT are 4 possible length 4 sequences).
  420.  
  421. I derived two generating functions g(x) and h(x) in order to enumerate
  422. P(n,k), they are compactly represented by the following matrix
  423. polynomial.
  424.  
  425.  
  426.             _    _      _     _           _   _
  427.        | g(x) |    | 1   1 | (n-3)   |  4  |
  428.        |      | =  |       |         |     | 
  429.        | h(x) |    | 1   x |         |2+2x |    
  430.            |_    _|    |_     _|         |_   _|
  431.  
  432. The above is expressed as matrix generating function.  It can be shown
  433. that P(n,k) is the coefficient of the x^k in the polynomial
  434. (g(x)+h(x)).
  435.  
  436. For example, if n=4 we get (g(x)+h(x)) from the matrix generating
  437. function as (10+4x+2x^2).  Clearly, P(4,1) (coefficient of x) is 4 and
  438. P(4,2)=2 ( There are two such sequences THTH, and HTHT).
  439.  
  440. We can show that
  441.  
  442.    mean(k) = (n-2)/4 and sd= square_root(5n-12)/4
  443.  
  444. We need to generate "n" samples. This can be done by using sequences of length
  445. (n+2).  Then our new statistics would be
  446.  
  447.    mean = n/4
  448.  
  449.    sd = square_root(5n-2)/4
  450.  
  451. Similar approach can be followed for higher dimensional cases.
  452.  
  453. ==> probability/flips/twice.in.run.p <==
  454. What is the probability in n flips of a fair coin that there will be two
  455. heads in a row?
  456.  
  457. ==> probability/flips/twice.in.run.s <==
  458. Well, the question is then how many strings of n h's and t's contain
  459. hh?  I would guess right off hand that its going to be easier to
  460. calculate the number of strings that _don't_ contain hh and then
  461. subtract that from the total number of strings.
  462.  
  463. So we want to count the strings of n h's and t's with no hh in them.
  464. How many h's and t's can there be? It is fairly clear that there must
  465. be from 0 to n/2 h's, inclusive. (If there were (n/2+1) then there
  466. would have to be two touching.)
  467.  
  468. How many strings are there with 0 h's? 1
  469.  
  470. How many strings are there with 1 h? Well, there are (n-1) t's, so
  471. there are a total of n places to put the one h. So the are nC1 such
  472. strings.  How many strings are there with 2 h's? Well, there are (n-1)
  473. places to put the two h's, so there are (n-1)C2 such strings.
  474.  
  475. Finally, with n/2 h's there are (n/2+1) places to put them, so there
  476. are (n/2+1)C(n/2) such strings.
  477.  
  478. Therefore the total number of strings is
  479. Sum (from i=0 to n/2) of (n-i+1)C(i)
  480.  
  481. Now, here's where it get's interesting. If we play around with Pascal's
  482. triangle for a while, we see that this sum equals none other than the
  483. (n+2)th Fibonacci number.
  484.  
  485. So the probability that n coin tosses will give a hh is:
  486.  
  487. 2^n-f(n+2) 
  488. ----------
  489. 2^n
  490.  
  491. (where f(x) is the xth Fibanocci number (so that f(1) is and f(2) are both 1))
  492.  
  493. ==> probability/flips/unfair.p <==
  494. Generate even odds from an unfair coin.  For example, if you
  495. thought a coin was biased toward heads, how could you get the
  496. equivalent of a fair coin with several tosses of the unfair coin?
  497.  
  498. ==> probability/flips/unfair.s <==
  499. Toss twice.  If both tosses give the same result, repeat this process
  500. (throw out the two tosses and start again).  Otherwise, take the first
  501. of the two results.
  502.  
  503. ==> probability/flips/waiting.time.p <==
  504. Compute the expected waiting time for a sequence of coin flips, or the
  505. probabilty that one sequence of coin flips will occur before another.
  506.  
  507.  
  508. ==> probability/flips/waiting.time.s <==
  509. Here's a C program I had lying around that is relevant to the
  510. current discussion of coin-flipping.  The algorithm is N^3 (for N flips)
  511. but it could certainly be improved.  Compile with 
  512.  
  513.     cc -o flip flip.c -lm
  514.  
  515. -- Guy
  516.  
  517. _________________ Cut here ___________________
  518.  
  519. #include <stdio.h>
  520. #include <math.h>
  521.  
  522. char *progname;              /* Program name */
  523.  
  524. #define NOT(c) (('H' + 'T') - (c))
  525.  
  526.  
  527. /* flip.c -- a program to compute the expected waiting time for a sequence
  528.              of coin flips, or the probabilty that one sequence
  529.              of coin flips will occur before another.
  530.  
  531.     Guy Jacobson, 11/1/90
  532. */
  533.  
  534. main (ac, av) int ac; char **av;
  535. {
  536.     char *f1, *f2, *parseflips ();
  537.     double compute ();
  538.  
  539.     progname = av[0];
  540.  
  541.     if (ac == 2) {
  542.     f1 = parseflips (av[1]);
  543.     printf ("Expected number of flips until %s = %.1f\n",
  544.         f1, compute (f1, NULL));
  545.     }
  546.     else if (ac == 3) {
  547.  
  548.     f1 = parseflips (av[1]);
  549.     f2 = parseflips (av[2]);
  550.  
  551.     if (strcmp (f1, f2) == 0) {
  552.         printf ("Can't use the same flip sequence.\n");
  553.         exit (1);
  554.     }
  555.     printf ("Probability of flipping %s before %s = %.1f%%\n",
  556.         av[1], av[2], compute (f1, f2) * 100.0);
  557.     }
  558.     else
  559.       usage ();
  560. }
  561.  
  562. char *parseflips (s) char *s;
  563. {
  564.     char *f = s;
  565.  
  566.     while (*s)
  567.       if (*s == 'H' || *s == 'h')
  568.     *s++ = 'H';
  569.       else if (*s == 'T' || *s == 't')
  570.     *s++ = 'T';
  571.       else
  572.     usage ();
  573.  
  574.     return f;
  575. }
  576.  
  577. usage ()
  578. {
  579.     printf ("usage: %s {HT}^n\n", progname);
  580.     printf ("\tto get the expected waiting time, or\n");
  581.     printf ("usage: %s s1 s2\n\t(where s1, s2 in {HT}^n for some fixed n)\n",
  582.         progname);
  583.     printf ("\tto get the probability that s1 will occur before s2\n");
  584.     exit (1);
  585. }
  586.  
  587. /*
  588.   compute --  if f2 is non-null, compute the probability that flip
  589.               sequence f1 will occur before f2.  With null f2, compute
  590.           the expected waiting time until f1 is flipped
  591.  
  592.   technique:
  593.  
  594.     Build a DFA to recognize (H+T)*f1 [or (H+T)*(f1+f2) when f2
  595.     is non-null].  Randomly flipping coins is a Markov process on the
  596.     graph of this DFA.  We can solve for the probability that f1 precedes
  597.     f2 or the expected waiting time for f1 by setting up a linear system
  598.     of equations relating the values of these unknowns starting from each
  599.     state of the DFA.  Solve this linear system by Gaussian Elimination.
  600. */
  601.  
  602. typedef struct state {
  603.     char *s;                /* pointer to substring string matched */
  604.     int len;                /* length of substring matched */
  605.     int backup;             /* number of one of the two next states */
  606. } state;
  607.  
  608. double compute (f1, f2) char *f1, *f2;
  609. {
  610.     double solvex0 ();
  611.     int i, j, n1, n;
  612.  
  613.     state *dfa;
  614.     int nstates;
  615.  
  616.     char *malloc ();
  617.  
  618.     n = n1 = strlen (f1);
  619.     if (f2)
  620.     n += strlen (f2); /* n + 1 states in the DFA */
  621.  
  622.     dfa = (state *) malloc ((unsigned) ((n + 1) * sizeof (state)));
  623.  
  624.     if (!dfa) {
  625.     printf ("Ouch, out of memory!\n");
  626.     exit (1);
  627.     }
  628.  
  629.     /* set up the backbone of the DFA */
  630.  
  631.     for (i = 0; i <= n; i++) {
  632.     dfa[i].s = (i <= n1) ? f1 : f2;
  633.     dfa[i].len = (i <= n1) ? i : i - n1;
  634.     }
  635.  
  636.     /* for i not a final state, one next state of i is simply
  637.        i+1 (this corresponds to another matching character of dfs[i].s
  638.        The other next state (the backup state) is now computed.
  639.        It is the state whose substring matches the longest suffix
  640.        with the last character changed */      
  641.        
  642.     for (i = 0; i <= n; i++) {
  643.     dfa[i].backup = 0;
  644.     for (j = 1; j <= n; j++)
  645.       if ((dfa[j].len > dfa[dfa[i].backup].len)
  646.           && dfa[i].s[dfa[i].len] == NOT (dfa[j].s[dfa[j].len - 1])
  647.           && strncmp (dfa[j].s, dfa[i].s + dfa[i].len - dfa[j].len + 1,
  648.               dfa[j].len - 1) == 0)
  649.         dfa[i].backup = j;
  650.     }
  651.  
  652.     /* our dfa has n + 1 states, so build a system n + 1 equations
  653.        in n + 1 unknowns */
  654.  
  655.     eqsystem (n + 1);
  656.  
  657.     for (i = 0; i < n; i++)
  658.       if (i == n1)
  659.     equation (1.0, n1, 0.0, 0, 0.0, 0, -1.0);
  660.       else
  661.     equation (1.0, i, -0.5, i + 1, -0.5, dfa[i].backup, f2 ? 0.0 : -1.0);
  662.     equation (1.0, n, 0.0, 0, 0.0, 0, 0.0);
  663.  
  664.     free (dfa);
  665.  
  666.     return solvex0 ();
  667. }
  668.  
  669.  
  670. /* a simple gaussian elimination equation solver */
  671.  
  672. double *m, **M;
  673. int rank;
  674. int neq = 0;
  675.  
  676. /* create an n by n system of linear equations.  allocate space
  677.    for the matrix m, filled with zeroes and the dope vector M */
  678.  
  679. eqsystem (n) int n;
  680. {
  681.     char *calloc ();
  682.     int i;
  683.  
  684.     m = (double *) calloc (n * (n + 1), sizeof (double));
  685.     M = (double **) calloc (n, sizeof (double *));
  686.  
  687.     if (!m || !M) {
  688.     printf ("Ouch, out of memory!\n");
  689.     exit (1);
  690.     }
  691.  
  692.     for (i = 0; i < n; i++)
  693.       M[i] = &m[i * (n + 1)];
  694.     rank = n;
  695.     neq = 0;
  696. }
  697.  
  698. /* add a new equation a * x_na + b * x_nb + c * x_nc + d = 0.0
  699.    (note that na, nb, and nc are not necessarily all distinct.) */
  700.  
  701. equation (a, na, b, nb, c, nc, d) double a, b, c, d; int na, nb, nc;
  702. {
  703.     double *eq = M[neq++]; /* each row is an equation */
  704.     eq[na + 1] += a;
  705.     eq[nb + 1] += b;
  706.     eq[nc + 1] += c;
  707.     eq[0] = d;             /* column zero holds the constant term */
  708. }
  709.  
  710. /* solve for the value of variable x_0.  This will go nuts if
  711.    therer are errors (for example, if m is singular) */
  712.  
  713. double solvex0 ()
  714. {
  715.     register i, j, jmax, k;
  716.     register double  max, val;
  717.     register double *maxrow, *row;
  718.  
  719.  
  720.     for (i = rank; i > 0; --i) {      /* for each variable */
  721.  
  722.         /* find pivot element--largest value in ith column*/
  723.     max = 0.0;
  724.     for (j = 0; j < i; j++)
  725.         if (fabs (M[j][i]) > fabs (max)) {
  726.         max = M[j][i];
  727.         jmax = j;
  728.         }
  729.     /* swap pivot row with last row using dope vectors */
  730.     maxrow = M[jmax];
  731.     M[jmax] = M[i - 1];
  732.     M[i - 1] = maxrow;
  733.  
  734.     /* normalize pivot row */
  735.     max = 1.0 / max;
  736.     for (k = 0; k <= i; k++)
  737.       maxrow[k] *= max;
  738.  
  739.     /* now eliminate variable i by subtracting multiples of pivot row */
  740.     for (j = 0; j < i - 1; j++) {
  741.         row = M[j];
  742.         if (val = row[i])              /* if variable i is in this eq */
  743.         for (k = 0; k <= i; k++)
  744.           row[k] -= maxrow[k] * val;
  745.     }
  746.     }
  747.  
  748.     /* the value of x0 is now in constant column of first row
  749.        we only need x0, so no need to back-substitute */
  750.  
  751.     val = -M[0][0];
  752.  
  753.     free (M);
  754.     free (m);
  755.  
  756.     return val;
  757. }
  758.  
  759. _________________________________________________________________
  760. Guy Jacobson   (201) 582-6558              AT&T Bell Laboratories
  761.         uucp:  {att,ucbvax}!ulysses!guy       600 Mountain Avenue
  762.     internet:  guy@ulysses.att.com         Murray Hill NJ, 07974
  763.  
  764.  
  765.  
  766. ==> probability/flush.p <==
  767. Which set contains proportionately more flushes than the set of all
  768. possible poker hands?
  769. (1) Hands whose first card is an ace
  770. (2) Hands whose first card is the ace of spades
  771. (3) Hands with at least one ace
  772. (4) Hands with the ace of spades
  773.  
  774. ==> probability/flush.s <==
  775. An arbitrary hand can have two aces but a flush hand can't.  The
  776. average number of aces that appear in flush hands is the same as the
  777. average number of aces in arbitrary hands, but the aces are spread out
  778. more evenly for the flush hands, so set #3 contains a higher fraction
  779. of flushes.
  780.  
  781. Aces of spades, on the other hand, are spread out the same way over
  782. possible hands as they are over flush hands, since there is only one of
  783. them in the deck.  Whether or not a hand is flush is based solely on a
  784. comparison between different cards in the hand, so looking at just one
  785. card is necessarily uninformative.  So the other sets contain the same
  786. fraction of flushes as the set of all possible hands.
  787.  
  788. ==> probability/hospital.p <==
  789. A town has two hospitals, one big and one small.  Every day the big
  790. hospital delivers 1000 babies and the small hospital delivers 100
  791. babies.  There's a 50/50 chance of male or female on each birth.
  792. Which hospital has a better chance of having the same number of boys
  793. as girls?
  794.  
  795. ==> probability/hospital.s <==
  796. The small one.  If there are 2N babies born, then the probability of an
  797. even split is
  798.  
  799. (2N choose N) / (2 ** 2N) ,
  800.  
  801. where (2N choose N) = (2N)! / (N! * N!) .
  802.  
  803. This is a DECREASING function.
  804.  
  805. If there are two babies born, then the probability of a split is 1/2
  806. (just have the second baby be different from the first).  With 2N
  807. babies, If there is a N,N-1 split in the first 2N-1, then there is a
  808. 1/2 chance of the last baby making it an even split.  Otherwise there
  809. can be no even split.  Therefore the probability is less than 1/2
  810. overall for an even split.
  811.  
  812. As N goes to infinity the probability of an even split approaches zero
  813. (although it is still the most likely event).
  814.  
  815. ==> probability/icos.p <==
  816. The "house" rolls two 20-sided dice and the "player" rolls one
  817. 20-sided die.  If the player rolls a number on his die between the
  818. two numbers the house rolled, then the player wins.  Otherwise, the
  819. house wins (including ties).  What are the probabilities of the player
  820. winning?
  821.  
  822. ==> probability/icos.s <==
  823. It is easily seen that if any two of the three dice agree that the
  824. house wins.  The probability that this does not happen is 19*18/(20*20).
  825. If the three numbers are different, the probability of winning is 1/3.
  826. So the chance of winning is 19*18/(20*20*3) = 3*19/200 = 57/200.
  827.  
  828. ==> probability/intervals.p <==
  829. Given two random points x and y on the interval 0..1, what is the average
  830. size of the smallest of the three resulting intervals?
  831.  
  832. ==> probability/intervals.s <==
  833. In between these positions the surface forms a series of planes.
  834. Thus the volume under it consists of 2 pyramids each with an
  835. altitude of 1/3 and an (isosceles triangular) base of area 1/2,
  836. yielding a total volume of 1/9.
  837.  
  838. ==> probability/killers.and.pacifists.p <==
  839. You enter a town that has K killers and P pacifists.  When a
  840. pacifist meets a pacifist, nothing happens.  When a pacifist meets a
  841. killer, the pacifist is killed.  When two killers meet, both die.
  842. Assume meetings always occur between exactly two persons and the pairs
  843. involved are completely random.  What are your odds of survival?
  844.  
  845. ==> probability/killers.and.pacifists.s <==
  846. Regardless of whether you are a pacifist or a killer, you may disregard
  847. all events in which a pacifist other than yourself is involved and
  848. consider only events in which you are killed or a pair of killers other
  849. than yourself is killed.
  850.  
  851. Thus we may assume that there are N killers and yourself.  If N is odd,
  852. your odds of surviving are zero.  If N is even, it doesn't matter to
  853. you whether you are a killer or not.  So WLOG assume you are.  Then your
  854. probability of survival is 1/(N+1).
  855.  
  856. ==> probability/leading.digit.p <==
  857. What is the probability that the ratio of two random reals starts with a 1?
  858. What about 9?
  859.  
  860. ==> probability/leading.digit.s <==
  861. What is the distribution of y/x for (x,y) chosen with uniform distribution
  862. from the unit square?
  863.  
  864. First you want y/x in one of the intervals ... [0.01,0.02), [0.1,0.2),
  865. [1,2), [10/20), ... .  This corresponds to (x,y) lying in one of
  866. several triangles with height 1 and bases on either the right or top
  867. edges of the square.  The bases along the right edge have lengths 0.1
  868. (for 0.1 <= y/x < 0.2), 0.01, 0.001, ... ; the sum of this series is
  869. 1/9.  The bases along the top edge have lengths 0.5 (for 0.5 < x/y <=
  870. 1), 0.05, 0.005, ... ; the sum of this series is 5/9.  So you have a
  871. total base length of 6/9 = 2/3, height 1, so the area is 1/3.  The
  872. total area of the square is 1/3, so the probability that y/x starts
  873. with a 1 is 1/3 ~= 0.333333.
  874.  
  875. In the second case you have the same lengths (but in different places)
  876. on the right edge, total 1/9.  But on the top edge, 9 <= y/x < 10 gives
  877. you 1/10 < x/y <= 1/9 gives you a base of length 1/9 - 1/10 = 1/90,
  878. and the series proceeds 1/900, 1/9000, ... ; the sum is 1/81.  Total
  879. base length then is 9/81 + 1/81 = 10/81, height 1, total area (and
  880. hence probability of a leading 9) is 5/81 ~= 0.061728.
  881.  
  882.  
  883. ==> probability/lights.p <==
  884. Waldo and Basil are exactly m blocks west and n blocks north from
  885. Central Park, and always go with the green light until they run out of
  886. options.  Assuming that the probability of the light being green is 1/2
  887. in each direction, that if the light is green in one direction it is
  888. red in the other, and that the lights are not synchronized, find the
  889. expected number of red lights that Waldo and Basil will encounter.
  890.  
  891. ==> probability/lights.s <==
  892. Let E(m,n) be this number, and let (x)C(y) = x!/(y! (x-y)!).  A model
  893. for this problem is the following nxm grid:
  894.  
  895.      ^         B---+---+---+ ... +---+---+---+ (m,0)
  896.      |         |   |   |   |     |   |   |   |
  897.      N         +---+---+---+ ... +---+---+---+ (m,1)
  898. <--W + E-->    :   :   :   :     :   :   :   :
  899.      S         +---+---+---+ ... +---+---+---+ (m,n-1)
  900.      |         |   |   |   |     |   |   |   |
  901.      v         +---+---+---+ ... +---+---+---E (m,n)
  902.  
  903. where each + represents a traffic light.  We can consider each
  904. traffic light to be a direction pointer, with an equal chance of
  905. pointing either east or south.
  906.  
  907. IMHO, the best way to approach this problem is to ask:  what is the
  908. probability that edge-light (x,y) will be the first red edge-light
  909. that the pedestrian encounters?  This is easy to answer; since the
  910. only way to reach (x,y) is by going south x times and east y times,
  911. in any order, we see that there are (x+y)C(x) possible paths from
  912. (0,0) to (x,y).  Since each of these has probability (1/2)^(x+y+1)
  913. of occuring, we see that the the probability we are looking for is
  914. (1/2)^(x+y+1)*(x+y)C(x).  Multiplying this by the expected number
  915. of red lights that will be encountered from that point, (n-k+1)/2,
  916. we see that
  917.  
  918.             m-1
  919.            -----
  920.            \
  921. E(m,n) =    >    ( 1/2 )^(n+k+1) * (n+k)C(n) * (m-k+1)/2
  922.            /
  923.            -----
  924.             k=0
  925.  
  926.             n-1
  927.            -----
  928.            \
  929.       +     >    ( 1/2 )^(m+k+1) * (m+k)C(m) * (n-k+1)/2 .
  930.            /
  931.            -----
  932.             k=0
  933.  
  934.  
  935. Are we done?  No!  Putting on our Captain Clever Cap, we define
  936.  
  937.             n-1
  938.            -----
  939.            \
  940. f(m,n) =    >    ( 1/2 )^k * (m+k)C(m) * k 
  941.            /
  942.            -----
  943.             k=0
  944.  
  945. and
  946.  
  947.             n-1
  948.            -----
  949.            \
  950. g(m,n) =    >    ( 1/2 )^k * (m+k)C(m) .
  951.            /
  952.            -----
  953.             k=0
  954.  
  955. Now, we know that
  956.  
  957.              n
  958.            -----
  959.            \
  960. f(m,n)/2 =  >    ( 1/2 )^k * (m+k-1)C(m) * (k-1) 
  961.            /
  962.            -----
  963.             k=1
  964.  
  965. and since f(m,n)/2 = f(m,n) - f(m,n)/2, we get that
  966.  
  967.             n-1
  968.            -----
  969.            \
  970. f(m,n)/2 =  >    ( 1/2 )^k * ( (m+k)C(m) * k - (m+k-1)C(m) * (k-1) )
  971.            /
  972.            -----
  973.             k=1
  974.  
  975. - (1/2)^n * (m+n-1)C(m) * (n-1)
  976.  
  977.             n-2
  978.            -----
  979.            \
  980.  =          >    ( 1/2 )^(k+1) * (m+k)C(m) * (m+1)
  981.            /
  982.            -----
  983.             k=0
  984.  
  985. - (1/2)^n * (m+n-1)C(m) * (n-1)
  986.  
  987. = (m+1)/2 * (g(m,n) - (1/2)^(n-1)*(m+n-1)C(m)) - (1/2)^n*(m+n-1)C(m)*(n-1)
  988.  
  989. therefore
  990.  
  991. f(m,n) = (m+1) * g(m,n) - (n+m) * (1/2)^(n-1) * (m+n-1)C(m) .
  992.  
  993.  
  994. Now, E(m,n) = (n+1) * (1/2)^(m+2) * g(m,n) - (1/2)^(m+2) * f(m,n)
  995.  
  996. + (m+1) * (1/2)^(n+2) * g(n,m) - (1/2)^(n+2) * f(n,m)
  997.  
  998. = (m+n) * (1/2)^(n+m+1) * (m+n)C(m) + (m-n) * (1/2)^(n+2) * g(n,m)
  999.  
  1000. + (n-m) * (1/2)^(m+2) * g(m,n) .
  1001.  
  1002.  
  1003. Setting m=n in this formula, we see that
  1004.  
  1005.            E(n,n) = n * (1/2)^(2n) * (2n)C(n),
  1006.  
  1007. and applying Stirling's theorem we get the beautiful asymptotic formula
  1008.  
  1009.                   E(n,n) ~ sqrt(n/pi).
  1010.  
  1011. ==> probability/lottery.p <==
  1012. There n tickets in the lottery, k winners and m allowing you to pick another
  1013. ticket. The problem is to determine the probability of winning the lottery
  1014. when you start by picking 1 (one) ticket.
  1015.  
  1016. A lottery has N balls in all, and you as a player can choose m numbers
  1017. on each card, and the lottery authorities then choose n balls, define
  1018. L(N,n,m,k) as the minimum number of cards you must purchase to ensure that
  1019. at least one of your cards will have at least k numbers in common with the
  1020. balls chosen in the lottery.
  1021.  
  1022. ==> probability/lottery.s <==
  1023. This relates to the problem of rolling a random number
  1024. from 1 to 17 given a 20 sided die.  You simply mark the numbers 18,
  1025. 19, and 20 as "roll again".
  1026.  
  1027. The probability of winning is, of course, k/(n-m) as for every k cases
  1028. in which you get x "draw again"'s before winning, you get n-m-k similar
  1029. cases where you get x "draw again"'s before losing.  The example using
  1030. dice makes it obvious, however.
  1031.  
  1032. L(N,k,n,k) >= Ceiling((N-choose-k)/(n-choose-k)*
  1033.                    (n-1-choose-k-1)/(N-1-choose-k-1)*L(N-1,k-1,n-1,k-1))
  1034.             = Ceiling(N/n*L(N-1,k-1,n-1,k-1))
  1035.  
  1036. The correct way to see this is as follows:  Suppose you have an
  1037. (N,k,n,k) system of cards.  Look at all of the cards that contain the
  1038. number 1.  They cover all ball sets that contain 1, and therefore these
  1039. cards become an (N-1,k-1,n-1,k-1) covering upon deletion of the number
  1040. 1.  Therefore the number 1 must appear at least L(N-1,k-1,n-1,k-1).
  1041. The same is true of all of the other numbers.  There are N of them but
  1042. n appear on each card.  Thus we obtain the bound.
  1043.  
  1044. ==> probability/oldest.girl.p <==
  1045. You meet a stranger on the street, and ask how many children he has.  He
  1046. truthfully says two.  You ask "Is the older one a girl?"  He truthfully
  1047. says yes.  What is the probability that both children are girls?  What
  1048. would the probability be if your second question had been "Is at least
  1049. one of them a girl?", with the other conditions unchanged?
  1050.  
  1051. ==> probability/oldest.girl.s <==
  1052. There are four possibilities:
  1053.  
  1054.     Oldest child    Youngest child
  1055. 1.  Girl            Girl
  1056. 2.  Girl            Boy
  1057. 3.  Boy             Girl
  1058. 4.  Boy             Boy
  1059.  
  1060. If your friend says "My oldest child is a girl," he has eliminated cases
  1061. 3 and 4, and in the remaining cases both are girls 1/2 of the time.  If
  1062. your friend says "At least one of my children is a girl," he has
  1063. eliminated case 4 only, and in the remaining cases both are girls 1/3
  1064. of the time.
  1065.  
  1066.  
  1067. ==> probability/particle.in.box.p <==
  1068. A particle is bouncing randomly in a two-dimensional box.  How far does it
  1069. travel between bounces, on average?
  1070.  
  1071. Suppose the particle is initially at some random position in the box and is
  1072. traveling in a straight line in a random direction and rebounds normally
  1073. at the edges.
  1074.  
  1075. ==> probability/particle.in.box.s <==
  1076. Let theta be the angle of the point's initial vector.  After traveling a
  1077. distance r, the point has moved r*cos(theta) horizontally and r*sin(theta)
  1078. vertically, and thus has struck r*(sin(theta)+cos(theta))+O(1) walls.  Hence
  1079. the average distance between walls will be 1/(sin(theta)+cos(theta)).  We now
  1080. average this over all angles theta:
  1081.     2/pi * intg from theta=0 to pi/2 (1/(sin(theta)+cos(theta))) dtheta
  1082. which (in a computation which is left as an exercise) reduces to
  1083.     2*sqrt(2)*ln(1+sqrt(2))/pi = 0.793515.
  1084.  
  1085. ==> probability/pi.p <==
  1086. Are the digits of pi random (i.e., can you make money betting on them)?
  1087.  
  1088. ==> probability/pi.s <==
  1089. No, the digits of pi are not truly random, therefore you can win money
  1090. playing against a supercomputer that can calculate the digits of pi far
  1091. beyond what we are currently capable of doing.  This computer selects a
  1092. position in the decimal expansion of pi -- say, at 10^100.  Your job is
  1093. to guess what the next digit or digit sequence is.  Specifically, you
  1094. have one dollar to bet.  A bet on the next digit, if correct, returns
  1095. 10 times the amount bet; a bet on the next two digits returns 100 times
  1096. the amount bet, and so on.  (The dollar may be divided in any fashion,
  1097. so we could bet 1/3 or 1/10000 of a dollar.)  You may place bets in any
  1098. combination.  The computer will tell you the position number, let you
  1099. examine the digits up to that point, and calculate statistics for you.
  1100.  
  1101. It is easy to set up strategies that might win, if the supercomputer
  1102. doesn't know your strategy.  For example, "Always bet on 7" might win,
  1103. if you are lucky.  Also, it is easy to set up bets that will always
  1104. return a dollar.  For example, if you bet a penny on every two-digit
  1105. sequence, you are sure to win back your dollar.  Also, there are
  1106. strategies that might be winning, but we can't prove it.  For example,
  1107. it may be that a certain sequence of digits never occurs in pi, but we
  1108. have no way of proving this.
  1109.  
  1110. The problem is to find a strategy that you can prove will always win
  1111. back more than a dollar.
  1112.  
  1113. The assumption that the position is beyond the reach of calculation
  1114. means that we must rely on general facts we know about the sequence of
  1115. digits of pi, which is practically nil.  But it is not completely nil,
  1116. and the challenge is to find a strategy that will always win money.
  1117.  
  1118. A theorem of Mahler (1953) states that for all integers p, q > 1,
  1119.  
  1120.         -42
  1121.   |pi - p/q| > q
  1122.  
  1123. This says that pi cannot have a rational approximation that is
  1124. extremely tight.
  1125.  
  1126. Now suppose that the computer picks position N.  I know that the next
  1127. 41 * N digits cannot be all zero.   For if they were, then the first N
  1128. digits, treated as a fraction with denominator 10^N, satisfies:
  1129.  
  1130.   |pi - p / 10^N|  <  10^(-42 N)
  1131.  
  1132. which contradicts Mahler's theorem.
  1133.  
  1134. So, I split my dollar into 10^(41N) - 1 equal parts, and bet on each of
  1135. the sequences of 41N digits, except the one with all zeroes.  One of
  1136. the bets is sure to win, so my total profit is about 10(^-41N) of a
  1137. dollar!
  1138.  
  1139. This strategy can be improved a number of ways, such as looking for
  1140. other repeating patterns, or improvements to the bound of 42 -- but the
  1141. earnings are so pathetic, it hardly seems worth the effort.
  1142.  
  1143. Are there other winning strategies, not based on Mahler's theorem?  I
  1144. believe there are algorithms that generate 2N binary digits of pi,
  1145. where the computations are separate for each block of N digits.  Maybe
  1146. from something like this, we can find a simple subsequence of the
  1147. binary digits of pi which is always zero, or which has some simple
  1148. pattern.
  1149.  
  1150. ==> probability/random.walk.p <==
  1151. Waldo has lost his car keys!  He's not using a very efficient search;
  1152. in fact, he's doing a random walk.  He starts at 0, and moves 1 unit
  1153. to the left or right, with equal probability.  On the next step, he
  1154. moves 2 units to the left or right, again with equal probability.  For
  1155. subsequent turns he follows the pattern 1, 2, 1, etc.
  1156.  
  1157. His keys, in truth, were right under his nose at point 0.  Assuming
  1158. that he'll spot them the next time he sees them, what is the
  1159. probability that poor Waldo will eventually return to 0?
  1160.  
  1161. ==> probability/random.walk.s <==
  1162. I can show the probability that Waldo returns to 0 is 1.  Waldo's
  1163. wanderings map to an integer grid in the plane as follows.  Let
  1164. (X_t,Y_t) be the cumulative sums of the length 1 and length 2 steps
  1165. respectively taken by Waldo through time t.  By looking only at even t,
  1166. we get the ordinary random walk in the plane, which returns to the
  1167. origin (0,0) with probability 1.  In fact, landing at (2n, n) for any n
  1168. will land Waldo on top of his keys too.  There's no need to look at odd
  1169. t.
  1170.  
  1171. Similar considerations apply for step sizes of arbitrary (fixed) size.
  1172.  
  1173. ==> probability/reactor.p <==
  1174. There is a reactor in which a reaction is to take place. This reaction
  1175. stops if an electron is present in the reactor. The reaction is started
  1176. with 18 positrons; the idea being that one of these positrons would
  1177. combine with any incoming electron (thus destroying both). Every second,
  1178. exactly one particle enters the reactor. The probablity that this particle   
  1179. is an electron is 0.49 and that it is a positron is 0.51.
  1180.  
  1181. What is the probability that the reaction would go on for ever?
  1182.  
  1183. Note: Once the reaction stops, it cannot restart.
  1184.  
  1185. ==> probability/reactor.s <==
  1186. Let P(n) be the probability that, starting with n positrons, the
  1187. reaction goes on forever.  Clearly P'(n+1)=P'(0)*P'(n), where the
  1188. ' indicates probabilistic complementation; also note that
  1189. P'(n) = .51*P'(n+1) + .49*P'(n-1).  Hence we have that P(1)=(P'(0))^2
  1190. and that P'(0) = .51*P'(1) ==> P'(0) equals 1 or 49/51.  We thus get
  1191. that either P'(18) = 1 or (49/51)^19 ==> P(18) = 0 or 1 - (49/51)^19.
  1192.  
  1193. The answer is indeed the latter.  A standard result in random walks
  1194. (which can be easily derived using Markov chains) yields that if p>1/2
  1195. then the probability of reaching the absorbing state at +infinity
  1196. as opposed to the absorbing state at -1 is 1-r^(-i), where r=p/(1-p)
  1197. (p is the probability of moving from state n to state n-1, in our
  1198. case .49) and i equals the starting location + 1.  Therefore we have
  1199. that P(18) = 1-(.49/.51)^19.
  1200.  
  1201. ==> probability/roulette.p <==
  1202. You are in a game of Russian roulette, but this time the gun (a 6
  1203. shooter revolver) has three bullets _in_a_row_ in three of the
  1204. chambers.  The barrel is spun only once.  Each player then points the
  1205. gun at his (her) head and pulls the trigger.  If he (she) is still
  1206. alive, the gun is passed to the other player who then points it at his
  1207. (her) own head and pulls the trigger.  The game stops when one player
  1208. dies.
  1209.  
  1210. Now to the point:  would you rather be first or second to shoot?
  1211.  
  1212. ==> probability/roulette.s <==
  1213. All you need to consider are the six possible bullet configurations
  1214.  
  1215.     B B B E E E         -> player 1 dies
  1216.     E B B B E E         -> player 2 dies
  1217.     E E B B B E         -> player 1 dies
  1218.     E E E B B B         -> player 2 dies
  1219.     B E E E B B         -> player 1 dies
  1220.     B B E E E B         -> player 1 dies
  1221.  
  1222. One therefore has a 2/3 probability of winning (and a 1/3 probability of
  1223. dying) by shooting second.  I for one would prefer this option.
  1224.  
  1225. ==> probability/transitivity.p <==
  1226. Can you number dice so that die A beats die B beats die C beats die A?
  1227. What is the largest probability p with which each event can occur?
  1228.  
  1229. ==> probability/transitivity.s <==
  1230. Yes.  The actual values on the dice faces don't matter, only their
  1231. ordering.  WLOG we may assume that no two faces of the same or
  1232. different dice are equal.  We can assume "generalised dice", where the
  1233. faces need not be equally probable.  These can be approximated by dice
  1234. with equi-probable faces by having enough faces and marking some of
  1235. them the same.
  1236.  
  1237. Take the case of three dice, called A, B, and C.  Picture the different
  1238. values on the faces of the A die.  Suppose there are three:
  1239.  
  1240.             A       A       A
  1241.  
  1242. The values on the B die must lie in between those of the A die:
  1243.  
  1244.         B   A   B   A   B   A   B
  1245.  
  1246. With three different A values, we need only four different B values.
  1247.  
  1248. Similarly, the C values must lie in between these:
  1249.  
  1250.       C B C A C B C A C B C A C B C
  1251.       
  1252. Assume we want A to beat B, B to beat C, and C to beat A.  Then the above
  1253. scheme for the ordering of values can be simplified to:
  1254.  
  1255.       B C A B C A B C A B C
  1256.  
  1257. since for example, the first C in the previous arrangement can be moved
  1258. to the second with the effect that the probability that B beats C is
  1259. increased, and the probabilities that C beats A or A beats B are
  1260. unchanged.  Similarly for the other omitted faces.
  1261.  
  1262. In general we obtain for n dice A...Z the arrangement
  1263.  
  1264.     B ... Z A B ... Z ...... A B ... Z
  1265.  
  1266. where there are k complete cycles of B..ZA followed by B...Z.  k must be
  1267. at least 1.
  1268.  
  1269. CONJECTURE:  The optimum can be obtained for k=1.
  1270.  
  1271. So the arrangement of face values is B ... Z A B ... Z.  For three dice
  1272. it is BCABC.  Thus one die has just one face, all the other dice have two
  1273. (with in general different probabilities).
  1274.  
  1275. CONJECTURE:  At the optimum, the probabilities that each die beats the
  1276. next can be equal.
  1277.  
  1278. Now put probabilities into the BCABC arrangement:
  1279.  
  1280.     B  C  A  B  C
  1281.     x  y  1  x' y'
  1282.  
  1283. Clearly x+x' = y+y' = 1.
  1284.  
  1285. Prob. that A beats B = x'
  1286.            B beats C = x + x'y'
  1287.            C beats A = y
  1288.  
  1289. Therefore x' = y = x + x'y'
  1290.  
  1291. Solving for these gives x = y' = 1-y, x' = y = (-1 + sqrt(5))/2 = prob.
  1292. of each die beating the next = 0.618...
  1293.  
  1294. For four dice one obtains the probabilities:
  1295.  
  1296.     B  C  D  A  B  C  D
  1297.     x  y  z  1  x' y' z'
  1298.  
  1299. A beats B:  x'
  1300. B beats C:  x + x'y'
  1301. C beats D:  y + y'z'
  1302. D beats A:  z
  1303.  
  1304. CONJECTURE: for any number of dice, at the optimum, the sequence of
  1305. probabilities abc...z1a'b'c...z' is palindromic.
  1306.  
  1307. We thus have the equalities:
  1308.  
  1309.     x+x' = 1
  1310.     y+y' = 1
  1311.     z+z' = 1
  1312.     x' = z = x + x'y' = x + x'y'
  1313.     y = y' (hence both = 1/2)
  1314.  
  1315. Solving this gives x = 1/3, z = 2/3 = prob. of each die beating the next.
  1316.  Since all the numbers are rational, the limit is attainable with
  1317. finitely many equiprobable faces.  E.g. A has one face, marked 0.  C has
  1318. two faces, marked 2 and -2.  B has three faces, marked 3, -1, -1.  D has
  1319. three faces, marked 1, 1, -3.  Or all four dice can be given six faces,
  1320. marked with numbers in the range 0 to 6.
  1321.  
  1322. Finding the solution for 5, 6, or n dice is left as an exercise.
  1323.  
  1324. --                                  ____
  1325. Richard Kennaway                  __\_ /    School of Information Systems
  1326. Internet:  jrk@sys.uea.ac.uk      \  X/     University of East Anglia
  1327. uucp:  ...mcsun!ukc!uea-sys!jrk    \/       Norwich NR4 7TJ, U.K.
  1328.  
  1329.  
  1330. Martin Gardner (of course!) wrote about notransitive dice, see the Oct '74
  1331. issue of Scientific American, or his book "Wheels, Life and Other Mathematical
  1332. Amusements", ISBN 0-7167-1588-0 or ISBN 0-7167-1589-9 (paperback).
  1333.  
  1334. In the book, Gardner cites Bradley Efron of Stanford U. as stating that
  1335. the maximum number for three dice is approx .618, requiring dice with more
  1336. than six sides.  He also mentions that .75 is the limit approached as the
  1337. number of dice increases.  The book shows three sets of 6-sided dice, where
  1338. each set has 2/3 as the advantage probability.
  1339.  
  1340.