Title: | Functions for normalisation of RT-qPCR data |
---|---|
Description: | Functions for the selection of optimal reference genes and the normalisation of real-time quantitative PCR data. |
Authors: | Matthias Kohl, James Perkins, Nor Izayu Abdul Rahman |
Maintainer: | James Perkins <[email protected]> |
License: | LGPL-3 |
Version: | 1.53.0 |
Built: | 2024-11-29 06:25:35 UTC |
Source: | https://github.com/bioc/NormqPCR |
Functions for normalisation of real-time quantitative PCR data.
Package: | NormqPCR |
Type: | Package |
Version: | 1.7.1 |
Date: | 2014-08-13 |
Depends: | R(>= 2.14.0), stats, RColorBrewer, Biobase, methods, ReadqPCR, qpcR |
Imports: | ReadqPCR |
biocViews: | MicrotitrePlateAssay, GeneExpression, qPCR |
License: | LGPL-3 |
LazyLoad: | yes |
LazyData: | yes |
require(NormqPCR)
Matthias Kohl, James Perkins, Nor Izayu Abdul Rahman
Maintainer: James Perkins [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
Jo Vandesompele, Katleen De Preter, Filip Pattyn et al. (2002). Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology 2002. 3(7):research0034.1-0034.11. http://genomebiology.com/2002/3/7/research/0034/
Claus Lindbjerg Andersen, Jens Ledet Jensen and Torben Falck Orntoft (2004). Normalization of Real-Time Quantitative Reverse Transcription-PCR Data: A Model-Based Variance Estimation Approach to Identify Genes Suited for Normalization, Applied to Bladder and Colon Cancer Data Sets. CANCER RESEARCH 64, 5245-5250, August 1, 2004. http://cancerres.aacrjournals.org/cgi/content/full/64/15/5245
Kenneth Livak, Thomase Schmittgen (2001). Analysis of Relative Gene Expression Data Using Real-Time Quantitative PCR and the 2^ddCt Method. Methods 25, 402-408, 2001 http://www.ncbi.nlm.nih.gov/pubmed/11846609
## some examples are given in the vignette ## Not run: library(NormqPCR) vignette("NormqPCR") ## End(Not run)
## some examples are given in the vignette ## Not run: library(NormqPCR) vignette("NormqPCR") ## End(Not run)
This dataset was used in Andersen et al (2004) to demonstrate normalization of real-time quantitative RT-PCR data by geometric averaging of housekeeping genes.
data(Bladder)
data(Bladder)
A qPCRBatch object which contains an expression matrix with the expression of
14 genes measured in 28 samples. The sample information is saved in the
phenoData
slot with variables
Sample.no.
sample number.
Grade
Grade of bladder cancer.
The following information on the measured genes is saved in the variables
Symbol
and Gene.name
of the featureData
slot.
ATP5B
ATP synthase, H+ transporting, mitochondrial F1 complex, beta polypetide.
HSPCB
Heat shock 90-kDa protein 1, beta.
S100A6
S100 calcium-binding protein A6 (calcylin).
FLOT2
Flotillin 2.
TEGT
Testis enhanced gene transcript (BAX inhibitor 1).
UBB
Ubiquitin B.
TPT1
Tumor protein, translationally controlled 1.
CFL1
Cofilin 1 (non-muscle).
ACTB
Actin, beta.
RPS13
Ribosomal protein S13.
RPS23
Ribosomal protein S23.
GAPD
Glyceraldehyde-3-phosphate dehydrogenase.
UBC
Ubiquitin C.
FLJ20030
Hypothetical protein FLJ20030.
For a detailed annotation see Table 1 in Anderson et al. (2004).
The genes included in this data set were selected by screening 99 bladder sample expression profiles.
The data set was obtained from http://www.mdl.dk/Publications_sup1.htm
Claus Lindbjerg Andersen, Jens Ledet Jensen and Torben Falck Orntoft (2004). Normalization of Real-Time Quantitative Reverse Transcription-PCR Data: A Model-Based Variance Estimation Approach to Identify Genes Suited for Normalization, Applied to Bladder and Colon Cancer Data Sets. CANCER RESEARCH 64, 5245-5250, August 1, 2004. http://cancerres.aacrjournals.org/cgi/content/full/64/15/5245
data(Bladder) Bladder head(exprs(Bladder)) pData(Bladder) fData(Bladder)
data(Bladder) Bladder head(exprs(Bladder)) pData(Bladder) fData(Bladder)
This dataset was used in Andersen et al (2004) to demonstrate normalization of real-time quantitative RT-PCR data by geometric averaging of housekeeping genes.
data(BladderRepro)
data(BladderRepro)
A qPCRBatch object which contains an expression matrix with the expression of
8 genes measured in 26 samples. The sample information is saved in the
phenoData
slot with variables
Sample.no.
sample number.
Grade
Grade of bladder cancer.
The following information on the measured genes is saved in the variables
Symbol
and Gene.name
of the featureData
slot.
CD14
CD14 antigen.
FCN1
Ficolin (collagen/fibrinogen domain containing) 1.
CCNG2
Cyclin G2.
NPAS2
Neuronal PAS domain protein 2.
UBC
Ubiquitin C.
CFL1
Cofilin 1 (non-muscle).
ACTB
Actin, beta.
GAPD
Glyceraldehyde-3-phosphate dehydrogenase.
For a detailed annotation see Table 1 and Supplementary table 1 in Anderson et al. (2004).
This data set was used to check the reproducibility of the results obtained in Andersen et al (2004).
The data set was obtained from http://www.mdl.dk/Publications_sup1.htm
Claus Lindbjerg Andersen, Jens Ledet Jensen and Torben Falck Orntoft (2004). Normalization of Real-Time Quantitative Reverse Transcription-PCR Data: A Model-Based Variance Estimation Approach to Identify Genes Suited for Normalization, Applied to Bladder and Colon Cancer Data Sets. CANCER RESEARCH 64, 5245-5250, August 1, 2004. http://cancerres.aacrjournals.org/cgi/content/full/64/15/5245
data(BladderRepro) BladderRepro head(exprs(BladderRepro)) pData(BladderRepro) fData(BladderRepro)
data(BladderRepro) BladderRepro head(exprs(BladderRepro)) pData(BladderRepro) fData(BladderRepro)
This dataset was used in Andersen et al (2004) to demonstrate normalization of real-time quantitative RT-PCR data by geometric averaging of housekeeping genes.
data(Colon)
data(Colon)
A qPCRBatch object which contains an expression matrix with the expression of
13 genes measured in 40 samples. The sample information is saved in the
phenoData
slot with variables
Sample.no.
sample number.
Classification
Classification of colon cancer.
The following information on the measured genes is saved in the variables
Symbol
and Gene.name
of the featureData
slot.
UBC
Ubiquitin C.
UBB
Ubiquitin B.
SUI1
Putative translation initiation factor.
NACA
Nascent-polypeptide-associated complex alpha polypeptide.
FLJ20030
Hypothetical protein FLJ20030.
CFL1
Cofilin 1 (non-muscle).
ACTB
Actin, beta.
CLTC
Clathrin, heavy polypeptide (Hc).
RPS13
Ribosomal protein S13.
RPS23
Ribosomal protein S23.
GAPD
Glyceraldehyde-3-phosphate dehydrogenase.
TPT1
Tumor protein, translationally controlled 1.
TUBA6
Tubulin alpha 6.
For a detailed annotation see Table 1 in Anderson et al. (2004).
The genes included in this data set were selected by screening 161 colon sample expression profiles.
The data set was obtained from http://www.mdl.dk/Publications_sup1.htm
Claus Lindbjerg Andersen, Jens Ledet Jensen and Torben Falck Orntoft (2004). Normalization of Real-Time Quantitative Reverse Transcription-PCR Data: A Model-Based Variance Estimation Approach to Identify Genes Suited for Normalization, Applied to Bladder and Colon Cancer Data Sets. CANCER RESEARCH 64, 5245-5250, August 1, 2004. http://cancerres.aacrjournals.org/cgi/content/full/64/15/5245
data(Colon) Colon head(exprs(Colon)) pData(Colon) fData(Colon)
data(Colon) Colon head(exprs(Colon)) pData(Colon) fData(Colon)
Takes expression set of qPCR values containing technical replicates and combines them.
combineTechReps(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' combineTechReps(qPCRBatch, calc="arith")
combineTechReps(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' combineTechReps(qPCRBatch, calc="arith")
qPCRBatch |
Expression set containing qPCR data, read in by a ReadqPCR function and containing technical reps, denoted by |
... |
Extra arguments, detailed below |
calc |
use median, arithmetic or geometric mean for combining the values |
Takes exprs
of qPCR values containing technical replicates and combines them using a specified centrality measure.
qPCRBatch
with same number of samples, but with less features, since all technical replicates are replaced with a single value of their means.
James Perkins [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
path <- system.file("exData", package = "NormqPCR") qPCR.example.techReps <- file.path(path, "qPCR.techReps.txt") qPCRBatch.qPCR.techReps <- read.qPCR(qPCR.example.techReps) rownames(exprs(qPCRBatch.qPCR.techReps)) combinedTechReps <- combineTechReps(qPCRBatch.qPCR.techReps) rownames(exprs(combinedTechReps))
path <- system.file("exData", package = "NormqPCR") qPCR.example.techReps <- file.path(path, "qPCR.techReps.txt") qPCRBatch.qPCR.techReps <- read.qPCR(qPCR.example.techReps) rownames(exprs(qPCRBatch.qPCR.techReps)) combinedTechReps <- combineTechReps(qPCRBatch.qPCR.techReps) rownames(exprs(combinedTechReps))
Takes expression set of qPCR values containing technical replicates and combines them. In addition the appropriate standard deviation (SD) is computed.
combineTechRepsWithSD(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' combineTechRepsWithSD(qPCRBatch, calc="arith")
combineTechRepsWithSD(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' combineTechRepsWithSD(qPCRBatch, calc="arith")
qPCRBatch |
Expression set containing qPCR data, read in by a ReadqPCR function and containing technical reps, denoted by |
... |
Extra arguments, detailed below |
calc |
use median, arithmetic or geometric mean for combining the values |
Takes exprs
of qPCR values containing technical replicates and combines them using a specified centrality measure.
The arithmetic mean (calc="arith"
) is combined with the classical standard deviation.
In case of the geometric mean (calc="geom"
) the classical standard deviation of the log-values is exponentiated.
The median (calc="median"
) is calculated in connection with the MAD.
qPCRBatch
with same number of samples, but with less features, since all technical replicates are replaced with a single value of their means.
In addition the slot assayData
includes a matrix with SD values which can be accessed via se.exprs
.
Matthias Kohl [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
path <- system.file("exData", package = "NormqPCR") qPCR.example.techReps <- file.path(path, "qPCR.techReps.txt") qPCRBatch.qPCR.techReps <- read.qPCR(qPCR.example.techReps) rownames(exprs(qPCRBatch.qPCR.techReps)) combinedTechReps <- combineTechRepsWithSD(qPCRBatch.qPCR.techReps) rownames(exprs(combinedTechReps)) exprs(combinedTechReps) se.exprs(combinedTechReps)
path <- system.file("exData", package = "NormqPCR") qPCR.example.techReps <- file.path(path, "qPCR.techReps.txt") qPCRBatch.qPCR.techReps <- read.qPCR(qPCR.example.techReps) rownames(exprs(qPCRBatch.qPCR.techReps)) combinedTechReps <- combineTechRepsWithSD(qPCRBatch.qPCR.techReps) rownames(exprs(combinedTechReps)) exprs(combinedTechReps) se.exprs(combinedTechReps)
This function computes normalized relative quantities (NRQs) for a qPCRBatch
.
ComputeNRQs(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' ComputeNRQs(qPCRBatch, hkgs)
ComputeNRQs(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' ComputeNRQs(qPCRBatch, hkgs)
qPCRBatch |
an object of class |
hkgs |
Names of reference/housekeeping genes. |
... |
other parameters to be passed to downstream methods. |
Allows the user to normalized relative quantities as defined in Hellemanns et al. (2007).
Object of class "qPCRBatch"
.
Nor Izayu Abdul Rahman, Matthias Kohl [email protected]
Jan Hellemans, Geert Mortier, Anne De Paepe, Frank Speleman and Jo Vandesompele (2007). qBase relative quantification framework and software for management and automated analysis of real-time quantitative PCR data. Genome Biology, 8:R19
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
## Example data path <- system.file("exData", package = "ReadqPCR") qPCR.example <- file.path(path, "qPCR.example.txt") Cq.data <- read.qPCR(qPCR.example) ## combine technichal replicates Cq.data1 <- combineTechRepsWithSD(Cq.data) ## add efficiencies Effs <- file.path(path, "Efficiencies.txt") Cq.effs <- read.table(file = Effs, row.names = 1, header = TRUE) rownames(Cq.effs) <- featureNames(Cq.data1) effs(Cq.data1) <- as.matrix(Cq.effs[,"efficiency",drop = FALSE]) se.effs(Cq.data1) <- as.matrix(Cq.effs[,"SD.efficiency",drop = FALSE]) ## res <- ComputeNRQs(Cq.data1, hkgs = c("gene_az", "gene_gx")) ## NRQs exprs(res) ## SD of NRQs se.exprs(res)
## Example data path <- system.file("exData", package = "ReadqPCR") qPCR.example <- file.path(path, "qPCR.example.txt") Cq.data <- read.qPCR(qPCR.example) ## combine technichal replicates Cq.data1 <- combineTechRepsWithSD(Cq.data) ## add efficiencies Effs <- file.path(path, "Efficiencies.txt") Cq.effs <- read.table(file = Effs, row.names = 1, header = TRUE) rownames(Cq.effs) <- featureNames(Cq.data1) effs(Cq.data1) <- as.matrix(Cq.effs[,"efficiency",drop = FALSE]) se.effs(Cq.data1) <- as.matrix(Cq.effs[,"SD.efficiency",drop = FALSE]) ## res <- ComputeNRQs(Cq.data1, hkgs = c("gene_az", "gene_gx")) ## NRQs exprs(res) ## SD of NRQs se.exprs(res)
This function calculates Cq value and amplification efficiency for a CyclesSet
.
It is based on function pcrbatch
of package qpcR.
CqValues(object, ...) ## S4 method for signature 'CyclesSet' CqValues(object, Effmethod = "expfit", group = NULL, model = l5, check = "uni2", checkPAR = parKOD(), remove = "none", exclude = NULL, type = "cpD2", labels = NULL, norm = FALSE, baseline = NULL, basefac = 1, smooth = NULL, smoothPAR = list(span = 0.1), factor = 1, opt = FALSE, optPAR = list(sig.level = 0.05, crit = "ftest"), plot = FALSE, verbose = FALSE, ...)
CqValues(object, ...) ## S4 method for signature 'CyclesSet' CqValues(object, Effmethod = "expfit", group = NULL, model = l5, check = "uni2", checkPAR = parKOD(), remove = "none", exclude = NULL, type = "cpD2", labels = NULL, norm = FALSE, baseline = NULL, basefac = 1, smooth = NULL, smoothPAR = list(span = 0.1), factor = 1, opt = FALSE, optPAR = list(sig.level = 0.05, crit = "ftest"), plot = FALSE, verbose = FALSE, ...)
object |
an object of class |
Effmethod |
a character vector defining the methods for computing amplification efficiency. |
group |
a vector containing the grouping for possible replicates. |
model |
the model to be used for all runs. Default model is |
check |
the method for kinetic outlier detection in |
checkPAR |
parameters to be supplied to the |
remove |
indicates which runs to be removed. Either |
exclude |
indicates samples to be excluded from calculation, either "" for samples with missing column names or a regular expression defining columns (samples);
see 'Details' and 'Examples' in |
type |
the point on the amplification curve which is used for efficiency estimation; see |
labels |
a vector containing labels which define replicate groups. See more details in |
norm |
a logical value which determines whether the raw data should be normalized within [0, 1] before model fitting or not. |
baseline |
type of baseline subtraction. More details in |
basefac |
a factor when using averaged baseline cycles, such as |
smooth |
the curve smoothing method. See more details in |
smoothPAR |
parameters to be supplied to smoothing method in |
factor |
a multiplication factor for the fluorescence response values. |
opt |
a logical value which determines whether model selection should be applied to each model or not. |
optPAR |
parameters to be supplied for model selection in |
plot |
a logical value. If |
verbose |
a logical value. If |
... |
other parameters to be passed to downstream methods. |
Allows the user to compute Cq value and amplification efficiency. In addition, all values generated during the computations are saved. This function has four choices of methods for computing amplification efficiency values which are the methods provided by package qpcR.
More details on technical replication and normalization is given in the vignette NormqPCR
.
Object of class "qPCRBatch"
.
Nor Izayu Abdul Rahman, Matthias Kohl [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
pcrbatch
, CyclesSet-class
, qPCRBatch-class
## Read in the raw qPCR data from file "LC480_Example.txt" path <- system.file("exData", package = "ReadqPCR") LC480.example <- file.path(path, "LC480_Example.txt") cycData <- read.LC480(file = LC480.example) ## Read in the sample information data from file "LC480_Example_SampleInfo.txt". LC480.SamInfo <- file.path(path, "LC480_Example_SampleInfo.txt") samInfo <- read.LC480SampleInfo(LC480.SamInfo) ## Merge information cycData1 <- merge(cycData, samInfo) ## Compute Cq values ## 1) use sigmoidal model res1 <- CqValues(cycData1, Effmethod = "sigfit") res1 effs(res1) se.effs(res1) ## 2) fit exponential model (default) res2 <- CqValues(cycData1, Effmethod = "expfit") res2 effs(res2) se.effs(res2) ## 3) use window of linearity res3 <- CqValues(cycData1, Effmethod = "sliwin") res3 effs(res3) se.effs(res3) ## 4) linear regression of efficiency res4 <- CqValues(cycData1, Effmethod = "LRE") res4 effs(res4) se.effs(res4)
## Read in the raw qPCR data from file "LC480_Example.txt" path <- system.file("exData", package = "ReadqPCR") LC480.example <- file.path(path, "LC480_Example.txt") cycData <- read.LC480(file = LC480.example) ## Read in the sample information data from file "LC480_Example_SampleInfo.txt". LC480.SamInfo <- file.path(path, "LC480_Example_SampleInfo.txt") samInfo <- read.LC480SampleInfo(LC480.SamInfo) ## Merge information cycData1 <- merge(cycData, samInfo) ## Compute Cq values ## 1) use sigmoidal model res1 <- CqValues(cycData1, Effmethod = "sigfit") res1 effs(res1) se.effs(res1) ## 2) fit exponential model (default) res2 <- CqValues(cycData1, Effmethod = "expfit") res2 effs(res2) se.effs(res2) ## 3) use window of linearity res3 <- CqValues(cycData1, Effmethod = "sliwin") res3 effs(res3) se.effs(res3) ## 4) linear regression of efficiency res4 <- CqValues(cycData1, Effmethod = "LRE") res4 effs(res4) se.effs(res4)
Normalise qPCR eset using a given housekeeping gene as control, then perform differential expression analysis using the delta delta Ct method
deltaCt(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' deltaCt(qPCRBatch, hkgs, combineHkgs=FALSE, calc="arith") deltaCq(qPCRBatch, hkgs, combineHkgs=FALSE, calc="arith")
deltaCt(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' deltaCt(qPCRBatch, hkgs, combineHkgs=FALSE, calc="arith") deltaCq(qPCRBatch, hkgs, combineHkgs=FALSE, calc="arith")
qPCRBatch |
qPCR-specific expression set, containing qPCR data. |
... |
Extra arguments, detailed below |
hkgs |
String containing the name of the name of the housekeeping gene which will be used to normalise the rest of the genes. |
combineHkgs |
Logical - if TRUE, then as long as more than one housekeeper given for argument hkgs, it will combine the housekeepers by finding the geometric mean. Housekeepers can be found using geNorm or NormFinder algorithms. |
calc |
use arithmetic or geometric mean. |
Takes expression set of qPCR values and normalises them using a housekeeping gene. Returns a qPCRBatch with exprs set of the same dimensions but with the given hkg value subtracted.
qPCRBatch with exprs set of the same dimensions but with the given hkg value subtracted.
James Perkins [email protected]
Kenneth Livak, Thomase Schmittgen (2001). Analysis of Relative Gene Expression Data Using Real-Time Quantitative PCR and the 2^DDCt Method. Methods 25, 402-408, 2001 http://www.ncbi.nlm.nih.gov/pubmed/11846609
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
selectHKs
, deltaDeltaCq
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) hkgs<-"Actb-Rn00667869_m1" qPCRBatch.norm <- deltaCq(qPCRBatch = qPCRBatch.taqman, hkgs = hkgs, calc="arith") head(exprs(qPCRBatch.norm))
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) hkgs<-"Actb-Rn00667869_m1" qPCRBatch.norm <- deltaCq(qPCRBatch = qPCRBatch.taqman, hkgs = hkgs, calc="arith") head(exprs(qPCRBatch.norm))
Normalise qPCRBatch RT-qPCR data using housekeeping genes as control, then perform differential expression analysis using the delta delta Cq method.
deltaDeltaCt(qPCRBatch,...) ## S4 method for signature 'qPCRBatch' deltaDeltaCt(qPCRBatch, maxNACase=0, maxNAControl=0, hkgs, contrastM, case, control, paired=TRUE, hkgCalc="arith", statCalc="arith") deltaDeltaCq(qPCRBatch, maxNACase=0, maxNAControl=0, hkgs, contrastM, case, control, paired=TRUE, hkgCalc="arith", statCalc="arith")
deltaDeltaCt(qPCRBatch,...) ## S4 method for signature 'qPCRBatch' deltaDeltaCt(qPCRBatch, maxNACase=0, maxNAControl=0, hkgs, contrastM, case, control, paired=TRUE, hkgCalc="arith", statCalc="arith") deltaDeltaCq(qPCRBatch, maxNACase=0, maxNAControl=0, hkgs, contrastM, case, control, paired=TRUE, hkgCalc="arith", statCalc="arith")
qPCRBatch |
qPCR-specific expression set, containing qPCR data. |
... |
Extra arguments, detailed below |
maxNACase |
Maximum number of NA values allowed before a detector's reading is discarded for samples designated as case. |
maxNAControl |
Maximum number of NA values allowed before a detector's reading is discarded for samples designated as control. |
hkgs |
String containing the name of th name of the housekeeping gene which will be used to normalise the rest of the genes. |
contrastM |
A binary matrix which designates case and control samples. |
case |
The name of the column in contrastM that corresponds to the case samples. |
control |
The name of the column in contrastM that corresponds to the control samples. |
paired |
Logical - if TRUE the detectors and housekeepers in the same sample will be paired for calculating standard deviation, effectively meaning we will be calculating standard deviation of the differences. If FALSE, there will be no pairing, and standard deviation will be pooled between the detector and housekeepers. |
hkgCalc |
String - either "arith" or "geom", details how the different housekeeper genes should be combined - either by using the arithmetic or geometric mean. |
statCalc |
String - either "arith" or "geom", details how genes should be combined - either by using the arithmetic or geometric mean. |
Takes expression set of qPCR values and normalises them using different housekeeping genes. Returns seperate sets of values for each housekeeping gene given.
matrix with columns containing the detector ids, 2^delta Cq values for the sample of interest and the callibrator sample, alongside their respective standard deviations, the 2^delta delta Cq values and the minimum and maximum values (ie the values that are 1 sd away )
James Perkins [email protected]
Kenneth Livak, Thomase Schmittgen (2001). Analysis of Relative Gene Expression Data Using Real-Time Quantitative PCR and the 2^DDCt Method. Methods 25, 402-408, 2001 http://www.ncbi.nlm.nih.gov/pubmed/11846609
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
selectHKs
, deltaCq
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) hkg <- "Actb-Rn00667869_m1" contM <- cbind(c(0,0,1,1,0,0,1,1),c(1,1,0,0,1,1,0,0)) colnames(contM) <- c("interestingPhenotype","wildTypePhenotype") rownames(contM) <- sampleNames(qPCRBatch.taqman) ddCq.taqman <- deltaDeltaCq(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkg=hkg, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", statCalc="geom", hkgCalc="arith") head(ddCq.taqman)
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) hkg <- "Actb-Rn00667869_m1" contM <- cbind(c(0,0,1,1,0,0,1,1),c(1,1,0,0,1,1,0,0)) colnames(contM) <- c("interestingPhenotype","wildTypePhenotype") rownames(contM) <- sampleNames(qPCRBatch.taqman) ddCq.taqman <- deltaDeltaCq(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkg=hkg, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", statCalc="geom", hkgCalc="arith") head(ddCq.taqman)
This data set was used in Vandesompele et al (2002) to demonstrate normalization of real-time quantitative RT-PCR data by geometric averaging of housekeeping genes.
data(geNorm)
data(geNorm)
A qPCRBatch object which contains an expression matrix with 85 observations on the following 10 variables which stand for expression data of ten potential reference/housekeeping genes
ACTB
actin, beta
B2M
beta-2-microglobulin
GAPD
glyceraldehyde-3-phosphate dehydrogenase
HMBS
hydroxymethylbilane synthase
HPRT1
hypoxanthine phosphoribosyltransferase 1
RPL13A
ribosomal protein L13a
SDHA
succinate dehydrogenase complex subunit A
TBP
TATA box binding protein
UBC
ubiquitin C
YWHAZ
tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein, zeta polypeptide
The row names of this data set indicate the various human tissues which were investigated.
9 normal bone-marrow samples
9 normal human tissues from pooled organs (heart, brain, fetal brain, lung, trachea, kidney, mammary gland, small intestine and uterus)
20 short-term cultured normal fibroblast samples from different individuals
13 normal leukocyte samples
34 neuroblastoma cell lines (independently prepared in different labs from different patients)
The data set was obtained from http://genomebiology.com/content/supplementary/gb-2002-3-7-research0034-s1.txt
Jo Vandesompele, Katleen De Preter, Filip Pattyn et al. (2002). Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology 2002. 3(7):research0034.1-0034.11. http://genomebiology.com/2002/3/7/research/0034/
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
data(geNorm) str(exprs(geNorm.qPCRBatch)) sampleNames(geNorm.qPCRBatch)
data(geNorm) str(exprs(geNorm.qPCRBatch)) sampleNames(geNorm.qPCRBatch)
Computation of the geometric mean.
geomMean(x, na.rm = TRUE)
geomMean(x, na.rm = TRUE)
x |
numeric vector of non-negative Reals |
na.rm |
a logical value indicating whether |
The computation of the geometric mean is done via prod(x)^(1/length(x))
.
geometric mean
A first version of this function appeared in package SLqPCR.
Matthias Kohl [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
x <- rlnorm(100) geomMean(x)
x <- rlnorm(100) geomMean(x)
Make all Cq values for a given detector NA when the number of NAs for that detector is above a given threshold
makeAllNAs(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' makeAllNAs(qPCRBatch, contrastM, sampleMaxM)
makeAllNAs(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' makeAllNAs(qPCRBatch, contrastM, sampleMaxM)
qPCRBatch |
Expression set containing qPCR data. |
... |
Extra arguments, detailed below |
contrastM |
Contrast Matrix like that used in |
sampleMaxM |
Sample Max Matrix. Columns represent the different sample types. There is one value per column, which represents the max number of NAs allowed for that sample type. |
Make all NAs when number of NAs above a given threshold
qPCRBatch
object with a new exprs slot, everything else equal
James Perkins [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
# read in the data path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) exprs(qPCRBatch.taqman)["Ccl20.Rn00570287_m1",] # values before # make contrastM a <- c(0,0,1,1,0,0,1,1) # one for each sample type, with 1 representing b <- c(1,1,0,0,1,1,0,0) # position of sample type in the samplenames vector contM <- cbind(a,b) colnames(contM) <- c("case","control") # then give the names of each sample type rownames(contM) <- sampleNames(qPCRBatch.taqman) # and the rows of the matrix contM # make sampleMaxM sMaxM <- t(as.matrix(c(3,3))) # now make the sample max matrix colnames(sMaxM) <- c("case","control") # make sure these line up with samples sMaxM # function qPCRBatch.taqman.replaced <- makeAllNAs(qPCRBatch.taqman, contM, sMaxM) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
# read in the data path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) exprs(qPCRBatch.taqman)["Ccl20.Rn00570287_m1",] # values before # make contrastM a <- c(0,0,1,1,0,0,1,1) # one for each sample type, with 1 representing b <- c(1,1,0,0,1,1,0,0) # position of sample type in the samplenames vector contM <- cbind(a,b) colnames(contM) <- c("case","control") # then give the names of each sample type rownames(contM) <- sampleNames(qPCRBatch.taqman) # and the rows of the matrix contM # make sampleMaxM sMaxM <- t(as.matrix(c(3,3))) # now make the sample max matrix colnames(sMaxM) <- c("case","control") # make sure these line up with samples sMaxM # function qPCRBatch.taqman.replaced <- makeAllNAs(qPCRBatch.taqman, contM, sMaxM) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
Make all Cq values for a given detector NA when the number of NAs for that detector is above a given threshold
makeAllNewVal(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' makeAllNewVal(qPCRBatch, contrastM, sampleMaxM, newVal)
makeAllNewVal(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' makeAllNewVal(qPCRBatch, contrastM, sampleMaxM, newVal)
qPCRBatch |
Expression set containing qPCR data. |
... |
Extra arguments, detailed below |
contrastM |
Contrast Matrix like that used in |
sampleMaxM |
Sample Max Matrix. Columns represent the different sample types. There is one value per column, which represents the max number of NAs allowed for that sample type. |
newVal |
New value to give the values in the group where the NAs are above the threshold. |
Make all a given value when number of NAs above a given threshold, with different thresholds for the different sample classes, using sMaxM and contM to provide this information, as detailed below.
qPCRBatch
object with a new exprs slot, everything else equal
James Perkins [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
# read in the data path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) exprs(qPCRBatch.taqman)["Ccl20.Rn00570287_m1",] # values before # make contrastM a <- c(0,0,1,1,0,0,1,1) # one for each sample type, with 1 representing b <- c(1,1,0,0,1,1,0,0) # position of sample type in the samplenames vector contM <- cbind(a,b) colnames(contM) <- c("case","control") # then give the names of each sample type rownames(contM) <- sampleNames(qPCRBatch.taqman) # and the rows of the matrix contM # make sampleMaxM sMaxM <- t(as.matrix(c(3,3))) # now make the sample max matrix colnames(sMaxM) <- c("case","control") # make sure these line up with samples sMaxM # function qPCRBatch.taqman.replaced <- makeAllNewVal(qPCRBatch.taqman, contM, sMaxM) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
# read in the data path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) exprs(qPCRBatch.taqman)["Ccl20.Rn00570287_m1",] # values before # make contrastM a <- c(0,0,1,1,0,0,1,1) # one for each sample type, with 1 representing b <- c(1,1,0,0,1,1,0,0) # position of sample type in the samplenames vector contM <- cbind(a,b) colnames(contM) <- c("case","control") # then give the names of each sample type rownames(contM) <- sampleNames(qPCRBatch.taqman) # and the rows of the matrix contM # make sampleMaxM sMaxM <- t(as.matrix(c(3,3))) # now make the sample max matrix colnames(sMaxM) <- c("case","control") # make sure these line up with samples sMaxM # function qPCRBatch.taqman.replaced <- makeAllNewVal(qPCRBatch.taqman, contM, sMaxM) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
Replace Cq values above a given threshold with a new value
replaceAboveCutOff(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' replaceAboveCutOff(qPCRBatch, newVal=NA, cutOff=38)
replaceAboveCutOff(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' replaceAboveCutOff(qPCRBatch, newVal=NA, cutOff=38)
qPCRBatch |
Expression set containing qPCR data. |
... |
Extra arguments, detailed below |
newVal |
The new value with which to replace the values above the cutoff |
cutOff |
the minimal threshold above which the values will be replaced |
Replaces values in the exprs slot of the qPCRBatch
object that are above a threshold value with a new number
qPCRBatch
object with a new exprs slot
James Perkins [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) exprs(qPCRBatch.taqman)["Ccl20.Rn00570287_m1",] qPCRBatch.taqman.replaced <- replaceAboveCutOff(qPCRBatch.taqman, newVal = NA, cutOff = 35) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) exprs(qPCRBatch.taqman)["Ccl20.Rn00570287_m1",] qPCRBatch.taqman.replaced <- replaceAboveCutOff(qPCRBatch.taqman, newVal = NA, cutOff = 35) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
Replace NAs with a given value
replaceNAs(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' replaceNAs(qPCRBatch, newNA)
replaceNAs(qPCRBatch, ...) ## S4 method for signature 'qPCRBatch' replaceNAs(qPCRBatch, newNA)
qPCRBatch |
Expression set containing qPCR data. |
... |
Extra arguments, detailed below |
newNA |
The new value to replace the NAs with |
Replaces NA values in the exprs slot of the qPCRBatch object with a given number
qPCRBatch object with a new exprs slot
James Perkins [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) qPCRBatch.taqman.replaced <- replaceNAs(qPCRBatch.taqman, newNA = 40) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
path <- system.file("exData", package = "NormqPCR") taqman.example <- file.path(path, "example.txt") qPCRBatch.taqman <- read.taqman(taqman.example) qPCRBatch.taqman.replaced <- replaceNAs(qPCRBatch.taqman, newNA = 40) exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
This function can be used to determine a set of reference/housekeeping (HK) genes for gene expression experiments
selectHKs(qPCRBatch, ...) ## S4 method for signature 'matrix' selectHKs(qPCRBatch, group, method = "geNorm", minNrHKs = 2, log = TRUE, Symbols, trace = TRUE, na.rm = TRUE) ## S4 method for signature 'qPCRBatch' selectHKs(qPCRBatch, group, method = "geNorm", minNrHKs = 2, log = TRUE, Symbols, trace = TRUE, na.rm = TRUE)
selectHKs(qPCRBatch, ...) ## S4 method for signature 'matrix' selectHKs(qPCRBatch, group, method = "geNorm", minNrHKs = 2, log = TRUE, Symbols, trace = TRUE, na.rm = TRUE) ## S4 method for signature 'qPCRBatch' selectHKs(qPCRBatch, group, method = "geNorm", minNrHKs = 2, log = TRUE, Symbols, trace = TRUE, na.rm = TRUE)
qPCRBatch |
matrix or qPCRBatch, containing the data (expression matrix) in the exprs slot |
... |
Extra arguments, detailed below |
group |
optional factor not used by all methods, hence may be missing |
method |
method to compute most stable genes |
minNrHKs |
minimum number of HK genes that should be considered |
log |
logical: is data on log-scale |
Symbols |
gene symbols |
trace |
logical, print additional information |
na.rm |
a logical value indicating whether |
This function can be used to determine a set of reference/housekeeping (HK) genes
for gene expression experiments. The default method "geNorm"
was proposed by Vandesompele et al. (2002).
Currently, the geNorm method by Vandesompele et al. (2002) and the NormFinder method of Andersen et al. (2004) are implemented.
Vandesompele et al. (2002) propose a cut-off value of 0.15 for the pairwise variation. Below this value the inclusion of an additional housekeeping gene is not required.
If method = "geNorm"
a list with the following components is
returned
ranking |
ranking of genes from best to worst where the two most stable genes cannot be ranked |
variation |
pairwise variation during stepwise selection |
meanM |
average expression stability M |
If method = "NormFinder"
a list with the following components is
returned
ranking |
ranking of genes from best to worst where the two most stable genes cannot be ranked |
rho |
stability measure rho of Andersen et al. (2004) |
Matthias Kohl [email protected]
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
Jo Vandesompele, Katleen De Preter, Filip Pattyn et al. (2002). Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology 2002. 3(7):research0034.1-0034.11. http://genomebiology.com/2002/3/7/research/0034/
Claus Lindbjerg Andersen, Jens Ledet Jensen and Torben Falck Orntoft (2004). Normalization of Real-Time Quantitative Reverse Transcription-PCR Data: A Model-Based Variance Estimation Approach to Identify Genes Suited for Normalization, Applied to Bladder and Colon Cancer Data Sets. CANCER RESEARCH 64, 5245-5250, August 1, 2004. http://cancerres.aacrjournals.org/cgi/content/full/64/15/5245
data(geNorm) tissue <- as.factor(c(rep("BM", 9), rep("FIB", 20), rep("LEU", 13), rep("NB", 34), rep("POOL", 9))) res.BM <- selectHKs(geNorm.qPCRBatch[,tissue == "BM"], method = "geNorm", Symbols = featureNames(geNorm.qPCRBatch), minNrHK = 2, log = FALSE)
data(geNorm) tissue <- as.factor(c(rep("BM", 9), rep("FIB", 20), rep("LEU", 13), rep("NB", 34), rep("POOL", 9))) res.BM <- selectHKs(geNorm.qPCRBatch[,tissue == "BM"], method = "geNorm", Symbols = featureNames(geNorm.qPCRBatch), minNrHK = 2, log = FALSE)
Computation of the gene expression stability value M for real-time quantitativ RT-PCR data. For more details we refer to Vandesompele et al. (2002).
stabMeasureM(x, log = TRUE, na.rm = TRUE)
stabMeasureM(x, log = TRUE, na.rm = TRUE)
x |
matrix or data.frame containing real-time quantitative RT-PCR data |
log |
logical: is data on log-scale |
na.rm |
a logical value indicating whether |
The gene expression stability value M is defined as the average pairwise normalization factor; i.e., one needs to specify data from at least two genes. For more details see Vandesompele et al. (2002). Note this dispatches on a transposed expression matrix, not a qPCRBatch object since it is only called from within the selectHKs method.
numeric vector with gene expression stability values
Matthias Kohl [email protected]
Jo Vandesompele, Katleen De Preter, Filip Pattyn et al. (2002). Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biology 2002. 3(7):research0034.1-0034.11. http://genomebiology.com/2002/3/7/research/0034/
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
selectHKs
data(geNorm) tissue <- as.factor(c(rep("BM", 9), rep("FIB", 20), rep("LEU", 13), rep("NB", 34), rep("POOL", 9))) res.BM <- selectHKs(geNorm.qPCRBatch[,tissue == "BM"], method = "geNorm", Symbols = featureNames(geNorm.qPCRBatch), minNrHK = 2, log = FALSE)
data(geNorm) tissue <- as.factor(c(rep("BM", 9), rep("FIB", 20), rep("LEU", 13), rep("NB", 34), rep("POOL", 9))) res.BM <- selectHKs(geNorm.qPCRBatch[,tissue == "BM"], method = "geNorm", Symbols = featureNames(geNorm.qPCRBatch), minNrHK = 2, log = FALSE)
Computation of the gene expression stability value rho for real-time quantitativ RT-PCR data. For more details we refer to Andersen et al. (2004).
stabMeasureRho(x,...) ## S4 method for signature 'x' stabMeasureRho(x, group, log = TRUE, na.rm = TRUE, returnAll = FALSE)
stabMeasureRho(x,...) ## S4 method for signature 'x' stabMeasureRho(x, group, log = TRUE, na.rm = TRUE, returnAll = FALSE)
x |
matrix containing real-time quantitative RT-PCR data, or qPCRBatch object |
... |
Extra arguments, detailed below |
group |
grouping factor, either a factor vector or a phenoData column called "Group" |
log |
logical: is data on log-scale |
na.rm |
a logical value indicating whether |
returnAll |
logical, return additional information. |
The gene expression stability value rho is computed. For more details see Andersen et al. (2004).
numeric vector with gene expression stability values
If returnAll == TRUE
a list with the following components is returned
rho |
stability measure rho of Andersen et al. (2004) |
d |
used by |
v |
used by |
Matthias Kohl [email protected]
Claus Lindbjerg Andersen, Jens Ledet Jensen and Torben Falck Orntoft (2004). Normalization of Real-Time Quantitative Reverse Transcription-PCR Data: A Model-Based Variance Estimation Approach to Identify Genes Suited for Normalization, Applied to Bladder and Colon Cancer Data Sets. CANCER RESEARCH 64, 5245-5250, August 1, 2004. http://cancerres.aacrjournals.org/cgi/content/full/64/15/5245
Perkins, JR, Dawes, JM, McMahon, SB, Bennett, DL, Orengo, C, Kohl, M (2012). ReadqPCR and NormqPCR: R packages for the reading, quality checking and normalisation of RT-qPCR quantification cycle (Cq) data. BMC Genomics, 13, 1:296.
selectHKs
data(Colon) Class <- pData(Colon)[,"Classification"] res.Colon <- stabMeasureRho(Colon, group = Class, log = FALSE) sort(res.Colon) # cf. Table 3 in Andersen et al (2004) data(Bladder) Grade <- pData(Bladder)[,"Grade"] res.Bladder <- stabMeasureRho(Bladder, group = Grade, log = FALSE) sort(res.Bladder) # cf. Table 3 in Andersen et al (2004)
data(Colon) Class <- pData(Colon)[,"Classification"] res.Colon <- stabMeasureRho(Colon, group = Class, log = FALSE) sort(res.Colon) # cf. Table 3 in Andersen et al (2004) data(Bladder) Grade <- pData(Bladder)[,"Grade"] res.Bladder <- stabMeasureRho(Bladder, group = Grade, log = FALSE) sort(res.Bladder) # cf. Table 3 in Andersen et al (2004)