Dear R'ers,
When I used S-plus i wrote a small program for a Mantel-Haenszel test
for trend (I think it worked). Unfortunately I can't get it working in
R.
It appears as if my use of 'el' is the problem but I can't sort it out.
Error in apply(array, c(, 2, 3), function(el) el * 1:s) :
argument is missing, with no default
Further down in the program I use 'el' as well
Can this be sorted out?
Fredrik
#######
"mantel.ext.test" <-
function(array)
{
# Mantel-Extension Test -- Testing for trend in the presence of
# confounding ref. Rosner : Fundamentals of Biostatistics, 5th
# ed, page 606-609
# array is an 3-dim array with as many strata (last dimension)
# you like case-control as second dimension and the dimension
# for which a trend should be tested as the first dimension
# Example:
# snoring <- array(c(196, 223, 103, 603, 486, 232, 188, 313,
# 232, 348, 383, 206), c(3,2,2))
# dimnames(snoring) <- list(age = c('30-39', '40-49', '50-60'), #
status = c('Case', 'Control'), sex = c('Women', 'Men'))
#
# mantel.ext.test(snoring)
# Data: snoring
# , , Women
# Case Control
# 30-39 196 603
# 40-49 223 486
# 50-60 103 232
# , , Men
# Case Control
# 30-39 188 348
# 40-49 313 383
# 50-60 232 206
# The trend of snoring with age , controlling for sex , has a
# Mantel-Extension-test chi-square = 35.06 (df = 1 ) p-value = 0
# Effect is " increasing " with increasing age
s <- dim(array)[1]
O <- sum(apply(apply(array, c(, 2, 3), function(el) el * 1:s)[,
1, ], 2, sum))
tot <- apply(array, c(3), sum)
sum.case <- apply(array, c(2, 3), sum)[1, ]
E <- sum((apply(apply(apply(array, c(1, 3), sum), 2, function(
el) el * 1:s), 2, sum) * sum.case)/tot)
s1 <- apply(apply(apply(array, c(1, 3), sum), 2, function(el)
el * 1:s), 2, sum, simplify = F)
s2 <- apply(apply(apply(array, c(1, 3), sum), 2, function(el)
el * (1:s)^2), 2, sum)
V <- sum((sum.case * (tot - sum.case) * (tot * s2 - s1^2))/(tot^
2 * (tot - 1)))
A <- abs(O - E) - 0.5
chi2 <- (abs(O - E) - 0.5)^2/V
p <- 1 - pchisq(chi2, 1)
incr <- ifelse(A > 0, "increasing", "decreasing")
cat("\n", "Data:", deparse(substitute(array)), "\n")
print(array)
cat("\n", "The trend of", deparse(substitute(array)), "with",
attr(dimnames(array), "names")[1], ", controlling for",
attr(dimnames(array), "names")[3], ",", "has a", "\n",
"Mantel-Extension-test chi-square =", round(chi2, 3),
"(df =", 1, ")", "p-value =", round(p, 6), "\n",
"Effect is \"", incr, "\" with increasing", attr(
dimnames(array), "names")[1], "\n")
}
######
########################
Fredrik Lundgren
fredrik.bg.lundg...@gmail.com
Engelbrektsgatan 31
582 21 Linköping
tel 013 - 47 30 117
mob 0706 - 86 39 29
Sommarhus: Ljungnäs 158
380 30 Rockneby
0480 - 650 98
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
