## ----style, eval=TRUE, echo=FALSE, results="asis"-------------------------- BiocStyle::latex() ## ----loading, fig.width=6, fig.height=6, echo=TRUE------------------------- library(Pigengene) data(aml) data(mds) d1 <- rbind(aml,mds) Labels <- c(rep("AML",nrow(aml)),rep("MDS",nrow(mds))) names(Labels) <- rownames(d1) Disease <- as.data.frame(Labels) h1 <- pheatmap.type(d1[,1:20],annRow=Disease,show_rownames=FALSE) ## ----oneStep, echo=TRUE---------------------------------------------------- p1 <- one.step.pigengene(Data=d1,saveDir='pigengene', bnNum=0, verbose=1, seed=1, Labels=Labels, toCompact=FALSE, doHeat=FALSE) ## ----tree, fig.width=5, fig.height=5, echo=TRUE---------------------------- plot(p1$c5treeRes$c5Trees[["34"]]) ## ----pigengene, fig.width=5, fig.height=5, echo=TRUE----------------------- dim(p1$pigengene$eigengenes) p1 <- pheatmap.type(p1$pigengene$eigengenes,annRow=Disease,show_rownames=FALSE) ## ----qc, echo=TRUE--------------------------------------------------------- ## QC: checked <- check.pigengene.input(Data=d1, Labels=Labels) DataI <- checked$Data LabelsI <- checked$Labels ## ----balance, echo=TRUE---------------------------------------------------- wData <- balance(Data=DataI, Labels=LabelsI, verbose=1)$balanced ## ----beta, echo=TRUE------------------------------------------------------- betaI <- calculate.beta(RsquaredCut=0.8, Data=wData, verbose=1)$power saveDir <- "steps" ## Results will be saved in this folder. dir.create(saveDir) ## ----wgcna, echo=TRUE------------------------------------------------------ ## WGCNA wgRes <- wgcna.one.step(Data=wData, seed=1, power=betaI, saveDir=saveDir, verbose=1) ## ----eigengenes, echo=TRUE------------------------------------------------- ## Eigengenes: pigengene <- compute.pigengene(Data=DataI, Labels=LabelsI, saveFile=file.path(saveDir, 'pigengene.RData'), modules=wgRes$modules, verbose=1) class(pigengene) print(dim(pigengene$eigengenes)) ##This is the eigengenes matrix. print(pigengene$eigengenes[1:3,1:4]) ## ----bn, echo=TRUE--------------------------------------------------------- ## Learning the BN structure: library(bnlearn) learnt <- learn.bn(pigengene=pigengene, bnPath=file.path(saveDir, "bn"), bnNum=10, ## In real applications, at least 100-1000. seed=1, verbose=1) BN <- learnt$consensus1$BN ## ----draw, echo=TRUE, eval=FALSE------------------------------------------- # drawn <- draw.bn(BN) # ## By construction, the Disease node can have no parents. ## ----bnParam, echo=TRUE---------------------------------------------------- ## Fit the parameters of the Bayesian network: fit <- bn.fit(x=BN, data=learnt$consensus1$Data, method="bayes", iss=10) ##where learnt$consensus1$Data is the discretized data matrix. ## The conditional probability table for one of the children of the Disease node: selectedFeatures <- children("Disease", x=BN) print(fit[[selectedFeatures[1]]]) ## ----bnPrediction, echo=TRUE----------------------------------------------- l2 <- predict(object=fit, node="Disease", data=learnt$consensus1$Data, method="bayes-lw") table(LabelsI, l2) ## ----trees, echo=TRUE------------------------------------------------------ ## Decision trees: treePath <- file.path(saveDir, 'C5Trees') dir.create(path=treePath) treeRes <- make.decision.tree(pigengene=pigengene, Data=DataI, selectedFeatures=selectedFeatures, saveDir=treePath, verbose=1, toCompact=FALSE) ## ----selectedFeatures, echo=TRUE------------------------------------------- ## Access tree results usedFeatures <- get.used.features(c5Tree=treeRes$c5Trees[["17"]]) moduleMembers <- setNames(paste0("ME", wgRes$modules), nm=sub("_[^_]+$", "", names(wgRes$modules))) modMembersUsed <- moduleMembers[moduleMembers %in% usedFeatures] moduList <- split(names(modMembersUsed), f=modMembersUsed) library(org.Hs.eg.db) pw1 <- get.enriched.pw(genes=moduList, idType="ENTREZID", pathwayDb="KEGG", Org=NULL, OrgDb=org.Hs.eg.db, outPath=saveDir, verbose=1) ## ----citation, results='asis', eval=TRUE----------------------------------- citation("Pigengene") ## ----sessionInfo, results='asis', eval=TRUE-------------------------------- toLatex(sessionInfo())