Проверка честности монеты

9

Мне задал следующий вопрос друг. Я не мог ей помочь, но я надеюсь, что кто-нибудь сможет мне это объяснить. Я не мог найти подобный пример. Спасибо за любую помощь и объяснение.

Q: Результаты 100 экспериментов с бросанием монет записываются как 0 = "Хвост" и 1 = "Голова". Выход x представляет собой строку из 0 и 1 длиной 100. И вычисляется количество раз, когда мы получаем 1-0-0 в x, и оно равно 20 (например, если x = (001001110100), 1-0-0 встречается 2 раза). Как вы думаете, это честная монета?

Джимми Дур
источник
1
Этот вопрос не звучит как настоящая реальная проблема. Это домашнее задание?
Секст Эмпирик
1
Я не уверена. Как я указал, это спросил меня друг. Я не мог помочь ей, но я все еще хочу узнать, как решить или ответить на этот вопрос. @ MartijnWeterings
Джимми Дур
2
stats.stackexchange.com/questions/158490/… содержит довольно полное описание ситуации. Более подробно см stats.stackexchange.com/...
whuber
1
@JimmyDur, ты не имеешь представления о толковании или значении вопроса. Например. Вы формулируете вопрос как «Вы думаете, что это честная монета?». Это звучит как вопрос с подвохом.
Секст Эмпирик
1
... Однако, с определенных точек зрения, это может быть неправильным способом решения этой проблемы, и может возникнуть желание использовать байесовский подход, например, если кто-то знает что-то вроде распределения вероятности справедливости монет до этого. Без каких-либо знаний о происхождении и обстоятельствах любое вычисление будет просто математическим упражнением, а не ответом на ваш явный вопрос «Как вы думаете, это честная монета?».
Секст Эмпирик

Ответы:

14

Решение проблемы с помощью симуляции

Моей первой попыткой было бы смоделировать это на компьютере, который может очень быстро перевернуть много монет. Ниже приведен пример с одним миллионом испытаний. Событие «то, что число раз Икс когда шаблон« 1-0-0 »происходит при Nзнак равно100 монетных бросках, составляет 20 или более», происходит примерно раз в три тысячи испытаний, поэтому то, что вы наблюдали, не очень вероятно (для справедливого монета).

Обратите внимание, что гистограмма предназначена для моделирования, а линия - это точное вычисление, объясненное ниже.

гистограмма

set.seed(1)

# number of trials
n <- 10^6

# flip coins
q <- matrix(rbinom(100*n, 1, 0.5),n)

# function to compute number of 100 patterns
npattern <- function(x) {
  sum((1-x[-c(99,100)])*(1-x[-c(1,100)])*x[-c(1,2)])
}

# apply function on data 
counts <- sapply(1:n, function(x) npattern(q[x,]))
hist(counts, freq = 0) 

# estimated probability
sum(counts>=20)/10^6
10^6/sum(counts>=20)

Решение проблемы с точным вычислением

Для аналитического подхода вы можете использовать тот факт, что «вероятность наблюдать 20 или более последовательностей« 1-0-0 »за 100 бросков монеты равна 1 минус вероятность того, что для получения 20 последовательностей потребуется более 100 бросков». , Это решается в следующих шагах:

Время ожидания вероятности переключения «1-0-0»

Распределение, еN,Иксзнак равно1(N) , количества раз, которое вам нужно перевернуть, пока вы не получите ровно одну последовательность '1-0-0', может быть вычислено следующим образом:

