Home > Mobile >  Simulation in R for probability puzzle
Simulation in R for probability puzzle

Time:09-25

Let two integers such that 1<a<7 and 2<b<5.We play the following game: We start with a euros and we roll a dice.If the result of the dice is at least b we earn 1 euro and we lose 1 euro otherwise. The game ends if we reach 8 euros (winners) or in 0 (losers). We prefer to be ?

A) Rich: we start with a equals 6 euros , but unlucky because b equals 5. B) Mediocre: we start with a equals 4 and b equals 4 C) Lucky: we start with a equals 2 but b equals 3.

How can I do it in R using simulation in order to decide If I want to be rich ,mediocre or lucky?

My effort


gamble <- function(a,n,p) {
  stake <- a
  while (stake > 0 & stake < n) {
    bet <- sample(c(-1,1),1,prob=c(1-p,p))
    stake <- stake   bet
  }
  if (stake == 0) return(1) else return(0)
}   

a <- 6
n <- 8
p <- 1/3
trials <- 100000
simlist <- replicate(trials, gamble(a, n, p))
mean(simlist) # Estimate of probability that gambler is ruined

a <- 4
n <- 8
p <- 1/2
trials <- 100000
simlist <- replicate(trials, gamble(a, n, p))
mean(simlist) 


a <- 2
n <- 8
p <- 2/3
trials <- 100000
simlist <- replicate(trials, gamble(a, n, p))
mean(simlist) 


CodePudding user response:

Here's an extremely over-engineered solution!

Your loop logic looks fine to me, the main issue is that you were originally rolling a die and just adding that many euros to your stake, and your function wasn't returning anything.

The key line below is euros <- euros (sample(1:6, 1) >= roll_to_win)*2 - 1 . This does the following:

  • Gets a random number between 1 and 6,
  • Checks to see if it's bigger than or equal to roll_to_win, which corresponds to b,
  • Converts a TRUE result to 1 and a FALSE result to -1 (with trickery! TRUE and FALSE behave like 1 and 0, so if you multiply TRUE by 2 and subtract 1 you get 1, and you get -1 for FALSE.
library(tidyverse)

# set up our gambling function. We've hard-coded that it's 8 euros to win.
gamble <- function(start_euros,
                   roll_to_win){
  
  euros <- start_euros
  win_euros <- 8
  
  while(euros > 0 & euros < win_euros){
    # for debugging, print euros to console
    #message(euros)
    
    # roll a 6-sided die, see if the result is
    # bigger than the value we need to win, and
    # add or subtract 1 from euros based on result
    euros <- euros   (sample(1:6, 1) >= roll_to_win)*2 - 1    
    
  }
  
  # return 0 for loss, 1 for win
  return (euros / win_euros)
}

# set up the number of reps we'll do
reps <- 10000

# build a nested tibble with our values for a and b and generate trials
results <- tribble(~a,~ b, 
                   6, 5,
                   4, 4,
                   2, 3) %>%
  mutate(trials = purrr::map2(a,b, function(x, y) replicate(reps, gamble(x,y))))

# get observed frequency by taking mean 
results %>%
  mutate(prob = purrr::map_dbl(trials, mean))

And here's what it looks like for me:

# A tibble: 3 x 4
      a     b trials          prob
  <dbl> <dbl> <list>         <dbl>
1     6     5 <dbl [10,000]> 0.240
2     4     4 <dbl [10,000]> 0.498
3     2     3 <dbl [10,000]> 0.751

No clue if that's right or not though :)

  • Related