Question

Blocked design and Surrogate Variables Analsysis(SVA)

2

Entering edit mode

kodream ▴ 20

@kodream-6952

Last seen 7.7 years ago

United States

I would like to try to use SVA with a block design. I've been using SVA in separate non-blocked experiments with similar data, with known artifacts which SVA removes quite well, however I am not quite sure how to use it with blocking. SVA doesn't seem to be aware of blocking, so I just run it first, then added the surrogate vectors to the design, before I run duplicate correlation. This is what I have so far, it seems like it will miss something:

svobj<-sva(exprs(eset),modMMR, modMMR0, method="irw", n.sv=n.sv)

modsMMR <-cbind(modMMR,svobj$sv)

#my version of duplicateCorrelation that uses the parallel foreach instead of for
corfit <- parallel_dupCor(eset,modsMMR,block=pData(eset)$Pair)
corfit$consensus

lfit <- lmFit(exprs(eset), design=modsMMR, block = pData(eset)$Pair, cor = corfit$consensus)

fit<-eBayes(lfit)

limma sva block • 2.2k views

ADD COMMENT • link updated 2.7 years ago by stuartthafass • 0 • written 10.1 years ago by kodream ▴ 20

0

Entering edit mode

I know this post is very old but could you post your code that runs duplicateCorrelation in parallel? I have a very big dataset and can't quite figure out how to update the rho variable in parallel... Thanks!

ADD REPLY • link 2.7 years ago stuartthafass • 0

score 1 · Answer 1 · 2014-11-10

1

Entering edit mode

Jeff Leek ▴ 650

@jeff-leek-5015

Last seen 3.8 years ago

United States

To handle the blocking, you need to include the blocking variables in both mod and mod0, then you should be able to calculate the surrogate variables and include them in the downstream analysis at the lmFit stage in the model matrix. If you show modMMR and modMMR0, I may be able to help you construct the exact model matrices.

ADD COMMENT • link 10.1 years ago Jeff Leek ▴ 650

1

Entering edit mode

Hi Jeff,

I have the same question regarding a blocked design. I am applying SVA to a large EPIC methylation array dataset with paired samples (patients/controls with enrollment and follow up). I am currently incorporating the block at the lmfit stage as the original poster, as I wasn't quite sure how to incorporate it in the mod correctly. My code is provided below. If you incorporate the blocking variable in the mod, then do you exclude corfit$consensus and block at the lmfit stage? Any help you could provide me with would be very greatly appreciated!

Best,

Katie

#Followup and Enrollment Pooled

############################

mval <- edata

mod <- model.matrix(~as.factor(SubjectType), data=pheno)

mod0 <- model.matrix(~1, data=pheno)

sva.results <- sva(mval, mod, mod0)

sv = sva.results$sv

modSv = cbind(mod,sv)

fit <- lmFit(edata,modSv, block=pheno$patientID,correlation=corfit$consensus)

contrast.matrix <- cbind("PDvsControl"=c(0,1,rep(0,sva.results$n.sv)))

fitContrasts = contrasts.fit(fit,contrast.matrix)

eb = eBayes(fitContrasts)

topTableF(eb, adjust="BH")

summary(decideTests(eb))

#Followup and enrollment separated by condition

########################################

mval <- edata

condition_visit <- paste(pheno$SubjectType, pheno$Visit_Name2, sep="_" ) #4 levels: PD_Enrollment, PD_Followup, HC_Enrollment, HC_Followup

mod <- model.matrix(~0+as.factor(condition_visit), data=pheno)

mod0 <- model.matrix(~1, data=pheno)

sva.results <- sva(mval, mod, mod0)

sv = sva.results$sv

modSv = cbind(mod,sv)

fit <- lmFit(edata,modSv, block=pheno$patientID,correlation=corfit$consensus)

contrast.matrix <- cbind("PDvsHC_Enrollment"=c(-1,0,1,0,rep(0,sva.results$n.sv)),

"PDvsHC_Followup"=c(0,-1,0,1,rep(0,sva.results$n.sv)),

"PD_FollowupvsEnrollment" = c(0,0,-1,1,rep(0,sva.results$n.sv)),

"HC_FollowupvsEnrollment" =c(-1,1,0,0,rep(0,sva.results$n.sv)))

fitContrasts = contrasts.fit(fit,contrast.matrix)

eb = eBayes(fitContrasts)

topTableF(eb, adjust="BH")

summary(decideTests(eb))

ADD REPLY • link 7.2 years ago kfisch ▴ 10