1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
|
\name{simulation}
\alias{simulation}
\title{Simulation studies}
\usage{simulation(replications, DataGen, listOfProcedures, discardProcInput=FALSE)}
\description{This function generates data according to a specified function and parameters
and then applies specified procedures to the generated data. The generated data is stored
in $data and the results are stored in $results. In order to recognize which
results and generated data belong together every result contains an element
$data.set.number. The data generating function must return a list which contains a list
$procInput. Every element in $procInput will be used as an input parameter
for the procedures provided by listOfProcedures.}
\value{A list with 2 elements. $data contains the objects generated by DataGen. $results
contains the objects generated by the procedures augmented by the data set number and
the parameter constellation.}
\author{MarselScheer}
\arguments{\item{replications}{The number of replications. This means how many simulation runs will be performed.}
\item{DataGen}{A list that contains the function and parameters for generating data
which will be analyzed be the procedures in listOfProcedures. It must contain the elements
$funName (character) $fun (function)}
\item{listOfProcedures}{A list of lists which contains the procedures with their parameters
for the simulation.}
\item{discardProcInput}{A list of lists which contains the procedures with their parameters}}
\examples{# this function generates pValues
myGen <- function(n, n0) {
list(procInput=list(pValues = c(runif(n-n0, 0, 0.01), runif(n0))),
groundTruth = c(rep(FALSE, times=n-n0), rep(TRUE, times=n0)))
}
sim <- simulation(replications = 10, list(funName="myGen", fun=myGen, n=200, n0=c(50,100)),
list(list(funName="BH",
fun=function(pValues, alpha) BH(pValues, alpha, silent=TRUE), alpha=c(0.25, 0.5)),
list(funName="holm", fun=holm, alpha=c(0.25, 0.5),silent=TRUE)))
# the following happend:
# Call myGen(200,50) and append result to sim$data
# Apply bonferroni and holm each with alpha 0.25 and 0.5 to this data set
# Append the results to sim$restults
# Repeat this 10 times.
# Call myGen(200, 100) and append restult to sim$data
# Apply bonferroni and holm each with alpha 0.25 and 0.5 to this data set
# Append the results to sim$restults
# Repeat this 10 times.
length(sim$data)
length(sim$results)
# we now reproduce the 6th item in results
print(sim$results[[6]]$data.set.number)
print(sim$results[[6]]$parameters)
all(BH(sim$data[[2]]$procInput$pValues, 0.5, silent=TRUE)$adjPValues == sim$results[[6]]$adjPValues)
#
# Just calculating some statistics and making some plots
NumberOfType1Error <- function(data, result) sum(data$groundTruth * result$rejected)
result.all <- gatherStatistics(sim, list(NumOfType1Err = NumberOfType1Error))
result <- gatherStatistics(sim, list(NumOfType1Err = NumberOfType1Error),
list(median=median, mean=mean, sd=sd))
print(result)
require(lattice)
histogram(~NumOfType1Err | method*alpha, data = result.all$statisticDF)
barchart(NumOfType1Err.median + NumOfType1Err.mean ~ method | alpha, data = result$statisticDF)}
|
