Problem 19

#The question if the binomial model suits this problem could be argued in both ways:
#Its suitable since we are looking for a binary outcome but the question asked could also be interpreted if the subjects are generally suffering from some sort of colour blindness.
phat_ml=8/95

print(phat_ml)
## [1] 0.08421053
binom.test(8 ,95 )
## 
##  Exact binomial test
## 
## data:  8 and 95
## number of successes = 8, number of trials = 95, p-value < 2.2e-16
## alternative hypothesis: true probability of success is not equal to 0.5
## 95 percent confidence interval:
##  0.03705905 0.15920028
## sample estimates:
## probability of success 
##             0.08421053
prop.test(8,95)
## 
##  1-sample proportions test with continuity correction
## 
## data:  8 out of 95, null probability 0.5
## X-squared = 64.042, df = 1, p-value = 1.218e-15
## alternative hypothesis: true p is not equal to 0.5
## 95 percent confidence interval:
##  0.0396893 0.1639164
## sample estimates:
##          p 
## 0.08421053
# A p-value (given a certain confidence level alpha) gives us the probability of gettin a certain value (or greater) which can be can help us to determine if we reject or not to reject a hypothesis. I basically gives us an "area" where we can reject a certain hypothesis if the p-value is as small or smaller than a predetermined probability value.

# if p-value < alpha --> reject H0. (or otherwise, depending on the hypothesis)
Wilson_CI <- function(x,n) {
  mid_wil <- (x + 1.96^2/2)/(n + 1.96^2)
  se_wil <- sqrt(n)/(n + 1.96^2) * sqrt(x/n * (1 - x/n) + 1.96^2/(4 * n))
  return(round(cbind(pmax(0, mid_wil - 1.96 * se_wil), pmin(1, mid_wil + 1.96 * se_wil)), 4))
}

Wilson_CI(8,95)
##        [,1]   [,2]
## [1,] 0.0433 0.1575

Problem 20

# Preparing Wald confidence intervall

Wald_CI <- function(x, n) {
se_wall <- sqrt(x * (n - x)/n^3)
return(round(cbind(pmax(0, x/n - 1.96 * se_wall), pmin(1, x/n + 1.96 * se_wall)), 4))
}

#########


# Wilson CI for Treatment A

Wilson_CI(9,27)
##        [,1]   [,2]
## [1,] 0.1864 0.5218
# Walden CI for Treatment A

Wald_CI(9,27)
##        [,1]   [,2]
## [1,] 0.1555 0.5111
# Wilson CI for Treatment B

Wilson_CI(5,27)
##        [,1]  [,2]
## [1,] 0.0818 0.367
# Walden CI for Treatment B

Wald_CI(5,27)
##        [,1]   [,2]
## [1,] 0.0387 0.3317
data_treat <- rbind(c(9,5),c(18,22))

data_treat
##      [,1] [,2]
## [1,]    9    5
## [2,]   18   22
prop.test(data_treat)
## 
##  2-sample test for equality of proportions with continuity correction
## 
## data:  data_treat
## X-squared = 0.86786, df = 1, p-value = 0.3515
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  -0.1499190  0.5356333
## sample estimates:
##    prop 1    prop 2 
## 0.6428571 0.4500000
# With the p-value of 0.3515 we fail to reject the Null-hypothesis which implied no difference between the two treatment groups. As a conclusion to this specific proportion test the two treatments do not differ in effectiveness.

p_CI95 <- function(p1,p2,n1,n2) {
  se <- sqrt((p1*(1-p1))/n1 + (p2*(1-p2))/n2)
  p_del <- p1-p2
  return(round(cbind(pmax(0, p_del - 1.96 * se), pmin(1, p_del + 1.96 * se)), 4))
}


p1 <- 9/18

p2 <- 5/22

n1 <- 27

n2 <- 27

p_CI95(p1=p1,p2=p2,n1=n1,n2=n2)
##        [,1]   [,2]
## [1,] 0.0266 0.5188
# With the approximate approach, there seems to be differing conclusion than with the proportion test. With the approximate approach we reject HO (A and B in this test indicate a significant difference) since O is not in the CI.
rr_ci95 <- function(x1, x2, n1, n2) {
    p1 <- x1/(n1)
    p2 <- x2/(n1)
    log_rr <- log(p1) - log(p2)
    se_log <- sqrt(1/x1 - 1/n1 + 1/x2 - 1/n2)
    return(round(c(exp(log_rr - 1.96 * se_log), exp(log_rr + 1.96 * se_log)), 4))
}
    
or_ci95 <- function(x1, x2, n1, n2) { 
    p1 <- x1/n1
    p2 <- x2/n2
    log_or <- log((p1/(1 - p1))/(p2/(1 - p2)))
    se_log <- sqrt(1/x1 + 1/x2 + 1/(n1 - x1) + 1/(n2 - x2))
    return(round(c(exp(log_or - 1.96 * se_log), exp(log_or + 1.96 * se_log)), 4))
}


rbind(c(rr_ci95(x1 = 9, x2 = 5, n1 = 27, n2 = 27)),c(or_ci95(x1 = 9, x2 = 5, n1 = 27, n2 = 27)))
##        [,1]   [,2]
## [1,] 0.6932 4.6741
## [2,] 0.6251 7.7424
# Both CI would fail to reject HO which implied no increased risk ( no risk association) since 1 is in the CI.