[Qpcr-commits] r126 - in pkg: NormqPCR NormqPCR/R NormqPCR/inst/doc NormqPCR/man QCqPCR/R QCqPCR/inst/doc
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri Dec 24 03:05:26 CET 2010
Author: jperkins
Date: 2010-12-24 03:05:25 +0100 (Fri, 24 Dec 2010)
New Revision: 126
Modified:
pkg/NormqPCR/NAMESPACE
pkg/NormqPCR/R/combineTechReps.R
pkg/NormqPCR/R/deltaCt.R
pkg/NormqPCR/R/deltaDeltaCt.R
pkg/NormqPCR/R/plotDCt.R
pkg/NormqPCR/R/plotDdCt.R
pkg/NormqPCR/R/selectHKs.R
pkg/NormqPCR/R/stabMeasureM.R
pkg/NormqPCR/R/stabMeasureRho.R
pkg/NormqPCR/inst/doc/NormqPCR.Rnw
pkg/NormqPCR/man/combineTechReps.Rd
pkg/NormqPCR/man/makeAllNAs.Rd
pkg/NormqPCR/man/replaceAboveCutOff.Rd
pkg/NormqPCR/man/replaceNAs.Rd
pkg/QCqPCR/R/PseudoPlot.R
pkg/QCqPCR/R/QCnormedData.R
pkg/QCqPCR/inst/doc/QCqPCR.Rnw
Log:
a few things..it's been a while
Modified: pkg/NormqPCR/NAMESPACE
===================================================================
--- pkg/NormqPCR/NAMESPACE 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/NAMESPACE 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,2 +1,2 @@
importClasses(qPCRSet)
-export(geomMean,stabMeasureM,stabMeasureRho,selectHKs,deltaCt,replaceNAs,deltaDeltaCt,makeAllNAs,combineTechReps,replaceAboveCutOff,makeAllNAs,replaceNAs,plotDCt)
+export(geomMean,stabMeasureM,stabMeasureRho,selectHKs,deltaCt,replaceNAs,deltaDeltaCt,makeAllNAs,combineTechReps,replaceAboveCutOff,makeAllNAs,replaceNAs,plotDCt,plotDdCt)
Modified: pkg/NormqPCR/R/combineTechReps.R
===================================================================
--- pkg/NormqPCR/R/combineTechReps.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/combineTechReps.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -9,14 +9,11 @@
if (FALSE %in% grepl("_TechReps", origDetectors)) stop("These are not tech reps")
newDetectors <- unique(gsub("_TechReps.\\d","", origDetectors))
NewExpM <- matrix(nrow = length(newDetectors), ncol = dim(expM)[2], dimnames = list(newDetectors,colnames(expM)))
- for(detector in newDetectors){
- cat(detector,"\t")
- cat(expM[gsub("_TechReps.\\d","", origDetectors) %in% origDetectors],"\t")
- dValues <- colMeans(expM[grepl(detector, origDetectors),],na.rm=TRUE)
- cat(dValues,"\n")
- NewExpM[detector,] <- dValues
+ for(detector in newDetectors){ # for each detector
+ dValues <- colMeans(expM[gsub("_TechReps.\\d","", origDetectors) %in% detector,],na.rm=TRUE) # find everything that begins with new detector, and find the mean
+ NewExpM[detector,] <- dValues # add values to the matrix
}
- NewExpM[is.na(NewExpM)] <- NA
+ NewExpM[is.na(NewExpM)] <- NA # make NAs real NAs
qPCRBatch <- new("qPCRBatch", exprs = NewExpM, phenoData = phenoData(qPCRBatch))
return(qPCRBatch)
}
Modified: pkg/NormqPCR/R/deltaCt.R
===================================================================
--- pkg/NormqPCR/R/deltaCt.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/deltaCt.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,34 +1,44 @@
setGeneric("deltaCt",
- function(qPCRBatch, hkgs, combineHkgs=FALSE)
+ function(qPCRBatch, hkgs, combineHkgs=FALSE, calc="arith")
standardGeneric("deltaCt")
)
setMethod("deltaCt", signature = "qPCRBatch", definition =
- function(qPCRBatch, hkgs, combineHkgs) {
+ function(qPCRBatch, hkgs, combineHkgs, calc) {
hkgs <- make.names(hkgs)
+# cat("hkgs",hkgs)
if(FALSE %in% (hkgs %in% featureNames(qPCRBatch))) stop ("given housekeeping gene, ", hkgs," not found in file. Ensure entered housekeeping genes appear in the file")
expM <- exprs(qPCRBatch)
hkgM <- expM[hkgs, ]
- print(hkgM)
+# print(hkgM)
if(length(hkgs) > 1) {
if (TRUE %in% apply(hkgM, 1, is.na)) {
- warning("NAs present in housekeeping genes readings")
- if (0 %in% apply(! apply(hkgM, 1, is.na),2,sum)) stop("Need at least 1 non NA for each housekeeper")
+ warning("NAs present in housekeeping genes readings")
+ if (0 %in% apply(! apply(hkgM, 1, is.na),2,sum)) stop("Need at least 1 non NA for each housekeeper")
}
# hkgV <- apply(hkgM,2,geomMean,na.rm=TRUE) CANT CURRENTLY DO THIS BECAUS EOF PROBLEMS WITH DOUBLE NAS
hkgV <- vector(length = dim(hkgM)[2])
- for(i in 1:dim(hkgM)[2]) {
+ if(calc == "arith") {
+ for(i in 1:dim(hkgM)[2]) {
if(! FALSE %in% is.na(hkgM[,i])) next # go to next sequence if we have only NAs
- hkgV[i] <- mean(hkgM[,i], na.rm=TRUE)
+ hkgV[i] <- mean(hkgM[,i], na.rm=TRUE)
+ }
}
- }
+ else {
+ for(i in 1:dim(hkgM)[2]) {
+ if(! FALSE %in% is.na(hkgM[,i])) next # go to next sequence if we have only NAs
+ hkgV[i] <- geomMean(hkgM[,i], na.rm=TRUE)
+ }
+ }
+# cat("hkgV",hkgV,"\n")
+ }
else {
if(TRUE %in% is.na(hkgM)) {
warning("NAs present in housekeeping gene readings")
if(! FALSE %in% is.na(hkgM)) stop("Need at least 1 non NA for the housekeeper")
}
hkgV <- hkgM # Because it's a vector really anyway
- cat(hkgV)
+# cat(hkgV)
}
exprs(qPCRBatch) <- t(t(exprs(qPCRBatch)) - hkgV)
return(qPCRBatch)
Modified: pkg/NormqPCR/R/deltaDeltaCt.R
===================================================================
--- pkg/NormqPCR/R/deltaDeltaCt.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/deltaDeltaCt.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,17 +1,14 @@
setGeneric("deltaDeltaCt",
- function(qPCRBatch, maxNACase=0, maxNAControl=0, hkgs, contrastM, case, control, paired=TRUE, combineHkgs=FALSE, calc="arith")
+ function(qPCRBatch, maxNACase=0, maxNAControl=0, hkgs, contrastM, case, control, paired=TRUE, hkgCalc="arith", statCalc="arith")
standardGeneric("deltaDeltaCt")
)
setMethod("deltaDeltaCt", signature = "qPCRBatch", definition =
- function(qPCRBatch, maxNACase, maxNAControl, hkgs, contrastM, case, control, paired, combineHkgs, calc) {
+ function(qPCRBatch, maxNACase, maxNAControl, hkgs, contrastM, case, control, paired, hkgCalc, statCalc) {
hkgs <- make.names(hkgs)
-#case <- as.character(case)
-#control <- as.character(control)
- if(combineHkgs == TRUE) {
- if(length(hkgs) == 1) stop("Not enough hkgs given")
- }
-# for(hkg in hkgs) {
-# if(! hkg %in% featureNames(qPCRBatch)) stop("invalid housekeeping gene")
+
+# if(combineHkgs == TRUE) {
+ if(length(hkgs) < 1) stop("Not enough hkgs given")
+# }
if(FALSE %in% (hkgs %in% featureNames(qPCRBatch))) stop("invalid housekeeping gene given")
for(hkg in hkgs){
if(sum(is.na(hkg)) > 0) warning(hkg, " May be a bad housekeeping gene to normalise with since it did not produce a reading ", sum(is.na(hkg)), "times out of", length(hkg))
@@ -22,21 +19,23 @@
caseM <- expM[,cases]
controlM <- expM[,controls]
-# hkgMCase <- caseM[hkgs, ]
-# hkgMControl <- controlM[hkgs, ]
- if(combineHkgs == TRUE) {
+ if(length(hkgs) > 1) {
hkgMCase <- caseM[hkgs, ]
hkgMControl <- controlM[hkgs, ]
- hkgVCase <- apply(hkgMCase, 2, geomMean, na.rm=TRUE)
- hkgVControl <- apply(hkgMControl, 2, geomMean, na.rm=TRUE)
- } else {
- hkg <- hkgs[1]
+ #hkgVCase <- apply(hkgMCase, 2, geomMean, na.rm=TRUE)
+ #hkgVControl <- apply(hkgMControl, 2, geomMean, na.rm=TRUE)
+ if(hkgCalc == "arith") {
+ hkgVCase <- apply(hkgMCase, 2, mean, na.rm=TRUE)
+ hkgVControl <- apply(hkgMControl, 2, mean, na.rm=TRUE)
+ } else {
+ hkgVCase <- apply(hkgMCase, 2, geomMean, na.rm=TRUE)
+ hkgVControl <- apply(hkgMControl, 2, geomMean, na.rm=TRUE)
+ }
+ } else { # Just use the first HKG
+ hkgVCase <- caseM[hkgs[1], ]
+ hkgVControl <- controlM[hkgs[1], ]
}
- hkgVCase <- caseM[hkg, ]
- hkgVControl <- controlM[hkg, ]
-
-
sdHkgCase <- sd(hkgVCase, na.rm=TRUE)
sdHkgControl <- sd(hkgVControl, na.rm=TRUE)
@@ -65,16 +64,16 @@
# dCtControl <- NA
sdCase <- NA
} else {
- if(calc == "arith") {
+ if(statCalc == "geom") {
dCtCase <- mean(2^-(VCase - hkgVCase), na.rm=TRUE)
sdCase <- sd(2^-(VCase - hkgVCase), na.rm=TRUE)
}
- if(calc == "geom") {
+ if(statCalc == "arith") {
dCtCase <- mean(VCase, na.rm=TRUE) - mean(hkgVCase, na.rm=TRUE)
if (paired == TRUE) {
sdCase <- sd(VCase - hkgVCase, na.rm=TRUE)
} else {
- sdCase <- sqrt(sd(VCase, na.rm=TRUE)^2 + sdHkgCase^2)
+ sdCase <- sqrt(sd(VCase, na.rm=TRUE)^2 + sdHkgCase^2)
}
}
}
@@ -89,11 +88,11 @@
# dCtCase <- NA
sdControl <- NA
} else {
- if(calc == "arith") {
+ if(statCalc == "geom") {
dCtControl <- mean(2^-(VControl - hkgVControl), na.rm=TRUE)
sdControl <- sd(2^-(VControl - hkgVControl), na.rm=TRUE)
}
- if(calc == "geom") {
+ if(statCalc == "arith") {
dCtControl <- mean(VControl, na.rm=TRUE) - mean(hkgVControl, na.rm=TRUE)
if (paired == TRUE) {
sdControl <- sd(VControl - hkgVControl, na.rm=TRUE)
@@ -108,10 +107,10 @@
if(sum(is.na(VControl)) > maxNAControl) {
dCtControl <- NA
}
- if(calc == "arith") {
+ if(statCalc == "geom") {
ddCt <- dCtCase / dCtControl
}
- if(calc == "geom") {
+ if(statCalc == "arith") {
ddCt <- dCtCase - dCtControl
}
if(is.na(ddCt)) {
@@ -128,12 +127,12 @@
maxddCts[i] <- NA
}
else {
- if(calc == "arith") {
- minddCts[i] <- (ddCt - sdCase)
- maxddCts[i] <- (ddCt + sdCase)
+ if(statCalc == "geom") {
+ minddCts[i] <- NA
+ maxddCts[i] <- NA
ddCts[i] <- ddCt
}
- if(calc == "geom") {
+ if(statCalc == "arith") {
minddCts[i] <- 2 ^ -(ddCt + sdCase)
maxddCts[i] <- 2 ^ -(ddCt - sdCase)
ddCts[i] <- 2^-ddCt
@@ -141,14 +140,14 @@
}
}
- if(calc == "arith") {
+ if(statCalc == "geom") {
# ddCts[i] <- (dCtCase/dCtControl)
dCtCases[i] <- dCtCase
sdCtCases[i] <- sdCase
dCtControls[i] <- dCtControl
sdCtControls[i] <- sdControl
}
- if(calc == "geom") {
+ if(statCalc == "arith") {
dCtCases[i] <- 2^-dCtCase
sdCtCases[i] <- sdCase
dCtControls[i] <- 2^-dCtControl
@@ -157,8 +156,9 @@
}
i <- i+1
}
- ddCtTable <- as.data.frame(cbind(featureNames(qPCRBatch),dCtCases,sdCtCases,dCtControls,sdCtControls,ddCts,minddCts,maxddCts))
- names(ddCtTable) <- c("ID", case, paste(case,"sd",sep="."), control, paste(control,"sd",sep="."),"ddCt","ddCt.min", "ddCt.max")
+
+ ddCtTable <- as.data.frame(cbind(featureNames(qPCRBatch),format(dCtCases, digits=4),format(sdCtCases, digits=4),format(dCtControls, digits=4),format(sdCtControls, digits=4),format(ddCts,digits=4),format(minddCts, digits=4),format(maxddCts, digits=4)))
+ names(ddCtTable) <- c("ID", paste("2^-dCt",case,sep="."), paste(case,"sd",sep="."), paste("2^-dCt",control,sep="."), paste(control,"sd",sep="."),"2^-ddCt","2^-ddCt.min", "2^-ddCt.max")
return(ddCtTable)
# stop(ddCtTable)
}
Modified: pkg/NormqPCR/R/plotDCt.R
===================================================================
--- pkg/NormqPCR/R/plotDCt.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/plotDCt.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,29 +1,26 @@
-plotDCt <- function(ddCtTable, detectors="", logFC = FALSE) {
+plotDCt <- function(..., ddCtTable, detectors="", statCalc="arith") {
if(detectors[1]!="") {
ddCtTable <- ddCtTable[ddCtTable$ID %in% detectors, ,drop=FALSE]
} else {
ddCtTable <- ddCtTable
}
plotNames <- ddCtTable$ID
- plotCts <- ddCtTable[,c("case","control")]
- plotSds <- ddCtTable[,c("case.sd","control.sd")]
+ plotCts <- ddCtTable[,c(2,4)]
+ plotSds <- ddCtTable[,c(3,5)]
plotCts <- sapply(plotCts, function(x) as.numeric(as.character(x)))
plotSds <- sapply(plotSds, function(x) as.numeric(as.character(x)))
-# plotCts <- plotTable[,c("case","control")]
-# plotSds <- plotTable[,c("case.sd","control.sd")]
- plotU <- plotCts + plotSds
-# plotL <- plotCts - plotTable[,c("case.sd","control.sd")]
-# plotU <- plotCts + plotSds
- plotL <- plotCts - plotSds
-print(plotCts)
+ if(statCalc == "arith") {
+ plotU <- 2^(log2(plotCts) + plotSds)
+ plotL <- 2^(log2(plotCts) - plotSds)
+ } else {
+ plotU <- plotCts + plotSds
+ plotL <- plotCts - plotSds
+ }
+# print(plotCts)
+
+# cat(plotU,"\t")
+# cat(plotL,"\n")
-#cat(plotU,"\n")
-#cat(plotL,"\n")
-# plotCts <- sapply(ddCtTable, function(x) as.numeric(as.character(x)))
-# plotCtsd
-# plotCtMax
-# plotCtMin
-
- barplot2(height=t(plotCts), plot.ci=TRUE, ci.u=t(plotU), ci.l=t(plotL), beside=T)
+ barplot2(..., height=t(plotCts), plot.ci=TRUE, ci.u=t(plotU), ci.l=t(plotL), beside=T)
}
Modified: pkg/NormqPCR/R/plotDdCt.R
===================================================================
--- pkg/NormqPCR/R/plotDdCt.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/plotDdCt.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,58 +1,29 @@
-#library(gplots)
-plotDdCt <- function(qPCRSet, ddCtTable, detectors="", logFC = FALSE) {
-# no.plots <- length(detectors)
-cat("start")
-#cat(eval(samples))
-#expM
+plotDdCt <- function(..., ddCtTable, detectors="", logFC=TRUE) {
+ if(detectors[1]!="") {
+ ddCtTable <- ddCtTable[ddCtTable$ID %in% detectors, ,drop=FALSE]
+ } else {
+ ddCtTable <- ddCtTable
+ }
+ plotNames <- ddCtTable$ID
+ds <- gsub("\\.\\w*$","",plotNames)
+ plotdDCt <- as.numeric(as.vector(ddCtTable[,6]))
+ plotMin <- as.numeric(as.vector(ddCtTable[,7]))
+ plotMax <- as.numeric(as.vector(ddCtTable[,8]))
-if(detectors!="") {
- bbb <- ddCtTable[ddCtTable$ID %in% detectors,]
- qPCRSet <- qPCRSet[detectors]
-# expM <- matrix(exprs(qPCRSet)[detectors,], ncol=8)
-# rownames(expM) <- detectors
-# expM <- exprs(qPCRSet)
-} else {
-#cat("try")
- bbb <- ddCtTable
-# expM <- exprs(qPCRSet)
-}
-expM <- exprs(qPCRSet)
-#cat("out")
-#print(bbb)
-#cat("YEAH")
-ddCts <- as.vector(bbb$ddCt)
-ddCtsNames <- as.vector(bbb$ID)
-ddCtsMax <- log2(as.numeric(as.vector(bbb$ddCt.max)))
-ddCtsMin <- log2(as.numeric(as.vector(bbb$ddCt.min)))
+plotdDCt[grep("\\+",as.vector(ddCtTable[,6]))] <- max(plotMax, na.rm=T)
+plotdDCt[grep("\\-",as.vector(ddCtTable[,6]))] <- min(plotMin, na.rm=T)
+bar.colours <- rep("blue",length(plotdDCt))
+bar.colours[grep("\\+",as.vector(ddCtTable[,6]))] <- "red"
+bar.colours[grep("\\-",as.vector(ddCtTable[,6]))] <- "red"
-maxo <- ceiling(max(ddCtsMax,na.rm=T))
-mino <- floor(min(ddCtsMin,na.rm=T))
-nonNAs.back <- log2(as.numeric(ddCts))
-nonNAs.back[ddCts == "+"] <- maxo
-nonNAs.back[ddCts == "-"] <- mino
-ddCts <- nonNAs.back
-op <- par(las=2)
-bp <- barplot2(ddCts, plot.ci=TRUE, ci.l=ddCtsMin, ci.u=ddCtsMax, names.arg=ddCtsNames, ylab=expression("Fold change relative to control ("~Delta~Delta~"Cti)"))
-par(op)
-#for(i in rownames(qq)) { t<-c(t,(ww[i,]));points(rep(j[m],4),ww[i,]);m=m+1 }
-#qq <- qPCRSet
-#cat("bbbbbbl")
-#return(expM)
-controlSubbed <- expM[,c(3:4,7:8)] - apply(expM[,c(1:2,5:6)],1,mean,na.rm=TRUE)
-#two <- 2^(-controlSubbed)
-#three <- log2(two)
-#controlSubbed <- log2((2^-((expM[,c(3:4,7:8)]) - (apply(expM[,c(1:2,5:6)],1,mean,na.rm=TRUE)))))
-#controlSubbed <- three
-#print(controlSubbed)
-print(expM)
-controlSubbed <- -(expM[,c(3:4,7:8)] - apply(expM[,c(1:2,5:6)],1,mean,na.rm=TRUE))
-controlSubebed <- log2(controlSubbed)
-#if(length(detectors) == 1) expM <- matrix(expM, ncol=4)
-#cat(featureNames(qPCRSet))
-#rownames(expM) <- featureNames(qPCRSet)
-#print(expM)
-#return(expM)
-#source("/home/bsm/jperkins/qpcr/pkg/NormqPCR/R/plotDdCt.R")
-#j<-1; for(i in 1:length(featureNames(qPCRSet))) { print(controlSubbed[i,]); points(rep(bp[j],4),controlSubbed[i,]);j=j+1 }
+
+ if(logFC == TRUE) {
+ plotdDCt <- log2(plotdDCt)
+ plotMin <- log2(plotMin)
+ plotMax <- log2(plotMax)
+ }
+ barplot2(..., names=ds, height=plotdDCt, plot.ci=TRUE, ci.u=plotMax, ci.l=plotMin, las=2, col=bar.colours, ylab="Log2 Fold Change")
}
+
+
Modified: pkg/NormqPCR/R/selectHKs.R
===================================================================
--- pkg/NormqPCR/R/selectHKs.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/selectHKs.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -12,7 +12,6 @@
log = TRUE, Symbols, trace = TRUE, na.rm = TRUE){
if(method == "geNorm") {
if(class(x) == "qPCRBatch") x <- t(exprs(x))
-# else x <- t(x)
}
if(method == "NormFinder") {
if(class(x) == "qPCRBatch"){
@@ -21,9 +20,6 @@
}
else stop("'x' must be of class qPCRBatch")
}
-# else x <- t(x)
-# print(x)
-# x <- t(exprs(x))
if(!is.matrix(x) & !is.data.frame(x))
stop("'x' needs to be of class matrix or data.frame")
if(is.data.frame(x)) x <- data.matrix(x)
@@ -57,7 +53,7 @@
R[1:minNrHKs] <- Symbols
else
R[i] <- Symbols[ind]
-
+
if(i > 2){
if(log){
NF.old <- rowMeans(x)
@@ -69,7 +65,7 @@
V[n-i+1] <- sd(log2(NF.new/NF.old), na.rm = na.rm)
}
}
-
+
if(trace){
cat("###############################################################\n")
cat("Step ", n-i+1, ":\n")
@@ -88,10 +84,6 @@
}
if(method == "NormFinder"){
-# if(missing(group))
-# group <- as.factor(pData(x)[,"Group"])
-
-
NF <- stabMeasureRho(x, group = group, log = log, na.rm = na.rm, returnAll = TRUE)
k <- length(NF$rho)
R <- character(minNrHKs)
@@ -101,7 +93,7 @@
b[1] <- which.min(rho)
rho.min <- numeric(minNrHKs)
rho.min[1] <- rho[b[1]]
-
+
if(trace){
cat("###############################################################\n")
cat("Step ", 1, ":\n")
Modified: pkg/NormqPCR/R/stabMeasureM.R
===================================================================
--- pkg/NormqPCR/R/stabMeasureM.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/stabMeasureM.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -25,7 +25,6 @@
N[N < 1] <- NA
Mean <- colMeans(A, na.rm = na.rm)
M[j] <- mean(sqrt(rowSums((t(A) - Mean)^2, na.rm = na.rm)/(N - 1)))
-# M[j] <- mean(apply(A, 2, sd, na.rm = na.rm))
}else
M[j] <- sd(A, na.rm = na.rm)
}
Modified: pkg/NormqPCR/R/stabMeasureRho.R
===================================================================
--- pkg/NormqPCR/R/stabMeasureRho.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/R/stabMeasureRho.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -4,7 +4,6 @@
## na.rm: remove NA values
stabMeasureRho <- function(x, group, log = TRUE, na.rm = TRUE, returnAll = FALSE){
if(class(x) == "qPCRBatch") {
-# stop("Must supply qPCRBatch to this function")
if(missing(group)) group <- pData(x)[,"Group"]
x <- t(exprs(x))
}
Modified: pkg/NormqPCR/inst/doc/NormqPCR.Rnw
===================================================================
--- pkg/NormqPCR/inst/doc/NormqPCR.Rnw 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/inst/doc/NormqPCR.Rnw 2010-12-24 02:05:25 UTC (rev 126)
@@ -16,7 +16,7 @@
pagebackref,bookmarks,colorlinks,linkcolor=darkblue,citecolor=darkblue,%
pagecolor=darkblue,raiselinks,plainpages,pdftex]{hyperref}
%
-\markboth{\sl Package ``{\tt NormqPCR}''}{\slinst/exData/ Package ``{\tt NormqPCR}''}
+\markboth{\sl Package ``{\tt NormqPCR}''}{\sl Package ``{\tt NormqPCR}''}
%
%------------------------------------------------------------------------------
\newcommand{\code}[1]{{\tt #1}}
@@ -362,7 +362,7 @@
<<dCt>>=
hkgs<-"Actb-Rn00667869_m1"
-qPCRBatch.norm <- deltaCt(qPCRBatch = qPCRBatch.taqman, hkgs = hkgs)
+qPCRBatch.norm <- deltaCt(qPCRBatch = qPCRBatch.taqman, hkgs = hkgs, calc="arith")
head(exprs(qPCRBatch.norm))
@
@@ -381,7 +381,7 @@
<<dCt many genes>>=
hkgs<-c("Actb-Rn00667869_m1", "B2m-Rn00560865_m1", "Gapdh-Rn99999916_s1")
-qPCRBatch.norm <- deltaCt(qPCRBatch = qPCRBatch.taqman, hkgs = hkgs)
+qPCRBatch.norm <- deltaCt(qPCRBatch = qPCRBatch.taqman, hkgs = hkgs, calc="arith")
head(exprs(qPCRBatch.norm))
@
@@ -446,7 +446,7 @@
<< ddCt >>=
hkg <- "Actb-Rn00667869_m1"
-ddCt.taqman <- deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkg=hkg, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype")
+ddCt.taqman <- deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkg=hkg, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", statCalc="geom", hkgCalc="arith")
head(ddCt.taqman)
@
@@ -457,7 +457,7 @@
<< ddCt Avg >>=
hkg <- "Actb-Rn00667869_m1"
-ddCtAvg.taqman <- deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkg=hkg, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", paired=FALSE)
+ddCtAvg.taqman <- deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkg=hkg, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", paired=FALSE, statCalc="geom", hkgCalc="arith")
head(ddCtAvg.taqman)
@
%-------------------------------------------------------------------------------
@@ -475,7 +475,7 @@
colnames(contM) <- c("interestingPhenotype","wildTypePhenotype")
rownames(contM) <- sampleNames(qPCRBatch.taqman)
hkgs<-c("Actb-Rn00667869_m1", "B2m-Rn00560865_m1", "Gapdh-Rn99999916_s1")
-ddCt.gM.taqman <- deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkgs=hkgs, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", combineHkgs=TRUE)
+ddCt.gM.taqman <- deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkgs=hkgs, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", statCalc="arith", hkgCalc="arith")
head(ddCt.gM.taqman)
@
@@ -488,7 +488,7 @@
colnames(contM) <- c("interestingPhenotype","wildTypePhenotype")
rownames(contM) <- sampleNames(qPCRBatch.taqman)
hkgs<-c("Actb-Rn00667869_m1", "B2m-Rn00560865_m1", "Gapdh-Rn99999916_s1")
-ddAvgCt.gM.taqman <-deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkgs=hkgs, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", paired=FALSE, combineHkgs=TRUE)
+ddAvgCt.gM.taqman <-deltaDeltaCt(qPCRBatch = qPCRBatch.taqman, maxNACase=1, maxNAControl=1, hkgs=hkgs, contrastM=contM, case="interestingPhenotype", control="wildTypePhenotype", paired=FALSE, statCalc="arith", hkgCalc="arith")
head(ddAvgCt.gM.taqman)
@
Modified: pkg/NormqPCR/man/combineTechReps.Rd
===================================================================
--- pkg/NormqPCR/man/combineTechReps.Rd 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/man/combineTechReps.Rd 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,7 +1,7 @@
\name{combineTechReps}
\alias{combineTechReps}
\alias{combineTechReps,qPCRBatch-method}
-\title{ Takes expression set of qPCR values containing technical replicates and combines them. }
+\title{ Combines Technical Replicates }
\description{
Takes expression set of qPCR values containing technical replicates and combines them.
}
@@ -9,19 +9,25 @@
combineTechReps(qPCRBatch)
}
\arguments{
- \item{qPCRBatch}{ Expression set containing qPCR data, read in by a ReadqPCR function and containing technical reps, denoted by _TechRep.n suffix.
+ \item{qPCRBatch}{ Expression set containing qPCR data, read in by a ReadqPCR function and containing technical reps, denoted by \code{_TechRep.n} suffix.
}
}
\details{
- Takes expression set of qPCR values containing technical replicates and combines them.
+ Takes \code{exprs} of qPCR values containing technical replicates and combines them.
}
\value{
- qPCRBatch with same number of samples, but with features reduced by the number of technical replicates that were present on the card.
+ \code{qPCRBatch} with same number of samples, but with less features, since all technical replicates are replaced with a single value of their means.
}
%\references{ }
\author{ James Perkins \email{jperkins at biochem.ucl.ac.uk}}
%\note{}
%\seealso{}
-%\examples{
-%}
+\examples{
+ path <- system.file("exData", package = "NormqPCR")
+ qPCR.example.techReps <- paste(path,"/qPCR.techReps.txt", sep = "")
+ qPCRBatch.qPCR.techReps <- read.qPCR(qPCR.example.techReps)
+ rownames(exprs(qPCRBatch.qPCR.techReps))
+ combinedTechReps <- combineTechReps(qPCRBatch.qPCR.techReps)
+ rownames(exprs(combinedTechReps))
+}
\keyword{data}
Modified: pkg/NormqPCR/man/makeAllNAs.Rd
===================================================================
--- pkg/NormqPCR/man/makeAllNAs.Rd 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/man/makeAllNAs.Rd 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,16 +1,16 @@
\name{makeAllNAs}
\alias{makeAllNAs}
\alias{makeAllNAs,qPCRBatch-method}
-\title{ Make all NAs when number of NAs above a given threshold }
-\description{ Make all NAs when number of NAs above a given threshold
+\title{ Make all Ct values NA }
+\description{ Make all Ct values for a given detector NA when the number of NAs for that detector is above a given threshold
}
\usage{
-makeAllNAs(qPCRBatch, contrastM, sampleMaxM)
+ makeAllNAs(qPCRBatch, contrastM, sampleMaxM)
}
\arguments{
\item{qPCRBatch}{ Expression set containing qPCR data.
}
- \item{contrastM}{ Contrast Matrix like that used in Limma. Columns represent the different samples types, rows are the different samples,
+ \item{contrastM}{ Contrast Matrix like that used in \code{limma}. Columns represent the different samples types, rows are the different samples,
with a 1 or 0 in the matrix indicating which sample types the different samples belong to.
}
\item{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.
@@ -20,12 +20,34 @@
Make all NAs when number of NAs above a given threshold
}
\value{
- qPCRBatch object with a new exprs slot, everything else equal
+ \code{qPCRBatch} object with a new exprs slot, everything else equal
}
%\references{ }
\author{ James Perkins \email{jperkins at biochem.ucl.ac.uk}}
%\note{}
%\seealso{}
-%\examples{
-%}
-\keyword{data}
+\examples{
+ # read in the data
+ path <- system.file("exData", package = "NormqPCR")
+ taqman.example <- paste(path, "/example.txt", sep="")
+ 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",]
+}
+\keyword{data}
\ No newline at end of file
Modified: pkg/NormqPCR/man/replaceAboveCutOff.Rd
===================================================================
--- pkg/NormqPCR/man/replaceAboveCutOff.Rd 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/man/replaceAboveCutOff.Rd 2010-12-24 02:05:25 UTC (rev 126)
@@ -1,8 +1,8 @@
\name{replaceAboveCutOff}
\alias{replaceAboveCutOff}
\alias{replaceAboveCutOff,qPCRBatch-method}
-\title{ Replace values above a given threshold with a new value }
-\description{ Replace values above a given threshold with a new value
+\title{ Replace Ct values with new value}
+\description{ Replace Ct values above a given threshold with a new value
}
\usage{
replaceAboveCutOff(qPCRBatch, newVal=NA, cutOff=38)
@@ -16,15 +16,21 @@
}
}
\details{
- Replaces values in the exprs slot of the qPCRBatch object that are above a threshold value with a new number
+ Replaces values in the exprs slot of the \code{qPCRBatch} object that are above a threshold value with a new number
}
\value{
- qPCRBatch object with a new exprs slot
+ \code{qPCRBatch} object with a new exprs slot
}
%\references{ }
\author{ James Perkins \email{jperkins at biochem.ucl.ac.uk}}
%\note{}
%\seealso{}
-%\examples{
-%}
+\examples{
+ path <- system.file("exData", package = "NormqPCR")
+ taqman.example <- paste(path, "/example.txt", sep="")
+ 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",]
+}
\keyword{data}
Modified: pkg/NormqPCR/man/replaceNAs.Rd
===================================================================
--- pkg/NormqPCR/man/replaceNAs.Rd 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/NormqPCR/man/replaceNAs.Rd 2010-12-24 02:05:25 UTC (rev 126)
@@ -23,6 +23,11 @@
\author{ James Perkins \email{jperkins at biochem.ucl.ac.uk}}
%\note{}
%\seealso{}
-%\examples{
-%}
+\examples{
+ path <- system.file("exData", package = "NormqPCR")
+ taqman.example <- paste(path, "/example.txt", sep="")
+ qPCRBatch.taqman <- read.taqman(taqman.example)
+ qPCRBatch.taqman.replaced <- replaceNAs(qPCRBatch.taqman, newNA = 40)
+ exprs(qPCRBatch.taqman.replaced)["Ccl20.Rn00570287_m1",]
+}
\keyword{data}
Modified: pkg/QCqPCR/R/PseudoPlot.R
===================================================================
--- pkg/QCqPCR/R/PseudoPlot.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/QCqPCR/R/PseudoPlot.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -4,7 +4,6 @@
)
setMethod("PseudoPlot", signature = "qPCRBatch", definition =
function(qPCRBatch, plotType, writeToFile, cutOff, statType, plateToPlot) {
-# CHECKING - IS THERE A CLEVERER WAY TO DO THIS?
if (statType == "parametric"
|| statType == "non-parametric") {
}
@@ -17,15 +16,14 @@
ctsMat[ctsMat > cutOff] <- cutOff
}
else {
- warning("No cutOff value given, if you are calculating residuals, the program it will crash out ungracefully")
+ warning("No cutOff value given, if you are calculating residuals, the program will crash out ungracefully")
}
orderMat <- exprs.well.order(qPCRBatch)
plateVec <- as.vector(gsub("-.*", "", orderMat))
- whichPlates <- unique(plateVec)
- whichPlates <- sort(plateVec)
+ whichPlates <- sort(unique(plateVec))
if(plateToPlot != "AllPlates") whichPlates <- plateToPlot
wellVec <- as.numeric(gsub(".*-", "", orderMat))
-
+cat("ZZZZ",whichPlates,"here\n")
if (plotType == "Cts.Values") {
minVal <- 0
maxVal <- round(max(ctsMat, na.rm=TRUE), 2)
@@ -135,7 +133,7 @@
par(mar = c(5.1, 4.1, 1, 2))
image(x.bar, 1, matrix(x.bar, length(x.bar), 1), axes = FALSE, xlab = "", ylab = "", col = myCol)
Labels <- c("Min", "Max")
- axis(1, at = c(0,40), labels = Labels, las = 1)
+ axis(1, at = c(minVal,maxVal), labels = Labels, las = 1)
if(writeToFile == TRUE) {
dev.off()
Modified: pkg/QCqPCR/R/QCnormedData.R
===================================================================
--- pkg/QCqPCR/R/QCnormedData.R 2010-10-12 15:32:24 UTC (rev 125)
+++ pkg/QCqPCR/R/QCnormedData.R 2010-12-24 02:05:25 UTC (rev 126)
@@ -2,29 +2,60 @@
function(qPCRBatch, hkgs, writeToFile=FALSE)
standardGeneric("plotVsHkg")
)
+
+
setMethod("plotVsHkg", signature = "qPCRBatch", definition =
function(qPCRBatch, hkgs, writeToFile)
{
cts <- exprs(qPCRBatch) # this refers to the actual data
-# hkgs <- # these are the housekeeping genes the data has been normalised to
+ hkgs <- make.names(hkgs)
+
if(FALSE %in% hkgs %in% featureNames(qPCRBatch)) stop("one or more housekeeper not found in exprs matrix")
- plotFrame <- row.names(cts)
+ plotFrame <- matrix(nrow=length(featureNames(qPCRBatch)),ncol=length(hkgs), dimnames = list(featureNames(qPCRBatch), hkgs))
for (hkg in hkgs) {
- cat("wiv da", hkg, "\n")
+ cat("housekeeping genes", hkg, "\n\n")
+
+ dCts <- deltaCt(qPCRBatch = qPCRBatch, hkgs = hkg, calc="arith")
+ for(detector in featureNames(qPCRBatch)) {
+#cat("\n\n\n")
+# print(dCts)
+# cat(exprs(qPCRBatch)[detector,],"\t")
+ plotFrame[detector,hkg] <- mean(exprs(dCts)[detector,],na.rm=TRUE)
+#cat(mean(exprs(qPCRBatch)[detector,],na.rm=TRUE),"\t")
+#cat(plotFrame[detector,hkg],"\n")
# hkg <- gsub("-.+$", "", hkg) # cut off the stuff from the detector's name after the -
# hkg2ddct <- paste(hkg, "_2^DDCt", sep = "")
# ddct <- cts[, hkg2ddct]
-# plotFrame <- data.frame(plotFrame, ddct)
-
+# plotFrame[,hkgs] <- <- data.frame(plotFrame, ddct)
+ }
+#cat("\n")
}
- plotFrame <- plotFrame[, -1] # take off first column of detector names
- for (i in seq(hkgs)) { # now order by each hkg and print
- plotFrame <- plotFrame[order(plotFrame[, i]), ] # order by 1st column
- if (writeToFile) jpeg(file = paste("OrderedBy", hkgs[i], ".jpg", sep = ""))
- matplot(plotFrame, type = "l", pch = seq(hkgs), lty = seq(hkgs), log = "y", main = paste("Ordered By ", hkgs[i], sep = ""),xlab=paste("rank order of ",hkgs[i]))
- legend(1, tail(seq(hkgs), n = 1), hkgs, lty = seq(hkgs), col = seq(hkgs))
- if (writeToFile) dev.off()
+
+# for (hkg in hkgs[1]) {
+# for(detectors in featureNames(qPCRBatch)) {
+# cat(plotFrame[detector,hkg])
+
+# }
+# }
+#stop()
+# plotFrame <- plotFrame[, -1] # take off first column of detector names
+ for (hkg in hkgs) { # now order by each hkg and print
+ ord.plotFrame <- plotFrame[order(plotFrame[, hkg]), ] # order by 1st column
+ if (writeToFile) jpeg(file = paste("mean.deltaCt.ordered.by.", hkg, ".jpg", sep = ""))
+
+
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/qpcr -r 126
More information about the Qpcr-commits
mailing list