Title: | The ddCt Algorithm for the Analysis of Quantitative Real-Time PCR (qRT-PCR) |
---|---|
Description: | The Delta-Delta-Ct (ddCt) Algorithm is an approximation method to determine relative gene expression with quantitative real-time PCR (qRT-PCR) experiments. Compared to other approaches, it requires no standard curve for each primer-target pair, therefore reducing the working load and yet returning accurate enough results as long as the assumptions of the amplification efficiency hold. The ddCt package implements a pipeline to collect, analyse and visualize qRT-PCR results, for example those from TaqMan SDM software, mainly using the ddCt method. The pipeline can be either invoked by a script in command-line or through the API consisting of S4-Classes, methods and functions. |
Authors: | Jitao David Zhang, Rudolf Biczok, and Markus Ruschhaupt |
Maintainer: | Jitao David Zhang <[email protected]> |
License: | LGPL-3 |
Version: | 1.63.0 |
Built: | 2024-11-29 05:39:18 UTC |
Source: | https://github.com/bioc/ddCt |
Barplot with error bars.
barploterrbar(y, yl, yh, barcol="orange", errcol="black", horiz=FALSE, w=0.2,theCut=NULL,columnForDiffBars=TRUE,cex.axis = par("cex.axis"),zeroForNA=TRUE,legend=FALSE,groups = NULL, order=FALSE, ...)
barploterrbar(y, yl, yh, barcol="orange", errcol="black", horiz=FALSE, w=0.2,theCut=NULL,columnForDiffBars=TRUE,cex.axis = par("cex.axis"),zeroForNA=TRUE,legend=FALSE,groups = NULL, order=FALSE, ...)
y |
Numeric vector. |
yl |
Numeric vector of same length as y. |
yh |
Numeric vector of same length as y. |
barcol |
Color of the bars. |
errcol |
Color of the error bars. |
horiz |
Logical. As in |
w |
Size of the error bar ticks. |
theCut |
The cut value |
columnForDiffBars |
Whether the matrix should be transposed (by default the rows are for diff bars) |
zeroForNA |
Draw 0 instead of NA |
cex.axis |
Axis font cex |
legend |
Sould a legend be plotted ? |
groups |
a factor - if specified the bars are collored according to the group they belong to |
order |
plot sample values in descending order |
... |
Further arguments that get passed on to
|
The function calls barplot
with
y
and decorates it with error bars according to yl
and
yh
.
The function is called for its side effect, producing a plot.
Markus Ruschhaupt, Florian Hahne
y <- matrix(runif(80), ncol=5) ym <- apply(y, 2, mean) dy <- apply(y, 2, sd)*2/sqrt(nrow(y)) barploterrbar(ym, ym-dy, ym+dy, barcol="#0000c0", errcol="orange")
y <- matrix(runif(80), ncol=5) ym <- apply(y, 2, mean) dy <- apply(y, 2, sd)*2/sqrt(nrow(y)) barploterrbar(ym, ym-dy, ym+dy, barcol="#0000c0", errcol="orange")
absolute quantification for Taqman data
ddCtAbsolute(raw.table, addData, type = "mean", ADD = -30.234, DIV = -1.6268, sampleInformation = NULL, toZero = FALSE, filename = "warning.output.txt")
ddCtAbsolute(raw.table, addData, type = "mean", ADD = -30.234, DIV = -1.6268, sampleInformation = NULL, toZero = FALSE, filename = "warning.output.txt")
raw.table |
data frame. It must contain columns with the following names:'Ct','Sample','Detector','Platename'. The column 'Ct' must contain numeric values. |
addData |
add data |
type |
character of length 1. 'mean' or 'median'- which method should be used for the aggregation of the repicates |
ADD |
Add constant |
DIV |
Div constant |
sampleInformation |
if specified it must be an object of class
|
toZero |
boolean - if there is only one replication should the error be treated as zero ? (only if 'type' is mean) |
filename |
character of length 1. The name of the file the warnings should be stored in. |
A an object of class
eSet
. The assayData
has the following components: exprs, error, Ct, Ct.error, Difference, number\_NA, number, Plate.
Markus Ruschhaupt mailto:[email protected]
~put references to the literature/web site here ~
This class is a subclass of ExpressionSet
and
represents objects which are produced
by the ddCt algorithm in the ddCtExpression
method
Class ExpressionSet
, directly.
Class eSet
, by class "ExpressionSet", distance 2.
Class VersionedBiobase
, by class "ExpressionSet", distance 3.
Class Versioned
, by class "ExpressionSet", distance 4.
signature(object = "ddCtExpression")
: returns the
Ct value of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the error number of the Ct value of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns the
dCt value of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the error number of the dCt value of this ddCtExpressionobject
signature(object = "ddCtExpression")
:returns the
ddCt value of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the error number of the ddCt value of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the levels in this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the error number of the levens in this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the Ct number of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns
the NA number of this ddCtExpressionobject
signature(object = "ddCtExpression")
: returns a
data frame which represents this expression object
signature(object = "ddCtExpression", file =
"character")
: writes ddCtExpression object into a file
Rudolf Biczok mailto:[email protected]
SDMFrame
: reader for SDM files
elist
, elistWrite
: utility functions for
ddCtExpression objects
ddCtExpression
: the method which invokes the ddCt algorithm
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## call ddCtExpression method to get a ddCt calculated expression result <- ddCtExpression(sampdat, calibrationSample="Sample1", housekeepingGenes=c("Gene1","Gene2")) ## use getter methods ddCt(result) ddCtErr(result)
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## call ddCtExpression method to get a ddCt calculated expression result <- ddCtExpression(sampdat, calibrationSample="Sample1", housekeepingGenes=c("Gene1","Gene2")) ## use getter methods ddCt(result) ddCtErr(result)
Apply the ddCt algorithm for a given data set
object |
SDMFrame Data object which holds a data set containing columns with the following names: 'Ct','Sample','Detector','Platename'. The column 'Ct' must contain numeric values. |
algorithm |
character. Name of the calibration samples. |
warningStream |
character of length 1. The name of the file the warnings should be stored in. |
calibrationSample |
character. Name of the calibration samples. |
housekeepingGenes |
character. Name of the housekeeping genes. |
type |
character of length 1. 'mean' or 'median'- which method should be used for the aggregation of the repicates |
sampleInformation |
if specified it must be an object of class |
toZero |
boolean - if there is only one replication should the error be treated as zero ? (only if 'type' is mean) |
efficiencies |
n.V. |
efficiencies.error |
n.V. |
A an object of class
ddCtExpression
.
ddCtExpression(object, warningStream = "warning.output.txt", algorithm="ddCt" calibrationSample, housekeepingGenes, type="mean", sampleInformation=NULL, toZero=TRUE, efficiencies = NULL, efficiencies.error = NULL)
An object of InputFrame
,
constructed with the method InputFrame
Rudolf Biczok mailto:[email protected]
Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta -Delta C(T)) Method. KJ Livak and TD Schmittgen, Methods, Vol. 25, No. 4. (December 2001), pp. 402-408
InputFrame
: reader for SDM files
ddCtExpression
: representation for ddCt
calculated expressions
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## call ddCtExpression method from class SDMFrame ## to get a ddCt calculated expression result <- ddCtExpression(sampdat, calibrationSample="Sample1", housekeepingGenes=c("Gene1","Gene2")) result
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## call ddCtExpression method from class SDMFrame ## to get a ddCt calculated expression result <- ddCtExpression(sampdat, calibrationSample="Sample1", housekeepingGenes=c("Gene1","Gene2")) result
ddCtExpression
object contains a list of matrices as
the results of ddCt
method. elist
combines these
lists into one data frame, and elistWrite
writes the data frame
into file.
summary
is a wrapper for the elist
method
elist(object,...) summary(object,...) elistWrite(object,file,...)
elist(object,...) summary(object,...) elistWrite(object,file,...)
object |
an |
file |
output file. |
... |
additional arguments passed to write.table. |
elist
is a wrapper to as(object, "data.frame")
function.
A data frame or output file.
Jitao David Zhang [email protected]
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## call ddCtExpression method from class SDMFrame ## to get a ddCt calculated expression result <- ddCtExpression(sampdat, calibrationSample="Sample1", housekeepingGenes=c("Gene1","Gene2")) ## call elist elistResult <- elist(result) elistResult
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## call ddCtExpression method from class SDMFrame ## to get a ddCt calculated expression result <- ddCtExpression(sampdat, calibrationSample="Sample1", housekeepingGenes=c("Gene1","Gene2")) ## call elist elistResult <- elist(result) elistResult
Draw barchart (with error-bars) of relative expression level
represented in ddCtExpression
object. The
barchart is implemented as grid plot by lattice
package, where
each panel represents one sample and the relative expression values of detectors (as
well as their standard errors) are depicted as bars.
Detectors which are not determined are marked by grey ND.
Two types of figures are supported: either condition on samples (by="Sample") or on detectors (by="Detector").
An object of ddCtExpression
,
constructed with the method ddCtExpression
Parameter object for errBarchart
Objects can be created by calls of the form
new("errBarchartParameter", ...)
. So far the object is only
internally used, but in the near future it will be exported.
exprsUndeterminedLabel
:Object of class
"character"
, specifying the text label when the expression
level is ‘Undetermined’
signature(object =
"errBarchartParameter")
: getting the text label when the
expression level is ‘Undetermined’
signature(object = "errBarchartParameter")
:
print method
So far it is only internallly used
Jitao David Zhang <[email protected]>
## Internally used ## param <- new("errBarchartParameter") ## exprsUndeterminedLabel(param)
## Internally used ## param <- new("errBarchartParameter") ## exprsUndeterminedLabel(param)
getDir
creates a directory in case it does not exist and
returns the directory name.
getDir(dir, ...)
getDir(dir, ...)
dir |
Directory name |
... |
Other parameters passed to |
Auxillary functions
getDir
returns the directory name
Jitao David Zhang <[email protected]>
getDir(tempdir())
getDir(tempdir())
Generally an InputFrame is built from a ReaderClass
(e.g. InputReader
), or a data.frame. See the example below
for building an object from a valid data.frame.
InputFrame(object)
InputFrame(object)
object |
A data.frame with three columns: Sample, Detector, and Ct |
A object of class InputFrame
Jitao David Zhang mailto: [email protected]
testDf <- data.frame(Sample=rep(paste("Sample", 1:3), each=2), Detector=rep(paste("Gene", 1:2), 3), Ct=30+rnorm(6)) testInputFrame <- InputFrame(testDf)
testDf <- data.frame(Sample=rep(paste("Sample", 1:3), each=2), Detector=rep(paste("Gene", 1:2), 3), Ct=30+rnorm(6)) testInputFrame <- InputFrame(testDf)
The class InputFrame provides core functionalities to read gene and sample information from SDM files and calculate them with a ddCt algorithm.
The function InputFrame
reads the data given in the colums
'Detector','Sample' and 'Ct' of the specified SDM output files and
stores them as a data.frame. An additional column including the
respective filename is added.
coreData
:Object of class "data.frame"
: Holds
all the required data extracted from the SDM file
files
:Object of class "character"
contains the
source SDM files
signature(x = "InputFrame")
: primitive
accessors. Returns an object of InputFrame-class
with the subset data.
signature(x = "InputFrame")
: returns the column
names in this SDM object
signature(object = "InputFrame")
: runs
a ddCt algorithm with this SDM object and returns a object of
class ddCtExpression
signature(object="InputFrame")
: returns the
source SDM file names.
signature(object = "InputFrame")
: returns
the detector names in this SDM object
signature(object = "InputFrame", value =
"character")
: replaces the detector names in this SDM object
signature(object = "InputFrame")
: returns
the sample names in this SDM object
signature(object = "InputFrame", value =
"character")
: replaces the sample names in this SDM object
signature(object = "InputFrame")
:
returns a vector of unique detector names in this SDM object
signature(object = "InputFrame",
target = "missing", value = "character")
: replaces all detector
names given by the 'names' attribute in 'value' with new detector
names
signature(object = "InputFrame",
target = "character", value = "character")
: replaces all detector
names given by 'target' with new detector names
signature(object = "InputFrame",
target = "missing", value = "character")
: replaces all sample
names given by the 'names' attribute in 'value' with new sample
names
signature(object = "InputFrame",
target = "character", value = "character")
: replaces all sample
names given by 'target' with new sample names
signature(object = "InputFrame")
:
returns a vector of unique sample names in this SDM object
signature(object = "InputFrame", sample="character")
:
removes the sample(s) specified from the InputFrame object
signature(object = "InputFrame",
target="character", value="character")
:
replace the detectors equal to the target with the value. Both
target
and value
can be vectors of the same length,
then the replace takes place iteratively.
signature(object = "InputFrame",
target="character", value="character")
:
replace the samples equal to the target with the value. Both
target
and value
can be vectors of the same length,
then the replace takes place iteratively.
signature(object="InputFrame")
: pretty print of the
InputFrame instance.
signature(object="InputFrame",
threshold="numeric")
: Right censoring the Ct value, which
targets the data points above a certain value
(threshold
). High Ct values (higher than 40 or 45 by the rule
of thumb) are often not accurate and may indicate too weak
expression. The function performs the right censoring on the data and
set the value above the threshold as NA
(by default) or a
given value. See the example.
signature(object="InputFrame")
: returns the
data frame read from SDM file.
signature(object="InputFrame")
: replace the
data frame read from SDM file.
signature(object="InputFrame")
: returns the Ct value
of the SDM file.
signature(object="InputFrame", value="numeric")
:
replace the Ct value in the object with the new values, and return the object.
Rudolf Biczok mailto:[email protected], Jitao David Zhang mailto:[email protected]
SDMFrame
function reads in data from SDM
files. Data from SDM files is used to construct
ddCtExpression
objects to analyze differetial expression.
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## you can also write ## sampdat <- new("SDMFrame",system.file("extdata", "Experiment1.txt", ## package="ddCt")) ## use the getter methods sampleNames(sampdat) ## or the overloaded primitive accessors sampdat[1:3,"Sample"] ## see all unique samples uniqueSampleNames(sampdat) ## replace all sample names 'Sample1' and 'Sample2' in sampdat ## with 'NewSample1' and 'NewSample2' uniqueSampleNames(sampdat,c("Sample1","Sample2")) <- c("NewSample1","NewSample2") uniqueSampleNames(sampdat) ## or use this syntax to replace the gene names uniqueDetectorNames(sampdat) <- c(Gene1="NewGene1", Gene2="NewGene2") uniqueDetectorNames(sampdat) ## remove sample or detector removeSample(sampdat, "Sample1") removeDetector(sampdat, "Gene1") ## replace sample or detector replaceSample(sampdat, "Sample1", "Sample0") replaceDetector(sampdat, "Gene1", "PLCG1") ## right censoring the data rightCensoring(sampdat, 35) rightCensoring(sampdat, 35, 35)
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## you can also write ## sampdat <- new("SDMFrame",system.file("extdata", "Experiment1.txt", ## package="ddCt")) ## use the getter methods sampleNames(sampdat) ## or the overloaded primitive accessors sampdat[1:3,"Sample"] ## see all unique samples uniqueSampleNames(sampdat) ## replace all sample names 'Sample1' and 'Sample2' in sampdat ## with 'NewSample1' and 'NewSample2' uniqueSampleNames(sampdat,c("Sample1","Sample2")) <- c("NewSample1","NewSample2") uniqueSampleNames(sampdat) ## or use this syntax to replace the gene names uniqueDetectorNames(sampdat) <- c(Gene1="NewGene1", Gene2="NewGene2") uniqueDetectorNames(sampdat) ## remove sample or detector removeSample(sampdat, "Sample1") removeDetector(sampdat, "Gene1") ## replace sample or detector replaceSample(sampdat, "Sample1", "Sample0") replaceDetector(sampdat, "Gene1", "PLCG1") ## right censoring the data rightCensoring(sampdat, 35) rightCensoring(sampdat, 35, 35)
Abstract factory for data input
A virtual Class: No objects may be created from it.
files
:Input files
colmap
:Column mapping
Rudolf Biczok and Jitao David Zhang
showClass("InputReader")
showClass("InputReader")
Read QuantStudio file(s)
QuantStudioFrame(file) readQuantStudio(file)
QuantStudioFrame(file) readQuantStudio(file)
file |
Character vector of filenames |
This function reads the data given in the QuantStudio output files
A object of class SDMFrame
Jitao David Zhang mailto:[email protected]
sampdat <- QuantStudioFrame(system.file("extdata", c("QuantStudio_File1.txt", "QuantStudio_File2.txt"),package="ddCt")) ## use the getter methods sampleNames(sampdat) ## or the overloaded primitive accessors sampdat[1:3,"Sample"] ## see all unique samples uniqueSampleNames(sampdat) ## replace all sample names 'Sample1' and 'Sample2' in sampdat ## with 'NewSample1' and 'NewSample2' uniqueSampleNames(sampdat,c("Sample1","Sample2")) <- c("NewSample1","NewSample2") uniqueSampleNames(sampdat) ## or use this syntax to replace the gene names uniqueDetectorNames(sampdat) <- c(Hs00559368_m1="Gene1", Hs02576168_g1="Gene2") uniqueDetectorNames(sampdat)
sampdat <- QuantStudioFrame(system.file("extdata", c("QuantStudio_File1.txt", "QuantStudio_File2.txt"),package="ddCt")) ## use the getter methods sampleNames(sampdat) ## or the overloaded primitive accessors sampdat[1:3,"Sample"] ## see all unique samples uniqueSampleNames(sampdat) ## replace all sample names 'Sample1' and 'Sample2' in sampdat ## with 'NewSample1' and 'NewSample2' uniqueSampleNames(sampdat,c("Sample1","Sample2")) <- c("NewSample1","NewSample2") uniqueSampleNames(sampdat) ## or use this syntax to replace the gene names uniqueDetectorNames(sampdat) <- c(Hs00559368_m1="Gene1", Hs02576168_g1="Gene2") uniqueDetectorNames(sampdat)
NTC stands for Non-template controls. This method remove the NTC samples from the input object.
signature(object = "ddCtExpression")
An object hat has been analyzed with the ddCt
method
signature(object = "InputFrame")
An input object
The function replces (or updates) the items of a given vector by
checking the equality with the target
parameter. If found, the
item will be replaced by the value
parameter. The length of
both target
and value
must be the same and could be
longer than 1, in which case the replace will be iterated.
replaceVectorByEquality(vector, target, value)
replaceVectorByEquality(vector, target, value)
vector |
A vector to be replaced. The items of the vector must be atom types, since the equality is checked by '=='. |
target |
targets to be replaced, could be either single or a vector |
value |
values to be replaced at the positions of targets, must
be of the same length of |
A warning will be prompted if any item in the target
cannot be found
A vector of the same length as the parameter vector
Jitao David Zhang
==
for checking equality.
vector <- c("java", "perl", "python", "c#") replaceVectorByEquality(vector, target="c#", value="c/c++") replaceVectorByEquality(vector, target=c("c#","perl"), value=c("c/c++","R"))
vector <- c("java", "perl", "python", "c#") replaceVectorByEquality(vector, target="c#", value="c/c++") replaceVectorByEquality(vector, target=c("c#","perl"), value=c("c/c++","R"))
Read an SDM file: Data Output File for SDS, Version 2.1
SDMFrame(file) readSDM(file)
SDMFrame(file) readSDM(file)
file |
Character vector of filenames |
This function reads the data given in the colums 'Detector','Sample' and 'Ct' of the specified SDM output file(s) and stores them as a data.frame. An additional column including the respective filename is added.
This function is a wrapper for the SDMFrame constructor
A object of class SDMFrame
Rudolf Biczok mailto:[email protected]
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## you can also write ## sampdat <- new("SDMFrame",system.file("extdata", "Experiment1.txt", ## package="ddCt")) ## or with ## sampdat <- readSDM(system.file("extdata", "Experiment1.txt", ## package="ddCt")) ## use the getter methods sampleNames(sampdat) ## or the overloaded primitive accessors sampdat[1:3,"Sample"] ## see all unique samples uniqueSampleNames(sampdat) ## replace all sample names 'Sample1' and 'Sample2' in sampdat ## with 'NewSample1' and 'NewSample2' uniqueSampleNames(sampdat,c("Sample1","Sample2")) <- c("NewSample1","NewSample2") uniqueSampleNames(sampdat) ## or use this syntax to replace the gene names uniqueDetectorNames(sampdat) <- c(Gene1="NewGene1", Gene2="NewGene2") uniqueDetectorNames(sampdat)
## read a SDM file sampdat <- SDMFrame(system.file("extdata", "Experiment1.txt", package="ddCt")) ## you can also write ## sampdat <- new("SDMFrame",system.file("extdata", "Experiment1.txt", ## package="ddCt")) ## or with ## sampdat <- readSDM(system.file("extdata", "Experiment1.txt", ## package="ddCt")) ## use the getter methods sampleNames(sampdat) ## or the overloaded primitive accessors sampdat[1:3,"Sample"] ## see all unique samples uniqueSampleNames(sampdat) ## replace all sample names 'Sample1' and 'Sample2' in sampdat ## with 'NewSample1' and 'NewSample2' uniqueSampleNames(sampdat,c("Sample1","Sample2")) <- c("NewSample1","NewSample2") uniqueSampleNames(sampdat) ## or use this syntax to replace the gene names uniqueDetectorNames(sampdat) <- c(Gene1="NewGene1", Gene2="NewGene2") uniqueDetectorNames(sampdat)
Write a 'data.frame' into an html table within a html page
write.htmltable(x, file, title = "", sortby = NULL, decreasing = TRUE, open = "wt")
write.htmltable(x, file, title = "", sortby = NULL, decreasing = TRUE, open = "wt")
x |
'data.frame' |
file |
character. File name. |
title |
character. Title of html page |
sortby |
character. Name of column by which to sort the table rows |
decreasing |
logical. Should the sort order be increasing or decreasing? |
open |
character. This argument is passed onto 'file' |
The funciton is called for its side effect: writing a file
Wolfgang Huber
Write a 'data.frame' into a tab delimited file (not quoted and no-row-name TSV file)
writeSimpleTabCsv(x, file, ...)
writeSimpleTabCsv(x, file, ...)
x |
'data.frame' |
file |
character. File name. |
... |
Additional arguments passed onto the function |
The function is called for its side effect: writing a file
Wolfgang Huber