CODEWITHSUNDEEP
rconstraintsrdavegantriplot
Steve Crittenden
Steve Crittenden

Reputation: 31

R vegan RDA not all levels of constraint are displayed in triplot

In my RDA triplot I would like to display 'sites', 'species' and their constraints which in my case are Field and Trt. The problem is that not all levels of the constraints are displayed in the plot. There are two levels of each factor.

My RDA code is:

Dummy.rda <- rda(species.rda ~ Field + Trt,RDA.env, scale=TRUE)

summary(Dummy.rda, scaling=3)  #Here I see only one level of each reported in:Biplot scores for constraining variables. However all levels appear in: Centroids for factor constraints

anova.cca(Dummy.rda, step=100, by='margin') # degrees of freedom are correct for both factors (df=1)

plot(Dummy.rda, scaling = 3) #This displays all levels of Field and Trt but only one of each has an arrow

plot(Dummy.rda, display = "species", xlim = xlims, ylim = ylims, 
       scaling = 3)
text(Dummy.rda, scaling = 3, display = "bp")  # I want to customize the RDA plot, but this 'text' only displays 1 level of each of Field and Trt.

Upvotes: 3

Views: 4154

Answers (1)

Gavin Simpson
Gavin Simpson

Reputation: 174813

The missing levels are because you are trying to view the factor variables as if they were continuous ones - strictly they should not be displayed as biplot arrows I guess. Anyway, just as in regression with dummy variables, one of the levels of the factor cannot be included because it is linearly dependent on the dummy variables for the remaining levels in the model matrix. Consider this example:

require("vegan")
data(dune)
data(dune.env)

mod <- rda(dune ~ Management, data = dune.env)

> model.matrix(mod)
   ManagementHF ManagementNM ManagementSF
2             0            0            0
13            0            0            1
4             0            0            1
16            0            0            1
6             1            0            0
1             0            0            1
8             1            0            0
5             1            0            0
....<truncated>

What you see in the output from model.matrix() are the variables that went into the ordination. Notice that there are three variables in the model matrix but there are four levels in the Management factor:

> with(dune.env, levels(Management))
[1] "BF" "HF" "NM" "SF"

The convention in R is to use the first level of a factor as the reference level. In a regression this would be included in the intercept, but we don't have one of those in RDA. Notice how in the first row of the model.matrix() output, all values are 0; this indicates that that row was in the BF management group. But as only three variables went into the model proper, we can only represent them by three biplot arrows - that is just the way the maths works.

What we can do is plot the group centroids and that is what is shown in the summary() output you refer to and which can be extracted using scores():

> scores(mod, display = "cn")
                   RDA1       RDA2
ManagementBF -1.2321245  1.9945903
ManagementHF -1.1847246  0.7128810
ManagementNM  2.1149031  0.4258579
ManagementSF -0.5115703 -2.0172205
attr(,"const")
[1] 6.322924

Hence to add the centroids to the existing plot do:

text(mod, scaling = 3, display = "cn")

No matter what you do, you can't add a biplot arrow for the reference group.

I hope that explains the behaviour you are seeing?

Upvotes: 4

Related Questions