r 中的条件概率

Conditional probability in r

问题:

一种影响 0.05% 男性人口的疾病的筛查测试能够在 90% 的个体实际患有该病的情况下识别出该病。然而,该测试会产生 1% 的假阳性(当个体没有患病时给出阳性读数)。求一个人在给定的测试呈阳性的情况下患有这种疾病的概率。然后,求一个人在检测结果为阴性的情况下患病的概率。

我的错误尝试:

我首先让: • 如果男性检测呈阳性 • 可能是一个人的测试结果呈阴性 • D 是一个人实际上患有这种疾病的事件 • Dc 是男性没有患病的事件

因此我们需要找到 P(D|T) 和 P(D|Tc)

然后我写了这段代码:

set.seed(110)
sims = 1000

D = rep(0, sims)
Dc = rep(0, sims)
T = rep(0, sims)
Tc = rep(0, sims)

# run the loop
for(i in 1:sims){
  
  # flip to see if we have the disease
  flip = runif(1)
  
  # if we got the disease, mark it
  if(flip <= .0005){
    D[i] = 1
  }
  
  # if we have the disease, we need to flip for T and Tc, 
  if(D[i] == 1){
    
    # flip for S1
    flip1 = runif(1)
    
    # see if we got S1
    if(flip1 < 1/9){
      T[i] = 1
    }
    
    # flip for S2
    flip2 = runif(1)
    
    # see if we got S1
    if(flip2 < 1/10){
      Tc[i] = 1
    }
  }
}


# P(D|T)
mean(D[T == 1])

# P(D|Tc)
mean(D[Tc == 1])

我真的很挣扎,所以任何帮助将不胜感激!

也许思考像这样的条件概率问题的最佳方式是举一个具体的例子。

假设我们测试了人口中的一百万个人。那么预计将有 500 人(一百万的 0.05%)患有这种疾病,其中 450 人预计检测呈阳性,50 人检测呈阴性(因为假阴性率为 10%)。

相反,999,500 人预计不会患病(100 万人减去 确实 患病的 500 人),但由于其中 1% 的人检测呈阳性,那么我们预计会有 9,995 人(999,500 人中的 1%)出现假阳性结果。

因此,鉴于随机取得的阳性检测结果,它要么属于检测呈阳性的 450 名患病者中的一位,要么属于检测呈阳性的 9,995 名未患病者中的一位 - 我们不知道

这是第一个问题的情况,因为我们有一个阳性检测结果,但不知道是真阳性还是假阳性。我们的受试者给定他们的阳性测试的概率是他们是 10,445 个阳性结果中的 450 个真阳性之一的概率(9995 个假阳性 + 450 个真阳性).这归结为简单计算450/10,445即0.043,也就是4.3%。

同样,随机进行的阴性测试要么属于 989505 (999500 - 9995) 人没有 测试呈阴性的人之一,要么属于 50 人with 测试阴性的疾病,因此患有该疾病的概率为 50/989505,即 0.005%。

我认为这个问题证明了在解释测试结果时需要考虑疾病流行率的重要性,而与编程或 R 关系不大。它只需要一个计算器(最多)。

如果你真的想 运行 在 R 中进行模拟,你可以这样做:

set.seed(1) # This makes the sample reproducible

sample_size <- 1000000 # This can be changed to get a larger or smaller sample

# Create a large sample of 1 million "people", using a 1 to denote disease and
# a 0 to denote no disease, with probabilities of 0.0005 (which is 0.05%) and
# 0.9995 (which is 99.95%) respectively.
disease <- sample(x = c(0, 1), 
                  size = sample_size, 
                  replace = TRUE, 
                  prob = c(0.9995, 0.0005))

# Create an empty vector to hold the test results for each person
test <- numeric(sample_size)

# Simulate the test results of people with the disease, using a 1 to denote
# a positive test and 0 to denote a negative test. This uses a probability of
# 0.9 (which is 90%) of having a positive test and 0.1 (which is 10%) of having
# a negative test. We draw as many samples as we have people with the disease
# and put them into the "test" vector at the locations corresponding to the
# people with the disease.
test[disease == 1] <- sample(x = c(0, 1), 
                             size = sum(disease), 
                             replace = TRUE, 
                             prob = c(0.1, 0.9))

# Now we do the same for people without the disease, simulating their test
# results, with a 1% probability of a positive test.
test[disease == 0] <- sample(x = c(0, 1), 
                             size = 1e6 - sum(disease), 
                             replace = TRUE, 
                             prob = c(0.99, 0.01))

现在我们有了 运行 我们的模拟,我们可以通过创建偶然事件 table

来计算真阳性、假阳性、真阴性和假阴性
contingency_table <- table(disease, test)

contingency_table
#>        test
#> disease      0      1
#>       0 989566   9976
#>       1     38    420

并得到这种疾病的大致概率,给出这样的阳性测试:

contingency_table[2, 2] / sum(contingency_table[,2])
#> [1] 0.04040015

以及像这样得到阴性测试的疾病概率:

contingency_table[2, 1] / sum(contingency_table[,1])
#> [1] 3.83992e-05

您会注意到,由于某些抽样概率非常小,抽样的概率估计并不那么准确。您可以模拟更大的样本,但您的计算机可能需要一段时间才能 运行 它。

reprex package (v2.0.0)

于 2021-08-19 创建

扩展 Allan 的答案,但如果您愿意,可以将其与贝叶斯定理联系起来:

从题目中,你知道(将百分比转化为概率):

插入: