inst/doc/CellNOptR-vignette.R

In CellNOptR: Training of boolean logic models of signalling networks using prior knowledge networks and perturbation data

### R code from vignette source 'CellNOptR-vignette.Rnw'

###################################################
### code chunk number 1: installRBGL (eval = FALSE)
###################################################
## if (!requireNamespace("BiocManager", quietly=TRUE))
##     install.packages("BiocManager")
## BiocManager::install("RBGL")


###################################################
### code chunk number 2: installPackage (eval = FALSE)
###################################################
## if (!requireNamespace("BiocManager", quietly=TRUE))
##     install.packages("BiocManager")
## BiocManager::install("CellNOptR")


###################################################
### code chunk number 3: installPackage2 (eval = FALSE)
###################################################
## install.packages("path_to_CellNOptR/CellNOptR_1.0.0.tar.gz",
##     repos=NULL, type="source")


###################################################
### code chunk number 4: Ropts
###################################################
options(width=70)


###################################################
### code chunk number 5: loadLib
###################################################
library(CellNOptR)


###################################################
### code chunk number 6: newDir (eval = FALSE)
###################################################
## dir.create("CNOR_analysis")


###################################################
### code chunk number 7: quickstart (eval = FALSE)
###################################################
## # ---------------------- load the library and get a SIF and MIDAS file
## library(CellNOptR)
## #
## # ---------------------- examples are provided in CellNOptR
## data("ToyModel", package="CellNOptR")
## data("CNOlistToy", package="CellNOptR")
## pknmodel = ToyModel
## cnolist = CNOlist(CNOlistToy)
## #
## # ---------------------- alternatively you can read your own files:
## # pknmodel = readSIF("ToyModel.sif")
## # cnolist = CNOlist("ToyDataMMB.csv")
## #
## # ---------------------- preprocess the network
## model = preprocessing(cnolist, pknmodel)
## #
## # ---------------------- perform the analysis
## res = gaBinaryT1(cnolist, model, verbose=FALSE)
## #
## # ---------------------- plot the results
## cutAndPlot(cnolist, model, list(res$bString))


###################################################
### code chunk number 8: directory
###################################################
cpfile<-dir(system.file("ToyModel",package="CellNOptR"),full=TRUE)
file.copy(from=cpfile,to=getwd(),overwrite=TRUE)


###################################################
### code chunk number 9: getData
###################################################
dataToy<-readMIDAS("ToyDataMMB.csv", verbose=FALSE)
CNOlistToy<-makeCNOlist(dataToy,subfield=FALSE, verbose=FALSE)


###################################################
### code chunk number 10: getData2
###################################################
data(CNOlistToy,package="CellNOptR", verbose=FALSE)


###################################################
### code chunk number 11: cnolistClass
###################################################
CNOlistToy = CNOlist("ToyDataMMB.csv")


###################################################
### code chunk number 12: cnolistClass2 (eval = FALSE)
###################################################
## data(CNOlistToy,package="CellNOptR")
## CNOlistToy = CNOlist(CNOlistToy)


###################################################
### code chunk number 13: showCNO
###################################################
CNOlistToy


###################################################
### code chunk number 14: plotCNO
###################################################
plot(CNOlistToy)


###################################################
### code chunk number 15: ploCNOPDF
###################################################
plotCNOlistPDF(CNOlist=CNOlistToy,filename="ToyModelGraph.pdf")


###################################################
### code chunk number 16: getModel
###################################################
pknmodel<-readSIF("ToyPKNMMB.sif")


###################################################
### code chunk number 17: getModel2
###################################################
data(ToyModel,package="CellNOptR")


###################################################
### code chunk number 18: indices
###################################################
checkSignals(CNOlistToy,pknmodel)


###################################################
### code chunk number 19: plotModel
###################################################
plotModel(pknmodel, CNOlistToy)


