PCA (rda) with presence/absence data

[LaFra7]

Hi all, I would like to perform a PCA with the function rda from the vegan package. I know that if you don't include any other variables such as environmental variables, it works as a PCA. Until now, I worked with relative abundances of trnL sequences, but now I would like to know if it is possibile to do it with presence/absence data. If so, should I include any additional argument to the function?

Thank you!

[gavinsimpson]

PCA isn't an appropriate method to use with presence absence data. Instead consider correspondence analysis (CA, via cca()), or convert your data to dissimilarities using a coefficient/method that is appropriate to 0/1 data and then use non-metric multidimensional scaling via metaMDS().

[LaFra7]

Thank you very much for your reply, but sorry for the question.. How do I convert my presence-absence data to dissimilarities?

[gavinsimpson]

I would use

dij <- vegdist(data, method = "jaccard", binary  = TRUE)

to get the binary form of the Jaccard distance

[LaFra7]

Ok, thanks a lot for your help!

[jarioksa]

I agree that CA or NMDS probably are more adequate methods, but there is nothing that makes PCA unsuitable for 0/1 data. In distance terms, PCA is based on Euclidean distances and with 0/1 data such distances are √[S₁ + S₂ − 2 × S_1∩2], where S are numbers of species in sampling units 1 and 2, and index 1∩2 refers to the number of species occurring in both sampling units (shared species). This is a legitimate way of assessing community distances, and inherently used in PCA of 0/1 data. However, it is not the most efficient and informative way, especially when the numbers of species vary a lot in the data, and the distances may rather reflect these differences in species numbers (α diversity) than differences among sampling units (β diversity). CA handles naturally differences in richness (α) and there are many dissimilarity indices that do even better work, and can be analysed with NMDS (metaMDS) or metric methods such as dbrda in vegan, cmdscale in stats or wcmdscale in vegan.

As to your technical question: PCA (including rda) handles numbers as numbers, and if your input is binary (0/1), it will handle zeros and ones. You need no special argument, but just to provide binary data.

[LaFra7]

Thank you very much for your reply. But now I have some doubts: I analyzed my first database with relative abundances of trnL. So, I just selected the matrix where I had the relative abundances of each plant found in each sample (data_rra) and I applied the function:

gen.pca <- vegan::rda(data_rra, scale = FALSE)

You said that with 0/1 data "the distances may rather reflect these differences in species numbers (α diversity) than differences among sampling units (β diversity)". Does it happen also with relative abundances?

Because I want to reproduce my analysis with the same dataset but in this case with 0/1 data and I dont understand what could be the difference if I still use the rda function.

Thank you for your help!

[gavinsimpson]

I was being too dogmatic :-) PCA is just the operation of certain mathematical operations so it is OK to apply to data that are 0/1.

[jarioksa]

PCA retains Euclidean distances of data: the distances of sampling unit points is related to their Euclidean distances. The Euclidean distance between two quantitative vectors x and y expands to √[Σx² + Σy² - 2 Σxy]. For 0/1 data these sums are numbers of species (numbers of ones) and for the last crossproduct term the number of species occurring in both sampling units (shared species). Quantitative PCA with scale = FALSE suffers from the same problem as PCA of 0/1 data: the variable sums of squares can be more important than the last term of sum of crossproducts (that measures similarity of sampling units).

One solution is to transform data to equalize the sampling unit totals. The most often recommended method is to use Helllinger standardization. In vegan you can achieve this with:

gen.pca <- rda(decostand(data_rra, "hellinger"))

After this standardization the retained distances are called Hellinger distances which can be directly written as √[2 - 2 Σ√xy / √(ΣxΣy)] which no longer is dependent on variable sum of squares. Similar transformation can also be used with 0/1 data.

[LaFra7]

Ok, it is clear for 0/1 data.

But I'm not sure I'm getting the point for the relative abundances: should I transform to hellinger even if my abundances are normalized (the sum of each row is 1)?

[jarioksa]

I hope you know how your data are standardized: normalization means that sum of squares of each row is 1. If your data are standardized to row unit total (divided by the row sum), then square root of its elements will be Hellinger standardization. Hellinger goes one step further from row total standardization. With Hellinger standardization the equation simplifies to more concise form (se previous messages), you get rid of squared terms, and completely rid of effect of varying totals of sampling units (without square root, but with unit total, you still get a variable term that is related to Simpson diversity).

[LaFra7]

Ok, thank you for your help!