Fredrik Lundgren
2010-Apr-04 10:51 UTC
[R] "mantel.haenszel.test for trend in S-plus doesn't work i R"
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.lundgren@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]]
David Winsemius
2010-Apr-04 13:17 UTC
[R] "mantel.haenszel.test for trend in S-plus doesn't work i R"
On Apr 4, 2010, at 6:51 AM, Fredrik Lundgren wrote:> 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 wellDon't think it has anything do do with "el". In this line: " O <- sum(apply(apply(array, c(, 2, 3), function(el) el * 1:s)[, + 1, ], 2, sum))" ... try removing the extraneous "," in the dimension spec of the inner apply call. I cannot be sure that the logic is correct but at least it allows the function to perform calculation and make a reprot. (I also removed the quotes around the function name since that seemed very un-R-like.) It might require some re-arrangements of the data and changing the age variable to a numeric score, but you considered using: glm( c(Case, Control) ~ age + sex, family="binomial") -- David> > 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.lundgren at 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 at 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.David Winsemius, MD West Hartford, CT
Erich Neuwirth
2010-Apr-04 14:40 UTC
[R] "mantel.haenszel.test for trend in S-plus doesn't work i R"
That is the problem> c(,2,3)Error: argument is missing, with no default On 4/4/2010 12:51 PM, Fredrik Lundgren wrote:> 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.lundgren at 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 at 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.-- Erich Neuwirth, University of Vienna Faculty of Computer Science Computer Supported Didactics Working Group Visit our SunSITE at http://sunsite.univie.ac.at Phone: +43-1-4277-39464 Fax: +43-1-4277-39459