Steven R. Dunbar
Department of Mathematics
203 Avery Hall
Lincoln, NE 68588-0130
http://www.math.unl.edu
Voice: 402-472-3731
Fax: 402-472-8466

Stochastic Processes and

__________________________________________________________________________

A Coin Tossing Experiment

_______________________________________________________________________

Note: These pages are prepared with MathJax. MathJax is an open source JavaScript display engine for mathematics that works in all browsers. See http://mathjax.org for details on supported browsers, accessibility, copy-and-paste, and other features.

_______________________________________________________________________________________________ ### Rating

Student: contains scenes of mild algebra or calculus that may require guidance.

_______________________________________________________________________________________________ ### Section Starter Question

Suppose you start with a fortune of $10, and you want to gamble to achieve$20 before you go broke. Your gambling game is to ﬂip a fair coin successively, and gain $1 if the coin comes up “Heads” and lose$1 if the coin comes up “Tails”. What do you estimate is the probability of achieving $20 before going broke? How long do you estimate it takes before one of the two outcomes occurs? How do you estimate each of these quantities? _______________________________________________________________________________________________ ### Key Concepts 1. Performing an experiment to gain intuition about coin-tossing games. __________________________________________________________________________ ### Vocabulary 1. We call victory the outcome of reaching a fortune goal by chance before going broke. 2. We call ruin the outcome of going broke by chance before reaching a fortune goal. __________________________________________________________________________ ### Mathematical Ideas #### Pure Risk Modeling We need a better understanding of the paths that risky securities take. We shall make and investigate a greatly simpliﬁed model of randomness and risk. For our model, we assume: 1. Time is discrete, occurring at ${t}_{0}=0,{t}_{1},{t}_{2},\dots$. 2. No risk-free investments are available, (i.e. no interest-bearing bonds). 3. No options, and no ﬁnancial derivatives are available. 4. The only investments are risk-only, that is, our net fortune at the $n$th time is a random variable: ${T}_{n+1}={T}_{n}+{Y}_{n+1}$ where ${T}_{0}=0$ is our given initial fortune, and for simplicity, Our model is commonly called “gambling” and we will investigate the probabilities of making a fortune by gambling. #### Some Humor Figure 1: Welcome to my casino! Figure 2: Where not to lose at the casino! #### An Experiment The reader should perform the following experiment as a “gambler” to gain intuition about the coin-tossing game. We call victory the outcome of reaching a fortune goal by chance before going broke. We call ruin the outcome of going broke by chance before reaching a fortune goal. 1. Each “gambler” has a chart for recording the outcomes of each game (see below) and a sheet of graph paper. 2. Each “gambler” has a fair coin to ﬂip, say a penny. 3. Each “gambler” ﬂips the coin, and records a $+1$ (gains$1) if the coin comes up “Heads” and records $-1$ (loses $1) if the coin comes up “Tails”. On the chart, the player records the outcome of each ﬂip with the ﬂip number, the outcome as “H” or “T” and keeps track of the cumulative fortune of the gambler so far. Keep these records in a neat chart, since some problems refer to them later. 4. Each “gambler” should record $100$ ﬂips, which takes about 10 to 20 minutes.  Toss $n$ 1 2 3 4 5 6 7 8 9 10 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 11 12 13 14 15 16 17 18 19 20 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 21 22 23 24 25 26 27 28 29 30 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 31 32 33 34 35 36 37 38 39 40 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 41 42 43 44 45 46 47 48 49 50 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 51 52 53 54 55 56 57 58 59 60 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 61 62 63 64 65 66 67 68 69 70 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 71 72 73 74 75 76 77 78 79 80 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 81 82 83 84 85 86 87 88 89 90 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$  Toss $n$ 91 92 93 94 95 96 97 98 99 100 H or T ${Y}_{n}=+1,-1$ ${T}_{n}={\sum }_{i=1}^{n}{Y}_{i}$ #### Some Outcomes A simulation with 30 “gamblers” obtained the following results: 1. 14 gamblers reached a net loss of $-10$ before reaching a net gain of $+10$, that is, were “ruined”. 2. 9 gamblers reached a net gain of $+10$ before reaching a net loss of $-10$, that is, achieved “victory”. 3. 7 gamblers still had not reached a net loss of $+10$ or $-10$ yet. This roughly matches the predicted outcomes of $1∕2$ the gamblers being ruined, not as close for the predicted proportion of $1∕2$ for victory. The mean duration until victory or ruin was $41.6$. Compare this to your results. #### Sources This section is adapted from ideas in William Feller’s classic text, An Introduction to Probability Theory and Its Applications, Volume I, Third Edition. The GeoGebra script is adapted from a demonstration and discussion among “carlosgomes” of Amarante, Portugal; Noel Lambert of Vosges, France; and Juan Vicente Snchez Gaitero of Cadiz, Spain on geogebra.org ____________________________________________________________ ### Algorithms, Scripts, Simulations #### Algorithm The probability of Heads is $p$. In matrix-oriented languages, use $k$ columns of the matrix for the iterations of an experiment each with $n$ trials. This is more eﬃcient than using for-loops or the equivalent. Then each $n×1$ column vector contains the outcomes of the experiment. Use the random number generator to create a random number in $\left(0,1\right)$. To simulate the ﬂip of a coin use logical functions to score a 1, also counted as a Head, if the value is less than $p$; a 0, also counted as a Tail, if the value is greater than $p$. Column sums then give a $1×k$ vector of total numbers of Heads. Statistical measures can then be applied to this vector. For wins and losses, linearly transform the coin ﬂips to $1$ if Heads, a $-1$ if Tails. Cumulatively sum the entries to create the running total. Apply statistical measures to the last row of the totals. Finally, set values for victory and ruin. Check when each column ﬁrst hits victory or ruin. #### Experiment Scripts Geogebra + R 1p <- 0.5 2n <- 100 3k <- 30 4coinFlips <- array( 0+(runif(n*k) <= p), dim=c(n,k)) 5 # 0+ coerces Boolean to numeric 6headsTotal <- colSums(coinFlips) # 0..n binomial rv sample, size k 7muHeads <- mean(headsTotal) # Expected value is n/2 8sigmaSquaredHeads <- var(headsTotal) # Theoretical value is np(1-p) 9cat(sprintf("Empirical Mean of Heads: %f \n", muHeads )) 10cat(sprintf("Empirical Variance of Heads: %f \n", sigmaSquaredHeads )) 11 12winLose <- 2*coinFlips - 1 # -1 for Tails, 1 for Heads 13totals <- apply( winLose, 2, cumsum) 14 # -n..n (every other integer) binomial rv sample 15 # the second argument ‘‘2’’ means column-wise 16muWinLose <- mean( totals[n,]) # Expected value is 0 17sigmaSquaredWinLose <- var( totals[n,]) # Theoretical value is 4np(1-p) 18cat(sprintf("Empirical Mean of Wins minus Losses: %f \n", muWinLose )) 19cat(sprintf("Empirical Variance of Wins minus Losses: %f \n", sigmaSquaredWinLose )) 20 Octave 1p = 0.5; 2n = 100; 3k = 30; 4 5coinFlips = rand(n,k) <= p; 6headsTotal = sum(coinFlips); # 0..n binomial rv sample, size k 7muHeads = mean(headsTotal); # Expected value is n/2 8sigmaSquaredHeads = var(headsTotal); # Theoretical value is np(1-p) 9disp("Empirical Mean of Heads:"), disp(muHeads) 10disp("Empirical Variance of Heads:"), disp(sigmaSquaredHeads) 11 12winLose = 2*coinFlips - 1; # -1 for Tails, 1 for Heads 13totals = cumsum(winLose); # -n..n (every other integer) binomial rv sample 14muWinLose = mean( totals(n,:) ); # Expected value is 0 15sigmaSquaredWinLose = var( totals(n,:) ); # Theoretical value is 4np(1-p) 16disp("Empirical Mean of Wins minus Losses:"), disp(muWinLose) 17disp("Empirical Variance of Wins minus Losses:"), \ 18 disp(sigmaSquaredWinLose) 19 Perl 1$p = 0.5;
2$n = 100; 3$k = 30;
4
5$coinFlips = random($k, $n ) <=$p;    #note order of dims!!
6$headsTotal = 7$coinFlips->transpose->sumover;     # 0..n binomial r.v. sample, size k
8
9#note transpose, PDL likes x (row) direction for implicitly threaded operations
10
11( $muHeads,$prms, $median,$min, $max,$adev, $rmsHeads ) = 12$headsTotal->stats;
13
14# easiest way to get the descriptive statistics
15
16printf "Empirical Mean of Heads: %6.3f \n", $muHeads; 17 18# Theoretical value is np 19$varHeads = $rmsHeads *$rmsHeads;
20printf "Emprical Variance of Heads: %6.3f \n", $varHeads; 21 22# Theoretical value is np(1-p) 23 24$winLose = 2 * $coinFlips - 1; # note use of threading 25$totals = ( $winLose->xchg( 0, 1 )->cumusumover )->transpose; 26 27# note xchg transpose, PDL likes x (row) direction for implicity threaded operations 28 29($muWinLose, $prms,$median, $min,$max, $adev,$rmsWinLose ) =
30    $totals ( 0 :$k - 1, $n - 1 )->stats; # Expected value is 0 31 32# pdl are indexed from 0 like C and Perl, so thats why$k-1, $n-1 33 34printf "Empirical Mean number of Wins over Losses: %6.3f \n",$muWinLose;
35
36# Theoretical value is 0
37$varWinLose =$rmsWinLose * $rmsWinLose; # Theoretical value is 4np(1-p) 38printf "Empirical Variance of Wins over Losses: %6.3f \n ",$varWinLose;
39
SciPy
1import scipy
2
3p = 0.5
4n = 100
5k = 30
6
7coinFlips = scipy.random.random((n,k))<= p
8# Note Booleans True for Heads and False for Tails
9headsTotal = scipy.sum(coinFlips, axis = 0)
10# Note how Booleans act as 0 (False) and 1 (True)
16
17winLose = 2*coinFlips - 1
18totals = scipy.cumsum(winLose, axis = 0)
19muWinLose = scipy.mean(totals[n-1,range(k)])
20stddevWinLose = scipy.std(totals[n-1, range(k)])
21sigmaSquaredWinLose = stddevWinLose * stddevWinLose
22print "Empirical Mean of Wins minus Losses:", muWinLose
23print "Empirical Variance of Wins minus Losses:",  sigmaSquaredWinLose

