The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

Validation

Jakob Schöpe

2021-03-19

Test Suite

The following fictitious example of a prospective cohort study will be used to validate the correct estimation of the BSW package in R.

Exposed Non-Exposed
Cases 200 50
Non-Cases 50 200 
library(testthat)
library(BSW)
Lade nötiges Paket: Matrix
Lade nötiges Paket: matrixStats
Lade nötiges Paket: quadprog
df <- data.frame(y = rep(c(0, 1), each = 250), 
                 x = rep(c(0, 1, 0, 1), times = c(200, 50, 50, 200))
                 )
RR <- (200 * 250) / (50 * 250)
SE <- sqrt((1/200 + 1/50) - (1/250 + 1/250))
fit <- bsw(y ~ x, df)
out <- summary(fit)
Call:
bsw(formula = y ~ x, data = df)

Convergence: TRUE
Coefficients:
             Estimate Std. Error   z value     Pr(>|z|)  RR      2.5%    97.5%
(Intercept) -1.609438  0.1264911 -12.72372 4.365487e-37 0.2 0.1560848 0.256271
x            1.386294  0.1303840  10.63239 2.106386e-26 4.0 3.0979677 5.164676

Iterations: 37

The relative risk for exposed individuals compared to non-exposed individuals can be calculated from

\(RR = \displaystyle\frac{200 * 250}{50*250} = 4\). 
test_that(desc = "Estimated relative risk is equal to 4",
          code = {
                  expect_equal(object = unname(exp(coef(fit)[2])),
                               expected = RR)
            }
          )
Test passed 🥇

The standard error of the natural logarithm of the relative risk can be calculated from

\(SE(ln(RR)) = \displaystyle\sqrt{\Big(\frac{1}{200} + \frac{1}{50}\Big) - \Big(\frac{1}{250}+\frac{1}{250}\Big)} = 0.130384\). 
test_that(desc = "Estimated standard error is equal to 0.1303840",
          code = {
                  expect_equal(object = unname(out$std.err[2]), 
                               expected = SE)
            }
          )
Test passed 🌈

The z-value can be calculated from

\(z = \displaystyle\frac{1.386294}{0.130384} = 10.63239\). 
test_that(desc = "Estimated z-value is equal to 10.63239",
          code = {
                  expect_equal(object = unname(out$z.value[2]), 
                               expected = log(RR) / SE)
            }
          )
Test passed 😀

The 95% confidence interval limits can be calculated from

\(exp(1.386294 \pm 1.959964 * 0.1303840) = [3.097968; 5.164676]\). 
test_that(desc = "Estimated 95% confidence interval limits are equal to 3.097968 and 5.164676",
          code = {
                  expect_equal(object = unname(exp(confint(fit)[2,])), 
                               expected = exp(log(RR) + SE * qnorm(c(0.025, 0.975))))
            }
          )
Test passed 😀

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.