Hi !
I need to compare the transcriptome of 5 cases and 5 controls using one-color Agilent microarrays. Weights were computed by the following function:
#My weight function
myFlagFun <- function(x) {
#Weight only strongly positive spots 1, everything else 0
present <- x$gIsPosAndSignif == 1
probe <- x$ControlType == 0
manual <- x$IsManualFlag == 0
strong <- x$gIsWellAboveBG == 1
y <- as.numeric(present & probe & manual & strong)
#Weight weak spots 0.5
weak <- strong == FALSE
weak <- (present & probe & manual & weak)
weak <- grep(TRUE,weak)
y[weak] <- 0.5
#Weight flagged spots 0.5
sat <- x$gIsSaturated == 0
xdr <- x$gIsLowPMTScaledUp == 0
featureOL1 <- x$gIsFeatNonUnifOL == 0
featureOL2 <- x$gIsFeatPopnOL == 0
flagged <- (sat & xdr & featureOL1 & featureOL2)
flagged <- grep(FALSE, flagged)
good <- grep(TRUE, y==1)
flagged <- intersect(flagged, good)
y[flagged] <- 0.5
y
}
Then I used the following script:
RG <- read.maimages(targets, source="agilent",green.only=TRUE, wt.fun=myFlagFun)
RGb <- backgroundCorrect(RG, method="none")
MA <- normalizeBetweenArrays(RGb, method="quantile")
MA.ave <- avereps(MA, ID=MA$genes$ProbeName)
f <- factor(targets$Condition, levels = unique(targets$Condition))
design <- model.matrix(~0 + f)
colnames(design) <- levels(f)
fit <- lmFit(MA.ave, design)
If I understood correctly, the weights are used in gene-wise weighted least squares regression to estimate the coefficients in the linear model. So they will affect the logFC values. My problem is that I need to get back to normalized individual measures so that their means are in line with coefficients calculated by the lmFit function.
Here is an example of what I get on a single gene:
In controls:
MA.ave values | 5.61618123133188 | 5.0746766862945 | 4.8972404255748 | 4.90207357931074 | 5,68229237143083 |
MA.ave weights | 1 | 0.5 | 0.5 | 0 | 1 |
The corresponding coefficient is
5.42814405
The corresponding coefficient is
Please let me know how I can get back to individual normalized measures taking into account the attributed weights ? Their means should be similar to coefficients reported in the fit object. Thank you very much for your help, Best Regards, Stephane |
Well, obviously the weights have some effect on the regression, otherwise we wouldn't need to use them. If TTCA does linear regression but doesn't accept weights, then there's no general solution for modifying the observations to mimic the effect of the weights in the linear model fit. Even if TTCA does accept weights, it may not be using the weights in the same way, e.g., if it's not doing linear regression. I would ignore the weights and accept the fact that the TTCA analysis may not be optimal. It's not like TTCA and limma would agree perfectly anyway, so you'll have to deal with different results regardless of what you do.
P.S. Respond to answers using "add comment" or "add reply" rather than adding a new answer.