[R-sig-dyn-mod] failed attempt to use FME's modFit

Hi,

did you test for simultanaeous identifiability of your parameters?

If not, then this may be an idea, and you may also try to fix one or 
several parameters and try to fit the others. In addition, you may 
consider to constrain both, lower and upper parameter limits.

Thomas

try this

venz = read.table(file.choose(), header = T)
venz

time      pob
1  1873  1732411
    ...   ...
13 2001 23232553
14 2011 27150095

attach(venz)
x = time - time[1]
x

[1]   0   8  18  47  53  63  68  77  88  98 108 117 128 138

y = (pob -pob[1])/1000000
y

[1]  0.000000  0.272728  0.489161  0.747114  1.081720  1.631936  2.118360 
3.302427  5.791588  8.989111
[11] 12.784324 16.372854 21.500142 25.417684

rhs <- function(x, b0, b1,b2) {b0 /(1+b1*exp(-b2*x))}
ds = data.frame(x,y)
modvenz <- nls(y ~ rhs(x, alpha, beta, expco), data = ds, start = 
list(alpha = 10

+  ,beta = 7, expco = 0.2), trace = T)
858.1073 :  10.0  7.0  0.2
686.5785 :  9.52205514 5.82514721 0.06032698
557.9303 :  10.39670917  5.39120514  0.03527819
 ...         ...           ...         ...
0.3278482 :   34.076566 530.985497   0.053214
0.3278482 :   34.07660934 530.98407321   0.05321395

summary(modvenz)

Formula: y ~ rhs(x, alpha, beta, expco)

Parameters:
       Estimate Std. Error t value Pr(>|t|)
alpha  34.07661    0.81994   41.56 1.90e-13 ***
beta  530.98407   48.27397   11.00 2.83e-07 ***
expco   0.05321    0.00117   45.46 7.13e-14 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.1726 on 11 degrees of freedom

Number of iterations to convergence: 17
Achieved convergence tolerance: 3.295e-07

----- Original Message ----- 
From: "jose romero" <jlaurentum at yahoo.com>
To: <r-sig-dynamic-models at r-project.org>
Sent: Tuesday, July 03, 2012 11:35 PM
Subject: [R-sig-dyn-mod] failed attempt to use FME's modFit


Hello list:

I'm attempting to fit some population growth models to venezuelan census 
data:

"time" "pob"
1873 1732411
1881 2005139
1891 2221572
1920 2479525
1926 2814131
1936 3364347
1941 3850771
1950 5034838
1961 7523999
1971 10721522
1981 14516735
1990 18105265
2001 23232553
2011 27150095

For example, the code for the generalized logistic growth model (as given in 
Petzoldt's "Konstruktion ?kologischer Mo delle mit der Open Source Software 
R", page 55) is thus:

#Modelos de crecimiento para la poblacion venezolana
#Primero, se lee la data de censo
censo <- read.table("censos.csv",header=TRUE)
library("simecol")
library("FME")
#Mdelo de crecimiento logistico generalizado
censo7 <- new("odeModel",
main=function(time,init,parms, ...) {
x <- init
p <- parms
dpob <- 
x["pob"]^p["alfa"]*p["c_anual"]*(1-(x["pob"]/p["K"])^p["beta"])^p["gamma"]
list(dpob)},
init=c(pob=1732411),
parms=c(c_anual=0.02014121,K=100000000,alfa=1,beta=1,gamma=1),
times=seq(from=1873,to=2011,by=1),
solver="rk")
#fitOdeModel
resultado <- fitOdeModel(censo7,obstime=censo$time,
yobs=censo[2],fn=ssqOdeModel,
method = "PORT", lower = c(c_anual=0,K=0,alfa=0,beta=0,gamma=0),
scale=c(1/0.1,1/100000000,1/0.1,1/0.1) )
resultado
#FME
parms(censo7) <- 
c(c_anual=0.0000444,K=27390000,alfa=1.4,beta=1.6,gamma=0.2982)
fCosto <- function(p) {
parms(censo7) <- p
out <- out(sim(censo7))
cost <- modCost(out, censo, weight="std")
cost
}
ajuste <- modFit(f = fCosto, p = parms(censo7))
summary(ajuste)
#Fin

So my questions are as follows:

1) Other attempts to use fitOdeModels with a different solver other than 
"rk" ("lsoda" for example) fail. Why is this?
2) There was no way to get FME's modFit to work. Why is it giving me an 
error indicating that I need at least 2 non NA values to interpolate with 
approx? This error also resulted with fitOdeModel on other models with less 
"fortunate" guesses of initial parameter values.

3}) An earlier attempt to use the Blumberg model, which is essentially a 
special case of the generalized logistic growth model minus one parameter 
was succesfull with bot modFit and fitOdeModel. I tried feeding the 
resulting parameter values of the Blumberg model to the generalized model, 
but no success. An attempt to use the "pseudo" method in modFit first and 
the retrofeed the parameter values back into modFit for further searching, 
as suggested in the FME vignette was not successful either.

What I am doing wrong in this?

Thanks in advance,

jose romero




--------------------------------------------------------------------------------

_______________________________________________
R-sig-dynamic-models mailing list
R-sig-dynamic-models at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-dynamic-models