Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Color RDA vectors by groups in R

Tags:

r

vegan

I'm plotting a series of RDAs in R, each of which has 10+ environmental vectors. The environmental variables each fall in to one of 5 categories. I'd like the vector colors to reflect these categories. I did it in a somewhat ghetto fashion by producing the raw black and white plot, then tracing it in Powerpoint (super time consuming), but there are ALOT of inherent issues with this but it looks good at least. Here's that plot, note the vector colors. enter image description here

I'd like to have my output be an R native plot. The plot I use to make the base plot is currently:

#download the data#
env<- read.csv("environmenta data.csv")
bio<- read.csv("bio data.csv")

#break out the groups of environmental vectors#
#these are purley for example#
group1<-as.matrix(subset(env,select=c(a,b,c))
group2<-as.matrix(subset(env,select=c(d,e,f,h))
group3<-as.matrix(subset(env,select=c(g))
group4<-as.matrix(subset(env,select=c(j,k,l,o))
group5<-as.matrix(subset(env,select=c(m,n,))

#run the RDA#
rda1<-rda(bio,env)

#Plot it#
plot(rda1,type="n",bty="n",main="",
     xlab="XX% variance explained",
     ylab="XX% variance explained",
     col.main="black",col.lab="black", col.axis="white",
     xaxt="n",yaxt="n")
#points(rda1,display="species",col="gray",pch=20)
#option to display species points
#text(rda2,display="species",col="gray") 
#option for species labels

points(rda1,display="cn",col="black",lwd=2)#<<< guessing this is the key statement for my issue, but not sure
text(rda1,display="cn",col="black",cex=0.5)

I'm assuming the answer falls in this points--> col="..." command, but I'm not sure how to tell it to subset by groups. Any and all help would be appreciated.

like image 258
Jesse001 Avatar asked Dec 05 '17 17:12

Jesse001


2 Answers

There are a couple of issues here:

  1. plotting the biplot scores, and
  2. scaling the biplot scores as per plot.cca

The latter doesn't really matter so much if you are only plotting the biplot scores, but a general solution to this kind of plotting requires that we consider that other scores may also be required to be draw.

Here's a reproducible example using data built into

library('vegan')
data(varespec, varechem)

Now we fit a silly ordination (don't do this at home folks)

ord <- rda(varespec ~ ., data = varechem)

Next, extract the biplot scores

bp <- scores(ord, display = 'bp')

and assume we have a variable that contains the groups we want to assign each biplot score to

f <- factor(sample(1:5, nrow(bp), replace = TRUE))

and the associated vector of colours we want to use

cols <- c('red','black','green','navy','purple')

and expand this vector of colours into one of length nrow(bp), i.e. one colour per biplot score/variable

cols <- cols[f]

Next start to plot, by preparing an empty plot region into which we can add the arrows

plot(ord, type = 'n')

In the above example, we leave set up the plot region so that species and site/sample scores are accommodated.

We now calculate a multiplier that allows the biplot scores to occupy a specified proportion (fill = 0.75) of the plot region

mul <- ordiArrowMul(bp, fill = 0.75)

Finally, add the biplot arrows, coloured as per their group

arrows(0, 0, mul * bp[,1], mul * bp[,2],
       length = 0.05, col = cols)

and add labels to the biplot arrows

labs <- rownames(bp)
text(ordiArrowTextXY(mul * bp, labs), labs, col = cols)

This produces

enter image description here

like image 146
Gavin Simpson Avatar answered Sep 28 '22 08:09

Gavin Simpson


library(vegan)

#DATA
data(varespec)
data(varechem)
env = varechem
bio = varespec
rm(list = c("varechem", "varespec"))

#Assign colors based on groups for the columns of env
groups = rep(c("red", "blue"), length.out = NCOL(env))

rda1 = rda(X = bio, Y = env)
plot(rda1, type = "n")
points(rda1)
text(rda1, display = "bp", col = groups)
like image 42
d.b Avatar answered Sep 28 '22 07:09

d.b