SEQanswers

Go Back   SEQanswers > Applications Forums > RNA Sequencing



Similar Threads
Thread Thread Starter Forum Replies Last Post
DESeq error clemen Bioinformatics 22 03-27-2014 04:59 AM
with or w/o pseudocounts in DESeq? ekimmike Bioinformatics 3 10-04-2012 01:19 PM
DESeq with SPIA billstevens Bioinformatics 7 09-23-2012 03:17 PM
Is more than two conditions possible in DESEQ? greener RNA Sequencing 5 05-09-2011 03:10 PM
DESeq question gfmgfm Bioinformatics 2 04-18-2011 03:15 AM

Reply
 
Thread Tools
Old 03-08-2013, 02:59 AM   #1
skoddo
Junior Member
 
Location: Germany, Freising

Join Date: Jul 2012
Posts: 1
Default DESeq plotPCA changes

Hi everyone,
I have a couple of questions concerning the plotPCA function in the DESeq package. I'm relatively new to R, so maybe they are relatively easy to answer.

1. Is it possible to extract the calculated coordinates from the PCA? I'm talking about the x-axis and y-axis values for the different samples.

2. How can I change the colors in the PCA? I can see from the function code

Code:
function (x, intgroup = "condition", ntop = 500) 
{
    rv = rowVars(exprs(x))
    select = order(rv, decreasing = TRUE)[seq_len(ntop)]
    pca = prcomp(t(exprs(x)[select, ]))
    fac = factor(apply(pData(x)[, intgroup, drop = FALSE], 1, 
        paste, collapse = " : "))
    if (length(fac) >= 3) 
        colours = brewer.pal(nlevels(fac), "Paired")
    else colours = c("green", "blue")
    xyplot(PC2 ~ PC1, groups = fac, data = as.data.frame(pca$x), 
        pch = 16, cex = 2, aspect = "iso", col = colours, main = draw.key(key = list(rect = list(col = colours), 
            text = list(levels(fac)), rep = FALSE)))
}
<environment: namespace:DESeq>
that it uses the "Paired" color palette from the RColorBrewer package. When I try to replicate the function with a different palette, I get an error that it can't find the "rowVars" function. I cannot find this function in R or in the DESeq manual. The default plotPCA is working fine though. I tried attaching my new function to DESeq namespace but I couldn't. I attached it to the package environment but it still didn't work. It seems to me that my problem is a fairly stupid one but even with the excellent manual and some intensive internet recherches I can't figure it out.

3. Is it possible to change the default colored circles to different shapes for different groups?

Thanks in advance for any help you can offer.
Kind Regards
Benedikt
skoddo is offline   Reply With Quote
Old 07-28-2014, 06:42 AM   #2
Nicki1984
Junior Member
 
Location: Scotland

Join Date: May 2013
Posts: 1
Default Same problem hacking plotPCA in DESeq

Hi there,

I have just encountered exactly the same problem with the rowVars function when trying to figure out how to extract the calculated co-ordinates. I was wondering if you ever received an answer to your question of managed to solve the problem yourself?

I would appreciate any help anyone can give!

Many thanks,
Nicki
Nicki1984 is offline   Reply With Quote
Old 08-14-2014, 12:55 AM   #3
rozitaa
Member
 
Location: Sweden

Join Date: Jun 2013
Posts: 51
Default

Hi Nicki and Benedikt,

I have the same problem. Have you got any answer for it?

Thanks,
Rozita
rozitaa is offline   Reply With Quote
Old 08-14-2014, 03:13 AM   #4
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

It's from the genefilter package, so just genefilter::rowVars() to use it. You could also just directly library(genefilter) and then not have to deal with the extra typing (the package is only loaded via a namespace, not attached).
dpryan is offline   Reply With Quote
Old 08-14-2014, 03:54 AM   #5
rozitaa
Member
 
Location: Sweden

Join Date: Jun 2013
Posts: 51
Default

