Hi Michael,

I'm trying to remove the batch effect from my samples doing like this:

APPOGGIO:

       Well Barcode day bioRep RAcond      TNx techRep
ATTACT   B4  ATTACT  D5     n1     RA NposTneg       1
TTATTG   D4  TTATTG  D5     n2     RA NposTneg       1
ATTCTA   B5  ATTCTA  D5     n1     RA NposTneg       2
TTGAAA   D5  TTGAAA  D5     n2     RA NposTneg       2
ATTTCA   B6  ATTTCA  D5     n1     RA NposTneg       3
TTTACA   D6  TTTACA  D5     n2     RA NposTneg       3
CAAATA   B7  CAAATA  D5     n1     RA NnegTpos       1
TTTCTT   D7  TTTCTT  D5     n2     RA NnegTpos       1
CATTAT   B8  CATTAT  D5     n1     RA NnegTpos       2
TTTTGA   D8  TTTTGA  D5     n2     RA NnegTpos       2
CTATAT   B9  CTATAT  D5     n1     RA NnegTpos       3
AGACCT   D9  AGACCT  D5     n2     RA NnegTpos       3

ddsIO <- DESeqDataSetFromMatrix(countData = CountMat_sel1,
                                colData = APPOGGIO,
                                design = ~bioRep + TNx )

ddsIO <-  DESeq(ddsIO, full=design(ddsIO), reduced = ~bioRep, test="LRT")

resLRT <- results(ddsIO, cooksCutoff=T, independentFiltering = T)

Doing that I got a ddsIO object and I want to go on and see the PCA plot of my samples... for this reason I do:

vsd <- vst(ddsIO)
plotPCA(vsd, "TNx")

This plot is showing me samples without the batch effetc removing...

Only if I use limma::removeBatchEffect (which is actually another funtion from another pkg) as below I can really appreciate the batch effect removing from the PCA plot

ddsIO$batch <- APPOGGIO$bioRep
assay(vsd) <- limma::removeBatchEffect(assay(vsd), vsd$batch)
plotPCA(vsd, "TNx")

How I can be sure about the batch effect removing step only with desing? There is a way to plot a vsd object without using limma::removeBatchEffect and get the "right" plot?

thanks in advance

G.

vst() does not remove batch effects. This is expected. It only transforms the data.

"How I can be sure about the batch effect removing step only with design?"

If you put the batch variable in the design, then you can be assured that the model is controlling for shifts in expression due to batch.

There isn't a way to visualize this process because it is being accounted for implicitly: the model does not run on modified data where batch variation has been removed in a preliminary step.