###################################################
### code chunk number 20: NONC
###################################################
indicesToy<-indexFinder(CNOlistToy,pknmodel,verbose=TRUE)
ToyNCNOindices<-findNONC(pknmodel,indicesToy,verbose=TRUE)
ToyNCNOcut<-cutNONC(pknmodel,ToyNCNOindices)
indicesToyNCNOcut<-indexFinder(CNOlistToy,ToyNCNOcut)


###################################################
### code chunk number 21: compress
###################################################
ToyNCNOcutComp<-compressModel(ToyNCNOcut,indicesToyNCNOcut)
indicesToyNCNOcutComp<-indexFinder(CNOlistToy,ToyNCNOcutComp)


###################################################
### code chunk number 22: expand
###################################################
model<-expandGates(ToyNCNOcutComp, maxInputsPerGate=3)


###################################################
### code chunk number 23: expand
###################################################
model <- preprocessing(CNOlistToy, pknmodel, expansion=TRUE,
    compression=TRUE, cutNONC=TRUE, verbose=FALSE)


###################################################
### code chunk number 24: resError
###################################################
resECNOlistToy<-residualError(CNOlistToy)


###################################################
### code chunk number 25: initbs
###################################################
initBstring<-rep(1,length(model$reacID))


###################################################
### code chunk number 26: optim
###################################################
ToyT1opt<-gaBinaryT1(CNOlist=CNOlistToy, model=model,
    initBstring=initBstring, verbose=FALSE)


###################################################
### code chunk number 27: resSim
###################################################
cutAndPlot(model=model, bStrings=list(ToyT1opt$bString),
    CNOlist=CNOlistToy,plotPDF=TRUE)


###################################################
### code chunk number 28: plotFit
###################################################
plotFit(optRes=ToyT1opt)


###################################################
### code chunk number 29: simFitPDF
###################################################
cutAndPlot(
    model=model,
    bStrings=list(ToyT1opt$bString),
    CNOlist=CNOlistToy,
    plotPDF=TRUE)
pdf("evolFitToyT1.pdf")
plotFit(optRes=ToyT1opt)
dev.off()


###################################################
### code chunk number 30: plotModelToy1
###################################################
plotModel(model, CNOlistToy, bString=ToyT1opt$bString)


###################################################
### code chunk number 31: plotModelToy1save
###################################################
plotModel(model, CNOlistToy, bString=ToyT1opt$bString, output="PDF", filename="plotModelToy1")


###################################################
### code chunk number 32: plotModelToy2
###################################################
bs = mapBack(model, pknmodel, ToyT1opt$bString)
plotModel(pknmodel, CNOlistToy, bs, compressed=model$speciesCompressed)


###################################################
### code chunk number 33: plotModelToy2save
###################################################
plotModel(pknmodel, CNOlistToy, bs, compressed=model$speciesCompressed,output="PDF", filename="plotModelToy2")


###################################################
### code chunk number 34: writeRes
###################################################
writeScaffold(
    modelComprExpanded=model,
    optimResT1=ToyT1opt,
    optimResT2=NA,
    modelOriginal=pknmodel,
    CNOlist=CNOlistToy)

writeNetwork(
    modelOriginal=pknmodel,
    modelComprExpanded=model,
    optimResT1=ToyT1opt,
    optimResT2=NA,
    CNOlist=CNOlistToy)

namesFilesToy<-list(
    dataPlot="ToyModelGraph.pdf",
    evolFitT1="evolFitToyT1.pdf",
    evolFitT2=NA,
    simResultsT1="SimResultsT1_1.pdf",
    simResultsT2=NA,
    scaffold="Scaffold.sif",
    scaffoldDot="Scaffold.dot",
    tscaffold="TimesScaffold.EA",
    wscaffold="weightsScaffold.EA",
    PKN="PKN.sif",
    PKNdot="PKN.dot",
    wPKN="TimesPKN.EA",
    nPKN="nodesPKN.NA")

writeReport(
    modelOriginal=pknmodel,
    modelOpt=model,
    optimResT1=ToyT1opt,
    optimResT2=NA,
    CNOlist=CNOlistToy,
    directory="testToy",
    namesFiles=namesFilesToy,
    namesData=list(CNOlist="Toy",model="ToyModel"))


