Browse code

1.9.3: plots, vignette, max warning and win. failure

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

Ramon Diaz-Uriarte authored on 19/06/2015 23:36:35
Showing 9 changed files

... ...
@@ -1,7 +1,7 @@
1 1
 Package: OncoSimulR
2 2
 Type: Package
3 3
 Title: Forward Genetic Simulation of Cancer Progresion with Epistasis 
4
-Version: 1.99.2
4
+Version: 1.99.3
5 5
 Date: 2015-06-19
6 6
 Author: Ramon Diaz-Uriarte. 
7 7
 Maintainer: Ramon Diaz-Uriarte <rdiaz02@gmail.com>
... ...
@@ -23,4 +23,4 @@ Imports: Rcpp (>= 0.11.1), parallel, data.table, graph, Rgraphviz, gtools, igrap
23 23
 Suggests: BiocStyle, knitr, Oncotree, testthat
24 24
 LinkingTo: Rcpp
25 25
 VignetteBuilder: knitr
26
-SystemRequirements: C++11
26
+# SystemRequirements: C++11
... ...
@@ -723,15 +723,16 @@ fitnessEffectsToIgraph <- function(rT, epistasis, orderEffects) {
723 723
 
724 724
 
725 725
 plot.fitnessEffects <- function(x, type = "graphNEL",
726
-                               layout = NULL,
727
-                               expandModules = FALSE,
728
-                               ...) {
726
+                                layout = NULL,
727
+                                expandModules = FALSE,
728
+                                autofit = FALSE,
729
+                                return_g = FALSE,
730
+                                ...) {
729 731
     ## some other layouts I find OK
730 732
     ## layout.circle
731 733
     ## layout.reingold.tilford if really a tree
732 734
     ## o.w. it will use the default
733 735
     g <- x$graph
734
-
735 736
     
736 737
     if(type == "igraph") {
737 738
         if(expandModules && (!x$gMOneToOne)) {
... ...
@@ -739,7 +740,18 @@ plot.fitnessEffects <- function(x, type = "graphNEL",
739 740
             vlabels <- x$geneToModule[V(g)$name]
740 741
             V(g)$label <- vlabels
741 742
         }
742
-        plot.igraph(g, layout = layout)
743
+        if(autofit) {
744
+            plot(0, type = "n", axes = FALSE, ann = FALSE)
745
+            ## ideas from http://stackoverflow.com/questions/14472079/match-vertex-size-to-label-size-in-igraph
746
+            ## vsize <- (strwidth(V(g)$label) + strwidth("oo")) * 200
747
+            ## but this is a kludge.
748
+            vsize <- (nchar(V(g)$label) + 3) * 4.5
749
+            plot.igraph(g, vertex.size = vsize, vertex.shape = "rectangle",
750
+                        layout = layout)
751
+        } else {
752
+            plot.igraph(g, layout = layout)
753
+        }
754
+        if(return_g) return(g)
743 755
     }
744 756
     else if (type == "graphNEL") {
745 757
         g1 <- igraph.to.graphNEL(g)
... ...
@@ -758,9 +770,19 @@ plot.fitnessEffects <- function(x, type = "graphNEL",
758 770
             names(nnodes) <- nodes(g1)
759 771
             nAttrs$label <- nnodes
760 772
         }
761
-        plot(g1, edgeAttrs = list(arrowsize = a1, lty = s1, lwd = lwd,
762
-                     color = c1),
763
-             nodeAttrs = nAttrs)
773
+        if(autofit) {
774
+            nAttrs$width <- (nchar(nAttrs$label) + 1)/10
775
+            names(nAttrs$width) <- names(nAttrs$label)
776
+            plot(g1, edgeAttrs = list(arrowsize = a1, lty = s1, lwd = lwd,
777
+                         color = c1), attrs=list(node=list(shape = "rectangle")),
778
+                 nodeAttrs = nAttrs)
779
+            
780
+        } else {
781
+            plot(g1, edgeAttrs = list(arrowsize = a1, lty = s1, lwd = lwd,
782
+                         color = c1),
783
+                 nodeAttrs = nAttrs)
784
+        }
785
+        if(return_g) return(g1)
764 786
     } else {
765 787
         stop("plot type not recognized")
766 788
     }
... ...
@@ -1,3 +1,9 @@
1
+Changes in version 1.99.3 (2015-06-19):
2
+	- More examples to vignette
3
+	- Using Makevars
4
+	- More functionality to plot.fitnessEffects
5
+	- Will Windoze work now?
6
+
1 7
 Changes in version 1.99.2 (2015-06-19):
2 8
 	- Fixed typos and other minor in vignett.
3 9
 
... ...
@@ -10,7 +10,7 @@
10 10
 }
11 11
 \usage{
12 12
 \method{plot}{fitnessEffects}(x, type = "graphNEL", layout = NULL,
13
-expandModules = FALSE, ...)
13
+expandModules = FALSE, autofit = FALSE, return_g = FALSE, ...)
14 14
 
15 15
 %%plotFitnessEffects(fe, type = "graphNEL", layout = NULL, expandModules = FALSE)
16 16
 }
... ...
@@ -31,6 +31,12 @@ expandModules = FALSE, ...)
31 31
    If there are modules with multiple genes, if you set this to TRUE