Давайте проанализируем способы получения «1-0-0» в виде цепочки Маркова. Мы следуем состояниям, описанным суффиксом строки переворотов: «1», «1-0» или «1-0-0». Например, если у вас есть следующие восемь сальто 10101100, то вы прошли по порядку следующие восемь состояний: «1», «1-0», «1», «1-0», «1», «1», «1-0», «1-0-0» и потребовалось восемь сальто, чтобы достичь «1-0-0». Обратите внимание , что вы не имеют равную вероятность достичь состояния «1-0-0» в каждом флип. Таким образом, вы не можете смоделировать это как биномиальное распределение . Вместо этого вы должны следовать дереву вероятностей. Состояние «1» может переходить в «1» и «1-0», состояние «1-0» может переходить в «1» и «1-0-0», и состояние «1-0-0» является поглощающим состоянием. Вы можете записать это как:

           number of flips
           1   2   3   4   5   6   7   8   9   ....   n

'1'        1   1   2   3   5   8  13  21  34   ....   F_n
'1-0'      0   1   1   2   3   5   8  13  21          F_{n-1}
'1-0-0'    0   0   1   2   4   7   12 20  33          sum_{x=1}^{n-2} F_{x}

и вероятность достижения паттерна «1-0-0» после броска первого «1» (вы начинаете с состояния «0», еще не перевернув головы), в течение N переворотов в полтора раза превышает вероятность быть в состоянии '1-0' в течение N-1 сальто:

еNс,Иксзнак равно1(N)знак равноFN-22N-1

где Fя - это я ое число Фибоначчи. Безусловная вероятность - это сумма

еN,Иксзнак равно1(N)знак равноΣКзнак равно1N-20,5КеNс,Иксзнак равно1(1+(N-К))знак равно0,5NΣКзнак равно1N-2FК

Время ожидания вероятности переворачивания К раз '1-0-0'

Это вы можете вычислить путем свертки.

еN,Иксзнак равноК(N)знак равноΣLзнак равно1NеN,Иксзнак равно1(L)еN,Иксзнак равно1(N-L)

Вы получите вероятность наблюдать 20 или более паттернов «1-0-0» (основываясь на гипотезе, что монета справедлива)

> # exact computation
> 1-Fx[20]
[1] 0.0003247105
> # estimated from simulation
> sum(counts>=20)/10^6
[1] 0.000337

Вот R-код для его вычисления:

# fibonacci numbers
fn <- c(1,1)
for (i in 3:99) {
  fn <- c(fn,fn[i-1]+fn[i-2])
}

# matrix to contain the probabilities
ps <- matrix(rep(0,101*33),33)

# waiting time probabilities to flip one pattern
ps[1,] <- c(0,0,cumsum(fn))/2^(c(1:101))

#convoluting to get the others
for (i in 2:33) {
  for (n in 3:101) {
     for (l in c(1:(n-2))) {
       ps[i,n] = ps[i,n] + ps[1,l]*ps[i-1,n-l]
     }  
  }
}

# cumulative probabilities to get x patterns in n flips
Fx <- 1-rowSums(ps[,1:100])

# probabilities to get x patterns in n flips
fx <- Fx[-1]-Fx[-33]

#plot in the previous histogram
lines(c(1:32)-0.5,fx)

Вычисление за нечестные монеты

ИксNп

Теперь мы используем обобщение чисел Фибоначчи:

