[R] constraining functions for optim()

[Guillaume Théroux Rancourt]

Hi Ravi,

> It makes no sense to set your constraints as the solution that you want.

The thing is that y.1 and y.2 are too far away from the measured values. These values (photosynthesis assimilation rates) should not be as far away from each other. With the data I sent where you proposed a very efficient solution, differences were very small in my complete data (under 1%). However, with the data here, y.1shows that there is a 30% difference between the measured and the estimated values. This is why I would like to reduce the difference and bring y.1 and y.2 closer to zero.

I have maybe used the wrong term whit "constraints".

> What was wrong with the approach that I already proposed?  I think it makes
> sense to minimize y1^2 + y2^2 for various reasons (e.g. roughly same scales,
> smoothness of objective function). 

Sorry! Nothing was wrong. I had just made changes to see what happened, but forgot to bring it back.

Thank you for your help, again.

Guillaume


> -----Original Message-----
> From: r-help-bounces at r-project.org [mailto:r-help-bounces at r-project.org] On
> Behalf Of Guillaume Théroux Rancourt
> Sent: Monday, February 08, 2010 2:40 PM
> To: r-help at r-project.org
> Subject: [R] constraining functions for optim()
> 
> Hi list!
> 
> I am optimizing one function with two data sets where the unknown values
> (Vcmax and gi) are the same in each dataset (see script below). The script
> works, but I would like to add other constraints to this function.
> 
> In this function, the optimization is carried out like:
> 
> ### code and functions before - see script below
> y.1 = (12.7313) - f.1
> y.2 = (10.52935) - f.2
> y = y.1^2 + y.2^2
> return(y)
> }
> res <- optim(par=c(0.1,130), fn=f, lower=c(0,0), upper=c(1,300),
> method="L-BFGS-B", control = list(maxit=10000))
> # I need lower and upper boundaries, so I use the L-BFGS-B method.
> 
> I would like to also optimize where y.1 and y.2 are equal or near zero. So
> the three conditions here are:
> y.1 = 0
> y.2 = 0
> y = 0
> 
> In the data set below, we get:
> 
>> res$value
> [1] 0.0004200628
> 
> but: 
>> y.1
> [1] 4.772095
>> y.2
> [1] -2.379089
> 
> Thank you for your help.
> 
> guillaume
> 
> 
> ################
> #### Complete script
> 
> f=function(x) { 
> 	
> Vcmax = x[2]
> gi = x[1]
> 		
> # First data set
> f.1=(-((Vcmax-0.89189)/gi+164.30765+272.38*(1+2*10/165.82))+sqrt(((Vcmax-0.8
> 9189)/gi+164.30765+272.38*(1+2*10/165.82))^2-4*(-1/gi)*(0.89189*(164.30765+2
> 72.38*(1+2*10/165.82))-Vcmax*(164.30765-(5*2/2.60546)))))/(-2/gi)
> 		
> if (is.nan(f.1)) f.1 = 1e30
> 		
> # Second data set
> f.2=(-((Vcmax-0.89189)/gi+161.61032+272.38*(1+10*10/165.82))+sqrt(((Vcmax-0.
> 89189)/gi+161.61032+272.38*(1+10*10/165.82))^2-4*(-1/gi)*(0.89189*(161.61032
> +272.38*(1+10*10/165.82))-Vcmax*(161.61032-(5*10/2.60546)))))/(-2/gi)
> 		
> if (is.nan(f.2)) f.2 = 1e30
> 		
> 		
> # Verification with measured values
> y.1 = (12.7313) - f.1
> y.2 = (10.52935) - f.2
> 		
> y = y.1^2 + y.2^2
> 		
> return(y)
> }
> 	
> res <- optim(par=c(0.1,130), fn=f, lower=c(0,0), upper=c(1,300),
> method="L-BFGS-B", control = list(maxit=10000))
> 
> 
> ## VERIFICATION 
> # First 
> f.1=function(gi, Vcmax){
> (-((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)}
> 	
> # Second 
> f.2=function(gi, Vcmax){
> (-((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)}
> 	
> 
> ### y.1 and y.2 should be equal or near zero
> 	y.1 = 12.7313 - f.1a(res$par[1],res$par[2])
> 	y.2 = 10.52935 - f.2a(res$par[1],res$par[2])
> 
> 
> 
> Guillaume Théroux Rancourt
> Ph.D. candidate --- Plant Biology
> Université Laval, Québec, QC, Canada
> guillaume.theroux-rancourt.1 at ulaval.ca
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>