#### Victory or Ruin Scripts

R
1    victory <- 10  # top boundary for random walk
2    ruin <- -10    # bottom boundary for random walk
3
4    victoryOrRuin <- array(0, dim = c(k))
5
6    hitVictory <- apply( (totals >= victory), 2, which)
7    hitRuin <- apply( (totals <= ruin), 2, which)
8
9    for (j in 1:k) {
10       if ( length(hitVictory[[j]]) == 0 && length(hitRuin[[j]]) == 0 ) {
11         # no victory, no ruin
12         # do nothing
13       }
14       else if ( length(hitVictory[[j]]) > 0 && length(hitRuin[[j]]) == 0 ) {
15         # victory, no ruin
16         victoryOrRuin[j] <- min(hitVictory[[j]])
17   }
18       else if ( length(hitVictory[[j]]) == 0 && length(hitRuin[[j]]) > 0 ) {
19         # no victory, ruin
20         victoryOrRuin[j] <- -min(hitRuin[[j]])
21         }
22       else # ( length(hitVictory[[j]]) > 0 && length(hitRuin[[j]]) > 0 )
23         # victory and ruin
24   if ( min(hitVictory[[j]]) < min(hitRuin[[j]]) ) { # victory first
25           victoryOrRuin[j] <- min(hitVictory[[j]])  # code hitting victory
26   }
27   else {  # ruin first
28           victoryOrRuin[j] <- -min(hitRuin[[j]]) # code hitting ruin as negative
29   }
30    }
31
32victoryBeforeRuin <- sum(0+(victoryOrRuin > 0))  # count exits through top
33ruinBeforeVictory <- sum(0+(victoryOrRuin < 0))  # count exits through bottom
34noRuinOrVictory   <- sum(0+(victoryOrRuin == 0))
35
36cat(sprintf("Victories: %i  Ruins: %i No Ruin or Victory: %i \n",
37    victoryBeforeRuin, ruinBeforeVictory, noRuinOrVictory))
38
39avgTimeVictoryOrRuin <- mean( abs(victoryOrRuin) )
40varTimeVictoryOrRuin <- var( abs(victoryOrRuin) )
41
42cat(sprintf("Average Time to Victory or Ruin: %f \n",
43    avgTimeVictoryOrRuin))
44cat(sprintf("Variance of Time to Victory or Ruin: %f \n",
45    varTimeVictoryOrRuin))
46
47hist( victoryOrRuin, nclass = 2*max(abs(victoryOrRuin))+1 ) # sample exit time distribution
48
Octave
1victory =  10;      # top boundary for random walk
2ruin    = -10;      # bottom boundary for random walk
3
4# Octave find treats a whole matrix as a single column,
5# but I need to treat experiment column by column, so here is a place
6# where for clarity, better to loop over columns, using find on each
7victoryOrRuin = zeros(1,k);
8for j = 1:k
9    hitVictory = find(totals(:,j) >= victory);
10    hitRuin  = find(totals(:,j) <= ruin);
11    if ( !rows(hitVictory) && !rows(hitRuin) )
12       # no victory, no ruin
13       # do nothing
14    elseif ( rows(hitVictory) && !rows(hitRuin) )
15       # victory, no ruin
16       victoryOrRuin(j) = hitVictory(1);
17    elseif ( !rows(hitVictory) && rows(hitRuin) )
18       # no victory, but hit ruin
19       victoryOrRuin(j) = -hitRuin(1);
20    else # ( rows(hitvictory) && rows(hitruin) )
21       # victory and ruin
22       if ( hitVictory(1) < hitRuin(1) )
23   victoryOrRuin(j) = hitVictory(1); # code hitting victory
24       else
25   victoryOrRuin(j) = -hitRuin(1); # code hitting ruin as negative
26       endif
27    endif
28endfor
29
30victoryOrRuin     # display vector of k hitting times
31
32victoryBeforeRuin = sum( victoryOrRuin > 0 ) # count exits through top
33ruinBeforeVictory = sum( victoryOrRuin < 0 ) # count exits through bottom
34noRuinOrVictory   = sum( victoryOrRuin == 0 )
35
36avgTimeVictoryOrRuin = mean( abs(victoryOrRuin) )
37varTimeVictoryOrRuin = var( abs(victoryOrRuin) )
38
39max(totals)
40hist( victoryOrRuin, max(max(victoryOrRuin))) # sample exit time distribution
41
Perl
1use PDL::NiceSlice;
2
3$victory = 10; # top boundary for random walk 4$ruin    = -10;    # bottom boundary for random walk
5
6$victoryOrRuin = zeros($k);
7
8$hitVictory = ($totals >= $victory ); 9$hitRuin    = ( $totals <=$ruin );
10
11$whenVictory =$hitVictory ( 0 : $k - 1, : )->xchg( 0, 1 )->maximum_ind; 12$whenRuin    = $hitRuin ( 0 :$k - 1, : )->xchg( 0, 1 )->maximum_ind;
13foreach $i ( 0 ..$k - 1 ) {
14    if ( $whenVictory ($i) == 0 and $whenRuin ($i) == 0 ) {
15
16        # no victory, no ruin, do nothing
17        $victoryOrRuin ($i) = 0;
18    }
19    elsif ( $whenVictory ($i) > 0 and $whenRuin ($i) == 0 ) {
20        $victoryOrRuin ($i) .= $whenVictory ($i);
21    }    # victory, no ruin
22    elsif ( $whenVictory ($i) == 0 and $whenRuin ($i) > 0 ) {
23        $victoryOrRuin ($i) .= -$whenRuin ($i);
24    }    # no victory, ruin
25    else {
26        if ( $whenVictory ($i) < $whenRuin ($i) ) {    # victory first
27            $victoryOrRuin ($i) .= $whenVictory ($i);   # code hitting victory
28        }
29        else {
30            $victoryOrRuin ($i)
31                .= -$whenRuin ($i);    # code hitting ruin as negative
32        }
33    }
34}
35
36$victoryBeforeRuin = 37 ($victoryOrRuin > 0 )->sumover;    # count exits through top
38$ruinBeforeVictory = 39 ($victoryOrRuin < 0 )->sumover;    # count exits through bottom
40
41(   $muTimeVictoryOrRuin,$prms, $median,$min, $max,$adev,
42    $rmsTimeVictoryOrRuin 43) =$victoryOrRuin->abs->stats;
44
45print "Mean time to victory or ruin", $muTimeVictoryOrRuin, "\n"; 46print "Standard deviation of time to victory or ruin",$rmsTimeVictoryOrRuin,
47    "\n";
48
49$histLimit =$victoryOrRuin->abs->max;
50print hist( $victoryOrRuin, -$histLimit, \$histLimit, 1 ), "\n";
51
SciPy
1victory = 10;                   # top boundary for random walk
2ruin = -10;                     # bottom boundary of random walk
3
4victoryOrRuin = scipy.zeros(k,int);
5hitVictory = scipy.argmax( (totals >= victory), axis=0)
6hitRuin = scipy.argmax( (totals <= ruin), axis=0)
7for i in range(k):
8    if hitVictory[i] == 0 and  hitRuin[i] == 0:
9        # no victory, no ruin, do nothing
10        victoryOrRuin[i] = 0
11    elif hitVictory[i] > 0 and  hitRuin[i] == 0:
12        victoryOrRuin[i] = hitVictory[i]
13    elif hitVictory[i] == 0 and  hitRuin[i] > 0:
14        victoryOrRuin[i] = -hitRuin[i]
15    else:
16        if hitVictory[i] < hitRuin[i]:
17            victoryOrRuin[i] = hitVictory[i]
18        else:
19            victoryOrRuin[i] = -hitRuin[i]
20
21victoryBeforeRuin = sum( victoryOrRuin > 0 ) # count exits through top
22ruinBeforeVictory = sum( victoryOrRuin < 0 ) # count exits through bottom
23noRuinOrVictory   = sum( victoryOrRuin == 0 )
24
25print "Victories:", victoryBeforeRuin, "Ruins:", ruinBeforeVictory,
26print "No Ruin or Victory:",  noRuinOrVictory
27
28avgTimeVictoryOrRuin = scipy.mean( abs(victoryOrRuin) )
29stdTimeVictoryOrRuin = scipy.std( abs(victoryOrRuin) )
30varTimeVictoryOrRuin = stdTimeVictoryOrRuin * stdTimeVictoryOrRuin
31
32print "Average Time to Victory or Ruin:", avgTimeVictoryOrRuin
33print "Variance of Time to Victory or Ruin:",  varTimeVictoryOrRuin

