SHrinkage covariance Incorporating Prior knowledge
SHIP.Rd
The SHIP-package implements the shrinkage estimator of a covariance matrix given any covariance target, such as described by Schaefer and Strimmer in 2005. In addition, it proposes several targets based on biological knowledge extracted from the public database KEGG.
Details
To use the shrinkage estimator, one should just have at hand a data set in the form of a \(n \times p\) matrix, and a covariance target.
If one wishes to use the proposed targets, the data set should be compatible with KEGG, i.e. it should be possible to extract for each gene the pathways it belongs to. This information, for example, can be found in libraries such as hgu133plus2.db.
References
J. Schaefer and K. Strimmer, 2005. A shrinkage approach to large-scale covariance matrix estimation and implications for functional genomics. Statist. Appl. Genet. Mol. Biol. 4:32.
M. Jelizarow, V. Guillemot, A. Tenenhaus, K. Strimmer, A.-L. Boulesteix, 2010. Over-optimism in bioinformatics: an illustration. Bioinformatics. Accepted.
Examples
# A short example on a toy dataset
# require(SHIP)
data(expl)
attach(expl)
sig1 <- shrink.estim(x,targetD(x))
#> Error in shrink.estim(x, targetD(x)): could not find function "shrink.estim"
sig2 <- shrink.estim(x,targetF(x))
#> Error in shrink.estim(x, targetF(x)): could not find function "shrink.estim"
sig3 <- shrink.estim(x,targetCor(x,genegroups))
#> Error in shrink.estim(x, targetCor(x, genegroups)): could not find function "shrink.estim"
sig4 <- shrink.estim(x,targetG(x,genegroups))
#> Error in shrink.estim(x, targetG(x, genegroups)): could not find function "shrink.estim"
paste(sig1[[2]],collapse=" ")
#> Error in eval(expr, envir, enclos): object 'sig1' not found
paste(sig2[[2]],collapse=" ")
#> Error in eval(expr, envir, enclos): object 'sig2' not found
paste(sig3[[2]],collapse=" ")
#> Error in eval(expr, envir, enclos): object 'sig3' not found
paste(sig4[[2]],collapse=" ")
#> Error in eval(expr, envir, enclos): object 'sig4' not found
if (FALSE) {
# Example on how to get the gene groups lists
require(hgu95av2.db)
# e.g. we have some interesting gene names :
vec <- c("MYC","ID2","PTGER4","ATF4","FGFR1","MET","HLA-DRB6")
# we then want to convert them into Probe Sets
symb <- as.list(hgu95av2SYMBOL)
pbsets <- names(symb[unlist(sapply(vec,function(x,l) which(l==x)[1],symb))])
# Probe Sets which are themselves converted into a gene groups list
genegroups <- as.list(hgu95av2PATH)[pbsets]
}