I am working with a set of RNASeq dataset. I have about 4000 observations (genes) on 20 samples and plotting a PCA I found the clustering doesn't vary much when I use different number of genes, but the % variance of PC1 and PC2 does. I have ~1300 DEG. But how do you select a sensible cutoff for the number of genes/observations to explain PCA variance? Would it be correct to select the number of genes that are DEG which account for >50% of the variance between PC1 and 2?

Few examples:

  • All genes: PC1(34%), PC2(14%)

  • Top1000: PC1(45%), PC2(17%)

  • Top500: PC1(49%), PC2(19%)

  • Top50: PC1(55%), PC2(24%)

  • Top 5: PC1(75%), PC2(24%)

  • 3
    $\begingroup$ Hi @Ecg ... this is a good question. What we need to understand is why not dump all the data into the PCA - it can certainly handle it. What you appear to be saying is 'how can I re-jig PCA to give a result when I get no result'. If you can edit your question to help here that would be great rather than a priori assume PC1 and PC2 comprise max resolution $\endgroup$ – M__ Nov 29 '20 at 18:17
  • 2
    $\begingroup$ Hi, yes I'd guessed this. The more data the fuzzier the signal. There's possibly a way foward, its not to reduce the amount of data .. in 'big data' its not cool. I can't remember its implemented though @StupidWolf might know. Hmmj looking at 'all genes' trapping a limited' variance that is weird. $\endgroup$ – M__ Nov 29 '20 at 21:02
  • 3
    $\begingroup$ Of course, you can subset however you want. You just have to adjust the message and be clear when communicating such a result. In this case, you have to specify that your PC1 explains 50% of the variance of all DE genes, or the top 1000 most variable genes.You however can no longer claim the PC explains 50% of the dataset. $\endgroup$ – Bastian Schiffthaler Nov 29 '20 at 21:23
  • 1
    $\begingroup$ If you are using DESeq2, you are performing a PCA on the top 1000 most variables genes and plotting that. Some of your DEGs might be among these genes but there might be other factors other than the condition of your interest. Agree with Bastian its a a matter of what you want to convey with this PCA $\endgroup$ – StupidWolf Nov 29 '20 at 22:45
  • 1
    $\begingroup$ This is not really what PCA is for, if I understand the question correctly. You don't select observations to "explain" PCA variance, that is more a feature of inferential tools like e.g. ANOVA. I would prefer to use some other measure to discard genes that are simply not informative or measured accurately. For example, first throw out low-expressed or non-expressed genes that are in the noise range. If what you want is to find a minimal set of genes that discriminate between two "clusters", then there are much better tools to do that than PCA, e.g. LDA or something like that. $\endgroup$ – Maximilian Press Nov 29 '20 at 23:31

The general goal of PCA in RNA-seq can be stated as, "I'd like a low-dimension representation of my data to allow easy assessment of the gross structure of my samples, specifically for assessing missing batch effects, samples swaps, etc." Ideally, this low-dimension representation can fit into a plot or two. In other words, we'd like the first 2 (maybe 3) PCs to explain a reasonable amount of variation. A few things to note from this alone:

  • The number of differentially expressed genes were not mentioned. They have no relevance here.
  • There was no mention made of the treatment groups clustering together. Seeing this is NOT a goal of PCA.

As you have noticed, the more input rows (genes, assuming you haven't already transposed it) the less variation explained by the first 2 principal components (PCs). On the flip side, the more rows you input the better PCA is describing the structure of your entire dataset.

Question: Do I really care about the gross structure of the entire dataset?

Answer: No, I just want to see problems that I might need to account for.

So you really don't have to use all of the data in the PCA, just the bits that will allow you to see any issues. In general that'll be the few hundred most-variable genes (after some reasonable transformation so you're not just looking at the most highly expressed genes!). The exact number shouldn't really matter very much, something in the range of 300-500 most variable genes should suffice. A common practice is to just randomly pick a nice round number somewhere in or around this range. This is not a number that needs optimization unless perhaps you have a very large number of samples and therefore PC1 and PC2 only account for a small (say 10-20%) percentage of the variation, in which case you're better off looking at either a number of difference principal components or using a different dimension reduction technique, such as UMAP.

  • $\begingroup$ Thanks, this was a very clear answer. Like you said, when using all the genes, both PC1 and 2 account already for almosth 50% of variance which is "strong" enough to my belief and to support my dataset. Thnaks again for the discussion above. $\endgroup$ – Ecg Nov 30 '20 at 11:54

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.