Есть ли способ определить точную долю из большой выборки модели Монте-Карло, используя R?

Question

В прошедшие выходные я участвовал в ежегодной MIT Puzzle Hunt. Эта головоломка имела 12 вероятностных проблем, которые нам пришлось решить. Хотя мы не успели завершить его вовремя, я решил, что хочу попытаться решить это сам (так что никаких спойлеров).

Задача 8 в настоящее время заставляет меня споткнуться. В нем говорится: Алиса получила 16 монет по советам, а Боб - 3. Они решили сыграть в следующую игру. Каждый игрок выбирает одну из своих монет. Они подбрасывают свои монеты, пока одна не станет головой, а другая - хвостом. Человек с головами держит обе монеты. Это повторяется до тех пор, пока один игрок не получит все монеты. Какова вероятность того, что Алиса потеряет свои советы Бобу?

Я не силен в статистике, и я уверен, что есть способ рассчитать это, используя формулу (любое руководство здесь также быть оцененным), но я понял, что потенциально могу рассчитать ответ, имитируя результаты в R.

Вот код, который я запустил для симуляции (предупреждение, я новичок в R, поэтому извиняюсь за грубость и плохое форматирование) :

coin_battle <- function(alice_points,bob_points){
Alice <- alice_points
Bob <- bob_points
    while(Alice !=0 & Bob !=0){
        result <- sample(1:2,1)
        if(result == 1){
            Alice <- Alice + 1
            Bob <- Bob - 1
        }
        else{
            Alice <- Alice - 1
            Bob <- Bob + 1
        }
    }
if(Alice==0){
    return("Bob")
}
if(Bob==0){
    return("Alice")
}
}
simulate <- function(sample_size){
n <- sample_size
res <- rep(0,n)
    for(i in 1:n){
        res[i] <- coin_battle(16,3)
    }
    return(res)
}
vector1 <- simulate(10000)
table(vector1)

Я только что запустил этот код с размером выборки 10 000 000, и в результате Алиса выиграла игру 8 420 661 раз, а Боб выиграл 1 579 339 раз, поэтому Боб выиграет 15,79339% от время. У меня был друг, подтвердивший, что этот десятичный ответ правильный. Однако для целей головоломки ответ должен быть в форме дроби. Я подозреваю, что как числитель, так и знаменатель уменьшенной формы будут составлять 26 или меньше, поскольку загадка, вероятно, будет попросите меня перевести дробь в буквы. Я мог бы, вероятно, угадать и проверить потенциальные варианты, используя числа от 1 до 26, но я хотел бы знать, есть ли более обобщенное решение для этого.

Спасибо!

Answer 1

Я не могу придумать прямой способ оценить точную долю. Но у вас есть большая симуляция, поэтому вы можете установить доверительный интервал вокруг истинной вероятности. Отсюда можно сделать предположение, что знаменатель меньше некоторого произвольного значения, и проверить все дроби, которые соответствуют этим критериям, чтобы увидеть, попадают ли они в доверительный интервал. Я написал код для этого, потому что этот вопрос заинтриговал меня. Будьте осторожны - это очень плохой код, потому что я немного болен и плохо соображаю.

ps <- prop.test(1579339,10000000)

# arbitrarily assume the denominator is <= 100
denominators <- 1:100
# numerators
num <- lapply(denominators, function(x){sequence(x)})
# the decimal values
vals <- lapply(denominators, function(x){sequence(x)/x})

# helper function to pluck the values that fall between the interval
between <- function(x, a, b){
  a <= x & x <= b
}

# find the values that fall in the interval
btw <- lapply(vals, between, a = ps$conf.int[1], b = ps$conf.int[2])
# not too proud of this.  it's ugly and there is a better way
output <- data.frame(percs = unlist(vals),
                     num = unlist(num), 
                     den = rep(denominators, denominators), 
                     between = unlist(btw))
# only grab the output that meets our criteria
possible <- output[output$between,]
possible

, что дает нам

> possible
         percs num den between
174  0.1578947   3  19    TRUE
709  0.1578947   6  38    TRUE
1605 0.1578947   9  57    TRUE
2862 0.1578947  12  76    TRUE
4480 0.1578947  15  95    TRUE

, и все эти цифры уменьшаются до 3 /19.