## ----global_options, include=FALSE-------------------------------------------- knitr::opts_chunk$set(fig.width = 6, fig.height = 6, fig.path = 'figures/') ## ----load, message = FALSE---------------------------------------------------- library(ropls) ## ----sacurine, message = FALSE------------------------------------------------ data(sacurine) names(sacurine) ## ----attach_code, message = FALSE--------------------------------------------- attach(sacurine) ## ----strF--------------------------------------------------------------------- view(dataMatrix) view(sampleMetadata) view(variableMetadata) ## ----pca_code, eval = FALSE--------------------------------------------------- # sacurine.pca <- opls(dataMatrix) ## ----pca_result, echo = FALSE------------------------------------------------- sacurine.pca <- opls(dataMatrix, fig.pdfC = "none") ## ----pca_figure, echo = FALSE, fig.show = 'hold'------------------------------ plot(sacurine.pca) ## ----pca-col------------------------------------------------------------------ genderFc <- sampleMetadata[, "gender"] plot(sacurine.pca, typeVc = "x-score", parAsColFcVn = genderFc) ## ----pca-col-personalized----------------------------------------------------- plot(sacurine.pca, typeVc = "x-score", parAsColFcVn = genderFc, parLabVc = as.character(sampleMetadata[, "age"]), parPaletteVc = c("green4", "magenta")) ## ----plsda-------------------------------------------------------------------- sacurine.plsda <- opls(dataMatrix, genderFc) ## ----oplsda------------------------------------------------------------------- sacurine.oplsda <- opls(dataMatrix, genderFc, predI = 1, orthoI = NA) ## ----oplsda_subset, warning=FALSE--------------------------------------------- sacurine.oplsda <- opls(dataMatrix, genderFc, predI = 1, orthoI = NA, subset = "odd") ## ----train-------------------------------------------------------------------- trainVi <- getSubsetVi(sacurine.oplsda) confusion_train.tb <- table(genderFc[trainVi], fitted(sacurine.oplsda)) confusion_train.tb ## ----test--------------------------------------------------------------------- confusion_test.tb <- table(genderFc[-trainVi], predict(sacurine.oplsda, dataMatrix[-trainVi, ])) confusion_test.tb ## ----overfit, echo = FALSE---------------------------------------------------- set.seed(123) obsI <- 20 featVi <- c(2, 20, 200) featMaxI <- max(featVi) xRandMN <- matrix(runif(obsI * featMaxI), nrow = obsI) yRandVn <- sample(c(rep(0, obsI / 2), rep(1, obsI / 2))) layout(matrix(1:4, nrow = 2, byrow = TRUE)) for (featI in featVi) { randPlsi <- opls(xRandMN[, 1:featI], yRandVn, predI = 2, permI = ifelse(featI == featMaxI, 100, 0), fig.pdfC = "none", info.txtC = "none") plot(randPlsi, typeVc = "x-score", parCexN = 1.3, parTitleL = FALSE, parCexMetricN = 0.5) mtext(featI/obsI, font = 2, line = 2) if (featI == featMaxI) plot(randPlsi, typeVc = "permutation", parCexN = 1.3) } mtext(" obs./feat. ratio:", adj = 0, at = 0, font = 2, line = -2, outer = TRUE) ## ----vip---------------------------------------------------------------------- ageVn <- sampleMetadata[, "age"] pvaVn <- apply(dataMatrix, 2, function(feaVn) cor.test(ageVn, feaVn)[["p.value"]]) vipVn <- getVipVn(opls(dataMatrix, ageVn, predI = 1, orthoI = NA, fig.pdfC = "none")) quantVn <- qnorm(1 - pvaVn / 2) rmsQuantN <- sqrt(mean(quantVn^2)) opar <- par(font = 2, font.axis = 2, font.lab = 2, las = 1, mar = c(5.1, 4.6, 4.1, 2.1), lwd = 2, pch = 16) plot(pvaVn, vipVn, col = "red", pch = 16, xlab = "p-value", ylab = "VIP", xaxs = "i", yaxs = "i") box(lwd = 2) curve(qnorm(1 - x / 2) / rmsQuantN, 0, 1, add = TRUE, col = "red", lwd = 3) abline(h = 1, col = "blue") abline(v = 0.05, col = "blue") par(opar) ## ----get_se------------------------------------------------------------------- data(sacurine) sac.se <- sacurine[["se"]] ## ----se_plsda, echo = TRUE, results = "hide"---------------------------------- sac.se <- opls(sac.se, "gender") ## ----se_updated, message = FALSE---------------------------------------------- library(SummarizedExperiment) SummarizedExperiment::colData(sac.se) SummarizedExperiment::rowData(sac.se) ## ----se_model----------------------------------------------------------------- sac_opls.ls <- getOpls(sac.se) names(sac_opls.ls) plot(sac_opls.ls[["gender_PLSDA"]], typeVc = "x-score") ## ----get_eset----------------------------------------------------------------- data("sacurine") sac.set <- sacurine[["eset"]] # viewing the ExpressionSet # ropls::view(sac.set) ## ----eset_plsda, echo = TRUE, results = "hide"-------------------------------- # performing the PLS-DA sac.plsda <- opls(sac.set, "gender") ## ----eset_updated------------------------------------------------------------- sac.set <- getEset(sac.plsda) library(Biobase) head(Biobase::pData(sac.set)) ## ----nci60_mae, message = FALSE----------------------------------------------- data("NCI60") nci.mae <- NCI60[["mae"]] ## ----mae_plsda---------------------------------------------------------------- nci.mae <- opls(nci.mae, "cancer", predI = 2, fig.pdfC = "none") ## ----mae_updated-------------------------------------------------------------- SummarizedExperiment::colData(nci.mae[["agilent"]]) ## ----mae_model---------------------------------------------------------------- mae_opls.ls <- getOpls(nci.mae) names(mae_opls.ls) plot(mae_opls.ls[["agilent"]][["cancer_PLSDA"]], typeVc = "x-score") ## ----get_mds------------------------------------------------------------------ data("NCI60") nci.mds <- NCI60[["mds"]] ## ----mds_plsda, echo = TRUE, results = "hide"--------------------------------- # Restricting to the "agilent" and "hgu95" datasets nci.mds <- nci.mds[, c("agilent", "hgu95")] # Restricting to the 'ME' and 'LE' cancer types library(Biobase) sample_names.vc <- Biobase::sampleNames(nci.mds[["agilent"]]) cancer_type.vc <- Biobase::pData(nci.mds[["agilent"]])[, "cancer"] nci.mds <- nci.mds[sample_names.vc[cancer_type.vc %in% c("ME", "LE")], ] # Building PLS-DA models for the cancer type nci.plsda <- ropls::opls(nci.mds, "cancer", predI = 2) ## ----mds_getmset-------------------------------------------------------------- nci.mds <- ropls::getMset(nci.plsda) ## ----phenomis, eval = FALSE--------------------------------------------------- # install.packages("devtools") # devtools::install_github("https://github.com/odisce/phenomis") # data(sacurine) # sacurine.se <- sacurine[["se"]] # library(phenomis) # phenomis::writing(sacurine.se, dir.c = getwd()) ## ----detach------------------------------------------------------------------- detach(sacurine) ## ----se_build----------------------------------------------------------------- # Preparing the data (matrix) and sample and variable metadata (data frames): data(sacurine, package = "ropls") data.mn <- sacurine[["dataMatrix"]] # matrix: samples x variables samp.df <- sacurine[["sampleMetadata"]] # data frame: samples x sample metadata feat.df <- sacurine[["variableMetadata"]] # data frame: features x feature metadata # Creating the SummarizedExperiment (package SummarizedExperiment) library(SummarizedExperiment) sac.se <- SummarizedExperiment(assays = list(sacurine = t(data.mn)), colData = samp.df, rowData = feat.df) # note that colData and rowData main format is DataFrame, but data frames are accepted when building the object stopifnot(validObject(sac.se)) # Viewing the SummarizedExperiment # ropls::view(sac.se) ## ----mae_build_load----------------------------------------------------------- data("NCI60_4arrays", package = "omicade4") ## ----mae_build, message = FALSE, warning=FALSE-------------------------------- library(MultiAssayExperiment) # Building the individual SummarizedExperiment instances experiment.ls <- list() sampleMap.ls <- list() for (set.c in names(NCI60_4arrays)) { # Getting the data and metadata assay.mn <- as.matrix(NCI60_4arrays[[set.c]]) coldata.df <- data.frame(row.names = colnames(assay.mn), .id = colnames(assay.mn), stringsAsFactors = FALSE) # the 'cancer' information will be stored in the main colData of the mae, not the individual SummarizedExperiments rowdata.df <- data.frame(row.names = rownames(assay.mn), name = rownames(assay.mn), stringsAsFactors = FALSE) # Building the SummarizedExperiment assay.ls <- list(se = assay.mn) names(assay.ls) <- set.c se <- SummarizedExperiment(assays = assay.ls, colData = coldata.df, rowData = rowdata.df) experiment.ls[[set.c]] <- se sampleMap.ls[[set.c]] <- data.frame(primary = colnames(se), colname = colnames(se)) # both datasets use identical sample names } sampleMap <- listToMap(sampleMap.ls) # The sample metadata are stored in the colData data frame (both datasets have the same samples) stopifnot(identical(colnames(NCI60_4arrays[[1]]), colnames(NCI60_4arrays[[2]]))) sample_names.vc <- colnames(NCI60_4arrays[[1]]) colData.df <- data.frame(row.names = sample_names.vc, cancer = substr(sample_names.vc, 1, 2)) nci.mae <- MultiAssayExperiment(experiments = experiment.ls, colData = colData.df, sampleMap = sampleMap) stopifnot(validObject(nci.mae)) ## ----eset_build, message = FALSE, warning = FALSE----------------------------- # Preparing the data (matrix) and sample and variable metadata (data frames): data(sacurine) data.mn <- sacurine[["dataMatrix"]] # matrix: samples x variables samp.df <- sacurine[["sampleMetadata"]] # data frame: samples x sample metadata feat.df <- sacurine[["variableMetadata"]] # data frame: features x feature metadata # Creating the ExpressionSet (package Biobase) sac.set <- Biobase::ExpressionSet(assayData = t(data.mn)) Biobase::pData(sac.set) <- samp.df Biobase::fData(sac.set) <- feat.df stopifnot(validObject(sac.set)) # Viewing the ExpressionSet # ropls::view(sac.set) ## ----mset_build_load---------------------------------------------------------- data("NCI60_4arrays", package = "omicade4") ## ----mset_build, message = FALSE, warning=FALSE------------------------------- library(MultiDataSet) # Creating the MultiDataSet instance nci.mds <- MultiDataSet::createMultiDataSet() # Adding the two datasets as ExpressionSet instances for (set.c in names(NCI60_4arrays)) { # Getting the data expr.mn <- as.matrix(NCI60_4arrays[[set.c]]) pdata.df <- data.frame(row.names = colnames(expr.mn), cancer = substr(colnames(expr.mn), 1, 2), stringsAsFactors = FALSE) fdata.df <- data.frame(row.names = rownames(expr.mn), name = rownames(expr.mn), stringsAsFactors = FALSE) # Building the ExpressionSet eset <- Biobase::ExpressionSet(assayData = expr.mn, phenoData = new("AnnotatedDataFrame", data = pdata.df), featureData = new("AnnotatedDataFrame", data = fdata.df), experimentData = new("MIAME", title = set.c)) # Adding to the MultiDataSet nci.mds <- MultiDataSet::add_eset(nci.mds, eset, dataset.type = set.c, GRanges = NA, warnings = FALSE) } stopifnot(validObject(nci.mds)) ## ----sessionInfo, echo=FALSE-------------------------------------------------- sessionInfo()