This is a test version of Biostars. For the public version, visit https://www.biostars.org.
boxplot issue

boxplot

I am working DESeq2 I wanted I proper boxplot but I don't understand what wrong I am doing and I wanted to know how to properly plot the x and y axis what parameters should be taken for boxplot . here's my code:

library("DESeq2")
library("ggplot2")
counts<-read.delim("PC_1.csv",header = TRUE, row.names = 1, sep = ",")
counts

#condition c=control P=patient
condition <- factor(c("C","C","C","C","C","C","C","C","P","P","P","P","P","P","P","P"))
condition

#creating data frame of condition and samples 
coldata <- data.frame(row.names = colnames(counts), condition)
coldata

#?DESeqDataSetFromMatrix
dds <- DESeqDataSetFromMatrix(countData = counts, colData = coldata, design = ~condition)
dds <- DESeq(dds)

#vst-variance stabilizing transformation
vsdata <- vst(dds, blind = FALSE)
#QC of data
plotPCA(vsdata, intgroup="condition")
#dispersion of data- ability between replicates as a function of normalized read counts(red line below 1)
plotDispEsts(dds)

#table of DEG compared p vs c
res <- results(dds, contrast = c("condition","P","C"))
res
#log2FC= +ve- upregulated ,-ve downregulated samples
#padj=p adjusted value

#new table for filtering gene which has padj value less than 0.05
sigs <- na.omit(res)
sigs <- sigs[sigs$padj < 0.05,]
sigs

library(ComplexHeatmap)
library(org.Hs.eg.db)
sigs.df <- as.data.frame(sigs)
#filtering
sigs.df <- sigs.df[(sigs.df$baseMean>100) & (abs(sigs.df$log2FoldChange)>1),]
sigs.df

#normalized counts
mat <- counts(dds, normalized = TRUE)[rownames(sigs.df),]

#z-score for each row
t(apply(mat, 1, scale))
colnames(mat) <- rownames(coldata)
mat

heatmap(mat,cluster_rows = T, cluster_columns = T, column_lables=colnames(mat), name="z-score")

boxplot(mat, 
        xlab="", 
        ylab="",
        las=2,
        col="red")
deseq2 boxplot

Please get used to make your question reproducible and precise, meaning provide example data and reduce code to a minimum.

1 answer

First of all, you're not outputting the result of your scale to anything...

t(apply(mat, 1, scale))
colnames(mat) <- rownames(coldata)

needs to be:

mat <- t(apply(mat, 1, scale))
colnames(mat) <- rownames(coldata)

So if you want to use the unscaled normalised counts:

options("scipen"=999) # set the default to scientific notation for large numbers off 
# plot the boxplot with the y axil log-transformed
boxplot(mat, 
        xlab="", 
        ylab="",
        las=2,
        col="red",
        log = "y", 
        pch = 16)
options("scipen"=0) # reset the scientific notation settings

Alternatively, consider transforming your dds object with either varianceStabilizingTransformation() or rlogTransformation() which will give you an object that has ~log2 transformed data.

data.rlog <- rlogTransformation(dds)
mat <- assay(data.rlog)

Log in to answer this question.