In situations where researchers are unaware which variables or combinations of variables influence a given outcome, and when computing the relevant models is very expensive, Monte Carlo simulations can be used. In this case computer-driven experimental designs can be constructed with the goal of detecting (or screening) which factors affect a given outcome variable through statistical inference. Montgomery (2017) describes this class of screening experiments as “… tests in which many factors are considered and the objective is to identify those factors that have significant effects”, which are typically relevant in early stages of a given project.
The most common designs implemented (and indeed, the default used by ) are those which fully-cross the factor variables, leading to a set of \(C_1\times C_2 \times \cdots \times C_n\) unique conditions to study. Fully-crossed designs are beneficial when detecting main effects, two-way interactions, and up to \(n\) -way interactions (provided that the number of replications per unique condition is greater than 1).
For example, consider the following fully-crossed experimental design that combines five simulation factors in the context of an independent samples \(t\) -test’s ability to detect true mean differences: sample size, the magnitude of the mean difference, group variances, distributions, whether the group sizes were equal or unequal, if the distributions were symmetric/asymmetric.
Fully-crossed design example
library (SimDesign)
Design <- createDesign (sample_size= c (100 ,200 ),
mean_diff= c (.25 , 2 ),
variance.ratio= c (1 ,4 ),
equal_size= c (TRUE , FALSE ),
dists= c ('norm' , 'skew' ),
same_dists= c (TRUE , FALSE ),
symmetric= c (TRUE , FALSE ),
# remove same-normal combo
subset = ! (symmetric & dists == 'norm' ))
Design
# A tibble: 96 × 7
sample_size mean_diff variance.ratio equal_size dists same_dists symmetric
<dbl> <dbl> <dbl> <lgl> <chr> <lgl> <lgl>
1 100 0.25 1 TRUE skew TRUE TRUE
2 200 0.25 1 TRUE skew TRUE TRUE
3 100 2 1 TRUE skew TRUE TRUE
4 200 2 1 TRUE skew TRUE TRUE
5 100 0.25 4 TRUE skew TRUE TRUE
6 200 0.25 4 TRUE skew TRUE TRUE
7 100 2 4 TRUE skew TRUE TRUE
8 200 2 4 TRUE skew TRUE TRUE
9 100 0.25 1 FALSE skew TRUE TRUE
10 200 0.25 1 FALSE skew TRUE TRUE
# ℹ 86 more rows
#-------------------------------------------------------------------
Generate <- function (condition, fixed_objects = NULL ) {
Attach (condition)
N1 <- ifelse (equal_size, sample_size/ 2 , sample_size/ 4 )
N2 <- sample_size - N1
DV1 <- if (dists == 'norm' ) rnorm (N1, sd= 1 ) else rchisq (N1, df= 5 ) - 5
if (dists == 'norm' && same_dists){
DV2 <- rnorm (N2, sd= sqrt (variance.ratio))
} else {
DV2 <- (rchisq (N2, df= 5 ) - 5 ) * sqrt (variance.ratio)
if (! symmetric) DV2 <- abs (0 - DV2)
}
dat <- data.frame (DV= c (DV1, DV2 + mean_diff),
group= c (rep ('g1' , N1), rep ('g2' , N2)))
dat
}
Analyse <- function (condition, dat, fixed_objects = NULL ) {
test <- t.test (DV ~ group, data= dat, var.equal= TRUE )
ret <- c (detect= test$ p.value < .05 )
ret
}
#-------------------------------------------------------------------
# note the omission of the summarise function
res <- runSimulation (design= Design, replications= 20 ,
generate= Generate, analyse= Analyse, verbose= FALSE )
res
# A tibble: 1,920 × 8
sample_size mean_diff variance.ratio equal_size dists same_dists symmetric
<dbl> <dbl> <dbl> <lgl> <chr> <lgl> <lgl>
1 100 0.25 1 TRUE skew TRUE TRUE
2 100 0.25 1 TRUE skew TRUE TRUE
3 100 0.25 1 TRUE skew TRUE TRUE
4 100 0.25 1 TRUE skew TRUE TRUE
5 100 0.25 1 TRUE skew TRUE TRUE
6 100 0.25 1 TRUE skew TRUE TRUE
7 100 0.25 1 TRUE skew TRUE TRUE
8 100 0.25 1 TRUE skew TRUE TRUE
9 100 0.25 1 TRUE skew TRUE TRUE
10 100 0.25 1 TRUE skew TRUE TRUE
# ℹ 1,910 more rows
# ℹ 1 more variable: detect <lgl>
Given these collected results, ANOVA-based approaches can be employed to investigate main and interaction effects for the recorded outcome variable.
# main effects only
mod <- aov (detect ~ ., data= res)
summary (mod)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 0.2 0.2 2.42 0.120
mean_diff 1 62.7 62.7 661.28 < 2e-16 ***
variance.ratio 1 3.3 3.3 34.28 5.6e-09 ***
equal_size 1 0.6 0.6 5.98 0.015 *
dists 1 4.6 4.6 48.57 4.4e-12 ***
same_dists 1 6.2 6.2 65.25 1.2e-15 ***
symmetric 1 107.5 107.5 1133.88 < 2e-16 ***
Residuals 1912 181.3 0.1
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# all two-way effects
mod2 <- aov (detect ~ (.)^ 2 , data= res)
summary (mod2)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 0.2 0.2 3.87 0.04929 *
mean_diff 1 62.7 62.7 1056.73 < 2e-16 ***
variance.ratio 1 3.3 3.3 54.77 2.0e-13 ***
equal_size 1 0.6 0.6 9.56 0.00202 **
dists 1 4.6 4.6 77.62 < 2e-16 ***
same_dists 1 6.2 6.2 104.27 < 2e-16 ***
symmetric 1 107.5 107.5 1811.94 < 2e-16 ***
sample_size:mean_diff 1 0.0 0.0 0.01 0.92537
sample_size:variance.ratio 1 0.0 0.0 0.22 0.63955
sample_size:equal_size 1 0.0 0.0 0.08 0.77871
sample_size:dists 1 0.0 0.0 0.36 0.55112
sample_size:same_dists 1 0.1 0.1 1.48 0.22343
sample_size:symmetric 1 0.1 0.1 1.07 0.30191
mean_diff:variance.ratio 1 0.0 0.0 0.43 0.51205
mean_diff:equal_size 1 0.6 0.6 9.56 0.00202 **
mean_diff:dists 1 0.0 0.0 0.74 0.38926
mean_diff:same_dists 1 7.9 7.9 132.77 < 2e-16 ***
mean_diff:symmetric 1 39.6 39.6 666.44 < 2e-16 ***
variance.ratio:equal_size 1 0.3 0.3 4.64 0.03131 *
variance.ratio:dists 1 0.1 0.1 1.58 0.20833
variance.ratio:same_dists 1 0.1 0.1 1.97 0.16012
variance.ratio:symmetric 1 2.9 2.9 48.98 3.6e-12 ***
equal_size:dists 1 0.2 0.2 3.20 0.07385 .
equal_size:same_dists 1 0.0 0.0 0.08 0.77871
equal_size:symmetric 1 0.7 0.7 11.07 0.00089 ***
dists:same_dists 1 16.4 16.4 276.45 < 2e-16 ***
same_dists:symmetric 1 0.1 0.1 2.22 0.13598
Residuals 1892 112.3 0.1
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# all three-way effects
mod3 <- aov (detect ~ (.)^ 3 , data= res)
summary (mod3)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 0.2 0.2 4.66 0.03094 *
mean_diff 1 62.7 62.7 1273.28 < 2e-16 ***
variance.ratio 1 3.3 3.3 66.00 8.1e-16 ***
equal_size 1 0.6 0.6 11.52 0.00070 ***
dists 1 4.6 4.6 93.53 < 2e-16 ***
same_dists 1 6.2 6.2 125.64 < 2e-16 ***
symmetric 1 107.5 107.5 2183.25 < 2e-16 ***
sample_size:mean_diff 1 0.0 0.0 0.01 0.91811
sample_size:variance.ratio 1 0.0 0.0 0.26 0.60720
sample_size:equal_size 1 0.0 0.0 0.10 0.75774
sample_size:dists 1 0.0 0.0 0.43 0.51292
sample_size:same_dists 1 0.1 0.1 1.79 0.18144
sample_size:symmetric 1 0.1 0.1 1.28 0.25715
mean_diff:variance.ratio 1 0.0 0.0 0.52 0.47172
mean_diff:equal_size 1 0.6 0.6 11.52 0.00070 ***
mean_diff:dists 1 0.0 0.0 0.89 0.34464
mean_diff:same_dists 1 7.9 7.9 159.98 < 2e-16 ***
mean_diff:symmetric 1 39.6 39.6 803.01 < 2e-16 ***
variance.ratio:equal_size 1 0.3 0.3 5.59 0.01812 *
variance.ratio:dists 1 0.1 0.1 1.91 0.16727
variance.ratio:same_dists 1 0.1 0.1 2.38 0.12312
variance.ratio:symmetric 1 2.9 2.9 59.02 2.5e-14 ***
equal_size:dists 1 0.2 0.2 3.85 0.04976 *
equal_size:same_dists 1 0.0 0.0 0.10 0.75774
equal_size:symmetric 1 0.7 0.7 13.34 0.00027 ***
dists:same_dists 1 16.4 16.4 333.11 < 2e-16 ***
same_dists:symmetric 1 0.1 0.1 2.68 0.10174
sample_size:mean_diff:variance.ratio 1 0.1 0.1 2.38 0.12312
sample_size:mean_diff:equal_size 1 0.0 0.0 0.01 0.91811
sample_size:mean_diff:dists 1 0.3 0.3 5.08 0.02430 *
sample_size:mean_diff:same_dists 1 0.1 0.1 2.38 0.12312
sample_size:mean_diff:symmetric 1 0.1 0.1 2.68 0.10174
sample_size:variance.ratio:equal_size 1 0.0 0.0 0.86 0.35483
sample_size:variance.ratio:dists 1 0.0 0.0 0.89 0.34464
sample_size:variance.ratio:same_dists 1 0.0 0.0 0.26 0.60720
sample_size:variance.ratio:symmetric 1 0.0 0.0 0.14 0.70560
sample_size:equal_size:dists 1 0.1 0.1 2.33 0.12693
sample_size:equal_size:same_dists 1 0.0 0.0 0.86 0.35483
sample_size:equal_size:symmetric 1 0.0 0.0 0.14 0.70560
sample_size:dists:same_dists 1 0.1 0.1 1.53 0.21656
sample_size:same_dists:symmetric 1 0.0 0.0 0.02 0.89979
mean_diff:variance.ratio:equal_size 1 0.0 0.0 0.52 0.47172
mean_diff:variance.ratio:dists 1 1.2 1.2 23.73 1.2e-06 ***
mean_diff:variance.ratio:same_dists 1 0.9 0.9 17.78 2.6e-05 ***
mean_diff:variance.ratio:symmetric 1 0.5 0.5 9.91 0.00167 **
mean_diff:equal_size:dists 1 0.4 0.4 8.04 0.00462 **
mean_diff:equal_size:same_dists 1 0.1 0.1 1.28 0.25813
mean_diff:equal_size:symmetric 1 1.1 1.1 21.71 3.4e-06 ***
mean_diff:dists:same_dists 1 14.6 14.6 296.98 < 2e-16 ***
mean_diff:same_dists:symmetric 1 0.0 0.0 0.02 0.89979
variance.ratio:equal_size:dists 1 0.0 0.0 0.13 0.71622
variance.ratio:equal_size:same_dists 1 0.0 0.0 0.26 0.60720
variance.ratio:equal_size:symmetric 1 0.3 0.3 5.73 0.01681 *
variance.ratio:dists:same_dists 1 0.5 0.5 10.71 0.00109 **
variance.ratio:same_dists:symmetric 1 0.0 0.0 0.14 0.70560
equal_size:dists:same_dists 1 0.0 0.0 0.43 0.51292
equal_size:same_dists:symmetric 1 0.0 0.0 0.78 0.37810
Residuals 1862 91.7 0.0
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
and so on.
Fractional factorial experiments
When interest is primarily in detecting main and two-way interaction effects only then substantial gains can be made by using fractional factorial designs. A fractional factorial design is a reduced version of the fully factorial design, allowing for a more efficient use of resources as it reduces the sample size and number of conditions studied by way of discarding higher-order factor combinations.
Note that to use the FrF2
package each factor must consist of exactly two-levels.
library (FrF2)
# help(FrF2)
# 7 factors with resolution 5
fr <- FrF2 (nfactors= 7 , resolution= 5 )
fr[1 : 6 ,]
A B C D E F G
1 1 -1 -1 1 1 1 1
2 1 -1 -1 1 -1 -1 1
3 1 -1 1 1 1 -1 1
4 1 1 1 1 1 1 1
5 -1 1 -1 -1 -1 -1 -1
6 1 1 1 -1 -1 -1 -1
# Create working simulation design given -1/1 combinations
fDesign <- createDesign (sample_size= c (100 ,200 ),
mean_diff= c (.25 , 2 ),
variance.ratio= c (1 ,4 ),
equal_size= c (TRUE , FALSE ),
dists= c ('norm' , 'skew' ),
same_dists= c (TRUE , FALSE ),
symmetric= c (TRUE , FALSE ),
# remove same-normal combo
subset = ! (symmetric & dists == 'norm' ),
fractional= fr)
fDesign
# A tibble: 48 × 7
sample_size mean_diff variance.ratio equal_size dists same_dists symmetric
<dbl> <dbl> <dbl> <lgl> <chr> <lgl> <lgl>
1 100 0.25 1 TRUE skew FALSE TRUE
2 200 2 1 FALSE skew FALSE FALSE
3 100 2 4 FALSE skew FALSE FALSE
4 100 2 1 TRUE skew TRUE TRUE
5 100 0.25 4 FALSE norm FALSE FALSE
6 200 0.25 1 TRUE norm TRUE FALSE
7 200 2 4 FALSE skew TRUE FALSE
8 200 2 4 TRUE skew TRUE TRUE
9 200 2 4 FALSE norm FALSE FALSE
10 100 0.25 4 FALSE skew TRUE FALSE
# ℹ 38 more rows
Notice that the fractional factor design has fewer condition than the fully-crossed version above. However, despite this reduced number of conditions valid main and two-way interaction effects are still estimable.
# note the omission of the summarise function
fres <- runSimulation (design= fDesign, replications= 20 ,
generate= Generate, analyse= Analyse, verbose= FALSE )
fres
# A tibble: 960 × 8
sample_size mean_diff variance.ratio equal_size dists same_dists symmetric
<dbl> <dbl> <dbl> <lgl> <chr> <lgl> <lgl>
1 100 0.25 1 TRUE skew FALSE TRUE
2 100 0.25 1 TRUE skew FALSE TRUE
3 100 0.25 1 TRUE skew FALSE TRUE
4 100 0.25 1 TRUE skew FALSE TRUE
5 100 0.25 1 TRUE skew FALSE TRUE
6 100 0.25 1 TRUE skew FALSE TRUE
7 100 0.25 1 TRUE skew FALSE TRUE
8 100 0.25 1 TRUE skew FALSE TRUE
9 100 0.25 1 TRUE skew FALSE TRUE
10 100 0.25 1 TRUE skew FALSE TRUE
# ℹ 950 more rows
# ℹ 1 more variable: detect <lgl>
# main effects only
mod <- aov (detect ~ ., data= fres)
summary (mod)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 1.4 1.4 14.88 0.00012 ***
mean_diff 1 35.3 35.3 388.81 < 2e-16 ***
variance.ratio 1 1.1 1.1 11.76 0.00063 ***
equal_size 1 0.3 0.3 3.72 0.05403 .
dists 1 1.9 1.9 20.67 6.2e-06 ***
same_dists 1 4.8 4.8 53.10 6.6e-13 ***
symmetric 1 57.6 57.6 635.03 < 2e-16 ***
Residuals 952 86.4 0.1
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# all two-way effects
mod2 <- aov (detect ~ (.)^ 2 , data= fres)
summary (mod2)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 1.4 1.4 27.23 2.2e-07 ***
mean_diff 1 35.3 35.3 711.31 < 2e-16 ***
variance.ratio 1 1.1 1.1 21.51 4.0e-06 ***
equal_size 1 0.3 0.3 6.81 0.0092 **
dists 1 1.9 1.9 37.82 1.2e-09 ***
same_dists 1 4.8 4.8 97.15 < 2e-16 ***
symmetric 1 57.6 57.6 1161.76 < 2e-16 ***
sample_size:mean_diff 1 0.3 0.3 6.81 0.0092 **
sample_size:variance.ratio 1 0.1 0.1 2.10 0.1475
sample_size:equal_size 1 0.1 0.1 1.34 0.2465
sample_size:dists 1 0.3 0.3 6.05 0.0141 *
sample_size:same_dists 1 0.0 0.0 0.76 0.3847
sample_size:symmetric 1 1.6 1.6 32.27 1.8e-08 ***
mean_diff:variance.ratio 1 0.2 0.2 4.12 0.0427 *
mean_diff:equal_size 1 0.0 0.0 0.34 0.5622
mean_diff:dists 1 0.1 0.1 2.06 0.1516
mean_diff:same_dists 1 5.1 5.1 102.95 < 2e-16 ***
mean_diff:symmetric 1 21.8 21.8 438.81 < 2e-16 ***
variance.ratio:equal_size 1 0.0 0.0 0.34 0.5622
variance.ratio:dists 1 0.0 0.0 0.04 0.8376
variance.ratio:same_dists 1 0.3 0.3 6.81 0.0092 **
variance.ratio:symmetric 1 0.8 0.8 15.25 0.0001 ***
equal_size:dists 1 0.0 0.0 0.38 0.5387
equal_size:same_dists 1 0.3 0.3 5.38 0.0206 *
equal_size:symmetric 1 0.3 0.3 6.18 0.0131 *
dists:same_dists 1 8.8 8.8 177.53 < 2e-16 ***
same_dists:symmetric 1 0.0 0.0 0.13 0.7226
Residuals 932 46.2 0.0
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Going even further, if only two-way interaction with respect to the first factor is of interest (sample_size
) then a resolution 4 input can be used, further reducing the total number of conditions to evaluate.
# 7 factors with resolution 4
fr <- FrF2 (nfactors= 7 , resolution= 4 )
# Create working simulation design given -1/1 combinations
fDesign <- createDesign (sample_size= c (100 ,200 ),
mean_diff= c (.25 , 2 ),
variance.ratio= c (1 ,4 ),
equal_size= c (TRUE , FALSE ),
dists= c ('norm' , 'skew' ),
same_dists= c (TRUE , FALSE ),
symmetric= c (TRUE , FALSE ),
# remove same-normal combo
subset = ! (symmetric & dists == 'norm' ),
fractional= fr)
fDesign
# A tibble: 12 × 7
sample_size mean_diff variance.ratio equal_size dists same_dists symmetric
<dbl> <dbl> <dbl> <lgl> <chr> <lgl> <lgl>
1 200 0.25 4 FALSE norm TRUE FALSE
2 100 0.25 4 TRUE skew TRUE FALSE
3 200 2 1 TRUE norm TRUE FALSE
4 100 2 4 TRUE norm FALSE FALSE
5 100 2 1 FALSE skew TRUE FALSE
6 100 0.25 1 FALSE norm FALSE FALSE
7 100 0.25 4 FALSE skew FALSE TRUE
8 200 2 4 FALSE skew FALSE FALSE
9 200 0.25 1 FALSE skew TRUE TRUE
10 100 2 1 TRUE skew FALSE TRUE
11 200 0.25 1 TRUE skew FALSE FALSE
12 200 2 4 TRUE skew TRUE TRUE
# note the omission of the summarise function
fres <- runSimulation (design= fDesign, replications= 20 ,
generate= Generate, analyse= Analyse, verbose= FALSE )
fres
# A tibble: 240 × 8
sample_size mean_diff variance.ratio equal_size dists same_dists symmetric
<dbl> <dbl> <dbl> <lgl> <chr> <lgl> <lgl>
1 200 0.25 4 FALSE norm TRUE FALSE
2 200 0.25 4 FALSE norm TRUE FALSE
3 200 0.25 4 FALSE norm TRUE FALSE
4 200 0.25 4 FALSE norm TRUE FALSE
5 200 0.25 4 FALSE norm TRUE FALSE
6 200 0.25 4 FALSE norm TRUE FALSE
7 200 0.25 4 FALSE norm TRUE FALSE
8 200 0.25 4 FALSE norm TRUE FALSE
9 200 0.25 4 FALSE norm TRUE FALSE
10 200 0.25 4 FALSE norm TRUE FALSE
# ℹ 230 more rows
# ℹ 1 more variable: detect <lgl>
# main effects
mod <- aov (detect ~ ., data= fres)
summary (mod)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 1.35 1.35 20.29 1.1e-05 ***
mean_diff 1 9.60 9.60 144.31 < 2e-16 ***
variance.ratio 1 2.02 2.02 30.32 9.7e-08 ***
equal_size 1 4.80 4.80 72.16 2.4e-15 ***
dists 1 0.17 0.17 2.54 0.1126
same_dists 1 0.67 0.67 10.15 0.0016 **
symmetric 1 13.81 13.81 207.54 < 2e-16 ***
Residuals 232 15.43 0.07
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# two-way interaction effects for sample_size only
mod2 <- aov (detect ~ (.)^ 2 , data= fres)
summary (mod2)
Df Sum Sq Mean Sq F value Pr(>F)
sample_size 1 1.35 1.35 31.73 5.2e-08 ***
mean_diff 1 9.60 9.60 225.65 < 2e-16 ***
variance.ratio 1 2.02 2.02 47.40 5.5e-11 ***
equal_size 1 4.80 4.80 112.82 < 2e-16 ***
dists 1 0.17 0.17 3.97 0.04761 *
same_dists 1 0.67 0.67 15.87 9.1e-05 ***
symmetric 1 13.81 13.81 324.52 < 2e-16 ***
sample_size:mean_diff 1 2.55 2.55 59.99 3.1e-13 ***
sample_size:variance.ratio 1 0.63 0.63 14.69 0.00016 ***
sample_size:equal_size 1 0.06 0.06 1.32 0.25141
sample_size:dists 1 2.50 2.50 58.76 5.1e-13 ***
Residuals 228 9.70 0.04
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
References
Montgomery, Douglas C. (2017). Design and Analysis of Experiments (9th ed.) . John Wiley & Sons, Inc.