qSig {signatureSearch} | R Documentation |
qSig
ObjectIt builds a 'qSig' object to store the query signature, reference database and GESS method used for GESS methods
qSig(query, gess_method, refdb)
query |
If 'gess_method' is 'CMAP' or 'LINCS', it should be a list with
two character vectors named If 'gess_method' is 'gCMAP', the query is a matrix with a single column representing gene ranks from a biological state of interest. The corresponding gene labels are stored in the row name slot of the matrix. Instead of ranks one can provide scores (e.g. z-scores). In such a case, the scores will be internally transformed to ranks. If 'gess_method' is 'Fisher', the query is expected to be a list with two
character vectors named If 'gess_method' is 'Cor', the query is a matrix with a single numeric column and the gene labels in the row name slot. The numeric column can contain z-scores, LFCs, (normalized) gene expression intensity values or read counts. |
gess_method |
one of 'CMAP', 'LINCS', 'gCMAP', 'Fisher' or 'Cor' |
refdb |
character(1), can be one of "cmap", "cmap_expr", "lincs", or
"lincs_expr" when using the CMAP/LINCS databases from the affiliated
To use a custom signature database, it should be the file path to the HDF5
file generated with the When the |
qSig
object
build_custom_db
,
signatureSearchData
,
gmt2h5
, qSig-class
db_path <- system.file("extdata", "sample_db.h5", package = "signatureSearch") ## Load sample_db as `SummarizedExperiment` object library(SummarizedExperiment); library(HDF5Array) sample_db <- SummarizedExperiment(HDF5Array(db_path, name="assay")) rownames(sample_db) <- HDF5Array(db_path, name="rownames") colnames(sample_db) <- HDF5Array(db_path, name="colnames") ## get "vorinostat__SKB__trt_cp" signature drawn from sample database query_mat <- as.matrix(assay(sample_db[,"vorinostat__SKB__trt_cp"])) query = as.numeric(query_mat); names(query) = rownames(query_mat) upset <- head(names(query[order(-query)]), 150) downset <- tail(names(query[order(-query)]), 150) qsig_lincs <- qSig(query=list(upset=upset, downset=downset), gess_method="LINCS", refdb=db_path) qsig_gcmap <- qSig(query=query_mat, gess_method="gCMAP", refdb=db_path)