[Biomod-commits] overlay response curve plot

[Wilfried Thuiller]

Dear Emilie, dear all,

I have patched the response curve function to overlay the response curves for all selected models. 
This is really cool actually! The results are quite nice. 
Please find attached the new function. Once BIOMOD is loaded, source this one.

source("PATH/response.plot.R")

Compared to the former function, there is one argument more: plot=TRUE/FALSE
If you put plot=TRUE, this is the same than the former one
If you put plot=FALSE, a 3-D array object is returned. This 3-D array contains:

First dimension: the same number of rows than the original data
Second dimension: 2 columns: The first one is the explanatory variable to plot, the second one, the probability of occurrence to plot
Third dimension: The set of environmental variables for which the response.plot was asked to run. 

I have tried with the following code: 

### Only three models are selected (GLM, FDA and RandomForest)
Initial.State(Response=Sp.Env[,12:13], Explanatory=Sp.Env[,4:10])
Models(GLM = TRUE, TypeGLM = "quad", Test = "AIC", FDA = TRUE, RF = TRUE, NbRunEval = 0, Roc=T)

### Load back all models for the first species. 
LoadModels(Sp=1, PA='all', rep=1, models='all')

### Save the object return by the new response.plot function
GLM.RC <- response.plot(Sp290_GLM_full, Sp.Env[,4:10], plot=F)
FDA.RC <- response.plot(Sp290_FDA_full, Sp.Env[,4:10], plot=F)
RF.RC <- response.plot(Sp290_RF_full, Sp.Env[,4:10], plot=F)

### plot the response for GLM for the first variable (the 1 for the third dimension)
plot(GLM.RC[,1,1], GLM.RC[,2,1],  type="l", ylim=c(0,1), xlab=colnames(GLM.RC[,1,])[1], ylab="Probability of occurrence")
### then add the lines for the two other models with different colors
lines(FDA.RC[,1,1], FDA.RC[,2,1], col='blue')
lines(RF.RC[,1,1], RF.RC[,2,1], col='red')
### load the gplots package to make nice and user-friendly legends.
library(plots)
### add a legend
smartlegend(x="left",y="top", inset = 0, c("GLM", "FDA", "RF"), fill = c("black", "blue", "red"))


#### the same but with all variables plotted in the same plot window. 
library(gplots)
par(mfrow=c(3,3))
for(i in 1:ncol(Sp.Env[,4:10])){
	
	plot(GLM.RC[,1,i], GLM.RC[,2,i],  type="l", ylim=c(0,1), xlab=colnames(GLM.RC[,1,])[i], ylab="Probability of occurrence")
	lines(FDA.RC[,1,i], FDA.RC[,2,i], col='blue')
	lines(RF.RC[,1,i], RF.RC[,2,i], col='red')
	smartlegend(x="left",y="top", inset = 0, c("GLM", "FDA", "RF"), fill = c("black", "blue", "red"),cex=0.6)


	
	
}


Hope it helps. This is very likely that we will replace the former response.plot function by this new one. Thanks Emilie for the question! ;-)

Best
Wilfried

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/biomod-commits/attachments/20120307/17ce3490/attachment-0002.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: response.plot.R
Type: application/octet-stream
Size: 4666 bytes
Desc: not available
URL: <http://lists.r-forge.r-project.org/pipermail/biomod-commits/attachments/20120307/17ce3490/attachment-0001.obj>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/biomod-commits/attachments/20120307/17ce3490/attachment-0003.html>