Question

I am using vegan package to do some analysis and plots, like this

require(vegan)
data(varechem)
data(varespec)

rdam <- rda(varechem,scale=T)


scal=2
ef <- envfit(rdam ~ Cal.vul, data = varespec)
sf <- ordisurf(rdam ~ Cal.vul, data = varespec, plot = FALSE, scaling = scal)
plot(rdam, type="po",scaling = scal)
spe.sc <- scores(rdam, choices=1:2, display="sp",scaling = scal)
arrows(0, 0, spe.sc[, 1], spe.sc[, 2], length=0, lty=1, col="red")
ordilabel(rdam,dis="sp",cex=0.8,col="red",scaling = scal)
plot(ef)
plot(sf, col = "darkgreen", add = TRUE)

But I have to repeat the same plot changing only the variable Cal.vul so I would like to make a function:

plot_spcSurface <- function(rdam,speFram, speName, scal=2 )
{
  ef <- envfit(rdam ~ speName, data = speFram)
  sf <- ordisurf(rdam ~ speName, data = speFram, plot = FALSE, scaling = scal)
  plot(rdam, type="po",scaling = scal)
  spe.sc <- scores(rdam, choices=1:2, display="sp",scaling = scal)
  arrows(0, 0, spe.sc[, 1], spe.sc[, 2], length=0, lty=1, col="red")
  ordilabel(rdam,dis="sp",cex=0.8,col="red",scaling = scal)
  plot(ef)
  plot(sf, col = "darkgreen", add = TRUE)

}  

plot_spcSurface(rdam,varespec,Cal.vul)

it gives

objeto 'Cal.vul' not found 

if I call it like

plot_spcSurface(rdam,varespec,varespec$Cal.vul)

it works, but the name of the species is incorrect labelled "speName" instead of "Cal.vul" and I can't get the name from inside the function, what is the correct way to do this?

Thanks

Was it helpful?

Solution 2

Try this:

plot_spcSurface <- function(rdam,speFram, speName, scal=2 )
{
  form <- as.formula(paste("rdam~",speName))
  ef <- envfit(form, data = speFram)
  sf <- ordisurf(form, data = speFram, plot = FALSE, scaling = scal)
  plot(rdam, type="po",scaling = scal)
  spe.sc <- scores(rdam, choices=1:2, display="sp",scaling = scal)
  arrows(0, 0, spe.sc[, 1], spe.sc[, 2], length=0, lty=1, col="red")
  ordilabel(rdam,dis="sp",cex=0.8,col="red",scaling = scal)
  plot(ef)
  plot(sf, col = "darkgreen", add = TRUE)

}  

plot_spcSurface(rdam,varespec,"Cal.vul")

This passes the name of the predictor variable as a string and builds the formula internally.

OTHER TIPS

You are making this more complicated than it really needs to be; messing with formulas and passing unevaluated arguments always causes, if not problems, headaches when you start writing functions.

Instead, use the default versions of the respective functions, which take an ordination object and a variable to plot as the arguments instead of a formula.

Here is a version that does just that:

plot_spcSurface <- function(ord, speFrame, speName, scal=2 ) {
  ef <- envfit(ord, speFrame[[speName]])
  sf <- ordisurf(ord, speFrame[[speName]], plot = FALSE, scaling = scal)
  plot(ord, type="po", scaling = scal)
  spe.sc <- scores(ord, choices=1:2, display="sp", scaling = scal)
  arrows(0, 0, spe.sc[, 1], spe.sc[, 2], length=0, lty=1, col="red")
  ordilabel(ord, dis="sp", cex=0.8, col="red", scaling = scal)
  plot(ef, labels = speName)
  plot(sf, col = "darkgreen", add = TRUE)
  invisible()
}

plot_spcSurface(rdam, varespec, "Cal.vul")

## or in the development version of vegan >= 2.1-41
plot_spcSurface(rdam, varespec, "Callvulg")

The only issue now is that due to my not accounting for the case where there are no labels stored on the internal object that is returned from envfit(), the function above doesn't label the arrow correctly. It should, but there is a corner case that I didn't account for in plot.envfit(), which I'll now go an fix on R-forge. This was fixed in r2862 on R-forge. (Note that in that version, Jari has changed the species codes so you'll need "Callvulg" now instead of "cal.vulg".)

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top