Entering edit mode
thu.vo
•
0
@313d375a
Last seen 4.1 years ago
I would like to create a PCA plot with ellipse for grouping. I have a dataset from DESeqFromHTSeqCount. I can change the color of point but could not adjust the color of ellipse. My code is below:
> #Set working directory for data
> setwd("C:/Users/t_v026/Dropbox/Atlantic salmon/New RNAseq/DEG Analysis/RawCount")
> directory <- "C:/Users/t_v026/Dropbox/Atlantic salmon/New RNAseq/DEG Analysis/RawCount/"
> sampleTableT <- read.delim("SampleTableTotal.txt", sep = "\t")
>
> ###Create DESeqDataSet
> ddsHTSeqTotal <- DESeqDataSetFromHTSeqCount(sampleTable = sampleTableT, directory = directory, design = ~ Tissue)
Warning message:
In DESeqDataSet(se, design = design, ignoreRank) :
some variables in design formula are characters, converting to factors
> ddsHTSeqTotal
class: DESeqDataSet
dim: 48436 40
metadata(1): version
assays(1): counts
rownames(48436): 1433b 1433z ... zyx zzef1
rowData names(0):
colnames(40): BC26B187 BC27N219 ... P26B180 P26N210
colData names(2): Color Tissue
>
>
> #Keep only read that have counts >10
> keep1 <- rowSums(counts(ddsHTSeqTotal)) >= 10
> ddsHTSeqTotal1 <- ddsHTSeqTotal[keep1,]
> ###Generation of PCA
> vsdT <- vst(ddsHTSeqTotal1,blind=FALSE)
> pltT <- ggplot(pT, aes(PC1, PC2, color = Tissue, shape = Tissue)) +
+ geom_point(aes(color = Tissue, shape=Tissue, size=Tissue))+
+ scale_shape_manual(values=c(16, 17, 18))+
+ scale_color_manual(values=c('#FFFF00','#FF6600', '#3399FF')) +
+ scale_size_manual(values=c(2,2,3)) +
+ theme_bw() +
+ theme(axis.line = element_line(colour = "black"),
+ panel.grid.major = element_blank(),
+ panel.grid.minor = element_blank(),
+ panel.border = element_blank(),
+ panel.background = element_blank()) +
+ xlab(paste0("PC1: ", percentVarT[1], "% variance")) +
+ ylab(paste0("PC2: ", percentVarT[2], "% variance")) +
+ coord_fixed() + ylim(-40,40) +
+ stat_ellipse(aes(x=PC1,y=PC2,fill = Tissue, color = Tissue), geom = "polygon", linetype = "blank", level=0.95, alpha=0.2)
> pltT
# include your problematic code here with any corresponding output
# please also include the results of running the following in an R session
sessionInfo( )
![Result of PCA plot][1]
Could anyone please give me an advice? Thank you so much.
Thu