32 32
   modules will be replaced by their genes.
33 33
 }
34
+
35
+\item{autofit}{If TRUE, we try to fit the edges to the labels. This is a
36
+  very experimental feature, likely to be not very robust.}
37
+
38
+\item{return_g}{It TRUE, the graph object (graphNEL or igrap) is returned.}
39
+
34 40
  \item{\dots}{
35 41
     Other arguments passed to \code{plot}. Not used for now.
36 42
   }
... ...
@@ -45,6 +51,8 @@ expandModules = FALSE, ...)
45 51
   type. Order relationships have an arrow from the earlier to the later
46 52
   event and have a different dotted line (lty 3).
47 53
 
54
+  If \code{return_g} is TRUE, you are returned also the graph object
55
+  (igraph or graphNEL) so that you can manipulate it further.
48 56
 }
49 57
 \author{
50 58
 Ramon Diaz-Uriarte
51 59
new file mode 100644
... ...
@@ -0,0 +1,3 @@
1
+## This is a C++11 package
2
+CXX_STD = CXX11
3
+  
0 4
new file mode 100644
... ...
@@ -0,0 +1,3 @@
1
+## This is a C++11 package
2
+CXX_STD = CXX11
3
+  
... ...
@@ -40,18 +40,19 @@ double prodFitness(std::vector<double> s) {
40 40
 		    [](double x, double y) {return (x * std::max(0.0, (1 + y)));});
41 41
 }
42 42
 
43
-// This is a better idea? If yes, change code in nr_fitness so that
44
-// birth 0 is not extinction.
45
-double prodFitnessNegInf(std::vector<double> s) {
46
-  double f = 1.0;
47
-  for(auto y : s) {
48
-    if( y == -std::numeric_limits<double>::infinity()) {
49
-      return -std::numeric_limits<double>::infinity();
50
-    } else {
51
-      f *= std::max(0.0, (1 + y));
52
-    }
53
-  }
54
-}
43
+// // This is a better idea? If yes, change code in nr_fitness so that
44
+// // birth 0 is not extinction.
45
+// double prodFitnessNegInf(std::vector<double> s) {
46
+//   double f = 1.0;
47
+//   for(auto y : s) {
48
+//     if( y == -std::numeric_limits<double>::infinity()) {
49
+//       return -std::numeric_limits<double>::infinity();
50
+//     } else {
51
+//       f *= std::max(0.0, (1 + y));
52
+//     }
53
+//   }
54
+//   return f;
55
+// }
55 56
 
