dsrTest - examples

Michael Nelson

2022-06-19

To demonstrate the functionality of dsrTest the various methods to construct confidence interval are demonstrated on the example data drawn from Table 14.4 in Fleiss (1981) (also used in Fay and Feuer (1997)).

The example data come from a study of Down Syndrome and maternal age (Stark and Mantell, 1966).

library(dsrTest)
# the data are in `downs.mi`
data("downs.mi", package = "dsrTest")
# Birth order 5 +
b5 <- downs.mi[downs.mi$BirthOrder == 5, ]
# Gamma Method
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "gamma"))
## 
##  Directly standardized rate: Gamma Method for Weighted Sum of Poissons
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##   67.70212 188.30017
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Gamma Mid-p
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "gamma",
                 control = list(midp = TRUE)))
## 
##  Directly standardized rate: Confidence Distribution method for
##  Weighted Sum of Poissons
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##   59.71632 173.08170
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Dobson (exact)
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "dobson"))
## 
##  Directly standardized rate: Dobson's method for Weighted Sum of
##  Poissons with Exact two-sided Poisson test (central method)
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  67.63284 83.86703
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Dobson (Mid-p)
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "dobson",
                 control = list(midp = TRUE)))
## 
##  Directly standardized rate: Dobson's method for Weighted Sum of
##  Poissons with Exact two-sided Poisson test (central method), mid-p
##  version
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  67.70418 83.79030
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Asymptotic (no transformation)
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "asymptotic"))
## 
##  Directly standardized rate: Asymptotic method for Weighted Sum of
##  Poissons normal approximation of MLE
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  67.48907 83.56895
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Asymptotic (log transformation)
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "asymptotic",
                 control = list(trans = "log")))
## 
##  Directly standardized rate: Asymptotic method for Weighted Sum of
##  Poissons normal approximation of the log-transformed MLE
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  67.90220 84.01246
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Beta Method
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "beta"))
## 
##  Directly standardized rate: Beta method for Weighted Sum of Poissons
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  67.70196 83.77768
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901
# Approximate Bootstrap Method
with(b5, dsrTest(Cases, Births, Standard, mult = 1e5, method = "bootstrap"))
## 
##  Directly standardized rate: Approximate bootstrap method for Weighted
##  Sum of Poissons
## 
## data:  x = Cases, time bases: n = Births, weights: w = Standard
## number of summands = 6, Multiplier = 1e+05, p-value = NA
## alternative hypothesis: two.sided
## 95 percent confidence interval:
##  68.35707 84.61566
## sample estimates:
## Directly standardized rate per 1e+05 
##                             75.52901

Summaries

# A list of methods to implement
methods_list <- list(
  gamma = list(
    list(wmtype = "max"),
    list(midp = TRUE),
    list(wmtype = "tcz"),
    list(wmtype = "mean"),
    list(wmtype = "minmaxavg")),
  asymptotic = list(
    list(trans = "none"),
    list(trans = "log"),
    list(trans = "loglog"),
    list(trans = "logit")),
  dobson = list(
    list(midp = FALSE),
    list(midp = TRUE)),
  beta = list(
    list(wmtype = "none"),
    list(wmtype = "tcz"),
    list(wmtype = "mean"),
    list(wmtype = "minmaxavg"),
    list(wmtype = "max")),
  bootstrap = list(list())
)
# split out to allow call to mapply
methods <-rep(names(methods_list), times = lengths(methods_list))
controls <- do.call(c, unname(methods_list))
all_methods <- mapply(dsrTest,
  method = methods, control = controls,
  MoreArgs = list(mult = 1e5, x = b5$Cases, n = b5$Births, w = b5$Standard),
  SIMPLIFY = FALSE)
# create some "short" names
control_types <- unlist(controls)
control_names <- c(gsub("midp=FALSE", "Exact CI",
  gsub("=TRUE", "",
       sprintf("[%s=%s]", names(control_types), control_types))), "")
names(all_methods) <- paste(methods, control_names)
# combine CI into one data.frame
results <- do.call(rbind,lapply(all_methods,
  function(data) data.frame(
    estimate = data$estimate,
    lci = data$conf.int[1], 
    uci = data$conf.int[2])))
# and display
knitr::kable(results, digits = 3)
estimate lci uci
gamma [wmtype=max] 75.529 67.702 188.300
gamma [midp] 75.529 59.716 173.082
gamma [wmtype=tcz] 75.529 67.702 112.858
gamma [wmtype=mean] 75.529 67.702 96.465
gamma [wmtype=minmaxavg] 75.529 67.702 130.933
asymptotic [trans=none] 75.529 67.489 83.569
asymptotic [trans=log] 75.529 67.902 84.012
asymptotic [trans=loglog] 75.529 67.848 83.947
asymptotic [trans=logit] 75.529 67.902 84.012
dobson [Exact CI] 75.529 67.633 83.867
dobson [midp] 75.529 67.704 83.790
beta [wmtype=none] 75.529 67.702 83.778
beta [wmtype=tcz] 75.529 55.287 112.856
beta [wmtype=mean] 75.529 67.891 96.465
beta [wmtype=minmaxavg] 75.529 61.280 130.930
beta [wmtype=max] 75.529 52.767 188.290
bootstrap 75.529 68.357 84.616

References

Fleiss, JL (1981) Statistical Methods for Rates and Proportions, Wiley, New York.

Stark CR and Mantel N (1966) ‘Effects of maternal age and birth order on the risk of mongolism and leukemia’ J Natl Cancer Inst 37 (5) 687–698. https://doi.org/10.1093/jnci/37.5.687

Fay MP & Feuer EJ (1997). ’Confidence intervals for directly standardized rates: a method based on the gamma distribution. Statistics in Medicine. 16: 791–801. https://doi.org/10.1002/(SICI)1097-0258(19970415)16:7<791::AID-SIM500>3.0.CO;2-%23