The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
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[@]metanet.ch.

Testing type-I error rates

A valid statistical test is one where the amount of type-I errors (rejecting the null when it should not) does not exceed the test threshold (often 5%). We say it in bold because it seems to be commonly forgotten. Logistic regressions applied to proportions deviates massively from this rule, becoming very liberal tests.

ANOPA respect very closely this rule when the corrected statistics are consulted.

The following code shows how to test the type-I error rate for ANOPA.It uses simulated data with no differences. The only limitations herein is that no correlations is added when repeated measures are used.

Note that the simulations are actually not run in this vignette, as they take times. We wished to provide code in case you wished to test type-I error rate by yourself. The present code is also not optimized for speed (and in particular, is not parallelized); we wished to keep the code as simple as possible for readability. In all the following, the number of simulated participants per group is small (20) but can be varied.

Simulations with a single factor

Simulation with one between factor

frm      <- s ~ grp         # the formula
BSDesign <- list(grp = c("ctrl","plcbo")) #one factor, two groups
thePs    <- c(0.3, 0.3)     # the true proportions, equal

# test type-I error rate when no effect as is the case for factor 2
set.seed(41)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, BSDesign )
    w   <- anopa(frm, smp[,2:3] )
    res <- c(res, if(summarize(w)[1,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design B,     testing B: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulation with one within factor

frm      <- cbind(s.early, s.middle, s.late) ~ .
WSDesign <- list(moment = c("early","middle","late"))
thePs    <- c(0.3, 0.3, 0.3)

# test type-I error rate when no effect as is the case for factor 2
set.seed(42)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, NULL, WSDesign )
    w   <- anopa(frm, smp[,2:4] , WSFactors = "M(3)" )
    res <- c(res, if(summarize(w)[1,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design W,     testing W: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulations with two factors

Simulation with two factors, between design

frm      <- s ~ grp * eta
WSDesign <- list()
BSDesign <- list(eta = c("repue","ajun"), grp = c("early","middle","late"))
thePs    <- c(0.3, 0.3, 0.5, 0.5, 0.7, 0.7)

# test type-I error rate when no effect as is the case for factor 2
set.seed(41)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, BSDesign )
    w   <- anopa(frm, smp )
    res <- c(res, if(summarize(w)[2,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design BxB,   testing B: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulation with two factors, within design

frm      <- cbind(s.repue.early, s.ajun.early, 
                  s.repue.middle, s.ajun.middle, 
                  s.repue.late, s.ajun.late) ~ .
BSDesign <- list()
WSDesign <- list(eta = c("repue","ajun"), moment = c("early","middle","late"))
thePs    <- c(0.3, 0.3, 0.5, 0.5, 0.7, 0.7)
# thePs    <- c(0.3, 0.7, 0.3, 0.7, 0.3, 0.7) # or no effect on factor 1

# test type-I error rate when no effect as is the case for factor 2
set.seed(41)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, NULL, WSDesign )
    w   <- anopa(frm, smp, WSFactors = c("e(2)", "m(3)") )
    res <- c(res, if(summarize(w)[2,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design WxW,   testing W: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulation with two factors, mixed design

frm      <- cbind(s.early, s.middle, s.late) ~ grp
BSDesign <- list(grp = c("ctrl","plcbo"))
WSDesign <- list(moment = c("early","middle","late"))
thePs    <- c(0.3, 0.3, 0.5, 0.5, 0.7, 0.7)
# thePs    <- c(0.3, 0.7, 0.3, 0.7, 0.3, 0.7) # or no effect on factor 1

# test type-I error rate when no effect as is the case for factor 2
set.seed(41)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, BSDesign, WSDesign )
    w   <- anopa(frm, smp[,2:5] , WSFactors = "M(3)")
    res <- c(res, if(summarize(w)[2,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design WxB,   testing B: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulations with three factors

Simulation with three factors, all between design

frm      <- s ~ grp * eta * a
BSDesign <- list(eta = c("repue","ajun"), 
                 grp = c("early","middle","late"), a = c("1","2","3","4"))
thePs    <- rep(0.3, 24)

# test type-I error rate when no effect as is the case for factor 2
set.seed(41)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, BSDesign )
    w   <- anopa(frm, smp )
    res <- c(res, if(summarize(w)[2,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design BxBxB, testing B: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulation with three factors, within design

frm      <- cbind(s.repue.early.1,  s.ajun.early.1,  s.repue.middle.1,
            s.ajun.middle.1, s.repue.late.1,   s.ajun.late.1,   s.repue.early.2, 
            s.ajun.early.2,  s.repue.middle.2, s.ajun.middle.2, s.repue.late.2,
            s.ajun.late.2,   s.repue.early.3,  s.ajun.early.3,  s.repue.middle.3,
            s.ajun.middle.3, s.repue.late.3,   s.ajun.late.3,   s.repue.early.4, 
            s.ajun.early.4,  s.repue.middle.4, s.ajun.middle.4, s.repue.late.4,
            s.ajun.late.4 ) ~ .
WSDesign <- list(eta = c("repue","ajun"), grp = c("early","middle","late"), a = c("1","2","3","4"))
thePs    <- rep(0.3, 24)

# test type-I error rate when no effect as is the case for factor 2
set.seed(43)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, NULL, WSDesign )
    w   <- anopa(frm, smp, WSFactors = c("e(2)","g(3)", "a(4)") )
    res <- c(res, if(summarize(w)[2,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design WxWxW, testing W: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulation with three factors, mixed design, testing within

frm      <- cbind(s.repue.early,  s.ajun.early,  s.repue.middle,
s.ajun.middle, s.repue.late,  s.ajun.late ) ~ a
BSDesign <- list( a = c("1","2","3","4") )
WSDesign <- list(eta = c("repue","ajun"), grp = c("early","middle","late") )
thePs    <- rep(0.3, 24)

# test type-I error rate when no effect as is the case for factor 2
set.seed(43)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, BSDesign, WSDesign )
    w   <- anopa(frm, smp, WSFactors = c("e(2)","g(3)") )
    res <- c(res, if(summarize(w)[1,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design BxWxW, testing W: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

Simulation with three factors, mixed design, testing between

frm      <- cbind(s.repue.early,  s.ajun.early,  s.repue.middle,
s.ajun.middle, s.repue.late,  s.ajun.late ) ~ a
BSDesign <- list( a = c("1","2","3","4") )
WSDesign <- list(eta = c("repue","ajun"), grp = c("early","middle","late") )
thePs    <- rep(0.3, 24)

# test type-I error rate when no effect as is the case for factor 2
set.seed(42)
res <- c()
for (i in 1:nsim) {
    smp <- GRP( thePs, theN, BSDesign, WSDesign )
    w   <- anopa(frm, smp, WSFactors = c("e(2)","g(3)") )
    res <- c(res, if(summarize(w)[3,4]<.05) 1 else 0)
}
typeI <- mean(res)
cat( "Design BxWxW, testing B: ", typeI, "\n")

# tolerance is large as the number of simulations is small
expect_equal( typeI, .05, tolerance = 0.035)

The end

Let’s restore the warnings and messages before leaving:

options("ANOPA.feedback" = 'all')

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.