computeRSquared {AffiXcan} | R Documentation |
Compute R and R^2 between rows of two SummarizedExperiment assays
computeRSquared( realExpr, imputedExpr, assay, testingSamples = NULL, BPPARAM = bpparam() )
realExpr |
A SummarizedExperiment object containing expression data |
imputedExpr |
The returning object of affiXcanImpute() |
assay |
A string with the name of the object in SummarizedExperiment::assays(realExpr) that contains expression values |
testingSamples |
A vector of strings. The identifiers of the samples that have to be considered by the function; default is NULL; if is.null(testingSamples)==TRUE then no filtering is performed |
BPPARAM |
A BiocParallelParam object. Default is bpparam(). For details on BiocParallelParam virtual base class see browseVignettes("BiocParallel") |
A list of lists; inner lists are named after the rows for which the correlation between realExpr and imputedExpr have been computed; inner lists contain two objects:
rho: the pearson's correlation coefficient (R) between the real expression values and the imputed GReX for the cross-validation i on the expressed gene y, computed with cor()
rho.sq: the coefficient of determination (R^2) between the real expression values and the imputed GReX for the cross-validation i on the expressed gene y, computed as pearson^2
cor.test.p.val: the p-value of the cor.test() between the real expression values and the imputed GReX for the cross-validation i on the expressed gene y
if (interactive()) { trainingTbaPaths <- system.file("extdata","training.tba.toydata.rds", package="AffiXcan") data(exprMatrix) data(regionAssoc) data(trainingCovariates) assay <- "values" training <- affiXcanTrain(exprMatrix=exprMatrix, assay=assay, tbaPaths=trainingTbaPaths, regionAssoc=regionAssoc, cov=trainingCovariates, varExplained=80, scale=TRUE) imputedExpr <- affiXcanImpute(tbaPaths=trainingTbaPaths, affiXcanTraining=training, scale=TRUE) realExpr <- exprMatrix correlation <- computeRSquared(realExpr, imputedExpr, assay) }