R code for the adjustment of the Wald confidence interval for a
difference of proportions, with matched pairs.
This is the interval called Wald+2 in Agresti and Min, Statistics
in Medicine, 2005, which adds
0.5 to each cell before constructing the Wald CI. The CI is truncated
when it overshoots the boundary.
--------------------------------------------------------------------------
diffpropci <- function(b,c,n,conflev)
{
z <- qnorm(1-(1-conflev)/2)
diff <- (c-b)/(n+2)
sd <- sqrt((b+c+1)-(c-b)^2/(n+2))/(n+2)
ll <- diff - z*sd
ul <- diff + z*sd
if(ll < -1) ll = -1
if(ul > 1) ul = 1
c(ll, ul)
}
# Adjusted Wald interval for difference of proportions with matched pairs
# "conflev"=confidence coefficient, n=sample size, b,c = off-diag counts
---------------------------------------------------------------------------
R code for Tango's score confidence interval for a
difference of proportions (Statistics in Medicine 1998),
with matched pairs, described in Agresti and Min (Statistics
in Medicine, 2005)
-----------------------------------------------------------------------------
scoreci <- function(b,c,n,conflev)
{
pa = 2*n
z = qnorm(1-(1-conflev)/2)
if(c == n) {ul = 1}
else{
proot = (c-b)/n
dp = 1-proot
niter = 1
while(niter <= 50){
dp = 0.5*dp
up2 = proot+dp
pb = - b - c + (2*n-c+b)*up2
pc = -b*up2*(1-up2)
q21 = (sqrt(pb^2-4*pa*pc)-pb)/(2*pa)
score = (c-b-n*up2)/sqrt(n*(2*q21+up2*(1-up2)))
if(abs(score)<z){ proot = up2 }
niter=niter+1
if((dp<0.0000001) || (abs(z-score)<.000001)){
niter=51
ul=up2
}
}
}
if(b == n) {ll = -1}
else{
proot = (c-b)/n
dp = 1+proot
niter = 1
while(niter <= 50){
dp = 0.5*dp
low2 = proot-dp
pb = - b - c + (2*n-c+b)*low2
pc = -b*low2*(1-low2)
q21 = (sqrt(pb^2-4*pa*pc)-pb)/(2*pa)
score = (c-b-n*low2)/sqrt(n*(2*q21+low2*(1-low2)))
if(abs(score) < z){proot = low2}
niter = niter+1
if((dp<0.0000001) || (abs(z-score)<.000001)){
ll = low2
niter = 51
}
}
}
c(ll,ul)
}
# Tango score interval for difference of proportions with matched
# pairs, confidence coefficient = conflev, n = sample size,
# b,c = off-diagonal counts
-----------------------------------------------------------------------------
R code for the adapted binomial score confidence interval for the
subject-specific odds ratio,
with matched pairs, described in Agresti and Min, Statistics in
Medicine, 2004.
This uses the Wilson score CI for a binomial parameter with the
off-diagonal counts.
-----------------------------------------------------------------------------
oddsratioci <- function(b,c,conflev)
{
z <- qchisq(conflev,1)
A <- b + c + z
B <- 2*c + z
C <- c^2/(b+c)
l <- (B - sqrt(B^2-4*A*C))/(2*A)
u <- (B + sqrt(B^2-4*A*C))/(2*A)
ll <- l/(1-l)
ul <- u/(1-u)
c(ll,ul)
}
# adapts Wilson binomial score CI to form CI for subject-specific
# odds ratio with matched pairs data
# "conflev"=confidence coefficient, b,c = off-diag counts
-----------------------------------------------------------------------------