In these probabilities, + represents a positive test, - represents a negative test, D=0 indicates no disease, and D=1 indicates the disease is present.
Probability of having the disease given a positive test: \(Pr(D=1∣+)\) 99% test accuracy when disease is present: \(Pr(+∣D=1)=0.99\) 99% test accuracy when disease is absent: \(Pr(−∣D=0)=0.99\) Rate of cystic fibrosis: \(Pr(D=1)=0.00025\) Bayes' theorem can be applied like this:
\[Pr(D=1∣+)=\frac{Pr(+∣D=1)⋅Pr(D=1)}{Pr(+)}\] \[Pr(D=1∣+)=\frac{Pr(+∣D=1)⋅Pr(D=1)}{Pr(+∣D=1)⋅Pr(D=1)+Pr(+∣D=0)⋅Pr(D=0)}\]
Substituting known values, we obtain: \[Pr(D=1∣+)=\frac{0.99⋅0.00025}{0.99⋅0.00025+0.01⋅0.99975}=0.02\]
Code: Monte Carlo simulation
prev <- 0.00025 # disease prevalence
N <- 100000 # number of tests
outcome <- sample(c("Disease", "Healthy"), N, replace = TRUE, prob = c(prev, 1-prev))
N_D <- sum(outcome == "Disease") # number with disease
N_H <- sum(outcome == "Healthy") # number healthy
# for each person, randomly determine if test is + or -
accuracy <- 0.99
test <- vector("character", N)
test[outcome == "Disease"] <- sample(c("+", "-"), N_D, replace=TRUE, prob = c(accuracy, 1-accuracy))
test[outcome == "Healthy"] <- sample(c("-", "+"), N_H, replace=TRUE, prob = c(accuracy, 1-accuracy))
table(outcome, test)
#> test
#> outcome - +
#> Disease 0 23
#> Healthy 99012 965
N <- 100000 # number of tests
outcome <- sample(c("Disease", "Healthy"), N, replace = TRUE, prob = c(prev, 1-prev))
N_D <- sum(outcome == "Disease") # number with disease
N_H <- sum(outcome == "Healthy") # number healthy
# for each person, randomly determine if test is + or -
accuracy <- 0.99
test <- vector("character", N)
test[outcome == "Disease"] <- sample(c("+", "-"), N_D, replace=TRUE, prob = c(accuracy, 1-accuracy))
test[outcome == "Healthy"] <- sample(c("-", "+"), N_H, replace=TRUE, prob = c(accuracy, 1-accuracy))
table(outcome, test)
#> test
#> outcome - +
#> Disease 0 23
#> Healthy 99012 965
Page last modified on February 14, 2021, at 07:51 AM
Powered by
PmWiki