Browse code

2.3.11;\n - evalAllGenotypes: order = FALSE by default; \n - Clarified difference plotFitnessEffects and plotFitnessLandscape.

git-svn-id: file:///home/git/hedgehog.fhcrc.org/bioconductor/trunk/madman/Rpacks/OncoSimulR@120020 bc3139a8-67e5-0310-9ffc-ced21a209358

Ramon Diaz-Uriarte authored on 10/08/2016 15:47:33
Showing 15 changed files

... ...
@@ -1,8 +1,8 @@
1 1
 Package: OncoSimulR
2 2
 Type: Package
3 3
 Title: Forward Genetic Simulation of Cancer Progresion with Epistasis 
4
-Version: 2.3.10
5
-Date: 2016-07-14
4
+Version: 2.3.11
5
+Date: 2016-08-10
6 6
 Authors@R: c(person("Ramon", "Diaz-Uriarte", role = c("aut", "cre"),
7 7
 		     email = "rdiaz02@gmail.com"),
8 8
 	      person("Mark", "Taylor", role = "ctb", email = "ningkiling@gmail.com"))
... ...
@@ -1032,7 +1032,7 @@ evalGenotypeFitAndMut <- function(genotype,
1032 1032
 ## I am here: simplify this
1033 1033
 
1034 1034
 evalAllGenotypesORMut <- function(fmEffects,
1035
-                                  order = TRUE, max = 256,
1035
+                                  order = FALSE, max = 256,
1036 1036
                              addwt = FALSE,
1037 1037
                              model = "",
1038 1038
                              calledBy_ = "") {
... ...
@@ -1129,7 +1129,7 @@ evalAllGenotypesORMut <- function(fmEffects,
1129 1129
     return(df)
1130 1130
 }
1131 1131
 
1132
-evalAllGenotypes <- function(fitnessEffects, order = TRUE, max = 256,
1132
+evalAllGenotypes <- function(fitnessEffects, order = FALSE, max = 256,
1133 1133
                              addwt = FALSE,
1134 1134
                              model = "") {
1135 1135
     evalAllGenotypesORMut(
... ...
@@ -1190,7 +1190,7 @@ generateAllGenotypes <- function(fitnessEffects, order = TRUE, max = 256) {
1190 1190
 }
1191 1191
 
1192 1192
 evalAllGenotypesFitAndMut <- function(fitnessEffects, mutatorEffects,
1193
-                                   order = TRUE, max = 256,
1193
+                                   order = FALSE, max = 256,
1194 1194
                                    addwt = FALSE,
1195 1195
                                    model = "" ){
1196 1196
 ##                                   minimal = FALSE) {
... ...
@@ -59,8 +59,12 @@ rfitness <- function(g, c= 0.5,
59 59
         }
60 60
         m <- cbind(m, Fitness = fi)
61 61
         if(min_accessible_genotypes > 0) {
62
-            if(count_accessible_g(m, accessible_th) >= min_accessible_genotypes) {
62
+            num_accessible_genotypes <- count_accessible_g(m, accessible_th)
63
+            if(num_accessible_genotypes >= min_accessible_genotypes) {
63 64
                 done <- TRUE
65
+                attributes(m) <- c(attributes(m),
66
+                                   accessible_genotypes = num_accessible_genotypes,
67
+                                   accessible_th = accessible_th)
64 68
             } else {
65 69
                 ## Cannot start again with a fitness column
66 70
                 m <- m[, -ncol(m), drop = FALSE]
... ...
@@ -1,4 +1,8 @@
1
-Changes in version 2.3.9 (2017-07-14):
1
+Changes in version 2.3.11 (2017-08-10):
2
+	- evalAllGenotypes: order = FALSE by default.
3
+	- Clarified difference plotFitnessEffects and plotFitnessLandscape.
4
+	
5
+Changes in version 2.3.10 (2017-07-14):
2 6
 	- The windows check failure with mc.cores.
3 7
 	
4 8
 Changes in version 2.3.9 (2017-07-09):
... ...
@@ -25,7 +25,7 @@ model = "")
25 25
 
26 26
 evalGenotypeMut(genotype, mutatorEffects, verbose = FALSE, echo = FALSE)
27 27
 
28
-evalAllGenotypes(fitnessEffects, order = TRUE, max = 256, addwt = FALSE,
28
+evalAllGenotypes(fitnessEffects, order = FALSE, max = 256, addwt = FALSE,
29 29
 model = "")
30 30
 
31 31
 evalAllGenotypesMut(mutatorEffects, max = 256, addwt = FALSE)
... ...
@@ -36,7 +36,7 @@ evalGenotypeFitAndMut(genotype, fitnessEffects,
36 36
                       model = "")
37 37
 
38 38
 evalAllGenotypesFitAndMut(fitnessEffects, mutatorEffects,
39
-                          order = TRUE, max = 256, addwt = FALSE,
39
+                          order = FALSE, max = 256, addwt = FALSE,
40 40
                           model = "")
41 41
 
42 42
 
... ...
@@ -64,8 +64,11 @@ evalAllGenotypesFitAndMut(fitnessEffects, mutatorEffects,
64 64
 
65 65
   
66 66
   \item{order}{
67
-    (For \code{evalAllGenotypes}). Does order matter? If it does, then generate not only all possible
68
-    combinations of the genes, but all possible permutations for each combination.
67
+
68
+    (For \code{evalAllGenotypes}). If TRUE, then order matters. If order
69
+    matters, then generate not only all possible combinations of the genes, but
70
+    all possible permutations for each combination.
71
+    
69 72
 }
70 73
 \item{max}{
71 74
   (For \code{evalAllGenotypes}). By default, no output is shown if the
... ...
@@ -186,7 +189,7 @@ ofe2 <- allFitnessEffects(orderEffects = c("F > D" = -0.3, "D > F" = 0.4),
186 189
                                 "F" = "f1, f2, f3",
187 190
                                 "D" = "d1, d2") )
188 191
 
189
-evalAllGenotypes(ofe2, max = 325)[1:15, ]
192
+evalAllGenotypes(ofe2, order = TRUE, max = 325)[1:15, ]
190 193
 
191 194
 ## Next two are identical
192 195
 evalGenotype("d1 > d2 > f3", ofe2, verbose = TRUE)
... ...
@@ -6,7 +6,8 @@
6 6
 }
7 7
 \description{
8 8
   Plot the restriction table/graph of restrictions, the epistasis, and
9
-  the order effects in a fitnessEffects object.
9
+  the order effects in a fitnessEffects object. This is not a plot of
10
+  the fitness landscape; for that, see \code{\link{plotFitnessLandscape}}.
10 11
 }
11 12
 \usage{
12 13
 \method{plot}{fitnessEffects}(x, type = "graphNEL", layout = NULL,
... ...
@@ -81,7 +82,9 @@ Genes without interactions are not shown.
81 82
 }
82 83
 
83 84
 \seealso{
84
-  \code{\link{allFitnessEffects}}
85
+  \code{\link{allFitnessEffects}},
86
+  \code{\link{plotFitnessLandscape}}
87
+  
85 88
 }
86 89
 \examples{
87 90
 
... ...
@@ -10,7 +10,10 @@
10 10
 
11 11
 \description{ Show a plot of a fitness landscape. The plot is modeled
12 12
   after (actually, mostly a blatant copy of) that of MAGELLAN,
13
-  \url{http://wwwabi.snv.jussieu.fr/public/Magellan/}. 
13
+  \url{http://wwwabi.snv.jussieu.fr/public/Magellan/}.
14
+
15
+  Note: this is not a plot of the fitnessEffects object; for that, see
16
+  \code{\link{plot.fitnessEffects}}.
14 17
 
15 18
 }
16 19
 
... ...
@@ -183,7 +186,7 @@ Ramon Diaz-Uriarte
183 186
   of a total of 5000 genotypes you'll have to wait a while for the plot
184 187
   to finish.
185 188
 
186
- 
189
+  
187 190
 
188 191
 }
189 192
 
... ...
@@ -191,7 +194,8 @@ Ramon Diaz-Uriarte
191 194
   \code{\link{allFitnessEffects}},
192 195
   \code{\link{evalAllGenotypes}},
193 196
   \code{\link{allFitnessEffects}},
194
-  \code{\link{rfitness}}
197
+  \code{\link{rfitness}},
198
+  \code{\link{plot.fitnessEffects}}
195 199
 }
196 200
 \examples{
197 201
 
... ...
@@ -23,7 +23,7 @@ rfitness(g, c = 0.5, sd = 1, reference = "random", scale = NULL,
23 23
   \item{g}{Number of genes.}
24 24
 
25 25
   \item{c}{The decrease in fitness of a genotype per each unit increase
26
-    in Hamming distance from the reference genotype (\code{reference}).}
26
+    in Hamming distance from the reference genotype (see \code{reference}).}
27 27
 
28 28
   \item{sd}{The standard deviation of the random component (a normal
29 29
   distribution of mean 0 and standard deviation \code{sd}).}
... ...
@@ -34,7 +34,7 @@ rfitness(g, c = 0.5, sd = 1, reference = "random", scale = NULL,
34 34
   distance from this reference. If "random" a genotype will be randomly
35 35
   chosen as the reference. If "max" the genotype with all positions
36 36
   mutated will be chosen as the reference. If you pass a vector (e.g.,
37
-  \code{fittest = c(1, 0, 1, 0)}) that will be the reference genotype.}
37
+  \code{reference = c(1, 0, 1, 0)}) that will be the reference genotype.}
38 38
 
39 39
 \item{scale}{Either NULL (nothing is done) or a two-element vector. If a
40 40
   two-element vector, fitness is re-scaled between \code{scale[1]} (the
... ...
@@ -101,7 +101,12 @@ rfitness(g, c = 0.5, sd = 1, reference = "random", scale = NULL,
101 101
   column denotes gene mutated/not-mutated. (For ease of use in other
102 102
   functions, this matrix has class  "genotype_fitness_matrix".) 
103 103
 
104
+  If you have specified \code{min_accessible_genotypes > 0}, the return
105
+  object has added attributes \code{accessible_genotypes} and
106
+  \code{accessible_th} that show the number of accessible
107
+  genotypes under the specified  threshold.
104 108
 }
109
+  
105 110
 \references{
106 111
 
107 112
   Szendro I.~G. et al. (2013). Quantitative analyses of empirical
... ...
@@ -4,7 +4,9 @@ CMD check". They can be run by doing "test_dir()".
4 4
 They are kept separate because they are long running, many are overkills,
5 5
 and in a few it is possible some might fail even when there are no bugs
6 6
 (since we are comparing expected outcomes from a distribution, and we
7
-could end up in unlikely cases).
7
+could end up in unlikely cases). In fact, you should expect some tests to
8
+fail if you run them enough times (if they were to never fail, this could
9
+mean that the tests are not actually sensitive enough).
8 10
 
9 11
 Note that we will rarely want to use something like try_again (available
10 12
 from the github version of testthat) or repeated tries (e.g.,
... ...
@@ -6,14 +6,14 @@
6 6
 
7 7
 cat(paste("\n Starting sample-prob-long", date(), "\n"))
8 8
 
9
-p.value.threshold <- 0.01
9
+p.value.threshold <- 1e-7
10 10
 
11 11
 
12 12
 test_that("Increasing cPDetect decreases time, Exp" , {
13 13
     gi <- rep(0.1, 10)
14 14
     names(gi) <- letters[1:10]
15 15
     oi <- allFitnessEffects(noIntGenes = gi)
16
-    n <- 100
16
+    n <- 800
17 17
     max.tries <- 4  
18 18
     for(tries in 1:max.tries) {
19 19
         sa <- oncoSimulPop(n,
... ...
@@ -35,7 +35,8 @@ test_that("Increasing cPDetect decreases time, Exp" , {
35 35
                            onlyCancer = FALSE,
36 36
                            detectionDrivers = 99, mc.cores = 2)
37 37
         ta <- unlist(lapply(sa, function(x) x$FinalTime))
38
-        tb <- unlist(lapply(sb, function(x) x$FinalTime))         
38
+        tb <- unlist(lapply(sb, function(x) x$FinalTime))
39
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
39 40
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
40 41
         if(T1) break;
41 42
     }
... ...
@@ -48,7 +49,7 @@ test_that("Increasing p2 decreases time, Exp" , {
48 49
     gi <- rep(0.1, 10)
49 50
     names(gi) <- letters[1:10]
50 51
     oi <- allFitnessEffects(noIntGenes = gi)
51
-    n <- 100
52
+    n <- 1200
52 53
     max.tries <- 4  
53 54
     for(tries in 1:max.tries) {
54 55
         sa <- oncoSimulPop(n,
... ...
@@ -70,7 +71,8 @@ test_that("Increasing p2 decreases time, Exp" , {
70 71
                            onlyCancer = FALSE,
71 72
                            detectionDrivers = 99, mc.cores = 2)
72 73
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
73
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
74
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
75
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
74 76
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
75 77
         if(T1) break;
76 78
     }
... ...
@@ -79,12 +81,11 @@ test_that("Increasing p2 decreases time, Exp" , {
79 81
 })
80 82
 
81 83
 
82
-
83 84
 test_that("Increasing n2 increases time, Exp" , {
84 85
     gi <- rep(0.1, 10)
85 86
     names(gi) <- letters[1:10]
86 87
     oi <- allFitnessEffects(noIntGenes = gi)
87
-    n <- 100
88
+    n <- 1200
88 89
     max.tries <- 4  
89 90
     for(tries in 1:max.tries) {
90 91
         sa <- oncoSimulPop(n,
... ...
@@ -92,7 +93,7 @@ test_that("Increasing n2 increases time, Exp" , {
92 93
                            model = "Exp",
93 94
                            initSize = 5000,
94 95
                            keepEvery = -9,
95
-                           detectionProb = c(p2 = .6, n2 = 7000, cPDetect = NULL),
96
+                           detectionProb = c(p2 = .6, n2 = 8000, cPDetect = NULL),
96 97
                            finalTime = NA,
97 98
                            onlyCancer = FALSE,
98 99
                            detectionDrivers = 99, mc.cores = 2)
... ...
@@ -101,12 +102,13 @@ test_that("Increasing n2 increases time, Exp" , {
101 102
                            model = "Exp",
102 103
                            initSize = 5000,
103 104
                            keepEvery = -9,
104
-                           detectionProb = c(p2 = .6, n2 = 5502, cPDetect = NULL),
105
+                           detectionProb = c(p2 = .6, n2 = 6002, cPDetect = NULL),
105 106
                            finalTime = NA,
106 107
                            onlyCancer = FALSE,
107 108
                            detectionDrivers = 99, mc.cores = 2)
108 109
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
109
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
110
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
111
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
110 112
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
111 113
         if(T1) break;
112 114
     }
... ...
@@ -121,7 +123,7 @@ test_that("Increasing checkSizePEvery increases time, Exp" , {
121 123
     gi <- rep(0.1, 10)
122 124
     names(gi) <- letters[1:10]
123 125
     oi <- allFitnessEffects(noIntGenes = gi)
124
-    n <- 100
126
+    n <- 1200
125 127
     max.tries <- 4  
126 128
     for(tries in 1:max.tries) {
127 129
         sa <- oncoSimulPop(n,
... ...
@@ -143,7 +145,8 @@ test_that("Increasing checkSizePEvery increases time, Exp" , {
143 145
                            onlyCancer = FALSE,
144 146
                            detectionDrivers = 99, mc.cores = 2)
145 147
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
146
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
148
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
149
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
147 150
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
148 151
         if(T1) break;
149 152
     }
... ...
@@ -14,7 +14,7 @@ test_that("Root name in module table or not", {
14 14
                                            "B" = 0.5),
15 15
                              geneToModule = c("A" = "a1, a2",
16 16
                                               "B" = "b1")))
17
-    expect_identical(evalAllGenotypes(fee), evalAllGenotypes(fee2))
17
+    expect_identical(evalAllGenotypes(fee, order = TRUE), evalAllGenotypes(fee2, order = TRUE))
18 18
     expect_error(allFitnessEffects(epistasis = c("A" = 0.3,
19 19
                                                  "B" = 0.5),
20 20
                                    geneToModule = c(
... ...
@@ -161,9 +161,11 @@ test_that("Bauer example: identical values fitness classes, diff. order", {
161 161
 
162 162
 test_that("Order effects, entry of labels and separation", {
163 163
     o1 <- evalAllGenotypes(allFitnessEffects(
164
-        orderEffects = c("d>f" = 0.4, "f > d" = -0.3) ))
164
+        orderEffects = c("d>f" = 0.4, "f > d" = -0.3) ),
165
+        order = TRUE)
165 166
     o2 <- evalAllGenotypes(allFitnessEffects(
166
-        orderEffects = c("f > d" = -0.3, "d > f" = 0.4) ))
167
+        orderEffects = c("f > d" = -0.3, "d > f" = 0.4) ),
168
+        order = TRUE)
167 169
     o1s <- o1[order(o1$Genotype),   ]
168 170
     o2s <- o2[order(o2$Genotype),   ]
169 171
     expect_equal(o1s, o2s)
... ...
@@ -175,7 +177,7 @@ test_that("Order effects, modules 1", {
175 177
                                   c("Root" = "Root",
176 178
                                     "F" = "f1, f2",
177 179
                                     "D" = "d1, d2") )
178
-    ag <- evalAllGenotypes(ofe1)
180
+    ag <- evalAllGenotypes(ofe1, order = TRUE)
179 181
     expect_true(all.equal(ag[c(17, 39, 19, 29), "Fitness"], c(1.4, 0.7, 1.4, 0.7)))
180 182
     expect_true(all.equal(ag[c(43, 44), "Fitness"], c(1.4, 1.4)))
181 183
     expect_true(all(ag[41:52, "Fitness"] == 1.4))
... ...
@@ -188,7 +190,7 @@ test_that("Order effects, modules 2", {
188 190
                                   c("Root" = "Root",
189 191
                                     "F" = "f1, f2, f3",
190 192
                                     "D" = "d1, d2") )
191
-    ag2 <- evalAllGenotypes(ofe2, max = 326)
193
+    ag2 <- evalAllGenotypes(ofe2, max = 326, order = TRUE)
192 194
     oe <- c(grep("^f.*d.*", ag2[, 1]), grep("^d.*f.*", ag2[, 1]))
193 195
     expect_true(all(ag2[grep("^d.*f.*", ag2[, 1]), "Fitness"] == 1.4))
194 196
     expect_true(all(ag2[grep("^f.*d.*", ag2[, 1]), "Fitness"] == 0.7))
... ...
@@ -200,10 +202,10 @@ test_that("Order effects, modules 2", {
200 202
 test_that("Order effects, twisted module names", {
201 203
     o1 <- evalAllGenotypes(allFitnessEffects(
202 204
       orderEffects = c("F > D" = -0.3, "D > F" = 0.4),  
203
-        geneToModule = c("Root" = "Root", "F" = "d", "D" = "f")))
205
+        geneToModule = c("Root" = "Root", "F" = "d", "D" = "f")), order = TRUE)
204 206
     o2 <- evalAllGenotypes(allFitnessEffects(
205 207
       orderEffects = c("D>F" = 0.4, "F >D" = -0.3),  
206
-        geneToModule = c("Root" = "Root", "F" = "d", "D" = "f")))
208
+        geneToModule = c("Root" = "Root", "F" = "d", "D" = "f")), order = TRUE)
207 209
     o1s <- o1[order(o1$Genotype),   ]
208 210
     o2s <- o2[order(o2$Genotype),   ]
209 211
     expect_equal(o1s, o2s)
... ...
@@ -225,7 +227,7 @@ test_that("Order effects, three-gene-orders and modules 1", {
225 227
                                     "M" = "m",
226 228
                                     "F" = "f",
227 229
                                     "D" = "d") )
228
-    ag <- evalAllGenotypes(o3)
230
+    ag <- evalAllGenotypes(o3, order = TRUE)
229 231
     expect_true(all.equal(ag[, "Fitness"],
230 232
                         c(rep(1, 4),
231 233
                           1.1,
... ...
@@ -318,7 +320,7 @@ test_that("No interaction genes and order effects, 1", {
318 320
     foi1 <- allFitnessEffects(
319 321
         orderEffects = c("D>B" = -0.3, "B > D" = 0.3),
320 322
         noIntGenes = c("A" = 0.05, "C" = -.2, "E" = .1))
321
-    agoi1 <- evalAllGenotypes(foi1,  max = 325)
323
+    agoi1 <- evalAllGenotypes(foi1,  max = 325, order = TRUE)
322 324
     rn <- 1:nrow(agoi1)
323 325
     names(rn) <- agoi1[, 1]
324 326
     expect_true(all.equal(agoi1[rn[LETTERS[1:5]], "Fitness"], c(1.05, 1, 0.8, 1, 1.1)))
... ...
@@ -349,7 +351,7 @@ test_that("No interaction genes and order effects, 2", {
349 351
     foi1 <- allFitnessEffects(
350 352
         orderEffects = c("D>B" = -0.3, "B > D" = 0.3),
351 353
         noIntGenes = c("M" = 0.05, "A" = -.2, "J" = .1))
352
-    agoi1 <- evalAllGenotypes(foi1,  max = 325)
354
+    agoi1 <- evalAllGenotypes(foi1,  max = 325, order = TRUE)
353 355
     rn <- 1:nrow(agoi1)
354 356
     names(rn) <- agoi1[, 1]
355 357
     expect_true(all.equal(agoi1[rn[c("B", "D", "M", "A", "J")], "Fitness"],
... ...
@@ -883,7 +885,7 @@ test_that("a silly epistasis example", {
883 885
     expect_silent(sv <- allFitnessEffects(
884 886
                      orderEffects = c("A > B" = 0.1),
885 887
                      epistasis = c("A" = 1)))
886
-    expect_output(print(evalAllGenotypes(sv)), "Genotype")
888
+    expect_output(print(evalAllGenotypes(sv, order = TRUE)), "Genotype")
887 889
 })
888 890
 
889 891
 test_that("can run without keeping input", {
... ...
@@ -255,7 +255,10 @@ test_that("mut and fitness both needed when needed", {
255 255
     expect_error(evalGenotypeFitAndMut("a, b", mutatorEffects = fm),
256 256
                  'argument "fitnessEffects" is missing',
257 257
                  fixed = TRUE)
258
-    expect_error(evalAllGenotypesFitAndMut(fe),
258
+    expect_error(evalAllGenotypesFitAndMut(fe, order = TRUE),
259
+                 'argument "mutatorEffects" is missing',
260
+                 fixed = TRUE)
261
+    expect_error(evalAllGenotypesFitAndMut(fe, order = FALSE),
259 262
                  'argument "mutatorEffects" is missing',
260 263
                  fixed = TRUE)
261 264
     expect_error(evalAllGenotypesFitAndMut(fm),
... ...
@@ -332,6 +335,9 @@ test_that("We cannot pass mutator/fitness objects to the wrong functions", {
332 335
     expect_error(evalAllGenotypesMut(fe2),
333 336
                  "You are trying to get the mutator effects of a fitness specification.",
334 337
                  fixed = TRUE)
338
+    expect_error(evalAllGenotypesMut(fe2),
339
+                 "You are trying to get the mutator effects of a fitness specification.",
340
+                 fixed = TRUE)
335 341
     expect_error(evalAllGenotypes(fm2),
336 342
                  "You are trying to get the fitness of a mutator specification.",
337 343
                  fixed = TRUE)
... ...
@@ -368,7 +374,7 @@ test_that("evaluating genotype and mutator", {
368 374
                     c(1.5 * 1.1, 5))
369 375
     expect_equivalent(dplyr::filter(ou, Genotype == "a, b, c, e")[, c(2, 3)],
370 376
                     c(1.3 * 1.5 * 1.1, 10 * 5))
371
-    oo <- evalAllGenotypesFitAndMut(fe, fm)
377
+    oo <- evalAllGenotypesFitAndMut(fe, fm, order = TRUE)
372 378
     expect_equivalent(dplyr::filter(oo, Genotype == "a")[, c(2, 3)],
373 379
                     c(1, 10))
374 380
     expect_equivalent(dplyr::filter(oo, Genotype == "b")[, c(2, 3)],
... ...
@@ -1,12 +1,12 @@
1 1
 cat(paste("\n Starting sample-prob", date(), "\n"))
2 2
 
3
-p.value.threshold <- 0.01
3
+p.value.threshold <- 1e-4
4 4
 
5 5
 test_that("Increasing cPDetect decreases time" , {
6 6
     gi <- rep(0.1, 10)
7 7
     names(gi) <- letters[1:10]
8 8
     oi <- allFitnessEffects(noIntGenes = gi)
9
-    n <- 15
9
+    n <- 35
10 10
     max.tries <- 4  
11 11
     for(tries in 1:max.tries) {
12 12
         sa <- oncoSimulPop(n,
... ...
@@ -28,7 +28,8 @@ test_that("Increasing cPDetect decreases time" , {
28 28
                            onlyCancer = TRUE,
29 29
                            detectionDrivers = NA, mc.cores = 2)
30 30
         ta <- unlist(lapply(sa, function(x) x$FinalTime))
31
-        tb <- unlist(lapply(sb, function(x) x$FinalTime))         
31
+        tb <- unlist(lapply(sb, function(x) x$FinalTime))
32
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
32 33
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
33 34
         if(T1) break;
34 35
     }
... ...
@@ -41,7 +42,7 @@ test_that("Increasing p2 decreases time" , {
41 42
     gi <- rep(0.1, 10)
42 43
     names(gi) <- letters[1:10]
43 44
     oi <- allFitnessEffects(noIntGenes = gi)
44
-    n <- 25
45
+    n <- 35
45 46
     max.tries <- 4  
46 47
     for(tries in 1:max.tries) {
47 48
         sa <- oncoSimulPop(n,
... ...
@@ -78,7 +79,7 @@ test_that("Increasing n2 increases time" , {
78 79
     gi <- rep(0.1, 10)
79 80
     names(gi) <- letters[1:10]
80 81
     oi <- allFitnessEffects(noIntGenes = gi)
81
-    n <- 30
82
+    n <- 40
82 83
     max.tries <- 4  
83 84
     for(tries in 1:max.tries) {
84 85
         sa <- oncoSimulPop(n,
... ...
@@ -116,7 +117,7 @@ test_that("Increasing checkSizePEvery increases time" , {
116 117
     gi <- rep(0.1, 10)
117 118
     names(gi) <- letters[1:10]
118 119
     oi <- allFitnessEffects(noIntGenes = gi)
119
-    n <- 30
120
+    n <- 40
120 121
     max.tries <- 4  
121 122
     for(tries in 1:max.tries) {
122 123
         sa <- oncoSimulPop(n,
... ...
@@ -138,7 +139,8 @@ test_that("Increasing checkSizePEvery increases time" , {
138 139
                            onlyCancer = TRUE,
139 140
                            detectionDrivers = NA, mc.cores = 2)
140 141
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
141
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
142
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
143
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
142 144
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
143 145
         if(T1) break;
144 146
     }
... ...
@@ -155,7 +157,7 @@ test_that("Increasing cPDetect decreases time, Exp" , {
155 157
     gi <- rep(0.1, 10)
156 158
     names(gi) <- letters[1:10]
157 159
     oi <- allFitnessEffects(noIntGenes = gi)
158
-    n <- 20
160
+    n <- 30
159 161
     max.tries <- 4  
160 162
     for(tries in 1:max.tries) {
161 163
         sa <- oncoSimulPop(n,
... ...
@@ -177,7 +179,8 @@ test_that("Increasing cPDetect decreases time, Exp" , {
177 179
                            onlyCancer = TRUE,
178 180
                            detectionDrivers = NA, mc.cores = 2)
179 181
         ta <- unlist(lapply(sa, function(x) x$FinalTime))
180
-        tb <- unlist(lapply(sb, function(x) x$FinalTime))         
182
+        tb <- unlist(lapply(sb, function(x) x$FinalTime))
183
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
181 184
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
182 185
         if(T1) break;
183 186
     }
... ...
@@ -190,7 +193,7 @@ test_that("Increasing p2 decreases time, Exp" , {
190 193
     gi <- rep(0.1, 10)
191 194
     names(gi) <- letters[1:10]
192 195
     oi <- allFitnessEffects(noIntGenes = gi)
193
-    n <- 30
196
+    n <- 40
194 197
     max.tries <- 4  
195 198
     for(tries in 1:max.tries) {
196 199
         sa <- oncoSimulPop(n,
... ...
@@ -213,6 +216,7 @@ test_that("Increasing p2 decreases time, Exp" , {
213 216
                            detectionDrivers = NA, mc.cores = 2)
214 217
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
215 218
         (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
219
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
216 220
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
217 221
         if(T1) break;
218 222
     }
... ...
@@ -226,7 +230,7 @@ test_that("Increasing n2 increases time, Exp" , {
226 230
     gi <- rep(0.1, 10)
227 231
     names(gi) <- letters[1:10]
228 232
     oi <- allFitnessEffects(noIntGenes = gi)
229
-    n <- 30
233
+    n <- 40
230 234
     max.tries <- 4  
231 235
     for(tries in 1:max.tries) {
232 236
         sa <- oncoSimulPop(n,
... ...
@@ -248,7 +252,8 @@ test_that("Increasing n2 increases time, Exp" , {
248 252
                            onlyCancer = TRUE,
249 253
                            detectionDrivers = NA, mc.cores = 2)
250 254
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
251
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
255
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
256
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
252 257
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
253 258
         if(T1) break;
254 259
     }
... ...
@@ -263,7 +268,7 @@ test_that("Increasing checkSizePEvery increases time, Exp" , {
263 268
     gi <- rep(0.1, 10)
264 269
     names(gi) <- letters[1:10]
265 270
     oi <- allFitnessEffects(noIntGenes = gi)
266
-    n <- 30
271
+    n <- 40
267 272
     max.tries <- 4  
268 273
     for(tries in 1:max.tries) {
269 274
         sa <- oncoSimulPop(n,
... ...
@@ -285,7 +290,8 @@ test_that("Increasing checkSizePEvery increases time, Exp" , {
285 290
                            onlyCancer = TRUE,
286 291
                            detectionDrivers = NA, mc.cores = 2)
287 292
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
288
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
293
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
294
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
289 295
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
290 296
         if(T1) break;
291 297
     }
... ...
@@ -306,7 +312,7 @@ test_that("Increasing cPDetect decreases time" , {
306 312
     gi <- rep(0.0,  10)
307 313
     names(gi) <- letters[1:10]
308 314
     oi <- allFitnessEffects(noIntGenes = gi)
309
-    n <- 15
315
+    n <- 25
310 316
     max.tries <- 4  
311 317
     for(tries in 1:max.tries) {
312 318
         sa <- oncoSimulPop(n,
... ...
@@ -328,7 +334,8 @@ test_that("Increasing cPDetect decreases time" , {
328 334
                            onlyCancer = FALSE,
329 335
                            detectionDrivers = NA, mc.cores = 2)
330 336
         ta <- unlist(lapply(sa, function(x) x$FinalTime))
331
-        tb <- unlist(lapply(sb, function(x) x$FinalTime))         
337
+        tb <- unlist(lapply(sb, function(x) x$FinalTime))
338
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
332 339
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
333 340
         if(T1) break;
334 341
     }
... ...
@@ -341,7 +348,7 @@ test_that("Increasing p2 decreases time" , {
341 348
     gi <- rep(0.0,  10)
342 349
     names(gi) <- letters[1:10]
343 350
     oi <- allFitnessEffects(noIntGenes = gi)
344
-    n <- 20
351
+    n <- 40
345 352
     max.tries <- 4  
346 353
     for(tries in 1:max.tries) {
347 354
         sa <- oncoSimulPop(n,
... ...
@@ -363,7 +370,8 @@ test_that("Increasing p2 decreases time" , {
363 370
                            onlyCancer = FALSE,
364 371
                            detectionDrivers = NA, mc.cores = 2)
365 372
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
366
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
373
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
374
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
367 375
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
368 376
         if(T1) break;
369 377
     }
... ...
@@ -377,7 +385,7 @@ test_that("Increasing n2 increases time" , {
377 385
     gi <- rep(0.0,  10)
378 386
     names(gi) <- letters[1:10]
379 387
     oi <- allFitnessEffects(noIntGenes = gi)
380
-    n <- 30
388
+    n <- 40
381 389
     max.tries <- 4  
382 390
     for(tries in 1:max.tries) {
383 391
         sa <- oncoSimulPop(n,
... ...
@@ -415,7 +423,7 @@ test_that("Increasing checkSizePEvery increases time" , {
415 423
     gi <- rep(0.0,  10)
416 424
     names(gi) <- letters[1:10]
417 425
     oi <- allFitnessEffects(noIntGenes = gi)
418
-    n <- 30
426
+    n <- 40
419 427
     max.tries <- 4  
420 428
     for(tries in 1:max.tries) {
421 429
         sa <- oncoSimulPop(n,
... ...
@@ -437,7 +445,8 @@ test_that("Increasing checkSizePEvery increases time" , {
437 445
                            onlyCancer = FALSE,
438 446
                            detectionDrivers = NA, mc.cores = 2)
439 447
         (ta <- unlist(lapply(sa, function(x) x$FinalTime)))
440
-        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))         
448
+        (tb <- unlist(lapply(sb, function(x) x$FinalTime)))
449
+        print(suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value))
441 450
         T1 <- suppressWarnings(wilcox.test(ta, tb, alternative = "greater")$p.value < p.value.threshold)
442 451
         if(T1) break;
443 452
     }
... ...
@@ -316,7 +316,7 @@ sv2 <- allFitnessEffects(epistasis = c("-A : B" = sb,
316 316
                              "B" = "b",
317 317
                              "C" = "c"),
318 318
                          drvNames = c("a1", "a2", "b", "c"))
319
-evalAllGenotypes(sv2, order = FALSE, addwt = TRUE)
319
+evalAllGenotypes(sv2, addwt = TRUE)
320 320
 
321 321
 ## 2. Simulate the data. Here we use the "McFL" model and set explicitly
322 322
 ##    parameters for mutation rate, initial size and size of the population
... ...
@@ -413,7 +413,7 @@ pancr <- allFitnessEffects(
413 413
                sh = -0.9,
414 414
                typeDep = "MN"),
415 415
     drvNames = c("KRAS", "SMAD4", "CDNK2A", "TP53", "MLL3", "TGFBR2", "PXDN"))
416
-## How does it look like?
416
+## How does the fitnessEffects object look like?
417 417
 plot(pancr)
418 418
 @ 
419 419
 <<fig.width=6.5, fig.height=10>>=
... ...
@@ -606,11 +606,16 @@ fem4 <- allFitnessEffects(genotFitness = m4)
606 606
 (The message is just telling you what the program guess you wanted.)
607 607
 
608 608
 
609
-That's it. You can plot that object
609
+That's it. You can plot that fitnessEffects object
610 610
 <<out.width="6cm", out.height = "6cm">>=
611 611
 plot(fem4)
612 612
 @ 
613
-or see what OncoSimulR thinks the fitnesses are, with the
613
+
614
+In this case, that plot is not very interesting (compare with the
615
+\texttt{plot(pancr)} we saw in \qref{quickexample} or the plots in
616
+\qref{posetslong}).  
617
+
618
+You can also check what OncoSimulR thinks the fitnesses are, with the
614 619
 \Rfunction{evalAllGenotypes} function that we will use repeatedly below
615 620
 (of course, here we should see the same fitnesses we entered):
616 621
 
... ...
@@ -619,6 +624,12 @@ evalAllGenotypes(fem4)
619 624
 @ 
620 625
 
621 626
 
627
+And you can plot the fitness landscape:
628
+
629
+<<>>=
630
+plotFitnessLandscape(evalAllGenotypes(fem4))
631
+@ 
632
+
622 633
 
623 634
 To specify the mapping you can also use a matrix (or data frame) with
624 635
 $g + 1$ columns; each of the first $g$ columns contains a 1 or a 0
... ...
@@ -629,8 +640,9 @@ genotypes (we will assume that the missing genotypes have fitness 1):
629 640
 <<out.width="6cm", out.height = "6cm">>=
630 641
 m6 <- cbind(c(1, 1), c(1, 0), c(2, 3))
631 642
 fem6 <- allFitnessEffects(genotFitness = m6)
632
-evalAllGenotypes(fem6, addwt = TRUE, order = FALSE)
643
+evalAllGenotypes(fem6, addwt = TRUE)
633 644
 plot(fem6)
645
+plotFitnessLandscape(evalAllGenotypes(fem6))
634 646
 @ 
635 647
 
636 648
 
... ...
@@ -1532,7 +1544,7 @@ o3 <- allFitnessEffects(orderEffects = c(
1532 1544
                               "D" = "d") )
1533 1545
 
1534 1546
 
1535
-(ag <- evalAllGenotypes(o3, addwt = TRUE))
1547
+(ag <- evalAllGenotypes(o3, addwt = TRUE, order = TRUE))
1536 1548
 @ 
1537 1549
 
1538 1550
 %% <<>>=
... ...
@@ -1576,7 +1588,7 @@ ofe1 <- allFitnessEffects(orderEffects = c("F > D" = -0.3, "D > F" = 0.4),
1576 1588
                               c("F" = "f1, f2",
1577 1589
                                 "D" = "d1, d2") )
1578 1590
 
1579
-ag <- evalAllGenotypes(ofe1)
1591
+ag <- evalAllGenotypes(ofe1, order = TRUE)
1580 1592
 
1581 1593
 @ 
1582 1594
 
... ...
@@ -1644,7 +1656,7 @@ ofe2 <- allFitnessEffects(orderEffects = c("F > D" = -0.3, "D > F" = 0.4),
1644 1656
                           geneToModule =
1645 1657
                               c("F" = "f1, f2, f3",
1646 1658
                                 "D" = "d1, d2") )
1647
-ag2 <- evalAllGenotypes(ofe2, max = 325)
1659
+ag2 <- evalAllGenotypes(ofe2, max = 325, order = TRUE)
1648 1660
 
1649 1661
 @ 
1650 1662
 
... ...
@@ -1689,7 +1701,7 @@ We can get the fitness of all genotypes (we set $max = 325$ because that
1689 1701
 is the number of possible genotypes):
1690 1702
 
1691 1703
 <<>>=
1692
-agoi1 <- evalAllGenotypes(foi1,  max = 325)
1704
+agoi1 <- evalAllGenotypes(foi1,  max = 325, order = TRUE)
1693 1705
 head(agoi1)
1694 1706
 @ 
1695 1707
 
... ...
@@ -4009,7 +4021,7 @@ very different for genotypes with the same number of mutations, and does
4009 4021
 not increase in a simple way with drivers:
4010 4022
 
4011 4023
 <<>>=
4012
-evalAllGenotypes(examplesFitnessEffects[["o3"]], addwt = TRUE)
4024
+evalAllGenotypes(examplesFitnessEffects[["o3"]], addwt = TRUE, order = TRUE)
4013 4025
 @ 
4014 4026
 
4015 4027
 A few figures could help:
... ...
@@ -1,15 +1,15 @@
1 1
 \usepackage[%
2
-		shash={66cf4f1},
3
-		lhash={66cf4f192d83ab3aa9ff159405ab08c1ffd05be3},
4
-		authname={Ramon Diaz-Uriarte},
5
-		authemail={rdiaz02@users.noreply.github.com},
6
-		authsdate={2016-07-11},
7
-		authidate={2016-07-11 08:27:58 +0200},
8
-		authudate={1468218478},
9
-		commname={Ramon Diaz-Uriarte},
10
-		commemail={rdiaz02@users.noreply.github.com},
11
-		commsdate={2016-07-11},
12
-		commidate={2016-07-11 08:27:58 +0200},
13
-		commudate={1468218478},
2
+		shash={55f3e48},
3
+		lhash={55f3e487f322e90ab69c172ab452f5b5401fa8a2},
4
+		authname={ramon diaz-uriarte (at Phelsuma)},
5
+		authemail={rdiaz02@gmail.com},
6
+		authsdate={2016-07-30},
7
+		authidate={2016-07-30 11:53:43 +0200},
8
+		authudate={1469872423},
9
+		commname={ramon diaz-uriarte (at Phelsuma)},
10
+		commemail={rdiaz02@gmail.com},
11
+		commsdate={2016-07-30},
12
+		commidate={2016-07-30 11:53:43 +0200},
13
+		commudate={1469872423},
14 14
 		refnames={ (HEAD -> master, origin/master, origin/HEAD)}
15 15
 	]{gitsetinfo}
16 16
\ No newline at end of file