__________________________________________________________________________ ### Problems to Work for Understanding

1. How many heads were obtained in your sequence of 100 ﬂips? What is the class average of the number of heads in 100 ﬂips? What is the variance of the number of heads obtained by the class members?
2. What is the net win-loss total obtained in your sequence of 100 ﬂips? What is the class average of the net win-loss total in 100 ﬂips? What is the variance of the net win-loss total obtained by the class members?
3. How many of the class reached a net gain of $+10$ (call it victory) before reaching a net loss of $-10$ (call it ruin)?
4. How many ﬂips did it take before you reached a net gain of $+10$ (victory) or a net loss of $-10$ (ruin)? What is the class average of the number of ﬂips before reaching victory before ruin at $+10$ or a net loss of $-10$?
5. What is the maximum net value achieved in your sequence of ﬂips? What is the class distribution of maximum values achieved in the sequence of ﬂips?
6. Perform some simulations of the coin-ﬂipping game, varying $p$. How does the value of $p$ aﬀect the experimental probability of victory and ruin?

__________________________________________________________________________ ### References

   William Feller. An Introduction to Probability Theory and Its Applications, Volume I, volume I. John Wiley and Sons, third edition, 1973. QA 273 F3712.

   Emmanuel Lesigne. Heads or Tails: An Introduction to Limit Theorems in Probability, volume 28 of Student Mathematical Library. American Mathematical Society, 2005.

__________________________________________________________________________ 1. Virtual Laboratories in Probability and Statistics Red and Black Game. Search for Red and Black Game, both the explanation and the simulation.
2. University of California, San Diego, Department of Mathematics, A.M. Garsia. A java applet that simulates how long it takes for a gambler to go broke. You can control how much money you and the casino start with, the house odds, and the maximum number of games. Results are a graph and a summary table. Submitted by Matt Odell, September 8, 2003.

__________________________________________________________________________

I check all the information on each page for correctness and typographical errors. Nevertheless, some errors may occur and I would be grateful if you would alert me to such errors. I make every reasonable eﬀort to present current and accurate information for public use, however I do not guarantee the accuracy or timeliness of information on this website. Your use of the information from this website is strictly voluntary and at your risk.

I have checked the links to external sites for usefulness. Links to external websites are provided as a convenience. I do not endorse, control, monitor, or guarantee the information contained in any external website. I don’t guarantee that the links are active at all times. Use the links here with the same caution as you would all information on the Internet. This website reﬂects the thoughts, interests and opinions of its author. They do not explicitly represent oﬃcial positions or policies of my employer.

Information on this website is subject to change without notice.