# Tests some of implementation in logistic_sample_size.R using examples 
# from Hsieh et al. 

alpha <- 0.05;   
beta <- 0.05;                     
Zalpha <- qnorm(1 - alpha/2);     
Zbeta <- qnorm(1 - beta);         

# Binary covariates: Table 1
Ps <- c(0.45, 0.35, 0.35, 0.075, 0.06, 0.09); 
Bs <- c(0.5, 0.5, 0.5, 0.5, 0.2, 0.8); 
P1s <- c(0.4, 0.5, 0.2, 0.05, 0.05, 0.05); 
P2s <- c(0.5, 0.2, 0.5, 0.1, 0.1, 0.1); 
expected <- c(1282, 126, 126, 1437, 2186, 2257); 
for (i in 1:length(Ps)) {
    P <- Ps[i]; 
    B <- Bs[i]; 
    P1 <- P1s[i]; 
    P2 <- P2s[i]; 
    n1 <- ( Zalpha * (P*(1 - P)/B)^(1/2) + 
	    Zbeta * (P1*(1 - P1) + P2*(1 - P2)*(1 - B)/B)^(1/2) )^2 / 
	    ( (P1 - P2)^2 * (1 - B) ); 
    cat(paste("Actual = ", n1, ", Expected = ", expected[i], "\n", sep="")); 
}

# Continuous covariates: Table 2
P1s <- c(0.5, 0.4, 0.1); 
betaStars <- c(0.405, 0.405, 0.405); 
expected <- c(317, 330, 880); 
for (i in 1:length(P1s)) {
    P1 <- P1s[i]; 
    betaStar <- betaStars[i]; 
    n1 <- ( (Zalpha + Zbeta)^2 ) / ( P1*(1 - P1)*betaStar^2 ); 
    cat(paste("Actual = ", n1, ", Expected = ", expected[i], "\n", sep="")); 
}