Calculate PCoA scores for dataframe "x", based on the distance matrix of dataframe "y"

Question

I'm trying to use multivariate techniques to compare two datasets (same structure) that were collected using different sampling techniques. I'd like to compute a PCoA for the first dataset (D1), and then see how the data from dataset 2 (D2) compares, relative to the distance matrix computed for D1. Essentially, I want to calculate PCoA scores for D2, based on the distance matrix of D1.

I'm using the "vegan" package in R.

Thanks in advance.

#calculate distance matrix and PCoA for "D2"
serb.d2 <- vegdist(serbbiomB2, "bray") #distance matrix for D2
serb.pcoa.2 <- cmdscale(serb.d2, k=4) #PCoA scores for D2
#how do I calculate PCoA scores for "D1", based on "serb.d2"?

Answer 1

If you are going to go down that route, I would suggest you look at capscale() in the vegan package that you are using.

Put all your x and y data in the same object, say XY. Then have another data frame (DF) with the indicator variable you need, say technique, which has levels awesome and cheap (say) --- essentially, this is just a single column with values awesome or cheap indicating which technique was used for each row of XY.

Now you can do

fit <- capscale(XY ~ technique, data = DF, distance = "bray")
fit

This will fit a constrained PCoA (constrained analysis of principal coordinates) and test whether the amount of variation in the dissimilarities of your data that can be explained by which technique was used to generate the data is larger than you'd expect to see if there were no difference between techniques.

Here's an example from vegan using the Dutch dune data set

data(dune, dune.env)
fit <- capscale(dune ~ Management, data = dune.env, distance = "bray")
fit
set.seed(23)
anova(fit)

This gives us

> fit
Call: capscale(formula = dune ~ Management, data = dune.env, distance =
"bray")

              Inertia Proportion Rank
Total          4.2990                
Real Total     4.5940     1.0000     
Constrained    1.5000     0.3265    3
Unconstrained  3.0940     0.6735   14
Imaginary     -0.2950               5
Inertia is squared Bray distance 

Eigenvalues for constrained axes:
  CAP1   CAP2   CAP3 
0.8998 0.4533 0.1471 

Eigenvalues for unconstrained axes:
  MDS1   MDS2   MDS3   MDS4   MDS5   MDS6   MDS7   MDS8   MDS9  MDS10  MDS11 
1.2730 0.4874 0.3784 0.3100 0.2100 0.1507 0.0855 0.0729 0.0601 0.0322 0.0172 
 MDS12  MDS13  MDS14 
0.0098 0.0044 0.0024

and

> set.seed(23)
> anova(fit)
Permutation test for capscale under reduced model
Permutation: free
Number of permutations: 999

Model: capscale(formula = dune ~ Management, data = dune.env, distance = "bray")
         Df Variance     F Pr(>F)   
Model     3   1.5001 2.586  0.003 **
Residual 16   3.0938                
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

The anova() is using a permutation test to compare the ratio of Model to Residual variance (as indicated by the pseudo $F$ entry in the table).

If you do this with the current CRAN version of vegan you'll see a slightly different display for anova() --- I've got the development version (2.2-0 to be) running here.

You could also do this analysis on the raw dissimilarity matrices without embedding them in an Euclidean space (what PCoA does), using the adonis() function in vegan. See its help page ?adonis for details.