[R-sig-eco] Simultaneous estimation of nonlinear models

Thomas Petzoldt Thomas.Petzoldt at tu-dresden.de
Mon May 4 09:58:34 CEST 2009


Hi Fred,

the new development package "FME" (Flexible Modelling Environment) 
contains (among others) functions to support fitting nonlinear models. 
It is focused on numerically solved differential equation models in 
ecology and related sciences, but it works also for other models.

The package was not yet uploaded on CRAN but can be installed from the 
R-Forge development server:

 >  install.packages("FME", repos="http://r-forge.r-project.org")

Note that the latest binary version requires to have R version 2.9.0 (or 
newer).

Look at the first example ("Fitted with analytical solution") for 
function "modFit" or read chapter "Fitting the model to data" of the 
package vignette:

http://r-forge.r-project.org/plugins/scmsvn/viewcvs.php/*checkout*/pkg/FME/inst/doc/FME.pdf?root=fme

Hope it helps

Thomas Petzoldt

Fred Takahashi wrote:
 > Hello
 >
 > I work with plant ecophysiology and use R for about 5 years. Now I
 > need to fit nonlinear data of response of photosynthetic rate to CO2
 > concentration  according well know mechanistic models. Basically,
 > there are one equation for the first part of the curve, and in the
 > second part, another equation (different factors limiting
 > photosynthesis). But some recent research pointed out that is better
 > to fit the data simultaneously to that two equations, than arbitrarily
 > cut the data in the half (simultaneous estimation, also know as
 > segmented regression). That is possible if the algorithm, like
 > Levenberg-Marquardt, minimizes some function of two equations,
 > informed in the conditional syntax (ex. If (eq1>eq2) then eq1 else
 > eq2).
 >
 > Dubois et al. (2007) showed the procedure using SAS. The key point of
 > that is to inform the equations in the conditional syntax like:
 > “
 >  if rubisco <  RuBP then A =  rubisco;
 >  else A =  RuBP;
 >  ...
 > or the equivalent:
 >  ...
 > A =  min(rubisco, RuBP);
 >  ...
 > “
 > when “rubisco” and “RuBP” are objects representing the two equations.
 >
 > Also in Dubois et al. (2007), there are the entire syntax in SAS:
 > “
 > proc sort data=ACI;
 > by factor1 factor2 days;
 > run;
 > ods output ParameterEstimates=estimates;
 > proc NLIN data=ACI method=Marquardt outest=status noitprint;
 > by factor1 factor2 days;
 > rubisco= ((Vcmax*(Ci-Gstar))/(Ci+(Kc*(1+(O/Ko)))));
 > RuBP= ((J*(Ci-Gstar))/((4*Ci)+(8*Gstar)));
 > control O=210; *fixes the value of O;
 > bounds -3<Rd<50;
 > parameters Vcmax= 5 to 330 by 25
 > J= 15 to 505 by 35
 > Rd= 1E-8 to 10.1 by .5;
 > model A=min(rubisco,RuBP)-Rd;
 > output out=ACiOUT parms=Vcmax J Rd predicted=Ahat student=student;
 > run;
 > ods listing;
 > data temporary;
 > set status;
 > where(_TYPE_ ne'GRID' and _TYPE_ ne'COVB');
 > run;
 > data estimates;
 > merge estimates temporary(keep=factor1 factor2 days _status_ _iter_);
 > by factor1 factor2 days;
 > run;
 > ”
 >
 > I try to do that in R in several ways, using the function nls, without
 > success. Does someone have a clue?
 >
 > Thank you for your help!
 >
 > Regards
 > Fred Takahashi
 >
 > PhD student
 > Universidade de Brasilia – Brazil
 >
 >
 >
 >
 > Reference:
 > Dubois, Fiscus, Booker, Flowers and Reid (2007) Optimizing the
 > statistical estimation of the parameters of the Farquhar–von
 > Caemmerer–Berry model of photosynthesis. New Phytologist 176: 402–414
 >
 > _______________________________________________
 > R-sig-ecology mailing list
 > R-sig-ecology at r-project.org
 > https://stat.ethz.ch/mailman/listinfo/r-sig-ecology
 >


-- 
Thomas Petzoldt
Technische Universitaet Dresden
Institut fuer Hydrobiologie        thomas.petzoldt at tu-dresden.de
01062 Dresden                      http://tu-dresden.de/hydrobiologie/
GERMANY

http://tu-dresden.de/Members/thomas.petzoldt



More information about the R-sig-ecology mailing list