###################################################
### code chunk number 35: eraseDir
###################################################
unlink("testToy",recursive=TRUE)


###################################################
### code chunk number 36: getToy
###################################################
dataToy<-readMIDAS("ToyDataMMB.csv")
CNOlistToy<-makeCNOlist(dataToy,subfield=FALSE)
pknmodel<-readSIF("ToyPKNMMB.sif")


###################################################
### code chunk number 37: wrap1 (eval = FALSE)
###################################################
## res <- CNORwrap(
##     paramsList=NA,
##     name="Toy",
##     namesData=list(CNOlist="ToyData",model="ToyModel"),
##     data=CNOlistToy,
##     model=pknmodel)


###################################################
### code chunk number 38: eraseToyDir
###################################################
unlink("Toy",recursive=TRUE)


###################################################
### code chunk number 39: wrap2
###################################################
pList<-defaultParameters(CNOlistToy, pknmodel)
#pList$data = CNOlistToy
#pList$model = ToyModel
res <- CNORwrap(
    paramsList=pList,
    name="Toy1Step",
    namesData=list(CNOlist="ToyData",model="ToyModel"))


###################################################
### code chunk number 40: eraseData
###################################################
unlink("ToyDataMMB.csv")
unlink("ToyPKNMMB.sif")
unlink("Toy1Step",recursive=TRUE)
unlink("ToyModelMMB2.sif")


###################################################
### code chunk number 41: getDREAM
###################################################
#Option 1: copy the SIF and MIDAS files (followed by readMIDAS, makeCNOlist and readSIF)
cpfile<-dir(system.file("DREAMModel",package="CellNOptR"),full=TRUE)
file.copy(from=cpfile,to=getwd(),overwrite=TRUE)
#Option 2: load the CNOlist and model objects
data(CNOlistDREAM,package="CellNOptR")
data(DreamModel,package="CellNOptR")


###################################################
### code chunk number 42: DREAMAnalysis (eval = FALSE)
###################################################
## model = preprocessing(CNOlistDREAM, DreamModel, verbose=FALSE)
## res = gaBinaryT1(CNOlistDREAM, model, verbose=FALSE, maxTime=10)
## cutAndPlot(CNOlistDREAM, model, bStrings=list(res$bString), plotPDF=TRUE,
##     plotParams=list(maxrow=25, margin=0.1, width=20, height=20))
## 
## 


###################################################
### code chunk number 43: t2load
###################################################
data(CNOlistToy2,package="CellNOptR")
data(ToyModel2,package="CellNOptR")
pknmodel = ToyModel2
cnolist = CNOlist(CNOlistToy2)


###################################################
### code chunk number 44: t2plotCNOlist
###################################################
plot(cnolist)
plotCNOlistPDF(cnolist,filename="ToyModelGraphT2.pdf")


###################################################
### code chunk number 45: t2Opt1
###################################################
model = preprocessing(cnolist, pknmodel, verbose=FALSE)
T1opt <- gaBinaryT1(cnolist, model, stallGenMax=10, maxTime=60, verbose=FALSE)


###################################################
### code chunk number 46: t2OptT1plot
###################################################
cutAndPlot(model=model, bStrings=list(T1opt$bString),
    CNOlist=cnolist, plotPDF=TRUE)
pdf("evolFitToy2T1.pdf")
plotFit(optRes=T1opt)
dev.off()
plotFit(optRes=T1opt)


###################################################
### code chunk number 47: t2OptT2
###################################################
T2opt<-gaBinaryTN(cnolist, model,  bStrings=list(T1opt$bString),
    stallGenMax=10, maxTime=60, verbose=FALSE)


###################################################
### code chunk number 48: resSimT2
###################################################
cutAndPlot(
    model=model,
    bStrings=list(T1opt$bString, T2opt$bString),
    CNOlist=cnolist,
    plotPDF=TRUE, plotParams=list(cex=0.8, cmap_scale=0.5, margin=0.2))


