BayesianTest

Overview and summary

Probability of texting:

• You are asked to compute the probability that the driver of a car is texting at a specific intersection.

• Nationally the cumulative probability that a driver is texting is:

P = 0.5, at x = 0.1

P = 0.75 at x = 0.3

• You observe cars at a location three times and note the number of texting drivers:

• 2 texting out of 20 drivers
• 4 texting out of 20 drivers

• 1 texting out of 20 drivers

– Given these data

• Compute the Beta prior, and report the coefficients

• Plot the prior, likelihood and posterior three times as you update your belief based on collecting more data

• Simulate the final posterior distribution and do the following:

1. Plot the posterior with the 90% HDI shown
2. Report the upper and lower limits of the 90% HDI
3. Of the next hundred drivers what are the number of texting drivers in the 90% HDI?
4. Are the drivers in this area better or worse that the national figures indicate?

Note: Following packages are required to run the below report.

• LearnBayes
• Knitr

Lets create a beta distribution, named beta.par.

beta.par <- beta.select(list(p=0.5, x=0.1), list(p=0.75, x=.3))

Q1)

• 2 texting out of 20 drivers
• 4 texting out of 20 drivers

• 1 texting out of 20 drivers

• Compute the Beta prior, and report the coefficients

• Plot the prior, likelihood and posterior three times as you update your belief based on collecting more data

par(mfrow = c(3,1))
beta.par + c(2, 18)

## [1]  2.41 19.73

triplot(beta.par, c(2, 18))

beta.par + c(4, 16)

## [1]  4.41 17.73

triplot(beta.par, c(4, 16))

beta.par + c(1, 19)

## [1]  1.41 20.73

triplot(beta.par, c(1, 19))

par(mfrow = c(1,1))

Q2)

Simulate the final posterior distribution and do the following:

1. Plot the posterior with the 90% HDI shown
2. Report the upper and lower limits of the 90% HDI
3. Of the next hundred drivers what are the number of texting drivers in the 90% HDI?

options(repr.plot.width=8, repr.plot.height=5)
beta.post.par <- beta.par + c(7, 53)
post.sample <- rbeta(100, beta.post.par[1], beta.post.par[2])
par(mfrow = c(1,2))
quants = quantile(post.sample, c(0.05, 0.95))
breaks = seq(min(post.sample), max(post.sample), length.out = 41)
hist(post.sample, breaks = breaks, 
     main = 'Distribution of samples \n with 90% HDI',
     xlab = 'Sample value',
     ylab = 'Density')
abline(v = quants[1], lty = 3, col = 'red', lwd = 3)
abline(v = quants[2], lty = 3, col = 'red', lwd = 3)
qqnorm(post.sample)

par(mfrow = c(1,1))

predplot(beta.post.par, 7, 53)

n <- 100
s <- 0:n
pred.probs <- pbetap(beta.post.par, n, s)
pred.probs.beta <- pbetap(beta.par, n, s)

plot(s, pred.probs, type="h", 
     main = paste('Probability distribution of successes in', as.character(n), 'trials'),
     xlab = 'Successes')

discint(cbind(s, pred.probs.beta) , 0.90)

## $prob
## [1] 0.9038973
## 
## $set
##  [1]  0  1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22
## [24] 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45
## [47] 46 47 48 49 50 51 52 53 54 55

discint(cbind(s, pred.probs), 0.90)

## $prob
## [1] 0.9098956
## 
## $set
##  [1]  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20

Visualization:

par(mfrow = c(2,1))


plot(s, pred.probs.beta, type="h", 
     main = paste('Probability distribution of Prior beleifs in', as.character(n), 'trials'),
     xlab = 'Successes')

plot(s, pred.probs, type="h", 
     main = paste('Probability distribution of posterier beleifs in', as.character(n), 'trials'),
     xlab = 'Successes')

par(mfrow = c(1,1))

Conclusion:

Based on the above results and visualizations, Cumulative probability of the city is better than the Nationally cumulative probability that a driver is texting when driving.