--- title: "Statistical test workflows" output: rmarkdown::pdf_document vignette: > %\VignetteIndexEntry{Statistical test workflows} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r, echo = FALSE} knitr::opts_chunk$set(collapse = TRUE, comment = "#>", eval = TRUE) ``` ## Introduction `testflow` provides statistical testing workflows organized by study design. ## One numerical variable ```{r} library(testflow) cardio <- make_cardio_data() test_one_sample(cardio, sbp_3m, mu = 140) ``` ## Two independent groups ```{r} test_two_groups(sbp_3m ~ sex, data = cardio) ``` ## Paired measurements ```{r} test_paired(sbp_3m ~ sbp_baseline, data = cardio) ``` ## More than two groups ```{r} test_groups(sbp_3m ~ treatment, data = cardio) ``` ## Factorial designs ```{r} test_factorial(sbp_3m ~ sex * treatment, data = cardio) ``` ## Repeated measurements ```{r} test_repeated(cardio, c(sbp_baseline, sbp_3m, sbp_6m), id = id) ``` The repeated numeric workflow chooses repeated-measures ANOVA when the within-time normality checks are acceptable and Friedman otherwise. Post-hoc comparisons are paired t-tests for the parametric branch and paired Wilcoxon tests for the non-parametric branch. ## Categorical outcomes ```{r} test_categorical(treatment ~ controlled_3m, data = cardio) ``` ## Repeated categorical outcomes ```{r} test_repeated_categorical(cardio, c(controlled_baseline, controlled_3m, controlled_6m)) ``` The repeated categorical workflow uses Cochran Q for binary repeated outcomes and pairwise McNemar tests for follow-up comparisons. ## References - Fisher, R. A. (1925). \emph{Statistical Methods for Research Workers}. - Gosset, W. S. (1908). The probable error of a mean. - Welch, B. L. (1947). Generalization of Student's problem with unequal variances. - Wilcoxon, F. (1945). Individual comparisons by ranking methods. - Mann, H. B., & Whitney, D. R. (1947). On a test of whether one of two random variables is stochastically larger than the other. - Levene, H. (1960). Robust tests for equality of variances. - Kruskal, W. H., & Wallis, W. A. (1952). Use of ranks in one-criterion variance analysis. - Tukey, J. W. (1949). Comparing individual means in the analysis of variance. - Dunn, O. J. (1964). Multiple comparisons using rank sums. - Friedman, M. (1937). The use of ranks to avoid the assumption of normality implicit in the analysis of variance. - Cochran, W. G. (1950). The comparison of percentages in matched samples. - McNemar, Q. (1947). Note on the sampling error of the difference between correlated proportions or percentages. - Pearson, K. (1895, 1900). - Spearman, C. (1904). The proof and measurement of association between two things. - Kendall, M. G. (1938). A new measure of rank correlation. - Cramer, H. (1946). \emph{Mathematical Methods of Statistics}. - Clopper, C. J., & Pearson, E. S. (1934). The use of confidence or fiducial limits illustrated in the case of the binomial. - Cohen, J. (1988). \emph{Statistical Power Analysis for the Behavioral Sciences}. - Draper, N. R., & Smith, H. (1998). \emph{Applied Regression Analysis} (3rd ed.). - Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). \emph{Applied Logistic Regression} (3rd ed.). - Kaplan, E. L., & Meier, P. (1958). Nonparametric estimation from incomplete observations. - Mantel, N. (1966). Evaluation of survival data and two new rank order statistics arising in its consideration. - Cox, D. R. (1972). Regression models and life-tables. - Grambsch, P. M., & Therneau, T. M. (1994). Proportional hazards tests and diagnostics based on weighted residuals. - Simel, D. L., Samsa, G. P., & Matchar, D. B. (1991). Likelihood ratios with confidence: sample size estimation for diagnostic test studies. - Altman, D. G., & Bland, J. M. (1994). Diagnostic tests 1: sensitivity and specificity. - Hanley, J. A., & McNeil, B. J. (1982). The meaning and use of the area under a receiver operating characteristic (ROC) curve. - Youden, W. J. (1950). Index for rating diagnostic tests. - Fleiss, J. L., Cohen, J., & Everitt, B. S. (1969). Large sample standard errors of kappa and weighted kappa. - Landis, J. R., & Koch, G. G. (1977). The measurement of observer agreement for categorical data. - Shrout, P. E., & Fleiss, J. L. (1979). Intraclass correlations: uses in assessing rater reliability. - McGraw, K. O., & Wong, S. P. (1996). Forming inferences about some intraclass correlation coefficients. - Koo, T. K., & Li, M. Y. (2016). A guideline of selecting and reporting intraclass correlation coefficients for reliability research. ## Correlation ```{r} test_correlation(sbp_3m ~ age, data = cardio) ``` ## Regression ```{r} test_linear_regression(sbp_3m ~ age + ldl, data = cardio) ``` Linear regression reports the overall model F-test as the primary result, the per-term coefficient table (with confidence intervals) in `alternative_tests$coefficients`, and R squared / adjusted R squared as the effect size. Residual normality (Shapiro-Wilk), homoscedasticity (Breusch-Pagan), and multicollinearity (variance inflation factor, when there is more than one predictor) are checked as assumptions. ```{r} test_logistic_regression(controlled_3m ~ age + ldl, data = cardio) ``` Logistic regression reports the likelihood-ratio test against the intercept-only model as the primary result, coefficients on both the log-odds and odds-ratio scale in `alternative_tests`, and McFadden's pseudo R squared as the effect size. Multicollinearity and influential observations (Cook's distance) are checked as assumptions. ## Survival analysis `make_cardio_data()` does not include time-to-event data, so this section uses a small synthetic example instead. ```{r} set.seed(1) n <- 100 survival_dat <- tibble::tibble( time = rexp(n, 0.1), status = rbinom(n, 1, 0.7), arm = rep(c("control", "treatment"), each = n / 2), age = rnorm(n, 60, 10) ) test_survival(Surv(time, status) ~ arm, data = survival_dat) ``` `test_survival()` reports the log-rank test as the primary result and a companion univariate Cox hazard ratio as the effect size; the hazard ratio is not itself part of the log-rank test and can disagree with it under strongly non-proportional hazards. `Surv()` is re-exported from the `survival` package, so it is available after `library(testflow)` alone. ```{r} test_cox(Surv(time, status) ~ age + arm, data = survival_dat) ``` `test_cox()` reports the overall likelihood-ratio test as the primary result, hazard ratios per term, the concordance index as the effect size, and checks the proportional-hazards assumption via the Schoenfeld residual test (`survival::cox.zph()`). ## Diagnostic test accuracy and agreement `make_cardio_data()` has no diagnostic-test or multi-rater columns, so this section uses small synthetic examples instead. ```{r} set.seed(1) diag_dat <- tibble::tibble( test = c(rep("positive", 55), rep("negative", 98)), reference = c(rep("positive", 45), rep("negative", 10), rep("positive", 8), rep("negative", 90)) ) test_diagnostic(diag_dat, test, reference) ``` `test_diagnostic()` reports sensitivity, specificity, predictive values, and likelihood ratios (each with an exact confidence interval) in `alternative_tests$diagnostic_table`, and tests overall accuracy against the no-information rate as the primary result, following the same convention as `caret::confusionMatrix()`. ```{r} roc_dat <- tibble::tibble( marker = c(rnorm(60, 2, 1), rnorm(50, 0, 1)), disease = c(rep("yes", 60), rep("no", 50)) ) test_roc(roc_dat, marker, disease) ``` `test_roc()` computes the AUC from the Mann-Whitney relationship, a closed-form Hanley-McNeil confidence interval, and the Youden's-J-optimal threshold. ```{r} agree_dat <- tibble::tibble( rater1 = sample(c("mild", "moderate", "severe"), 100, replace = TRUE), rater2 = sample(c("mild", "moderate", "severe"), 100, replace = TRUE) ) test_agreement(agree_dat, rater1, rater2) ``` `test_agreement()` reports Cohen's kappa with the Fleiss-Cohen-Everitt (1969) large-sample confidence interval. ```{r} icc_dat <- tibble::tibble( rater1 = rnorm(30, 50, 10), rater2 = rnorm(30, 50, 10), rater3 = rnorm(30, 50, 10) ) test_icc(icc_dat, c(rater1, rater2, rater3)) ``` `test_icc()` reports ICC(2,1) (two-way random effects, absolute agreement, single measurement) as the primary result, following the reliability-study recommendation of Koo & Li (2016), alongside ICC(1,1) and ICC(3,1) for comparison in `alternative_tests$icc_table`. ## Outliers ```{r} test_outliers(c(sbp_3m, ldl, crp), data = cardio) ``` ## Reporting and plotting Every workflow returns a `testflow` object. Use `report(x)`, `plot(x)`, and `as_tibble(x)`. See `effect-size-formulas.Rmd` for the exact formulas used by the reported effect-size estimates.