Package MKinfer
Matthias Kohl
2022-09-06
Introduction
Package MKinfer includes a collection of functions for the computation of various confidence intervals (Altman et al. 2000; Hedderich and Sachs 2018) including bootstrapped versions (Davison and Hinkley 1997) as well as Hsu (Hedderich and Sachs 2018), permutation (Janssen 1997), bootstrap (Davison and Hinkley 1997) and multiple imputation (Barnard and Rubin 1999) t-test.
We first load the package.
Confidence Intervals
Binomial Proportion
There are several functions for computing confidence intervals. We can compute 12 different confidence intervals for binomial proportions (DasGupta, Cai, and Brown 2001); e.g.
##
## wilson confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## prob 0.1429739 0.3741268
##
## sample estimate:
## prob
## 0.24
##
## additional information:
## standard error of prob
## 0.05896867
##
## clopper-pearson confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## prob 0.1306099 0.3816907
##
## sample estimate:
## prob
## 0.24
## [1] 0.1306099 0.3816907
## attr(,"conf.level")
## [1] 0.95
For all intervals implemented see the help page of function binomCI. One can also compute bootstrap intervals via function boot.ci of package boot (Davison and Hinkley 1997) as well as one-sided intervals.
Difference of Two Binomial Proportions
There are several functions for computing confidence intervals. We can compute different confidence intervals for the difference of two binomial proportions (independent (Newcombe 1998a) and paired case (Newcombe 1998b)); e.g.
##
## wilson confidence interval (independent proportions)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## difference of independent proportions -0.03813715 0.19256
##
## sample estimate:
## difference of proportions
## 0.08928571
##
## additional information:
## proportion of group 1 proportion of group 2
## 0.08928571 0.00000000
##
## wilson-cc confidence interval (paired data)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## difference of proportions (paired data) -0.1141915 -0.05543347
##
## sample estimate:
## difference of proportions
## -0.08491281
##
## additional information:
## proportion of group 1 proportion of group 2
## 0.2699014 0.3548143
For all intervals implemented see the help page of function binomDiffCI. One can also compute boostrap intervals via function boot.ci of package boot (Davison and Hinkley 1997) as well as one-sided intervals.
Mean and SD
We can compute confidence intervals for mean and SD (Altman et al. 2000, @Davison1997).
##
## Exact confidence interval(s)
##
## 95 percent confidence intervals:
## 2.5 % 97.5 %
## mean 1.013980 2.601463
## sd 2.333029 3.480362
##
## sample estimates:
## mean sd
## 1.807722 2.792929
##
## additional information:
## SE of mean
## 0.3949798
##
## Exact confidence interval(s)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## mean 1.01398 2.601463
##
## sample estimates:
## mean sd
## 1.807722 2.792929
##
## additional information:
## SE of mean
## 0.3949798
##
## Exact confidence interval(s)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## sd 2.333029 3.480362
##
## sample estimates:
## mean sd
## 1.807722 2.792929
##
## additional information:
## SE of mean
## 0.3949798
##
## Exact confidence interval(s)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## mean 0.9761795 2.639264
##
## sample estimate:
## mean
## 1.807722
##
## additional information:
## SE of mean
## 0.4242641
##
## Exact confidence interval(s)
##
## 95 percent confidence interval:
## 0 % 95 %
## sd 0 3.356325
##
## sample estimate:
## sd
## 2.792929
##
## Bootstrap confidence interval(s)
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level Normal Basic Studentized
## 95% ( 1.043, 2.578 ) ( 1.034, 2.565 ) ( 1.030, 2.607 )
##
## Level Percentile BCa
## 95% ( 1.050, 2.582 ) ( 1.056, 2.586 )
## Calculations and Intervals on Original Scale
##
## sample estimates:
## mean sd
## 1.807722 2.792929
Difference in Means
We can compute confidence interval for the difference of means (Altman et al. 2000; Hedderich and Sachs 2018; Davison and Hinkley 1997).
##
## Confidence interval (paired)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## mean of differences -0.9182056 1.000461
##
## sample estimates:
## mean of differences sd of differences
## 0.0411276 2.0497941
##
## additional information:
## SE of mean of differences
## 0.4583479
##
## Exact confidence interval(s)
##
## 95 percent confidence intervals:
## 2.5 % 97.5 %
## mean -0.9182056 1.000461
## sd 1.5588488 2.993871
##
## sample estimates:
## mean sd
## 0.0411276 2.0497941
##
## additional information:
## SE of mean
## 0.4583479
##
## Bootstrap confidence interval (paired)
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level Normal Basic Studentized
## 95% (-0.8325, 0.9082 ) (-0.8388, 0.8927 ) (-0.9203, 0.9866 )
##
## Level Percentile BCa
## 95% (-0.8105, 0.9211 ) (-0.8068, 0.9230 )
## Calculations and Intervals on Original Scale
##
## sample estimates:
## mean of differences sd of differences
## 0.0411276 2.0497941
##
## Classical confidence interval (unpaired)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## difference in means -2.723246 0.008406416
##
## sample estimate:
## difference in means
## -1.35742
##
## additional information:
## SE of difference in means Cohen's d (SMD)
## 0.6667748 -0.7884618
##
## mean of x SD of x mean of y SD of y
## 0.4030681 1.1849233 1.7604879 2.5014008
##
## Welch confidence interval (unpaired)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## difference in means -3.192146 0.4773068
##
## sample estimate:
## difference in means
## -1.35742
##
## additional information:
## SE of difference in means Cohen's d (SMD)
## 0.8342079 -0.4904221
##
## mean of x SD of x mean of y SD of y
## 0.4030681 1.1849233 1.7604879 2.5014008
##
## Hsu confidence interval (unpaired)
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## difference in means -3.244529 0.5296895
##
## sample estimate:
## difference in means
## -1.35742
##
## additional information:
## SE of difference in means Cohen's d (SMD)
## 0.8342079 -0.4904221
##
## mean of x SD of x mean of y SD of y
## 0.4030681 1.1849233 1.7604879 2.5014008
##
## Bootstrap confidence interval (equal variances, unpaired)
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level BCa
## 95% (-3.015, 0.106 )
## Calculations and Intervals on Original Scale
##
## sample estimate:
## difference in means
## -1.35742
##
## additional information:
## SE of difference in means Cohen's d (SMD)
## 0.6667748 -0.7884618
##
## mean of x SD of x mean of y SD of y
## 0.4030681 1.1849233 1.7604879 2.5014008
##
## Bootstrap confidence interval (unequal variances, unpaired)
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level BCa
## 95% (-2.999, 0.089 )
## Calculations and Intervals on Original Scale
##
## sample estimate:
## difference in means
## -1.35742
##
## additional information:
## SE of difference in means Cohen's d (SMD)
## 0.8342079 -0.4904221
##
## mean of x SD of x mean of y SD of y
## 0.4030681 1.1849233 1.7604879 2.5014008
In case of unequal variances and unequal sample sizes per group the classical confidence interval may have a bad coverage (too long or too short), as is indicated by the small Monte-Carlo simulation study below.
M <- 100
CIhsu <- CIwelch <- CIclass <- matrix(NA, nrow = M, ncol = 2)
for(i in 1:M){
x <- rnorm(10)
y <- rnorm(30, sd = 0.1)
CIclass[i,] <- normDiffCI(x, y, method = "classical")$conf.int
CIwelch[i,] <- normDiffCI(x, y, method = "welch")$conf.int
CIhsu[i,] <- normDiffCI(x, y, method = "hsu")$conf.int
}
## coverage probabilies
## classical
sum(CIclass[,1] < 0 & 0 < CIclass[,2])/M
## [1] 0.68
## [1] 0.94
## [1] 0.94
Coefficient of Variation
We provide 12 different confidence intervals for the (classical) coefficient of variation (Gulhar, Kibria, and Ahmed 2012; Davison and Hinkley 1997); e.g.
##
## Miller (1991) confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## CV 0.1784019 0.2390245
##
## sample estimate:
## CV
## 0.2087132
##
## additional information:
## standard error of CV
## 0.01546524
##
## Gulhar et al (2012) confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## CV 0.1832516 0.242457
##
## sample estimate:
## CV
## 0.2087132
##
## Bootstrap confidence interval
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level Normal Basic
## 95% ( 0.1710, 0.2507 ) ( 0.1696, 0.2498 )
##
## Level Percentile BCa
## 95% ( 0.1676, 0.2479 ) ( 0.1753, 0.2586 )
## Calculations and Intervals on Original Scale
##
## sample estimate:
## CV
## 0.2087132
For all intervals implemented see the help page of function cvCI.
Quantiles, Median and MAD
We start with the computation of confidence intervals for quantiles (Hedderich and Sachs 2018; Davison and Hinkley 1997).
##
## exact confidence interval
##
## 95.14464 percent confidence interval:
## lower upper
## 95 % quantile 5.53907 8.231155
##
## sample estimate:
## 95 % quantile
## 7.33884
##
## asymptotic confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## 95 % quantile 5.53907 8.550966
##
## sample estimate:
## 95 % quantile
## 7.33884
##
## bootstrap confidence interval
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level Normal Basic
## 95% ( 5.960, 9.557 ) ( 6.526, 9.312 )
##
## Level Percentile BCa
## 95% ( 5.366, 8.152 ) ( 5.543, 8.231 )
## Calculations and Intervals on Original Scale
##
## sample estimate:
## 95 % quantile
## 7.33884
Next, we consider the median.
##
## exact confidence interval
##
## 95.23684 percent confidence intervals:
## lower upper
## median 1.067520 2.505523
## median 1.351246 2.674241
##
## sample estimate:
## median
## 1.800282
##
## asymptotic confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## median 1.33018 2.548754
##
## sample estimate:
## median
## 1.800282
##
## bootstrap confidence interval
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level Normal Basic
## 95% ( 1.038, 2.438 ) ( 1.052, 2.260 )
##
## Level Percentile BCa
## 95% ( 1.341, 2.549 ) ( 1.333, 2.534 )
## Calculations and Intervals on Original Scale
##
## sample estimate:
## median
## 1.800282
It often happens that quantile confidence intervals are not unique. Here the minimum length interval might be of interest.
##
## minimum length exact confidence interval
##
## 95.23684 percent confidence interval:
## lower upper
## median 1.351246 2.674241
##
## sample estimate:
## median
## 1.800282
Finally, we take a look at MAD (median absolute deviation) where by default the standardized MAD is used (see function mad).
##
## exact confidence interval
##
## 95.23684 percent confidence intervals:
## lower upper
## MAD 1.565245 2.240502
## MAD 1.658633 2.427747
##
## sample estimate:
## MAD
## 1.907953
##
## asymptotic confidence interval
##
## 95 percent confidence interval:
## 2.5 % 97.5 %
## MAD 1.612338 2.280246
##
## sample estimate:
## MAD
## 1.907953
##
## bootstrap confidence interval
##
## BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
## Based on 9999 bootstrap replicates
##
## CALL :
## boot.ci(boot.out = boot.out, conf = 1 - alpha, type = bootci.type)
##
## Intervals :
## Level Normal Basic
## 95% ( 1.576, 2.203 ) ( 1.536, 2.182 )
##
## Level Percentile BCa
## 95% ( 1.634, 2.280 ) ( 1.623, 2.269 )
## Calculations and Intervals on Original Scale
##
## sample estimate:
## MAD
## 1.907953
##
## exact confidence interval
##
## 95.23684 percent confidence intervals:
## lower upper
## MAD 1.055743 1.511198
## MAD 1.118733 1.637493
##
## sample estimate:
## MAD
## 1.286897
Hsu Two-Sample t-Test
The Hsu two-sample t-test is an alternative to the Welch two-sample t-test using a different formula for computing the degrees of freedom of the respective t-distribution (Hedderich and Sachs 2018). The following code is taken and adapted from the help page of the t.test function.
##
## Welch Two Sample t-test
##
## data: 1:10 and c(7:20)
## t = -5.4349, df = 21.982, p-value = 1.855e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -11.052802 -4.947198
## sample estimates:
## mean of x mean of y
## 5.5 13.5
##
## Welch Two Sample t-test
##
## data: 1:10 and c(7:20, 200)
## t = -1.6329, df = 14.165, p-value = 0.1245
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -47.242900 6.376233
## sample estimates:
## mean of x mean of y
## 5.50000 25.93333
##
## Hsu Two Sample t-test
##
## data: 1:10 and c(7:20)
## t = -5.4349, df = 9, p-value = 0.0004137
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -11.329805 -4.670195
## sample estimates:
## mean of x mean of y SD of x SD of y
## 5.50000 13.50000 3.02765 4.18330
##
## Hsu Two Sample t-test
##
## data: 1:10 and c(7:20, 200)
## t = -1.6329, df = 9, p-value = 0.1369
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -48.740846 7.874179
## sample estimates:
## mean of x mean of y SD of x SD of y
## 5.50000 25.93333 3.02765 48.32253
##
## Welch Two Sample t-test
##
## data: extra[group == 1] and extra[group == 2]
## t = -1.8608, df = 17.776, p-value = 0.07939
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.3654832 0.2054832
## sample estimates:
## mean of x mean of y
## 0.75 2.33
##
## Hsu Two Sample t-test
##
## data: extra[group == 1] and extra[group == 2]
## t = -1.8608, df = 9, p-value = 0.09569
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.5007773 0.3407773
## sample estimates:
## mean of x mean of y SD of x SD of y
## 0.750000 2.330000 1.789010 2.002249
##
## Welch Two Sample t-test
##
## data: extra by group
## t = -1.8608, df = 17.776, p-value = 0.07939
## alternative hypothesis: true difference in means between group 1 and group 2 is not equal to 0
## 95 percent confidence interval:
## -3.3654832 0.2054832
## sample estimates:
## mean in group 1 mean in group 2
## 0.75 2.33
##
## Hsu Two Sample t-test
##
## data: extra by group
## t = -1.8608, df = 9, p-value = 0.09569
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.5007773 0.3407773
## sample estimates:
## mean of x mean of y SD of x SD of y
## 0.750000 2.330000 1.789010 2.002249
Bootstrap t-Test
One and two sample bootstrap t-tests with equal or unequal variances in the two sample case (Efron and Tibshirani 1993).
##
## Bootstrapped Welch Two Sample t-test
##
## data: 1:10 and c(7:20)
## bootstrapped p-value = 2e-04
## 95 percent bootstrap percentile confidence interval:
## -10.742857 -5.270714
##
## Results without bootstrap:
## t = -5.4349, df = 21.982, p-value = 1.855e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -11.052802 -4.947198
## sample estimates:
## mean of x mean of y
## 5.5 13.5
##
## Bootstrapped Welch Two Sample t-test
##
## data: 1:10 and c(7:20, 200)
## bootstrapped p-value = 0.0234
## 95 percent bootstrap percentile confidence interval:
## -46.635000 -5.966667
##
## Results without bootstrap:
## t = -1.6329, df = 14.165, p-value = 0.1245
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -47.242900 6.376233
## sample estimates:
## mean of x mean of y
## 5.50000 25.93333
##
## Bootstrapped Welch Two Sample t-test
##
## data: extra[group == 1] and extra[group == 2]
## bootstrapped p-value = 0.08541
## 95 percent bootstrap percentile confidence interval:
## -3.17 0.01
##
## Results without bootstrap:
## t = -1.8608, df = 17.776, p-value = 0.07939
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.3654832 0.2054832
## sample estimates:
## mean of x mean of y
## 0.75 2.33
##
## Bootstrapped Welch Two Sample t-test
##
## data: extra by group
## bootstrapped p-value = 0.07861
## 95 percent bootstrap percentile confidence interval:
## -3.1605 0.0300
##
## Results without bootstrap:
## t = -1.8608, df = 17.776, p-value = 0.07939
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.3654832 0.2054832
## sample estimates:
## mean in group 1 mean in group 2
## 0.75 2.33
Permutation t-Test
One and two sample permutation t-tests with equal (Efron and Tibshirani 1993) or unequal variances (Janssen 1997) in the two sample case.
##
## Permutation Welch Two Sample t-test
##
## data: 1:10 and c(7:20)
## (Monte-Carlo) permutation p-value = 1e-04
## 95 percent (Monte-Carlo) permutation percentile confidence interval:
## -12.400000 -3.657143
##
## Results without permutation:
## t = -5.4349, df = 21.982, p-value = 1.855e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -11.052802 -4.947198
## sample estimates:
## mean of x mean of y
## 5.5 13.5
##
## Permutation Welch Two Sample t-test
##
## data: 1:10 and c(7:20, 200)
## (Monte-Carlo) permutation p-value < 2.2e-16
## 95 percent (Monte-Carlo) permutation percentile confidence interval:
## -36.700000 1.966667
##
## Results without permutation:
## t = -1.6329, df = 14.165, p-value = 0.1245
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -47.242900 6.376233
## sample estimates:
## mean of x mean of y
## 5.50000 25.93333
##
## Permutation Welch Two Sample t-test
##
## data: extra[group == 1] and extra[group == 2]
## (Monte-Carlo) permutation p-value = 0.07961
## 95 percent (Monte-Carlo) permutation percentile confidence interval:
## -3.36 0.16
##
## Results without permutation:
## t = -1.8608, df = 17.776, p-value = 0.07939
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.3654832 0.2054832
## sample estimates:
## mean of x mean of y
## 0.75 2.33
##
## Permutation Welch Two Sample t-test
##
## data: extra by group
## (Monte-Carlo) permutation p-value = 0.07551
## 95 percent (Monte-Carlo) permutation percentile confidence interval:
## -3.34 0.16
##
## Results without permutation:
## t = -1.8608, df = 17.776, p-value = 0.07939
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.3654832 0.2054832
## sample estimates:
## mean in group 1 mean in group 2
## 0.75 2.33
In case of skewed distributions, one may use function p2ses to compute an alternative standardized effect size (SES) as proposed by Botta-Dukat (Botta-Dukát 2018).
## [1] 1.754212
Multiple Imputation t-Test
Function mi.t.test can be used to compute a multiple imputation t-test by applying the approch of Rubin (Rubin 1987) in combination with the adjustment of Barnard and Rubin (Barnard and Rubin 1999).
## Generate some data
set.seed(123)
x <- rnorm(25, mean = 1)
x[sample(1:25, 5)] <- NA
y <- rnorm(20, mean = -1)
y[sample(1:20, 4)] <- NA
pair <- c(rnorm(25, mean = 1), rnorm(20, mean = -1))
g <- factor(c(rep("yes", 25), rep("no", 20)))
D <- data.frame(ID = 1:45, variable = c(x, y), pair = pair, group = g)
## Use Amelia to impute missing values
library(Amelia)
## Lade nötiges Paket: Rcpp
## ##
## ## Amelia II: Multiple Imputation
## ## (Version 1.8.0, built: 2021-05-26)
## ## Copyright (C) 2005-2022 James Honaker, Gary King and Matthew Blackwell
## ## Refer to http://gking.harvard.edu/amelia/ for more information
## ##
##
## Welch Two Sample t-test
##
## data: variable by group
## t = -6.9214, df = 33.69, p-value = 5.903e-08
## alternative hypothesis: true difference in means between group no and group yes is not equal to 0
## 95 percent confidence interval:
## -2.628749 -1.435125
## sample estimates:
## mean in group no mean in group yes
## -1.0862469 0.9456901
##
## Multiple Imputation Welch Two Sample t-test
##
## data: Variable variable: group no vs group yes
## t = -6.6457, df = 25.142, p-value = 5.63e-07
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.668011 -1.405862
## sample estimates:
## mean (no) SD (no) mean (yes) SD (yes)
## -1.0809584 0.9379070 0.9559786 1.0203219
##
## Two Sample t-test
##
## data: variable by group
## t = -6.8168, df = 34, p-value = 7.643e-08
## alternative hypothesis: true difference in means between group no and group yes is not equal to 0
## 95 percent confidence interval:
## -2.637703 -1.426171
## sample estimates:
## mean in group no mean in group yes
## -1.0862469 0.9456901
##
## Multiple Imputation Two Sample t-test
##
## data: Variable variable: group no vs group yes
## t = -6.5736, df = 26.042, p-value = 5.656e-07
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.673828 -1.400046
## sample estimates:
## mean (no) SD (no) mean (yes) SD (yes)
## -1.0809584 0.9379070 0.9559786 1.0203219
##
## Multiple Imputation Welch Two Sample t-test
##
## data: Variable variable: group no vs group yes
## t = -6.6457, df = 25.142, p-value = 2.815e-07
## alternative hypothesis: true difference in means is less than 0
## 95 percent confidence interval:
## -Inf -1.5135
## sample estimates:
## mean (no) SD (no) mean (yes) SD (yes)
## -1.0809584 0.9379070 0.9559786 1.0203219
##
## Multiple Imputation Welch Two Sample t-test
##
## data: Variable variable: group no vs group yes
## t = -6.6457, df = 25.142, p-value = 1
## alternative hypothesis: true difference in means is greater than 0
## 95 percent confidence interval:
## -2.560374 Inf
## sample estimates:
## mean (no) SD (no) mean (yes) SD (yes)
## -1.0809584 0.9379070 0.9559786 1.0203219
##
## One Sample t-test
##
## data: D$variable[D$group == "yes"]
## t = 4.5054, df = 19, p-value = 0.0002422
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 0.5063618 1.3850184
## sample estimates:
## mean of x
## 0.9456901
##
## Multiple Imputation One Sample t-test
##
## data: Variable variable
## t = 4.6847, df = 18.494, p-value = 0.0001725
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
## 0.5280752 1.3838820
## sample estimates:
## mean SD
## 0.9559786 1.0203219
##
## Multiple Imputation One Sample t-test
##
## data: Variable variable
## t = 9.5851, df = 18.494, p-value = 1
## alternative hypothesis: true mean is less than -1
## 95 percent confidence interval:
## -Inf 1.309328
## sample estimates:
## mean SD
## 0.9559786 1.0203219
##
## Multiple Imputation One Sample t-test
##
## data: Variable variable
## t = 9.5851, df = 18.494, p-value = 6.655e-09
## alternative hypothesis: true mean is greater than -1
## 95 percent confidence interval:
## 0.6026297 Inf
## sample estimates:
## mean SD
## 0.9559786 1.0203219
##
## Paired t-test
##
## data: D$variable and D$pair
## t = -1.3532, df = 35, p-value = 0.1847
## alternative hypothesis: true mean difference is not equal to 0
## 95 percent confidence interval:
## -0.6976921 0.1395993
## sample estimates:
## mean difference
## -0.2790464
##
## Multiple Imputation Paired t-test
##
## data: Variables variable and pair
## t = -1.0455, df = 39.974, p-value = 0.3021
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.5680741 0.1807237
## sample estimates:
## mean of difference SD of difference
## -0.1936752 1.2426515
Function mi.t.test also works with package mice by applying function mids2datlist of package miceadds.
##
## Attache Paket: 'mice'
## Das folgende Objekt ist maskiert 'package:stats':
##
## filter
## Die folgenden Objekte sind maskiert von 'package:base':
##
## cbind, rbind
## * miceadds 3.13-12 (2022-05-30 15:14:07)
##
## Multiple Imputation Welch Two Sample t-test
##
## data: Variable variable: group no vs group yes
## t = -7.2111, df = 29.735, p-value = 5.289e-08
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.6351 -1.4716
## sample estimates:
## mean (no) SD (no) mean (yes) SD (yes)
## -1.0771373 0.8585436 0.9762127 1.0261920
Repeated Measures One-Way ANOVA
We provide a simple wrapper function that allows to compute a classical repeated measures one-way ANOVA as well as a respective mixed effects model. In addition, the non-parametric Friedman and Quade tests can be computed.
##
## Repeated measures 1-way ANOVA
##
## data: outcome , timepoints and patients
## F = 6.5898, num df = 2, denom df = 18, p-value = 0.007122
##
## Mixed-effects 1-way ANOVA
##
## data: outcome , timepoints and patients
## F = 7.3674, num df = 2, denom df = 18, p-value = 0.004596
##
## Friedman rank sum test
##
## data: x, g and id
## Friedman chi-squared = 12.8, df = 2, p-value = 0.001662
##
## Quade test
##
## data: x, g and id
## Quade F = 7.2906, num df = 2, denom df = 18, p-value = 0.004795
Volcano Plots
Volcano plots can be used to visualize the results when many tests have been applied. They show a measure of effect size in combination with (adjusted) p values.
## Generate some data
x <- matrix(rnorm(1000, mean = 10), nrow = 10)
g1 <- rep("control", 10)
y1 <- matrix(rnorm(500, mean = 11.75), nrow = 10)
y2 <- matrix(rnorm(500, mean = 9.75, sd = 3), nrow = 10)
g2 <- rep("treatment", 10)
group <- factor(c(g1, g2))
Data <- rbind(x, cbind(y1, y2))
## compute Hsu t-test
pvals <- apply(Data, 2, function(x, group) hsu.t.test(x ~ group)$p.value,
group = group)
## compute log-fold change
logfc <- function(x, group){
res <- tapply(x, group, mean)
log2(res[1]/res[2])
}
lfcs <- apply(Data, 2, logfc, group = group)
volcano(lfcs, p.adjust(pvals, method = "fdr"),
effect.low = -0.25, effect.high = 0.25,
xlab = "log-fold change", ylab = "-log10(adj. p value)")
Imputation of Standard Deviations for Changes from Baseline
The function imputeSD can be used to impute standard deviations for changes from baseline adopting the approach of the Cochrane handbook (Higgins et al. 2019, Section 6.5.2.8).
SD1 <- c(0.149, 0.022, 0.036, 0.085, 0.125, NA, 0.139, 0.124, 0.038)
SD2 <- c(NA, 0.039, 0.038, 0.087, 0.125, NA, 0.135, 0.126, 0.038)
SDchange <- c(NA, NA, NA, 0.026, 0.058, NA, NA, NA, NA)
imputeSD(SD1, SD2, SDchange)
## SD1 SD2 SDchange Cor SDchange.min SDchange.mean SDchange.max
## 1 0.149 0.149 NA NA 0.04491608 0.05829769 0.06913600
## 2 0.022 0.039 NA NA 0.01915642 0.02050235 0.02176520
## 3 0.036 0.038 NA NA 0.01132754 0.01460887 0.01727787
## 4 0.085 0.087 0.026 0.9545639 0.02600000 0.02600000 0.02600000
## 5 0.125 0.125 0.058 0.8923520 0.05800000 0.05800000 0.05800000
## 6 NA NA NA NA NA NA NA
## 7 0.139 0.135 NA NA 0.04148755 0.05374591 0.06368696
## 8 0.124 0.126 NA NA 0.03773311 0.04894677 0.05803262
## 9 0.038 0.038 NA NA 0.01145511 0.01486787 0.01763200
Correlations can also be provided via argument corr. This option may particularly be useful, if no complete data is available.
## SD1 SD2 SDchange Cor SDchange.min SDchange.mean SDchange.max
## 1 0.149 0.149 NA NA 0.04711794 0.06663483 0.08161066
## 2 0.022 0.039 NA NA 0.01935975 0.02146159 0.02337520
## 3 0.036 0.038 NA NA 0.01186592 0.01666133 0.02035682
## 4 0.085 0.087 NA NA 0.02726720 0.03850974 0.04714340
## 5 0.125 0.125 NA NA 0.03952847 0.05590170 0.06846532
## 6 NA NA NA NA NA NA NA
## 7 0.139 0.135 NA NA 0.04350287 0.06139218 0.07513654
## 8 0.124 0.126 NA NA 0.03957777 0.05593568 0.06849234
## 9 0.038 0.038 NA NA 0.01201666 0.01699412 0.02081346
Pairwise Comparisons
Function pairwise.fun enables the application of arbitrary functions for pairwise comparisons.
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
## Warning in wilcox.test.default(xi, xj, paired = paired, ...): kann bei Bindungen
## keinen exakten p-Wert Berechnen
##
## Pairwise comparisons using Wilcoxon rank sum test with continuity correction
##
## data: airquality$Ozone and airquality$Month
##
## 5 6 7 8
## 6 0.19250 - - -
## 7 3e-05 0.01414 - -
## 8 0.00012 0.02591 0.86195 -
## 9 0.11859 0.95887 0.00074 0.00325
##
## P value adjustment method: none
## Package 'exactRankTests' is no longer under development.
## Please consider using package 'coin' instead.
## 5 vs 6 5 vs 7 5 vs 8 5 vs 9 6 vs 7 6 vs 8
## 1.897186e-01 1.087205e-05 6.108735e-05 1.171796e-01 1.183753e-02 2.333564e-02
## 6 vs 9 7 vs 8 7 vs 9 8 vs 9
## 9.528659e-01 8.595683e-01 5.281796e-04 2.714694e-03
sessionInfo
## R version 4.2.1 Patched (2022-08-19 r82735)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Linux Mint 20.3
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/libf77blas.so.3.10.3
## LAPACK: /home/kohlm/RTOP/Rbranch/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=de_DE.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=de_DE.UTF-8 LC_COLLATE=C
## [5] LC_MONETARY=de_DE.UTF-8 LC_MESSAGES=de_DE.UTF-8
## [7] LC_PAPER=de_DE.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=de_DE.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] exactRankTests_0.8-35 miceadds_3.13-12 mice_3.14.0
## [4] Amelia_1.8.0 Rcpp_1.0.9 MKinfer_0.7
##
## loaded via a namespace (and not attached):
## [1] tidyselect_1.1.2 xfun_0.31 bslib_0.4.0 mitools_2.4
## [5] purrr_0.3.4 lattice_0.20-45 colorspace_2.0-3 vctrs_0.4.1
## [9] generics_0.1.3 htmltools_0.5.3 yaml_2.3.5 gmp_0.6-5
## [13] utf8_1.2.2 rlang_1.0.4 jquerylib_0.1.4 pillar_1.8.0
## [17] foreign_0.8-82 glue_1.6.2 withr_2.5.0 DBI_1.1.3
## [21] arrangements_1.1.9 lifecycle_1.0.1 stringr_1.4.0 munsell_0.5.0
## [25] gtable_0.3.0 evaluate_0.15 labeling_0.4.2 knitr_1.39
## [29] fastmap_1.1.0 fansi_1.0.3 highr_0.9 broom_1.0.0
## [33] scales_1.2.0 backports_1.4.1 cachem_1.0.6 jsonlite_1.8.0
## [37] farver_2.1.1 ggplot2_3.3.6 digest_0.6.29 stringi_1.7.8
## [41] dplyr_1.0.9 grid_4.2.1 cli_3.3.0 tools_4.2.1
## [45] magrittr_2.0.3 sass_0.4.2 tibble_3.1.8 MKdescr_0.7
## [49] tidyr_1.2.0 pkgconfig_2.0.3 assertthat_0.2.1 rmarkdown_2.14
## [53] R6_2.5.1 boot_1.3-28 nlme_3.1-158 compiler_4.2.1
References
Altman, Douglas G., David Machin, Trevor N. Bryant, and Martin J. Gardner, eds. 2000. Statistics with Confidence: Confidence Intervals and Statistical Guidelines. BMJ Books.
Efron, Bradley, and Robert J Tibshirani. 1993. An Introduction to the Bootstrap (Chapman & Hall/Crc Monographs on Statistics and Applied Probability). Chapman; Hall/CRC.
Gulhar, Monika, B M Golam Kibria, and Nasar Ahmed. 2012. “A Comparison of Some Confidence Intervals for Estimating the Population Coefficient of Variation: A Simulation Study.” SORT 36 (January): 45–68.
Higgins, Julian P.T., James Thomas, Jacqueline Chandler, Miranda Cumpston, Tianjing Li, Matthew J. Page, and Vivian A. Welch, eds. 2019. Cochrane Handbook for Systematic Reviews of Interventions. Wiley. https://doi.org/10.1002/9781119536604.
———. 1998b. “Improved confidence intervals for the difference between binomial proportions based on paired data.” Stat Med 17 (22): 2635–50.
Rubin, Donald B. 1987. Multiple Imputation for Nonresponse in Surveys (Wiley Series in Probability and Statistics). Wiley.