自举大小和绘制不同观察的概率

机器算法验证可能性引导程序

2022-03-12 10:38:40

假设我们从 n 个编号的球中抽取 k 个样本而，并得到个（唯一）球的样本。然后，我们将罐子里的所有球都放回原处，并再次重复绘制以再次获得个（唯一）球。 $k$ $n$ $N_1$ $N_2$

我们对最终只被抽出一次的球感兴趣（在两个实验过程中）。很明显，对于其中一个球来说，事件“恰好被抽到一次”的概率是

2 k / n (1 - k / n)

$2k/n(1-k/n)$

在处最大化。 $k^*=n/2$

现在考虑相同的实验，但有替换，我们再次寻找的值（现在表示为），它使事件“恰好绘制一次”的频率最大化。现在，每次平局中不同球的预期数量是，所以我预计。运行一些模拟，我发现的实际位置似乎要低得多，可能在和之间。我很难理解我的错误在哪里。 $k$ $\bar{k}^*$ $1-e^{-1}\approx0.632n$ $\bar{k}^*=k^*/(1-e^{-1})\approx0.79n$ $\bar{k}^*$ $0.65n$ $0.7n$

更一般地说（第二个问题）我想知道在第二个实验（有替换的那个）中恰好绘制一次的球的比例作为的函数的分布函数是什么？通过数值实验，我得到以下曲线（）： $k$ $n=100$

在此处输入图像描述

编辑

这是重现上述示例的 R 代码：

fx02<-function(ll,n,k){
    a1<-matrix(0,n,2)
    a1[sample(1:n,k,replace=TRUE),1]<-1
    a1[sample(1:n,k,replace=TRUE),2]<-1
    sum(rowSums(a1)==1)/n
}

ss<-(1:60)*5   #The grid of values of k for which we'll compute the probability.
a4<-matrix(NA,length(ss),2)
for(i in 1:length(ss)){
    a2<-ss[i]
    a3<-c(lapply(1:1000,fx02,n=100,k=a2),recursive=TRUE);
    a4[i,]<-c(a2,mean(a3))
}
plot(a4,xlab="k",ylab="frequency of distinct draw")

1个回答

在看到您的模拟代码后进行编辑

我认为您没有考虑到某个数量也可以在一个样本中多次出现的事实。考虑到这一点，我们得到相同的结果。这是您修改后的代码：

fx02<-function(ll,n,k){
  a1<-matrix(0,n,2)
  samp1 <- sample(1:n, k, replace=TRUE)
  samp2 <- sample(1:n, k, replace=TRUE)

  a1[sort(unique(samp1)), 1] <- as.numeric(table(samp1))
  a1[sort(unique(samp2)), 2] <- as.numeric(table(samp2))

  sum(rowSums(a1)==1)/n
}

ss<-(1:60)*5   #The grid of values of k for which we'll compute the probability.
a4<-matrix(NA,length(ss),2)
for(i in 1:length(ss)){
  a2<-ss[i]
  a3<-c(lapply(1:1000,fx02,n=100,k=a2),recursive=TRUE);
  a4[i,]<-c(a2,mean(a3))
}
plot(a4,xlab="k",ylab="frequency of distinct draw", pch=16, las=1)
abline(h=0.3697296, v=50)

max(a4[,2])
[1] 0.36971

user603 代码图

原始答案

对于实验2（有放回），我认为某个数字恰好画一次的概率是：或

P_{once} = 2 k (n - 1)^{(2 k - 1)} (1 / n)^{2 k}

$P_{\text{once}}=2k(n-1)^{(2k-1)}(1/n)^{2k}$

P_{once} = 2 k (n - 1)^{(2 k - 1)} n^{- 2 k}

$P_{\text{once}}=2k(n-1)^{(2k-1)}n^{-2k}$

这可以通过一个简单的例子来检查，其中和。在这种情况下，恰好抽取某个数字一次的概率是。 $n=3$ $k=2$ $2(4/9)^{2}\approx0.395$

上述公式的最大值出现在：或其中被四舍五入到最接近的整数。正如@user603 在评论中已经提到的，对于较大的来说，这大约是对于，所以 50。和的最大概率将在左右（正如评论中已经计算的那样@user603）。

k_{max} = ⌊ - \frac{1}{[2 (\log (n + 1) + \log (1 / n))]} ⌉

$k_{\text{max}}=\lfloor-\frac{1}{[2(\log(n + 1) + \log(1/n))]}\rceil$

k_{max} = ⌊ - \frac{1}{[2 (\log (n + 1) - \log (n))]} ⌉

$k_{\text{max}}=\lfloor-\frac{1}{[2(\log(n + 1) - \log(n))]}\rceil$

k_{max}

$k_{\text{max}}$

n / 2

$n/2$

n

$n$

n = 100

$n=100$

k_{max} \approx 49.75

$k_{\text{max}}\approx 49.75$

n = 100

$n=100$

k_{max} = 50

$k_{\text{max}}=50$

0.3697

$0.3697$

我设置了一个模拟来检查这个结果R：

prob.once <- vector()

draw.once <- function(n, k, sim=10000, repl=TRUE){

  for ( i in 1:sim ) {

    samp1 <- sample(1:n, size=k, replace=repl)
    samp2 <- sample(1:n, size=k, replace=repl)

    if ((is.element(1, samp1) & !is.element(1, samp2) & !is.element(1, samp1[duplicated(samp1)])) |
          (!is.element(1, samp1) & is.element(1, samp2) & !is.element(1, samp2[duplicated(samp2)])) ){
      prob.once[i] <- 1

    } else {
      prob.once[i] <- 0
    }    
  }      
  mean(prob.once)   
}

krepl <- 1:300
probs.repl <- sapply(krepl, FUN=draw.once, n=100, sim=20000, repl=TRUE)

plot(probs.repl~krepl, pch=16, type="p", lwd=2, las=1, ylab="Probability", xlab="k", col="steelblue")

abline(h=0.3697296)
abline(v=50)

模拟图

模拟结果似乎证实了上述考虑。

其它你可能感兴趣的问题

上一篇如何组合多个估算数据集？下一篇在社会科学中何时以及如何使用权重进行序列分析？