###################################################
### code chunk number 49: quickstart (eval = FALSE)
###################################################
## # ---------------------- load the library and get a SIF and MIDAS file
## library(CellNOptR)
## library(stringr)
## #
## # ---------------------- examples are provided in CellNOptR
## data("ToyModel", package="CellNOptR")
## data("CNOlistToy", package="CellNOptR")
## pknmodel = ToyModel
## cnolist = CNOlist(CNOlistToy)
## #
## # ---------------------- alternatively you can read your own files:
## # pknmodel = readSIF("ToyModel.sif")
## # cnolist = CNOlist("ToyDataMMB.csv")
## #
## # ---------------------- preprocess the network
## model = preprocessing(cnolist, pknmodel)
## #
## # ---------------------- perform the analysis
## cplexPath = "path/to/cplex"
## resILP = ilpBinaryT1(cnolist = cnolist, model = model, 
##                      numSolutions = 3, relGap = 0.05,
##                      cplexPath = cplexPath) # asking to retrieve 3 
##                                             # equivalent solutions with 
##                                             # tolerance gap relGap=0.05
## #
## # ---------------------- plot the results (optimized models + fits)
## cutAndPlot(CNOlist = cnolist, model = model, bStrings = list(resILP$bitstringILP[[1]]))
## plotModel(model, cnolist, bString=resILP$bitstringILP[[1]])
## 
## cutAndPlot(CNOlist = cnolist, model = model, bStrings = list(resILP$bitstringILP[[2]]))
## plotModel(model, cnolist, bString=resILP$bitstringILP[[2]])
## 
## cutAndPlot(CNOlist = cnolist, model = model, bStrings = list(resILP$bitstringILP[[3]]))
## plotModel(model, cnolist, bString=resILP$bitstringILP[[3]])


###################################################
### code chunk number 50: quickstart (eval = FALSE)
###################################################
## # ---------------------- load the library and get a SIF and MIDAS file
## library(CellNOptR)
## #
## # ---------------------- examples are provided in CellNOptR
## data(PKN_ToyPB, package="CellNOptR")
## data(CNOlist_ToyPB, package="CellNOptR")
## #
## # ---------------------- preprocess the network
## model = preprocessing(data = cnodata,model = pknmodel,compression = T,expansion = T)
## plotModel(model,cnodata)
## #
## #------------------------ original CNOlist contains many timepoints, we use only a subset
## plot(cnodata)
## selectedTime = c(0,10)
## cnodata_prep = cutCNOlist(cnodata ,model = model,
##                           cutTimeIndices = which(!cnodata@timepoints %in% selectedTime))
## plot(cnodata_prep)
## #
## # ---------------------- perform the analysis
## opt = gaBinaryT1(CNOlist = cnodata_prep,model = model)
## #
## # ---------------------- 10-fold cross-validation procedure using T1 data
## library(doParallel)
## doParallel::registerDoParallel(cores=3)
## system.time({R1=crossvalidateBoolean(CNOlist = cnodata_prep, model = model, 
##                                      type="datapoint", nfolds=10, parallel = TRUE)})
## system.time({R2=crossvalidateBoolean(CNOlist = cnodata_prep, model=model, 
##                                      type="experiment", nfolds=10, parallel = TRUE)})
## system.time({R3=crossvalidateBoolean(CNOlist = cnodata_prep, model=model, 
##                                      type="observable", nfolds=10, parallel = TRUE)})
## system.time({R4=crossvalidateBoolean(CNOlist = cnodata_prep, model=model, 
##                                      type="datapoint", nfolds=10, parallel = FALSE)})
## system.time({R5=crossvalidateBoolean(CNOlist = cnodata_prep, model=model, 
##                                      type="experiment", nfolds=10, parallel = FALSE)})
## system.time({R6=crossvalidateBoolean(CNOlist = cnodata_prep, model=model, 
##                                      type="observable", nfolds=10, parallel = FALSE)})