Optimisation function

[Guillaume Théroux Rancourt]

Hi DRUG users,

I'm trying to find the fit to measured data 2 variables involved in predicting photosynthesis based on light response curves. Right now, my code gives me this error:

Error in optim(c(200, 0.5), .fun.abs_fixed) : objective function in optim evaluates to length 7 not 1

I get it that the error probably comes from Cc having a length of 7, while I'm trying to fit Jmax and curvature each of length. One way to maybe tackle this would be to have Jmax and curvature each of length=7, but I want a unique value.

Any thoughts? I feel lost. Any hints to tidy up this code is welcomed!

To explain the code of the final optimization:
- J is computed from the Jmax and curvature values;
- Cc is then fitted with the optimize function (which I know works because I use it somewhere else in my full script);
- An.pred (predicted values) are then computed from;
- Then, the difference between the measured and the predicted values is squared.

Thanks!

Guillaume

Vcmax=120; Kc=259; Ko=179
gm.total=0.4; O=200; Ci=300  
Sco=2.5; Rd.total=1
G.star=0.5 * (O / Sco)
Kco = Kc * (1 + O/Ko)
I = c(349.571984317016, 180.556171115112, 89.8671780143736, 89.8149127418517,
      45.8250868240358, 27.1946174978256, 9.08917854537962)
An.meas = c(21.8704644017378, 14.3670999006506, 7.88797540677535, 7.58768883616424,
            3.86654753565164, 2.0512337565626, -0.029012968206055)

## Optimisation to compute Cc in final optimization
.fun = function(Cc, Ci, gm, Vcmax, Kc, Ko, O, J, G.star, Rd) {
  ## squared diff between two computed An values
  Vc.c  = Cc * Vcmax  /  (Cc + Kc  * (1 + O/Ko) )
  Vc.j  = J  / (4 + 8 * G.star/Cc)
  .diff = ((gm * (Ci - Cc)) - ((1 - G.star / Cc)   *  min(Vc.c , Vc.j) - Rd ))^2
}

## Final optimization
.fun.abs_fixed = function(x) {
  Jmax = x[1]
  curvature = x[2]
  J = (I + Jmax - sqrt((I +  Jmax)^2  - 4 * curvature * I * Jmax))    /  (2 * curvature)

## Part that might cause the error (?)  
  for(j in 1:length(I)){
    Cc[j] = optimize(f=.fun, interval=c(0,500), Ci=Ci, gm=gm.total, Vcmax=Vcmax, Kc=Kc, Ko=Ko, O=O, J=J[j], G.star=G.star, Rd=Rd.total)$minimum
  }
 
  Vc.c  = Cc * Vcmax  /  (Cc + Kc  * (1 + O/Ko) )
  Vc.j  = J  / (4 + 8 * G.star/Cc)
  An.pred = ((1 - G.star / Cc)  *  min(Vc.c , Vc.j) - Rd.total )
 
  .diff = (An.meas - An.pred)^2
}

optim(c(200, 0.5), .fun.abs_fixed)

[Jaime Ashander]

Hi Guillaume,

The following changes to your final optimization eliminate this error. Note this makes your objective to minimize the sum of squared deviations

- Jaime

## Final optimization
.fun.abs_fixed = function(x) {
  Jmax = x[1]
  curvature = x[2]
  J = (I + Jmax - sqrt((I +  Jmax)^2  - 4 * curvature * I * Jmax))    /  (2 * curvature)
 

  Cc <- numeric(length(I)) ## needed to add this to avoid:
  ##Error in Cc[j] = optimize(f = .fun, interval = c(0, 500), Ci = Ci, gm = gm.total,  :
  ##                            object 'Cc' not found
  ##                          Calls: optim -> <Anonymous> -> fn

  ## Part that might cause the error (?) 
  for(j in 1:length(I)){
    Cc[j] = optimize(f=.fun, interval=c(0,500), Ci=Ci, gm=gm.total, Vcmax=Vcmax, Kc=Kc, Ko=Ko, O=O, J=J[j], G.star=G.star, Rd=Rd.total)$minimum
  }
 
  Vc.c  = Cc * Vcmax  /  (Cc + Kc  * (1 + O/Ko) )
  Vc.j  = J  / (4 + 8 * G.star/Cc)
  An.pred = ((1 - G.star / Cc)  *  min(Vc.c , Vc.j) - Rd.total )
 
  .diff = (An.meas - An.pred)^2

  sum(.diff)## optimize base on sum of squared differences
}

optim(c(200, 0.5), .fun.abs_fixed)

--
Check out our R resources at http://www.noamross.net/davis-r-users-group.html
---
You received this message because you are subscribed to the Google Groups "Davis R Users' Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to davis-rug+...@googlegroups.com.
Visit this group at http://groups.google.com/group/davis-rug.
For more options, visit https://groups.google.com/d/optout.