## Probability of a girl birth given placenta previa (BDA3 p. 37).

Simulate samples from Beta(438,544), draw a histogram with quantiles, and do the same for a transformed variable.

ggplot2 is used for plotting, tidyr for manipulating data frames

library(ggplot2)
theme_set(theme_minimal())
library(tidyr)

Sample from posterior Beta(438,544). Obtain all draws at once and store them in vector ‘theta’

a <- 438
b <- 544
theta <- rbeta(10000, a, b)

Compute odds ratio for all draws

phi <- (1 - theta) / theta

Compute 2.5% and 97.5% quantile approximation using samples

quantiles <- c(0.025, 0.975)
thetaq <- quantile(theta, quantiles)
phiq <- quantile(phi, quantiles)

Histogram plots with 30 bins for theta and phi

# merge the data into one data frame for plotting
df1 <- data.frame(phi,theta) %>% pivot_longer(everything())
# merge quantiles into one data frame for plotting
df2 <- data.frame(phi = phiq, theta = thetaq) %>% pivot_longer(everything())
ggplot() +
geom_histogram(data = df1, aes(value), bins = 30) +
geom_vline(data = df2, aes(xintercept = value), linetype = 'dotted') +
facet_wrap(~name, ncol = 1, scales = 'free_x')  +
labs(x = '', y = '') +
scale_y_continuous(breaks = NULL)

LS0tCnRpdGxlOiAiQmF5ZXNpYW4gZGF0YSBhbmFseXNpcyBkZW1vIDIuMyIKYXV0aG9yOiAiQWtpIFZlaHRhcmksIE1hcmt1cyBQYWFzaW5pZW1pIgpkYXRlOiAiYHIgZm9ybWF0KFN5cy5EYXRlKCkpYCIKb3V0cHV0OgogIGh0bWxfZG9jdW1lbnQ6CiAgICB0aGVtZTogcmVhZGFibGUKICAgIGNvZGVfZG93bmxvYWQ6IHRydWUKLS0tCiMjIFByb2JhYmlsaXR5IG9mIGEgZ2lybCBiaXJ0aCBnaXZlbiBwbGFjZW50YSBwcmV2aWEgKEJEQTMgcC4gMzcpLgoKU2ltdWxhdGUgc2FtcGxlcyBmcm9tIEJldGEoNDM4LDU0NCksIGRyYXcgYSBoaXN0b2dyYW0gd2l0aApxdWFudGlsZXMsIGFuZCBkbyB0aGUgc2FtZSBmb3IgYSB0cmFuc2Zvcm1lZCB2YXJpYWJsZS4KCmdncGxvdDIgaXMgdXNlZCBmb3IgcGxvdHRpbmcsIHRpZHlyIGZvciBtYW5pcHVsYXRpbmcgZGF0YSBmcmFtZXMKCmBgYHtyIHNldHVwLCBtZXNzYWdlPUZBTFNFLCBlcnJvcj1GQUxTRSwgd2FybmluZz1GQUxTRX0KbGlicmFyeShnZ3Bsb3QyKQp0aGVtZV9zZXQodGhlbWVfbWluaW1hbCgpKQpsaWJyYXJ5KHRpZHlyKQpgYGAKClNhbXBsZSBmcm9tIHBvc3RlcmlvciBCZXRhKDQzOCw1NDQpLgpPYnRhaW4gYWxsIGRyYXdzIGF0IG9uY2UgYW5kIHN0b3JlIHRoZW0gaW4gdmVjdG9yICd0aGV0YScKCmBgYHtyIH0KYSA8LSA0MzgKYiA8LSA1NDQKdGhldGEgPC0gcmJldGEoMTAwMDAsIGEsIGIpCmBgYAoKQ29tcHV0ZSBvZGRzIHJhdGlvIGZvciBhbGwgZHJhd3MKCmBgYHtyIH0KcGhpIDwtICgxIC0gdGhldGEpIC8gdGhldGEKYGBgCgpDb21wdXRlIDIuNSUgYW5kIDk3LjUlIHF1YW50aWxlIGFwcHJveGltYXRpb24gdXNpbmcgc2FtcGxlcwoKYGBge3IgfQpxdWFudGlsZXMgPC0gYygwLjAyNSwgMC45NzUpCnRoZXRhcSA8LSBxdWFudGlsZSh0aGV0YSwgcXVhbnRpbGVzKQpwaGlxIDwtIHF1YW50aWxlKHBoaSwgcXVhbnRpbGVzKQpgYGAKCkhpc3RvZ3JhbSBwbG90cyB3aXRoIDMwIGJpbnMgZm9yIHRoZXRhIGFuZCBwaGkKCmBgYHtyIH0KIyBtZXJnZSB0aGUgZGF0YSBpbnRvIG9uZSBkYXRhIGZyYW1lIGZvciBwbG90dGluZwpkZjEgPC0gZGF0YS5mcmFtZShwaGksdGhldGEpICU+JSBwaXZvdF9sb25nZXIoZXZlcnl0aGluZygpKQojIG1lcmdlIHF1YW50aWxlcyBpbnRvIG9uZSBkYXRhIGZyYW1lIGZvciBwbG90dGluZwpkZjIgPC0gZGF0YS5mcmFtZShwaGkgPSBwaGlxLCB0aGV0YSA9IHRoZXRhcSkgJT4lIHBpdm90X2xvbmdlcihldmVyeXRoaW5nKCkpCmdncGxvdCgpICsKICBnZW9tX2hpc3RvZ3JhbShkYXRhID0gZGYxLCBhZXModmFsdWUpLCBiaW5zID0gMzApICsKICBnZW9tX3ZsaW5lKGRhdGEgPSBkZjIsIGFlcyh4aW50ZXJjZXB0ID0gdmFsdWUpLCBsaW5ldHlwZSA9ICdkb3R0ZWQnKSArCiAgZmFjZXRfd3JhcCh+bmFtZSwgbmNvbCA9IDEsIHNjYWxlcyA9ICdmcmVlX3gnKSAgKwogIGxhYnMoeCA9ICcnLCB5ID0gJycpICsKICBzY2FsZV95X2NvbnRpbnVvdXMoYnJlYWtzID0gTlVMTCkKYGBgCgo=