R/parseformula.R

In HMSC: Hierarchical modelling of species community

parseformula <-
function(formula,param,data){
#### F. Guillaume Blanchet - Avril 2013, Juin 2013
##########################################################################################
	#### Extract the response variable (2L) and the model (3L) in the formula
	respFormula<-formula[[2L]]
	explanaFormula<-formula[[3L]]
	
	#### Check if the names in the formula and are associated to variables in the data
	varNames<-colnames(data)
	paramNames<-names(param)
	formulaNames<-all.vars(formula)
	
	varInFormula<-which(formulaNames%in%varNames)
	nvarInFormula<-length(varInFormula)
	
	paramInFormula<-which(formulaNames%in%paramNames)
	nparamInFormula<-length(paramInFormula)
	
	respInFormula<-which(formulaNames%in%as.character(respFormula))
	formulaAllNames<-formulaNames[-c(respInFormula,varInFormula,paramInFormula)]
	
	if(length(formulaAllNames)!=0){
		stop("Some variables (or regression parameters) used in the formula are not found in 'data' (or 'param')")
	}
	
	formulaNames<-all.vars(formula)
	#==================================================
	#### Construct a list including data and parameters
	#==================================================
	nVarParam<-length(formulaNames)
	listData<-vector("list",length=nVarParam)
	
	names(listData)[1:nvarInFormula]<-formulaNames[varInFormula]
	names(listData)[(nvarInFormula+1):(nparamInFormula+nvarInFormula)]<-formulaNames[paramInFormula]
	
	counter<-1
	for(i in which(formulaNames[varInFormula]%in%formulaNames)){
		listData[[counter]]<-data[,i]
		counter<-counter+1
	}
	
	for(i in which(formulaNames[paramInFormula]%in%paramNames)){
		listData[[counter]]<-param[i]
		counter<-counter+1
	}
	
	#=======================
	#### Calculate the model
	#=======================
	model<-eval(explanaFormula,listData)
	
	if(any(is.na(model))){
		stop("NAs were generated from in the evaluation of the model... There is a problem")
	}
	
	#### Extract the response variable
	y<-data[,which(varNames==respFormula)]
	
	return(list(y=y,model=model))
}