FN(Икс)знак равно{1если Nзнак равно1Иксесли Nзнак равно2Икс(FN-1+FN-2)если N>2

вероятности теперь как:

еNс,Иксзнак равно1,п(N)знак равно(1-п)N-1FN-2((1-п)-1-1)

а также

еN,Иксзнак равно1,п(N)знак равноΣКзнак равно1N-2п(1-п)К-1еNс,Иксзнак равно1,п(1+N-К)знак равноп(1-п)N-1ΣКзнак равно1N-2FК((1-п)-1-1)

Когда мы строим это, вы получаете:

разные р

пзнак равно0,5пзнак равно0,33

Таким образом, если вы думали перед экспериментом, что монета вряд ли несправедлива, то теперь вы все равно должны думать, что монета вряд ли несправедлива.


пчасеadsзнак равнопTaяLs

Можно было бы намного проще проверить вероятность получения справедливой монеты, подсчитав количество голов и хвостов и используя биномиальное распределение, чтобы смоделировать эти наблюдения и проверить, является ли наблюдение конкретным или нет.

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

п

коррелированная монета

пзнак равно0,45пп

# number of trials
set.seed(1)
n <- 10^6

p <- seq(0.3,0.6,0.02)
np <- length(p)
mcounts <- matrix(rep(0,33*np),33)

pb <- txtProgressBar(title = "progress bar", min = 0,
                     max = np, style=3)
for (i in 1:np) {
  # flip first coins
  qfirst <- matrix(rbinom(n, 1, 0.5),n)*2-1
  # flip the changes of the sign of the coin
  qrest <- matrix(rbinom(99*n, 1, p[i]),n)*2-1
  # determining the sign of the coins
  qprod <- t(sapply(1:n, function(x) qfirst[x]*cumprod(qrest[x,])))
  # representing in terms of 1s and 0s
  qcoins <- cbind(qfirst,qprod)*0.5+0.5
  counts <- sapply(1:n, function(x) npattern(qcoins[x,]))

  mcounts[,i] <- sapply(1:33, function(x) sum(counts==x))
  setTxtProgressBar(pb, i)
}
close(pb)

plot(p,colSums(mcounts[c(20:33),]),
     type="l", xlab="p same flip", ylab="counts/million trials", 
     main="observation of 20 or more times '1-0-0' pattern \n for coin with correlated flips")
points(p,colSums(mcounts[c(20:33),]))

Использование математики в статистике

Выше все хорошо, но это не прямой ответ на вопрос

"Как вы думаете, это честная монета?"

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

Один открытый вопрос: почему и как мы ищем модель «1-0-0».

  • Например, возможно, эта схема не была целью, о чем было решено до начала расследования. Может быть, это просто что-то, что «выделяется» в данных, и это было то, что привлекло внимание после эксперимента. В этом случае нужно учитывать, что эффективно проводить множественные сравнения .
  • Другая проблема заключается в том, что вычисленная выше вероятность представляет собой p-значение. Значение p-значения должно быть тщательно продумано. Это не вероятность, что монета справедлива. Это, скорее, вероятность наблюдать конкретный результат, если монета справедлива. Если у кого-то есть среда, в которой он знает какое-то распределение справедливости монет, или можно сделать разумное предположение, то можно принять это во внимание и использовать байесовское выражение .
  • Что справедливо, что несправедливо. В конце концов, при достаточном количестве испытаний можно обнаружить небольшую несправедливость. Но уместно ли это и не является ли такой поиск необъективным? Когда мы придерживаемся подхода, основанного на частоте, то следует описать нечто вроде границы, над которой мы считаем монету справедливой (некоторый соответствующий размер эффекта). Тогда можно было бы использовать что-то похожее на двухсторонний t-критерий , чтобы решить, является ли монета честной или нет (относительно шаблона «1-0-0»).
Секст Эмпирик
источник
1
Почему бы не искать решение в закрытой форме через MLE?
Digio
п^MLзнак равноaргмaИксп[п(Икс1,,,,,ИксN|п,N)]Икс~ВяNомяaL(N,п)п^ML
Digio
пп
Иксзнак равно{0,0,1,0,0,1,1,1,0,1,0,0}п^знак равноИксNзнак равно512знак равно0,421,96п^(1-п^)Nзнак равно0,280,42+0,28<0,50
3
@FrankHarrell также нет ничего в настройке проблемы, что заставляет нас думать, что сальто некоррелированы. Проблема установки относительно плоха с информацией. Мой пример, который касается корреляции сальто, просто для того, чтобы охватить суть вопроса. Я не о том , что это есть способ ответить на него. Но скажем, что кто-то (возможно, ОП) проводит исследования последовательностей ДНК или других проблемных ситуаций, где возможность корреляции имеет больше смысла, тогда у них есть пример того, как это может получиться.
Секст Эмпирик