56 57
 double prodDeathFitness(std::vector<double> s) {
57 58
   return accumulate(s.begin(), s.end(), 1.0,
... ...
@@ -239,10 +240,10 @@ std::vector<Poset_struct> rTable_to_Poset(Rcpp::List rt) {
239 240
     if(! is_sorted(Poset[i].parentsNumID.begin(), Poset[i].parentsNumID.end()) )
240 241
       throw std::logic_error("ParentsNumID not sorted");
241 242
     
242
-    if(isinf(Poset[i].s))
243
+    if(std::isinf(Poset[i].s))
243 244
       Rcpp::Rcout << "WARNING: at least one s is infinite" 
244 245
 		  << std::endl;
245
-    if(isinf(Poset[i].sh) && (Poset[i].sh > 0))
246
+    if(std::isinf(Poset[i].sh) && (Poset[i].sh > 0))
246 247
       Rcpp::Rcout << "WARNING: at least one sh is positive infinite" 
247 248
 		  << std::endl;
248 249
   }
... ...
@@ -938,9 +939,9 @@ void printPoset(const std::vector<Poset_struct>& Poset) {
938 939
       Rcpp::Rcout <<"\t\t typeDep = " << depToString(Poset[i].typeDep) << ' ' ;
939 940
       Rcpp::Rcout <<"\t s = " << Poset[i].s << " ";
940 941
       Rcpp::Rcout <<"\t sh = " << Poset[i].sh << std::endl;
941
-      if(isinf(Poset[i].sh))
942
+      if(std::isinf(Poset[i].sh))
942 943
 	++counterInfs;
943
-      if(isinf(Poset[i].sh) && (Poset[i].sh < 0))
944
+      if(std::isinf(Poset[i].sh) && (Poset[i].sh < 0))
944 945
 	++counterNegInfs;
945 946
       Rcpp::Rcout << "\t\t Number of parent modules or genes = " << 
946 947
 	Poset[i].parents.size() << std::endl;
... ...
@@ -1105,436 +1106,3 @@ double evalRGenotype(Rcpp::IntegerVector rG, Rcpp::List rFE,
1105 1106
 }
1106 1107
 
1107 1108
 
1108
-
1109
-// // [[Rcpp::export]]
1110
-// void ovalRGenotype(Rcpp::IntegerVector rG, Rcpp::List rFE, bool verbose) {
1111
-//   const Rcpp::List rF(rFE);
1112
-//   fitnessEffectsAll F = convertFitnessEffects(rF);
1113
-//   Genotype g = convertGenotypeFromR(rG, F);
1114
-//   vector<double> s = evalGenotypeFitness(g, F);
1115
-//   // if(verbose) {
1116
-//   //   Rcpp::Rcout << "\n Individual s terms are :";
1117
-//   //   for(auto i : s) Rcpp::Rcout << " " << i;
1118
-//   //   Rcpp::Rcout << std::endl;
1119
-//   // }
1120
-// }
1121
-
1122
-// // [[Rcpp::export]]
1123
-// void dummyconvertFitnessEffects(Rcpp::List rFE) {
1124
-//   // Yes, some of the things below are data.frames in R, but for
1125
-//   // us that is used just as a list.
1126
-
1127
-//   fitnessEffectsAll fe;
1128
-  
1129
-//   Rcpp::List rrt = rFE["long.rt"];
1130
-//   Rcpp::List re = rFE["long.epistasis"];
1131
-//   Rcpp::List ro = rFE["long.orderEffects"];
1132
-//   Rcpp::List rgi = rFE["long.geneNoInt"];
1133
-//   Rcpp::List rgm = rFE["geneModule"];
1134
-//   bool rone = rFE["gMOneToOne"];
1135
-  
1136
-
1137
-//   // if(rrt.size()) {
1138
-//   //   fe.Poset = rTable_to_Poset(rrt);
1139
-//   // } 
1140
-//   // if(re.size()) {
1141
-//   //   fe.Epistasis = convertEpiOrderEff(re);
1142
-//   // } 
1143
-//   // if(ro.size()) {
1144
-//   //   fe.orderE = convertEpiOrderEff(ro);
1145
-//   // } 
1146
-//   if(rgi.size()) {
1147
-//     fe.genesNoInt = convertNoInts(rgi);
1148
-//   } else {
1149
-//     fe.genesNoInt.shift = -9L;
1150
-//   }
1151
-//   // fe.Gene_Module_tabl = R_GeneModuleToGeneModule(rgm);
1152
-//   // fe.allOrderG = sortedAllOrder(fe.orderE);
1153
-//   // fe.allPosetG = sortedAllPoset(fe.Poset);
1154
-//   // fe.gMOneToOne = rone;
1155
-
1156
-//   //return fe;
1157
-// }
1158
-
1159
-
1160
-
1161
-// // [[Rcpp::export]]
1162
-// void dummy(Rcpp::List rt) { 
1163
-
1164
-//   // The restriction table, or Poset, has a first element
1165
-//   // with nothing, so that all references by mutated gene
1166
-//   // are simply accessing the Poset[mutated gene] without
1167
-//   // having to remember to add 1, etc.
1168
-
1169
-//   std::vector<Poset_struct> Poset;
1170
-  
1171
-//   Poset.resize(rt.size() + 1);
1172
-//   Poset[0].child = "0"; //should this be Root?? I don't think so.
1173
-//   Poset[0].childNumID = 0;
1174
-//   Poset[0].typeDep = Dependency::NA;
1175
-//   Poset[0].s = std::numeric_limits<double>::quiet_NaN();
1176
-//   Poset[0].sh = std::numeric_limits<double>::quiet_NaN();
1177
-//   Poset[0].parents.resize(0);
1178
-//   Poset[0].parentsNumID.resize(0);
1179
-
1180
-//   Rcpp::List rt_element;
1181
-//   // std::string tmpname;
1182
-//   Rcpp::IntegerVector parentsid;
1183
-//   Rcpp::CharacterVector parents;
1184
-
1185
-//   for(int i = 1; i != (rt.size() + 1); ++i) {
1186
-//     rt_element = rt[i - 1];
1187
-//     Poset[i].child = Rcpp::as<std::string>(rt_element["child"]);
1188
-//     Poset[i].childNumID = as<int>(rt_element["childNumID"]);
1189
-//     Poset[i].typeDep = stringToDep(as<std::string>(rt_element["typeDep"]));
1190
-//     Poset[i].s = as<double>(rt_element["s"]);
1191
-//     Poset[i].sh = as<double>(rt_element["sh"]);
1192
-
1193
-//     // if(i != Poset[i].childNumID) {
1194
-//     // Nope: this assumes we only deal with posets.
1195
-//     //   // Rcpp::Rcout << "\n childNumID, original = " << as<int>(rt_element["childNumID"]);
1196
-//     //   // Rcpp::Rcout << "\n childNumID, Poset = " << Poset[i].childNumID;
1197
-//     //   // Rcpp::Rcout << "\n i = " << i << std::endl;
1198
-//     //   throw std::logic_error("childNumID != index");
1199
-//     // }
1200
-//     // Rcpp::IntegerVector parentsid = as<Rcpp::IntegerVector>(rt_element["parentsNumID"]);
1201
-//     // Rcpp::CharacterVector parents = as<Rcpp::CharacterVector>(rt_element["parents"]);
1202
-
1203
-//     parentsid = as<Rcpp::IntegerVector>(rt_element["parentsNumID"]);
1204
-//     parents = as<Rcpp::CharacterVector>(rt_element["parents"]);
1205
-
1206
-//     if( parentsid.size() != parents.size() ) {
1207
-//       throw std::logic_error("parents size != parentsNumID size");
1208
-//     }
1209
-
1210
-//     for(int j = 0; j != parentsid.size(); ++j) {
1211
-//       Poset[i].parentsNumID.push_back(parentsid[j]);
1212
-//       Poset[i].parents.push_back( (Rcpp::as< std::string >(parents[j])) );
1213
-//       // tmpname = Rcpp::as< std::string >(parents[j]);
1214
-//       // Poset[i].parents.push_back(tmpname);
1215
-//     }
1216
-
1217
-//     // Should not be needed if R always does what is should. Disable later?
1218
-//     if(! is_sorted(Poset[i].parentsNumID.begin(), Poset[i].parentsNumID.end()) )
1219
-//       throw std::logic_error("ParentsNumID not sorted");
1220
-    
1221
-//     if(isinf(Poset[i].s))
1222
-//       Rcpp::Rcout << "WARNING: at least one s is infinite" 
1223
-// 		  << std::endl;
1224
-//     if(isinf(Poset[i].sh) && (Poset[i].sh > 0))
1225
-//       Rcpp::Rcout << "WARNING: at least one sh is positive infinite" 
1226
-// 		  << std::endl;
1227
-//   }
1228
-// }
1229
-
1230
-// // [[Rcpp::export]]
1231
-// void dummy2(Rcpp::List rFE) {
1232
-//   // Yes, some of the things below are data.frames in R, but for
1233
-//   // us that is used just as a list.
1234
-
1235
-//   fitnessEffectsAll fe;
1236
-  
1237
-//   Rcpp::List rrt = rFE["long.rt"];
1238
-//   Rcpp::List re = rFE["long.epistasis"];
1239
-//   Rcpp::List ro = rFE["long.orderEffects"];
1240
-//   Rcpp::List rgi = rFE["long.geneNoInt"];
1241
-//   Rcpp::List rgm = rFE["geneModule"];
1242
-//   bool rone = rFE["gMOneToOne"];
1243
-  
1244
-
1245
-//   if(rrt.size()) {
1246
-//     fe.Poset = rTable_to_Poset(rrt);
1247
-//   } 
1248
-//   if(re.size()) {
1249
-//     fe.Epistasis = convertEpiOrderEff(re);
1250
-//   } 
1251
-//   if(ro.size()) {
1252
-//     fe.orderE = convertEpiOrderEff(ro);
1253
-//   } 
1254
-//   // if(rgi.size()) {
1255
-//   //   fe.genesNoInt = convertNoInts(rgi);
1256
-//   // } else {
1257
-//   //   fe.genesNoInt.shift = -9L;
1258
-//   // }
1259
-  
1260
-//   // fe.Gene_Module_tabl = R_GeneModuleToGeneModule(rgm);
1261
-//   // fe.allOrderG = sortedAllOrder(fe.orderE);
1262
-//   // fe.allPosetG = sortedAllPoset(fe.Poset);
1263
-//   fe.gMOneToOne = rone;
1264
-// }
1265
-
1266
-
1267
-// // [[Rcpp::export]]
1268
-// void wrap_test_rt(Rcpp::List rtR, Rcpp::DataFrame rGM) {
1269
-//   std::vector<Poset_struct> Poset;
1270
-//   std::vector<geneToModule> geneModules;
1271
-//   std::vector<geneToModuleLong> geneModules;
1272
-
1273
-//   R_GeneModuleToGeneModule(rGM, geneModules, geneModules);
1274
-//   printGeneToModule(geneModules);
1275
-
1276
-//   rTable_to_Poset(rtR, Poset);
1277
-//   printPoset(Poset);
1278
-// }
1279
-
1280
-// // stretching vector of mutatedDrv unlikely to speed up;
1281
-// // access (seeing if a module gene in there) is O(1) but I do
1282
-// // that for every gene and then sum over the vector (linear in
1283
-// // number of positions?)
1284
-
1285
-// // Could be made much faster if we could assume lower numbered
1286
-// // positions cannot depend on higher numbered ones. Nope, not that much.
1287
-
1288
-// SEXP wrap_test_checkRestriction(Rcpp::List rtR, 
1289
-// 				Rcpp::DataFrame rGM,
1290
-// 				Rcpp::IntegerVector genotype) {
1291
-
1292
-
1293
-// Turn the following into a function called from R naturally.  Allows for
1294
-// testing AND allows users to understand the consequences of an rt and a
1295
-// genotype. The R function is evalGenotype
1296
-
1297
-
1298
-// // [[Rcpp::export]]
1299
-// SEXP eval_Genotype(Rcpp::List rtR, 
1300
-// 		   Rcpp::DataFrame rGM,
1301
-// 		   Rcpp::IntegerVector genotype) {
1302
-
1303
-  
1304
-//   // std::vector<geneToModule> geneModules;
1305
-//   std::vector<geneToModuleLong> geneModules;
1306
-
1307
-//   std::vector<int> Genotype;
1308
-//   std::vector<double> s_vector;
1309
-//   std::vector<double> sh_vector;
1310
-
1311
-//   std::vector<Poset_struct> Poset = rTable_to_Poset(rtR);
1312
-//   geneModules = R_GeneModuleToGeneModule(rGM, geneModules, geneModules);
1313
-
1314
-
1315
-//   Genotype = Rcpp::as< std::vector<int> >(genotype);
1316
-
1317
-//   evalPosetConstraints(Genotype, Poset, geneModules,  s_vector, sh_vector);
1318
-
1319
-//   // Rcpp::Rcout << "s_vector:\n ";
1320
-//   // for (auto c : s_vector)
1321
-//   //   Rcpp::Rcout << c << ' ';
1322
-//   // Rcpp::Rcout << std::endl;
1323
-//   // Rcpp::Rcout << "sh_vector:\n ";
1324
-//   // for (auto c : sh_vector)
1325
-//   //   Rcpp::Rcout << c << ' ';
1326
-//   // Rcpp::Rcout << std::endl;
1327
-
1328
-//   return
1329
-//     List::create(Named("s_vector") = s_vector,
1330
-// 		 Named("sh_vector") = sh_vector);
1331
-// }
1332
-
1333
-
1334
-// when mutation happens, if a passenger, insert in mutatedPass and give p_vector.
1335
-  // if(mutatedPos >= numDrivers) { //the new mutation is a passenger
1336
-  //   // do something: iterate through the passenger part of genotype
1337
-  //   // keeping the rest of the parent stuff.
1338
-  //   return;
1339
-  // } else {
1340
-
1341
-
1342
-// if driver, insert in Drv.
1343
-// and check restrictions
1344
-
1345
-
1346
-// static void printOrderEffects(const std::vector<epistasis>& orderE) {
1347
-//   // do something
1348
-// }
1349
-
1350
-
1351
-
1352
-
1353
-
1354
-
1355
-
1356
-
1357
-// // [[Rcpp::export]]
1358
-// void wrap_FitnessEffects(Rcpp::List rtR,
1359
-// 			Rcpp::List longEpistasis,
1360
-// 			Rcpp::List longOrderE,
1361
-// 			Rcpp::DataFrame rGM,
1362
-// 			Rcpp::DataFrame geneNoInt,
1363
-// 			Rcpp::IntegerVector gMOneToOne) {
1364
-
1365
-   
1366
-//   std::vector<Poset_struct> Poset;
1367
-//   std::vector<geneToModule> geneModules;
1368
-//   std::vector<geneToModuleLong> geneModules;
1369
-
1370
-//   convertFitnessEffects(rtR, longEpistasis, longOrderE,
1371
-// 			rGM, geneNoInt, gMOneToOne,
1372
-// 			Poset, geneModules, geneModules);
1373
-//   printGeneToModule(geneModules);
1374
-//   printPoset(Poset);
1375
-// }
1376
-
1377
-
1378
-
1379
-// void readFitnessEffects(Rcpp::List rtR,
1380
-// 			       Rcpp::List longEpistasis,
1381
-// 			       Rcpp::List longOrderE,
1382
-// 			       Rcpp::DataFrame rGM,
1383
-// 			       Rcpp::DataFrame geneNoInt,
1384
-// 			       Rcpp::IntegerVector gMOneToOne,
1385
-// 			       std::vector<Poset_struct>&  Poset,
1386
-// 			       std::vector<geneToModule>& geneModules,
1387
-// 			       std::vector<geneToModuleLong>& geneModules
1388
-// 			       // return objects
1389
-// 			       // something else
1390
-// 			       ) {
1391
-//    rTable_to_Poset(rtR, Poset);
1392
-//    R_GeneModuleToGeneModule(rGM, geneModules, geneModules);
1393
- 
1394
-// }
1395
-
1396
-
1397
-// epistasis and order can also be of modules, as we map modules to genes
1398
-// for that.
1399
-
1400
-
1401
-
1402
-// How we store a given genotype and how we store the restrictions need to
1403
-// have a 1-to-1 mapping.
1404
-
1405
-
1406
-// checkEpistasis: for every epistatic set, check if there?
1407
-
1408
-// checkOrderEffects: we find each member in the vector, but the indices
1409
-// must be correlative.
1410
-
1411
-
1412
-
1413
-
1414
-
1415
-
1416
-
1417
-
1418
-
1419
-
1420
-
1421
-
1422
-
1423
-
1424
-
1425
-
1426
-
1427
-
1428
-
1429
-// static void evalPosetConstraints(const std::vector<int>& Drv,
1430
-// 			     const std::vector<Poset_struct>& Poset,
1431
-// 			     const std::vector<geneToModule>& geneToModules,
1432
-// 			     std::vector<double>& s_vector,
1433
-// 			     std::vector<double>& sh_vector) {
1434
-//   size_t numDeps;
1435
-//   size_t sumDepsMet = 0;
1436
-//   int parent_module_mutated = 0;
1437
-//   std::vector<int> mutatedModules;
1438
-  
1439
-
1440
-//   // Map mutated genes to modules (mutatedModules) and then examine, for
1441
-//   // each mutatedModule, if its dependencies (in terms of modules) are met.
1442
-
1443
-//   mutatedModules = DrvToModule(Drv, geneToModules);
1444
-
1445
-//   for(auto it_mutatedModule = mutatedModules.begin();
1446
-//       it_mutatedModule != mutatedModules.end(); ++it_mutatedModule) {
1447
-//     if( (Poset[(*it_mutatedModule)].parentsNumID.size() == 1) &&
1448
-// 	(Poset[(*it_mutatedModule)].parentsNumID[0] == 0) ) { //Depends only on root.
1449
-//       // FIXME: isn't it enough to check the second condition?
1450
-//       s_vector.push_back(Poset[(*it_mutatedModule)].s);
1451
-//     } else {
1452
-//       sumDepsMet = 0;
1453
-//       numDeps = Poset[(*it_mutatedModule)].parentsNumID.size();
1454
-//       for(auto it_Parents = Poset[(*it_mutatedModule)].parentsNumID.begin();
1455
-// 	  it_Parents != Poset[(*it_mutatedModule)].parentsNumID.end();
1456
-// 	  ++it_Parents) {
1457
-// 	// if sorted, could use binary search
1458
-// 	// FIXME: try a set or sort mutatedModules?
1459
-
1460
-// 	//parent_module_mutated =
1461
-// 	//std::binary_search(mutatedModules.begin(), mutatedModules.end(),
1462
-// 	//(*it_Parents))
1463
-
1464
-// 	parent_module_mutated = 
1465
-// 	  (std::find(mutatedModules.begin(), 
1466
-// 		     mutatedModules.end(), 
1467
-// 		     (*it_Parents)) != mutatedModules.end());
1468
-// 	if(parent_module_mutated)  {
1469
-// 	  ++sumDepsMet;
1470
-// 	  if( Poset[(*it_mutatedModule)].typeDep == Dependency::semimonotone)
1471
-// 	    break;
1472
-// 	}
1473
-//       }
1474
-//       if( ((Poset[(*it_mutatedModule)].typeDep == Dependency::semimonotone) && 
1475
-// 	   (sumDepsMet)) ||
1476
-// 	  ((Poset[(*it_mutatedModule)].typeDep == Dependency::monotone) && 
1477
-// 	   (sumDepsMet == numDeps)) ) {
1478
-// 	s_vector.push_back(Poset[(*it_mutatedModule)].s);
1479
-//       } else {
1480
-// 	sh_vector.push_back(Poset[(*it_mutatedModule)].sh);
1481
-//       }
1482
-//     }
1483
-//   }
1484
-  
1485
-// }
1486
-
1487
-
1488
-
1489
-// std::vector<double> evalPosetConstraints0(const std::vector<int>& mutatedModules,
1490
-// 					 const std::vector<Poset_struct>& Poset) {
1491
-
1492
-//   size_t numDeps;
1493
-//   size_t sumDepsMet = 0;
1494
-//   int parent_module_mutated = 0;
1495
-//   // std::vector<int> mutatedModules;
1496
-  
1497
-//   std::vector<double> s;
1498
-  
1499
-//   // Map mutated genes to modules (mutatedModules) and then examine, for
1500
-//   // each mutatedModule, if its dependencies (in terms of modules) are met.
1501
-
1502
-  
1503
-
1504
-//   for(auto it_mutatedModule = mutatedModules.begin();
1505
-//       it_mutatedModule != mutatedModules.end(); ++it_mutatedModule) {
1506
-//     if( (Poset[(*it_mutatedModule)].parentsNumID.size() == 1) &&
1507
-// 	(Poset[(*it_mutatedModule)].parentsNumID[0] == 0) ) { //Depends only on root.
1508
-//       // FIXME: isn't it enough to check the second condition?
1509
-//       s.push_back(Poset[(*it_mutatedModule)].s);
1510
-//     } else {
1511
-//       sumDepsMet = 0;
1512
-//       numDeps = Poset[(*it_mutatedModule)].parentsNumID.size();
1513
-//       for(auto it_Parents = Poset[(*it_mutatedModule)].parentsNumID.begin();
1514
-// 	  it_Parents != Poset[(*it_mutatedModule)].parentsNumID.end();
1515
-// 	  ++it_Parents) {
1516
-
1517
-// 	parent_module_mutated = binary_search(mutatedModules.begin(), 
1518
-// 					      mutatedModules.end(),
1519
-// 					      (*it_Parents));
1520
-// 	  // (std::find(mutatedModules.begin(), 
1521
-// 	  // 	     mutatedModules.end(), 
1522
-// 	  // 	     (*it_Parents)) != mutatedModules.end());
1523
-// 	if(parent_module_mutated)  {
1524
-// 	  ++sumDepsMet;
1525
-// 	  if( Poset[(*it_mutatedModule)].typeDep == Dependency::semimonotone)
1526
-// 	    break;
1527
-// 	}
1528
-//       }
1529
-//       if( ((Poset[(*it_mutatedModule)].typeDep == Dependency::semimonotone) && 
1530
-// 	   (sumDepsMet)) ||
1531
-// 	  ((Poset[(*it_mutatedModule)].typeDep == Dependency::monotone) && 
1532
-// 	   (sumDepsMet == numDeps)) ) {
1533
-// 	s.push_back(Poset[(*it_mutatedModule)].s);
1534
-//       } else {
1535
-// 	s.push_back(Poset[(*it_mutatedModule)].sh);
1536
-//       }
1537
-//     }
1538
-//   }
1539
-//   return s;
1540
-// }
... ...
@@ -242,7 +242,7 @@ Before anything else, let us load the package. We also explicitly load
242 242
 \Biocpkg{graph} and \CRANpkg{igraph} for the vignette to work (you do not
243 243
 need that for your usual interactive work).
244 244
 
245
-<<echo=FALSE>>=
245
+<<results="hide">>=
246 246
 library(OncoSimulR)
247 247
 library(graph)
248 248
 library(igraph)
... ...
@@ -2314,6 +2314,53 @@ plot(pancr)
2314 2314
 Of course the ``s'' and ``sh'' are set arbitrarily here.
2315 2315
 
2316 2316
 
2317
+\subsection{Raphael and Vandin's modules}\label{raphael-ex}
2318
+
2319
+In \cite{Raphael2014a}, Raphael and Vandin show several progression models
2320
+in terms of modules. The extended poset for the colorectal cancer model in
2321
+their Figure 4.a is (s and sh are arbitrary):
2322
+
2323
+
2324
+<<fig.height = 4>>=
2325
+rv1 <- allFitnessEffects(data.frame(parent = c("Root", "A", "KRAS"),
2326
+                                    child = c("A", "KRAS", "FBXW7"),
2327
+                                    s = 0.1,
2328
+                                    sh = -0.01,
2329
+                                    typeDep = "MN"),
2330
+                         geneToModule = c("Root" = "Root",
2331
+                             "A" = "EVC2, PIK3CA, TP53",
2332
+                             "KRAS" = "KRAS",
2333
+                             "FBXW7" = "FBXW7"))
2334
+
2335
+plot(rv1, expandModules = TRUE, autofit = TRUE)
2336
+
2337
+@ 
2338
+
2339
+We use the (experimental) \Rcode{autofit} option to fit the labels to the
2340
+edges.
2341
+
2342
+
2343
+Their Figure 5b is
2344
+
2345
+<<fig.height=12>>=
2346
+
2347
+rv2 <- allFitnessEffects(data.frame(parent = c("Root", "1", "2", "3", "4"),
2348
+                                    child = c("1", "2", "3", "4", "ELF3"),
2349
+                                    s = 0.1,
2350
+                                    sh = -0.01,
2351
+                                    typeDep = "MN"),
2352
+                         geneToModule = c("Root" = "Root",
2353
+                             "1" = "APC, FBXW7",
2354
+                             "2" = "ATM, FAM123B, PIK3CA, TP53",
2355
+                             "3" = "BRAF, KRAS, NRAS",
2356
+                             "4" = "SMAD2, SMAD4, SOX9",
2357
+                             "ELF3" = "ELF3"))
2358
+plot(rv2, "igraph", expandModules = TRUE, 
2359
+     layout = layout.reingold.tilford,
2360
+     autofit = TRUE)
2361
+
2362
+@ 
2363
+
2317 2364
 \section{Running the simulations}\label{simul}
2318 2365
 
2319 2366
 
... ...
@@ -2629,7 +2676,7 @@ Often, you will want to simulate multiple runs of the same scenario, and
2629 2676
 do something with them. Conceptually, the first step is running multiple
2630 2677
 simulations and, then, sampling them.
2631 2678
 
2632
-We will use the ``pancreas'' example, above \section{pancreas}
2679
+We will use the ``pancreas'' example, above section \ref{pancreas}.
2633 2680
 <<>>=
2634 2681
 
2635 2682
 
... ...
@@ -1,15 +1,15 @@
1 1
 \usepackage[%
2
-		shash={4d8bb43},
3
-		lhash={4d8bb439bd91e1f2a877322a8d59485a34678e58},
4
-		authname={ramon diaz-uriarte (at Bufo)},
2
+		shash={1db4b3a},
3
+		lhash={1db4b3a4511805f1d9914d24b269f5e10f4753a1},
4
+		authname={Ramon Diaz-Uriarte (at Coleonyx)},
5 5
 		authemail={rdiaz02@gmail.com},
6 6
 		authsdate={2015-06-19},
7
-		authidate={2015-06-19 14:00:52 +0200},
8
-		authudate={1434715252},
9
-		commname={ramon diaz-uriarte (at Bufo)},
7
+		authidate={2015-06-19 22:51:22 +0200},
8
+		authudate={1434747082},
9
+		commname={Ramon Diaz-Uriarte (at Coleonyx)},
10 10
 		commemail={rdiaz02@gmail.com},
11 11
 		commsdate={2015-06-19},
12
-		commidate={2015-06-19 14:00:52 +0200},
13
-		commudate={1434715252},
14
-		refnames={ (HEAD, master)}
12
+		commidate={2015-06-19 22:51:22 +0200},
13
+		commudate={1434747082},
14
+		refnames={ (HEAD, origin/master, origin/HEAD, master)}
15 15
 	]{gitsetinfo}
16 16
\ No newline at end of file