Guillaume Théroux Rancourt
2010-Feb-05  23:48 UTC
[R] Append multiple optim result with a loop (or apply?)
Hello R-helpers,
Thanks to Ravi Varadhan, I have improve the function I am working on top
optimize two equations. Now, my objective is to do a series of optimization from
a data table, where each row is one data serie (i.e data from one plant) to be
optimized.
The function below works up to the for loops. Before that, it outputs the
desired values. However, what I want to do is to loop the function (fun() in
this case) where each of the loop should optimize the functions using the ith
elements. I am working with a data table, so each of the elements that are used
in the testF (below) are vectors, since each of them are columns of the data
table.
I am new at writing function and looping, and ?for and ?apply did not help me
understand my error, as well as R help etc. I was wondering if the apply family
functions could be of any help here.
Also, I might not be handling my data correctly (is the data frame correct or
should I convert the data into a matrix?). Any improvement is warmly welcome.
Thank you for your help.
Guillaume
Guillaume Th?roux Rancourt
Ph.D. candidate --- Plant Biology
Universit? Laval, Qu?bec, QC, Canada
guillaume.theroux-rancourt.1 at ulaval.ca
#####
Here is my complete code:
testF = function(A.a, A.b, Ci.a, Ci.b, O.a, O.b, Rd, Sco, id)
{
		
fun=function(x) {
	
f=function(x) { 
		
Vcmax = x[2]
gi = x[1]
			
# First 
f.1=(-((Vcmax-Rd)/gi+Ci.a+272.38*(1+O.a*10/165.82))+sqrt(((Vcmax-Rd)/gi+Ci.a+272.38*(1+O.a*10/165.82))^2-4*(-1/gi)*(Rd*(Ci.a+272.38*(1+O.a*10/165.82))-Vcmax*(Ci.a-(5*O.a/Sco)))))/(-2/gi)
			
if (is.nan(f.1)) f.1 = 1e30
			
# Second 
f.2=(-((Vcmax-Rd)/gi+Ci.b+272.38*(1+O.b*10/165.82))+sqrt(((Vcmax-Rd)/gi+Ci.b+272.38*(1+O.b*10/165.82))^2-4*(-1/gi)*(Rd*(Ci.b+272.38*(1+O.b*10/165.82))-Vcmax*(Ci.b-(5*O.b/Sco)))))/(-2/gi)
			
if (is.nan(f.2)) f.2 = 1e30
			
			
# Verification with measured A values
y.1 = A.a - f.1
y.2 = A.b - f.2
			
y = y.1^2 + y.2^2
			
return(y)
}
		
res <- optim(par=c(0.15,50), fn=f, lower=c(0,0), upper=c(1,250),
method="L-BFGS-B")
		
output = data.frame(res$par[2],res$par[1])
colnames(output) = c("Vcmax", "gi")
		
return(output)
		
}
	
# Looping
n = length(A.a)
	
for (i in 1:n) {
A.a = A.a[i]
A.b = A.b[i]
Ci.a = Ci.a[i]
Ci.b = Ci.b[i]
O.a = O.a[i]
O.b = O.b[i]
Rd = Rd[i]
Sco = Sco[i]
		
fn = fun(i)
fn
}
}
#### SAMPLE DATA
	O.21	A.21	Ci.21	O.2	A.2	Ci.2	O.10	A.10	Ci.10	Vcmax	Rd	gi	Gamma.s	Sco	Ci.s
1	21	7.478326633	164.6573484	2	12.73133946	164.3076501	10	10.52935048
161.6103191	64.48	0.89189	0.2835	40.3	2.605459057	37.15400353
2	21	5.797166726	162.4292412	2	11.12480648	165.8296147	10	8.704447278
165.0450986	53.232	0.9703	0.22335	40.3	2.605459057	35.95569734
pop=read.table("file.name", header=T)
testF(pop$A.21, pop$A.2, pop$Ci.21, pop$Ci.2, pop$O.21, pop$O.2, pop$Rd,
pop$Sco)