The rankGenes function is a generic function that can deal with mutilple types of inputs. Given a matrix of gene expression that has samples in columns, genes in rows, and values being gene expression intensity,rankGenes ranks gene expression intensities in each sample.

It can also work with S4 objects that have gene expression matrix as a component (i.e ExpressionSet, DGEList,SummarizedExperiment). It calls the rank function in the base package which ranks the gene expression matrix by its absolute expression level. If the input is S4 object of DGEList, ExpressionSet, or SummarizedExperiment, it will extract the gene expression matrix from the object and rank the genes. The default 'tiesMethod' is set to 'min'.

rankGenes(expreMatrix, tiesMethod = "min", stableGenes = NULL)

# S4 method for matrix
rankGenes(expreMatrix, tiesMethod = "min", stableGenes = NULL)

# S4 method for data.frame
rankGenes(expreMatrix, tiesMethod = "min", stableGenes = NULL)

# S4 method for DGEList
rankGenes(expreMatrix, tiesMethod = "min", stableGenes = NULL)

# S4 method for ExpressionSet
rankGenes(expreMatrix, tiesMethod = "min", stableGenes = NULL)

# S4 method for SummarizedExperiment
rankGenes(expreMatrix, tiesMethod = "min", stableGenes = NULL)

Arguments

expreMatrix

matrix, data.frame, ExpressionSet, DGEList or SummarizedExperiment storing gene expression measurements

tiesMethod

character, indicating what method to use when dealing with ties

stableGenes

character, containing a list of stable genes to be used to rank genes using expression of stable genes. This is required when using the stable genes dependent version of singscore (see details in simpleScore). Stable genes for solid cancers (carcinomas) and blood transcriptomes can be obtained using the getStableGenes function

Value

The ranked gene expression matrix that has samples in columns and genes in rows. Unit normalised ranks are returned if data is ranked using stable genes

See also

getStableGenes, simpleScore, rank, "ExpressionSet", "SummarizedExperiment", "DGEList"

Examples

rankGenes(toy_expr_se) # toy_expr_se is a gene expression dataset
#>    D_Ctrl_R1 D_TGFb_R1
#> 2          2         2
#> 9         13        12
#> 10         4         5
#> 12        10         8
#> 13         8         7
#> 14        18        18
#> 15         3         1
#> 16        20        20
#> 18         6        10
#> 19        14        14
#> 20         9         9
#> 21        11        11
#> 22        17        16
#> 23        19        19
#> 24         5         4
#> 25        16        17
#> 26         1         3
#> 27        12        13
#> 28         7         6
#> 29        15        15
#> attr(,"stable")
#> [1] FALSE

# ExpressionSet object
emat <- SummarizedExperiment::assay(toy_expr_se)
e <- Biobase::ExpressionSet(assayData = as.matrix(emat))
rankGenes(e)
#>    D_Ctrl_R1 D_TGFb_R1
#> 2          2         2
#> 9         13        12
#> 10         4         5
#> 12        10         8
#> 13         8         7
#> 14        18        18
#> 15         3         1
#> 16        20        20
#> 18         6        10
#> 19        14        14
#> 20         9         9
#> 21        11        11
#> 22        17        16
#> 23        19        19
#> 24         5         4
#> 25        16        17
#> 26         1         3
#> 27        12        13
#> 28         7         6
#> 29        15        15
#> attr(,"stable")
#> [1] FALSE

#scoring using the stable version of singscore
rankGenes(e, stableGenes = c('2', '20', '25'))
#>    D_Ctrl_R1 D_TGFb_R1
#> 2       0.25      0.25
#> 9       0.75      0.75
#> 10      0.50      0.50
#> 12      0.75      0.50
#> 13      0.50      0.50
#> 14      1.00      1.00
#> 15      0.50      0.25
#> 16      1.00      1.00
#> 18      0.50      0.75
#> 19      0.75      0.75
#> 20      0.50      0.50
#> 21      0.75      0.75
#> 22      1.00      0.75
#> 23      1.00      1.00
#> 24      0.50      0.50
#> 25      0.75      0.75
#> 26      0.25      0.50
#> 27      0.75      0.75
#> 28      0.50      0.50
#> 29      0.75      0.75
#> attr(,"stable")
#> [1] TRUE

if (FALSE) {
#for real cancer or blood datasets, use getStableGenes()
rankGenes(cancer_expr, stableGenes = getStableGenes(5))
rankGenes(blood_expr, stableGenes = getStableGenes(5, type = 'blood'))
}