B Código da função paired_proportions
Abaixo é exibido o código da função paired_proportions, utilizada no capítulo 17, seção 17.3.3.
library(dplyr)
<- function(data, id, row, col, row_ref, row_trt, col_ref,
paired_proportions case_control = FALSE, alpha = 0.05, ...) {
col_out, = id
id_subj99 = mutate(data,
summary cell = case_when(
== row_ref & data[, col] == col_ref ~ 1,
data[, row] == row_trt & data[, col] == col_ref ~ 2,
data[, row] == row_ref & data[, col] == col_out ~ 4,
data[, row] == row_trt & data[, col] == col_out ~ 8,
data[, row] TRUE ~ 0)) %>% group_by(data[, id_subj99]) %>%
summarize(cell_class = sum(cell))
= table(summary$cell_class)
tab if (case_control == TRUE) {
#
# tab["5"] -> (row_ref in both),
# tab["6"] -> (row_ref in col_out and row_trt in col_ref),
# tab["9"] -> (row_trt in col_out and row_ref in col_ref),
# tab["10"] -> (row_trt in both),
#
if (is.na(tab["5"])) tab["5"] = 0
if (is.na(tab["6"])) tab["6"] = 0
if (is.na(tab["9"])) tab["9"] = 0
if (is.na(tab["10"])) tab["10"] = 0
= matrix(c(tab["10"], tab["6"], tab["9"], tab["5"]), nrow = 2)
mat colnames(mat) = c(row_trt, row_ref)
rownames(mat) = c(row_trt, row_ref)
names(dimnames(mat)) = c(col_out,col_ref)
else {
} #
# tab["3"] -> (col_ref in both),
# tab["6"] -> (col_ref in row_trt and col_out in row_ref),
# tab["9"] -> (col_out in row_trt and col_ref in row_ref),
# tab["12"] -> (col_out in both),
#
if (is.na(tab["3"])) tab["3"] = 0
if (is.na(tab["6"])) tab["6"] = 0
if (is.na(tab["9"])) tab["9"] = 0
if (is.na(tab["12"])) tab["12"] = 0
= matrix(c(tab["12"], tab["6"], tab["9"], tab["3"]), nrow = 2)
mat colnames(mat) = c(col_out, col_ref)
rownames(mat) = c(col_out, col_ref)
names(dimnames(mat)) = c(row_trt, row_ref)
}#
# IC proportion differences - Wald with Bonett–Price Laplace adjustment
#
= sum(mat)
n = mat[1,1]
r = mat[2,1]
t = mat[1,2]
s = r+s
n1. .1 = r+t
n= n1. / n
p2 = n.1 / n
p1 = (s + 1) / (n + 2)
p2_ = (t + 1) / (n + 2)
p1_ = qnorm(1-alpha/2)
z = p2 - p1
D = max(-1, (p2_-p1_) - z * sqrt((p2_+p1_-(p2_-p1_)^2)/(n+2)))
Dinf = min((p2_-p1_) + z * sqrt((p2_+p1_-(p2_-p1_)^2)/(n+2)), 1)
Dsup = c(D, Dinf, Dsup)
prop_diff #
# IC odds ratio - Transformed Wilson score
#
= s/t
or_est = t+s
nd = 2*(nd+z^2)
div = (2*s + z^2)
num1 = z*sqrt(z^2 + 4*s*(1-s/nd))
num2 = (num1 + num2) / div
pu = (num1 - num2) / div
pi = c(or_est, pi/(1-pi), pu/(1-pu))
or #
# IC relative risk - Bonett–Price hybrid Wilson score
#
= p2/p1
rr_est = r + t + s
nn = sqrt((nd + 2) / ((n1.+1)*(n.1+1)))
A = sqrt((1-(n1.+1)/(nn+2))/(n1.+1))
B = sqrt((1-(n.1+1)/(nn+2))/(n.1+1))
C = A/(B+C) * z
z =2*(nn+z^2)
den = 2*n1.+z^2 - z*sqrt(z^2+4*n1.*(1-n1./nn))
l1 = 2*n1.+z^2 + z*sqrt(z^2+4*n1.*(1-n1./nn))
u1 = 2*n.1+z^2 - z*sqrt(z^2+4*n.1*(1-n.1/nn))
l2 = 2*n.1+z^2 + z*sqrt(z^2+4*n.1*(1-n.1/nn))
u2 = c(rr_est, l1/u2, u1/l2)
rr print(mat)
print(mcnemar.test(mat, ...))
= paste("Lower ", sprintf("%.0f", (1-alpha)*100),"% CI",
text_low sep = '')
= paste("Upper ", sprintf("%.0f", (1-alpha)*100),"% CI",
text_upper sep = '')
if (case_control == FALSE) {
names(prop_diff) = c("proportion differences", text_low, text_upper)
print(prop_diff)
cat("", sep = "\n")
names(rr) = c("relative risk", text_low, text_upper)
print(rr)
cat("", sep = "\n")
}names(or) = c("odds ratio", text_low,text_upper)
print(or)
cat("", sep="\n")
}