Quote:
Originally Posted by dpryan View Post
It's from the genefilter package, so just genefilter::rowVars() to use it. You could also just directly library(genefilter) and then not have to deal with the extra typing (the package is only loaded via a namespace, not attached).
Hej dpryan,

I didn't understand! This is DeSeq package we are talking about then how I can include gene filter package? There is no rowVar() even there! Can you please explain a bit more.

Thanks
rozitaa is offline   Reply With Quote
Old 08-14-2014, 04:01 AM   #6
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Yes, I know that you're talking about DESeq. DESeq doesn't provide that function, it uses genefilter. If you "library(genefilter)", you'll find that the rowVars() (not rowVar()) function is then useable.
dpryan is offline   Reply With Quote
Old 08-14-2014, 04:31 AM   #7
rozitaa
Member
 
Location: Sweden

Join Date: Jun 2013
Posts: 51
Default

Thanks! Do you also know if it is possible to change PCA plot colors and shapes in Deseq?
rozitaa is offline   Reply With Quote
Old 08-14-2014, 04:35 AM   #8
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

You just need to modify the plotPCA function. Colors are assigned by the "col=colours" parameter in the xyplot function and shapes are the "pch=16". You can provide a vector of shapes if you want and then different groups can have different shapes.
dpryan is offline   Reply With Quote
Old 11-23-2014, 01:56 AM   #9
bpb9
Member
 
Location: NYC

Join Date: Aug 2012
Posts: 24
Default

I am able to change the colors by modifying the plotPCA function (check out ?brewer.pal to see the options for colors--there are several color schemes if you don't like "paired"). Or you can just remove these lines of the plotPCA code altogether:
if (length(fac) >= 3)
colours = brewer.pal(nlevels(fac), "Paired")
else colours = c("green", "blue")

and instead do any of the R colors you like:
colours = c("red","blue","green","darkgoldenrod","hotpink") etc

The second part of the question--pulling the actual values for PC1 and PC2 that are plotted--I haven't quite figured out yet. If I enter parts of the plotPCA function line by line, I am able to print all principal components to the screen, but I cannot save the object (or any subset of columns from the object) as a table, even after trying to data.frame(pca).
bpb9 is offline   Reply With Quote
Old 11-23-2014, 02:48 AM   #10
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

I'm not by a computer with R at the moment, but I recall that it uses prcomp to compute the principal components. I recall that that returns a list, so converting it to a data.frame probably won't work well (what would reasonable dimensions be when the list elements include matrices?). Just read the help page for prcomp and you should be able to tell which list element has the values you're after.
dpryan is offline   Reply With Quote
Old 01-20-2015, 05:31 PM   #11
bpb9
Member
 
Location: NYC

Join Date: Aug 2012
Posts: 24
Default Pull Principal Components from DESeq plotPCA

#Load libraries
library(gplots)
library(RColorBrewer)
library(lattice)
library(genefilter)

#Define the plotPCA function
plotPCA4<-function (x, intgroup = "condition", ntop = 500)
{
rv = rowVars(exprs(x))
select = order(rv, decreasing = TRUE)[seq_len(ntop)]
pca = prcomp(t(exprs(x)[select, ]))
fac = factor(apply(pData(x)[, intgroup, drop = FALSE], 1,
paste, collapse = " : "))
if (length(fac) >= 3)
#colours = brewer.pal(nlevels(fac), "YlOrRd")
colours = brewer.pal(nlevels(fac), "Paired")
else colours = c("darkred", "darkblue")
xyplot(PC2 ~ PC1, groups = fac, data = as.data.frame(pca$x),
pch = 16, cex = 2, aspect = "iso", col = colours, main = draw.key(key = list(rect = list(col = colours),
text = list(levels(fac)), rep = FALSE)))
}

#To make a pca object using vsd transformed data, for example:
pca = prcomp(t(exprs(vsd)[select, ]))

#To view the actual values of each principal component in each individual:
pca$x[,1:10]

#This let me pull, for each of my samples, the first ten principal components.
bpb9 is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 07:40 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO