时间戳

首页 > 解决方案 > r中的条件概率

r - r中的条件概率

问题描述

问题：

一项影响 0.05% 男性人口的疾病筛查测试能够在 90% 的个体实际患有该疾病的病例中识别该疾病。然而，该测试会产生 1% 的假阳性（当个体没有患病时给出阳性读数）。找出一个人患上检测呈阳性的疾病的概率。然后，在一个人的测试结果为阴性的情况下，找出他患有该疾病的概率。

我的错误尝试：

我首先让： • T 是男性检测呈阳性的事件 • Tc 是男性检测呈阴性的事件 • 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])

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

标签： rprobability

解决方案

思考这样一个条件概率问题的最佳方法可能是举一个具体的例子。

假设我们测试了人口中的 100 万人。然后预计将有 500 人（100 万人的 0.05%）患有这种疾病，其中 450 人预计检测呈阳性，50 人检测呈阴性（因为假阴性率为 10%）。

相反，预计将有 999,500 人没有这种疾病（100 万人减去 500 人确实患有这种疾病），但由于其中 1% 的人检测呈阳性，因此我们预计有 9,995 人（999,500 人的 1%）出现假阳性结果.

因此，鉴于随机抽取的阳性检测结果，它要么属于检测呈阳性的 450 名患病者之一，要么属于检测呈阳性的 9,995 名未患病者之一——我们不知道是哪一个。

这是第一个问题的情况，因为我们有一个阳性测试结果，但不知道它是真阳性还是假阳性。鉴于他们的阳性测试，我们的受试者患有该疾病的概率是他们是 10,445 名阳性结果（9995 名假阳性 + 450 名真阳性）中的 450 名真阳性之一的概率。这归结为简单的计算 450/10,445 或 0.043，即 4.3%。

类似地，随机进行的阴性测试要么属于 989505 (999500 - 9995) 名无病检测呈阴性的人之一，要么属于 50 名检测呈阴性的有病患者之一，因此患该病的概率为 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))

现在我们已经运行了模拟，我们可以通过创建列联表来计算真阳性、假阳性、真阴性和假阴性

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 包于 2021-08-19 创建 (v2.0.0 )

推荐阅读