两枚神奇的硬币

KeyTo9_Fans

题目来自正在进行的比赛：
http://acm.sgu.ru/problem.php?contest=35&problem=B

开始时间：12月20日17点0分，结束时间：12月20日22点0分。

题目大意：

有两枚硬币，正面朝上的概率分别是p和q。

我们都不知道p和q的具体值，所以各抛了n1次和n2次。

它们正面朝上的次数分别是m1和m2。

根据上述结果求p>q的概率。

此题是KeyTo9_Fans喜欢的类型，所以发到论坛上，大家一起做做看。

mathe

这个结果应该同对p,q的先验分布有关系吧。当然这里我们可以假设p,q的先验分布都是[0,1]上均匀分布，于是就可以有答案了。但是我觉得这个题目不好，实际上按照我们的实际经验，p,q的先验分布信息应该都是非常集中到1/2附近的

KeyTo9_Fans

mathe说得对，先验分布确实很关键。

由于发帖太急躁，我疏忽了这个条件：

Random values p and q are both uniformly distributed on [0,1] and are independent.

p和q是独立的，且均匀分布在[0,1]区间。

看来一道题目出得好不好，背景描述很关键。

如果描述成这样的场景：

-----------------------------------------------
甲乙两人投篮，命中率分别是p和q。

投篮之前，我们完全不知道p和q的具体值，所以认为p和q是独立的，且均匀分布在[0,1]区间。

后来甲$n_1$投$m_1$中，乙$n_2$投$m_2$中。

根据此结果求p>q的概率。
-----------------------------------------------

这样就容易理解了。因为这个描述更贴切，而硬币再怎么神奇也不会随心所欲地以任意的概率抛到正面向上。

shshsh_0510

对于$0<q<p$,在$[0,p]$对$C_{n_2}^{m_2}q^{m_2}(1-q)^{n_2-m_2}$以$q$为自变量积分，再乘$C_{n_1}^{m_1}p^{m_1}(1-p)^{n_1-m_1}$，然后以$p$为变量在$[0,1]$积分

KeyTo9_Fans

按照你的原式帮你写成Tex格式了。因为我觉得这样好看一点。

方法是对的。但为什么前面的系数是$C_{n_2}^{m_2}$和$C_{n_1}^{m_1}$呢？

mathe

这是组合数呀，这是古典概率问题。$n_2$次共有$m_2$次朝上的不同方案就是$C_{n_2}^{m_2}$,其中每次的概率为$q^{m_2}(1-q)^{n_2-m_2}$,
所以就可以直接计算积分
$\int_{0<q<p<1}\int C_{n_2}^{m_2}q^{m_2}(1-q)^{n_2-m_2}C_{n_1}^{m_1}p^{m_1}(1-p)^{n_1-m_1}dpdq$

shshsh_0510是直接分成两次积分。
不过结果还应该除以p和q都在0,1上的积分，也就是$\beta(m_2+1,n_2-m_2+1)\beta(m_1+1,n_1-m_1+1)C_{n_2}^{m_2}C_{n_1}^{m_1}=1/{(n_2+1)(n_1+1)}$

KeyTo9_Fans

按照楼上的方法，套了几个简单的数据，算出来的结果与测试结果完全吻合。

为什么不仅要乘$C_{n_1}^{m_1}$和$C_{n_2}^{m_2}$，还要乘上$(n_1+1)$和$(n_2+1)$才是最终结果呢？

wiley

因为是条件概率, 用了Bayes theorem.
需要计算的是Pr[ p>q | p的记录n1/m1, q的记录n2/m2], 所以分母上需要: Pr[p的记录n1/m1, q的记录n2/m2]对于任意的0<=p,q<=1

mathe

这个问题的关键是如何比较有效的计算
$I(n_1,m_1,n_2,m_2,p,q)=\int_{0<q<p<1}\int q^{m_2}(1-q)^{n_2-m_2} p^{m_1}(1-p)^{n_1-m_1}dpdq$
由于
$\int_0^p q^{m_2}(1-q)^{n_2-m_2}dq={m_2}/{n_2-m_2+1}\int_0^p q^{m_2-1}(1-q)^{n_2-m_2+1}dq -{p^{m_2}(1-p)^{n_2-m_2+1}}/{n_2-m_2+1}$
我们可以得到
$I(n_1,m_1,n_2,m_2,p,q)={m_2}/{n_2-m_2+1}I(n_1,m_1,n_2,m_2-1,p,q)-{beta(m_2+m_1+1,n_1+n_2-m_1-m_2+2)}/{n_2-m_2+1}$
通过这个方法，我们可以对单个变量$m_1$或$m_2$进行递归计算

KeyTo9_Fans

我昨天干坏事了，把这张帖子每分钟刷新1次，连续刷了半天，结果这张帖子成了一周内最热门的帖子……

还有那张northwolves的帖子也是。

这道题Fans在12月20日23点30分才完美解决，比比赛结束时间晚了1小时30分。

所用的方法极其暴力。

首先令$n_2=0$，$m_2=0$，考虑$n_1$和$m_1$的变化，把测试结果列成一张表。



规律是显而易见的。

然后令$n_2=1$，$m_2=0$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第2张表。



规律也是显而易见的。

然后令$n_2=1$，$m_2=1$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第3张表。



规律就更显而易见了。

然后令$n_2=2$，$m_2=0$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第4张表。



规律不怎么显而易见，不过肯定还是能找到的。

然后令$n_2=2$，$m_2=1$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第5张表。



规律也是不怎么显而易见，不过肯定还是能找到的，继续努力。

接着令$n_2=2$，$m_2=2$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第6张表。



这个规律显而易见多了。

接下来令$n_2=3$，$m_2=0$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第7张表。



如果前面的规律都找到了，现在这些表格的规律应该一眼就看出来了。

接下来令$n_2=3$，$m_2=1$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第8张表。



如果前面的规律都找到了，现在这些表格的规律应该一眼就看出来了。

继续令$n_2=3$，$m_2=2$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第9张表。



如果前面的规律都找到了，现在这些表格的规律应该一眼就看出来了。

最后令$n_2=3$，$m_2=3$，同样考虑$n_1$和$m_1$的变化，把测试结果列成第10张表。



如果前面的规律都找到了，现在这些表格的规律应该一眼就看出来了。

好了，看完这10张表，表里的规律和表与表之间的变化规律是不是都找出来了呢？

如果都找出来了，那么就可以用含有$n_1$、$m_1$、$n_2$、$m_2$的式子来表示p>q的概率了：

p>q的概率为$\frac{\sum_{i=0}^{m_1}C_{m_2+i}^{m_2}C_{n_2-m_2+n_1+1-i}^{n_2-m_2}}{\sum_{i=0}^{n_1+1}C_{m_2+i}^{m_2}C_{n_2-m_2+n_1+1-i}^{n_2-m_2}}$

如果预先打好了$C_n^m$的表，上述表达式的结果是很好求的。

代码如下：

1. 
   
2. 
   
3. #include<cstdio>
   
4. #include<memory>
   
5. #include<math.h>
   
6. 
   
7. double c[1024][1024],x,y,s,p,q;
   
8. int i,j,k,t,n1,n2,p2,m1,m2;
   
9. 
   
10. int main()
    
11. {
    
12.         for(i=0;i<1024;i++)
    
13.         {
    
14.                 c[i][0]=1e-153;
    
15.                 c[i][i]=1e-153;
    
16.                 for(j=1;j<i;j++)
    
17.                         c[i][j]=c[i-1][j]+c[i-1][j-1];
    
18.         }
    
19.         while(scanf("%d%d%d%d",&n1,&m1,&n2,&m2)>3)
    
20.         {
    
21.                 p2=n2-m2;
    
22.                 p=0;
    
23.                 q=0;
    
24.                 for(i=0;i<n1+2;i++)
    
25.                 {
    
26.                         s=c[m2+i][m2]*c[p2+n1+1-i][p2];
    
27.                         q+=s;
    
28.                         if(i<=m1)p+=s;
    
29.                 }
    
30.                 printf("%.9lf\n",p/q);
    
31.         }
    
32.         return 0;
    
33. }
    
34. 
    

复制代码

这是该代码编译出来的可执行程序，没有编程工具也可以运行。

Coins2_2.rar (38.9 KB, 下载次数: 2)

2009-12-24 10:48 上传

点击文件名下载附件

由于Fans大牛水平不足，不知道如何从积分式出发，推导出上述表达式。

所以这个工作只好交给更高级的大牛们来完成了。