[R] Append multiple optim result with a loop (or apply?)
Guillaume Théroux Rancourt
Guillaume.Theroux-Rancourt at fsaa.ulaval.ca
Sat Feb 6 00:48:56 CET 2010
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)
More information about the R-help
mailing list