[Pomp-commits] r501 - in pkg: . R tests
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Fri May 27 21:43:00 CEST 2011
Author: kingaa
Date: 2011-05-27 21:42:59 +0200 (Fri, 27 May 2011)
New Revision: 501
Modified:
pkg/DESCRIPTION
pkg/R/nlf.R
pkg/tests/ou2-nlf.R
pkg/tests/ou2-nlf.Rout.save
Log:
- fix bug in 'nlf': nconverge was essentially ignored
Modified: pkg/DESCRIPTION
===================================================================
--- pkg/DESCRIPTION 2011-05-22 21:01:49 UTC (rev 500)
+++ pkg/DESCRIPTION 2011-05-27 19:42:59 UTC (rev 501)
@@ -1,8 +1,8 @@
Package: pomp
Type: Package
Title: Statistical inference for partially observed Markov processes
-Version: 0.38-1
-Date: 2011-05-22
+Version: 0.38-2
+Date: 2011-05-27
Author: Aaron A. King, Edward L. Ionides, Carles Breto, Steve Ellner, Bruce Kendall, Helen Wearing, Matthew J. Ferrari, Michael Lavine, Daniel C. Reuman
Maintainer: Aaron A. King <kingaa at umich.edu>
URL: http://pomp.r-forge.r-project.org
Modified: pkg/R/nlf.R
===================================================================
--- pkg/R/nlf.R 2011-05-22 21:01:49 UTC (rev 500)
+++ pkg/R/nlf.R 2011-05-27 19:42:59 UTC (rev 501)
@@ -62,8 +62,8 @@
## Vector of times to output the simulation
times <- seq(
- from=t0,
- length=nconverge+nasymp+1,
+ from=t0+nconverge*dt,
+ length=nasymp,
by=dt
)
@@ -88,50 +88,33 @@
return(-val)
}
- if (method == 'subplex') {
- opt <- subplex::subplex(
- par=guess,
- fn=nlf.objfun,
- object=object,
- params=params,
- par.index=par.index,
- times=times,
- t0=t0,
- lags=lags,
- period=period,
- tensor=tensor,
- seed=seed,
- transform=transform,
- nrbf=nrbf,
- verbose=verbose,
- bootstrap=bootstrap,
- bootsamp=bootsamp,
- control=list(...)
- )
- } else {
- opt <- optim(
- par=guess,
- fn=nlf.objfun,
- gr=gr,
- method=method,
- object=object,
- params=params,
- par.index=par.index,
- times=times,
- t0=t0,
- lags=lags,
- period=period,
- tensor=tensor,
- seed=seed,
- transform=transform,
- nrbf=nrbf,
- verbose=verbose,
- bootstrap=bootstrap,
- bootsamp=bootsamp,
- control=list(...)
- )
- }
+ optimizer <- switch(
+ method,
+ subplex=subplex::subplex,
+ stats::optim
+ )
+ opt <- optimizer(
+ par=guess,
+ fn=nlf.objfun,
+ object=object,
+ params=params,
+ par.index=par.index,
+ times=times,
+ t0=t0,
+ lags=lags,
+ period=period,
+ tensor=tensor,
+ seed=seed,
+ transform=transform,
+ nrbf=nrbf,
+ verbose=verbose,
+ bootstrap=bootstrap,
+ bootsamp=bootsamp,
+ control=list(...)
+ )
+
+ opt$est <- est
opt$value <- -opt$value
params[par.index] <- opt$par
opt$params <- params
Modified: pkg/tests/ou2-nlf.R
===================================================================
--- pkg/tests/ou2-nlf.R 2011-05-22 21:01:49 UTC (rev 500)
+++ pkg/tests/ou2-nlf.R 2011-05-27 19:42:59 UTC (rev 501)
@@ -3,12 +3,14 @@
set.seed(594861940L)
data(ou2)
+estnames=c("alpha.2","alpha.3","tau")
+theta.truth <- coef(ou2)
+theta.guess <- theta.truth
+theta.guess[estnames] <- theta.guess[estnames]*1.5
-p.truth <- coef(ou2)
-
m1 <- nlf(
object=ou2,
- start=p.truth,
+ start=theta.truth,
lags=c(4,6),
nconverge=100,
nasymp=2000,
@@ -18,12 +20,13 @@
seed=384886L,
lql.frac = 0.025
)
+se1 <- rep(NA,length(estnames))
print(m1)
m2 <- nlf(
object=ou2,
- start=p.truth,
- est=c("alpha.2","alpha.3","tau"),
+ start=theta.truth,
+ est=estnames,
lags=c(4,6),
nconverge=100,
nasymp=2000,
@@ -35,5 +38,46 @@
lql.frac = 0.025
)
-se <- m2$se
-print(cbind(truth=p.truth[names(se)],fit=m2$params[names(se)],se=se),digits=3)
+se2 <- m2$se
+
+m3 <- nlf(
+ object=ou2,
+ start=theta.guess,
+ lags=c(4,6),
+ nconverge=100,
+ nasymp=2000,
+ method="Nelder-Mead",
+ maxit=500,
+ trace=1,
+ verbose=TRUE,
+ eval.only=TRUE,
+ seed=384886L,
+ lql.frac = 0.025
+ )
+se3 <- rep(NA,length(estnames))
+
+m4 <- nlf(
+ object=ou2,
+ start=theta.guess,
+ est=estnames,
+ lags=c(4,6),
+ nconverge=100,
+ nasymp=2000,
+ method="Nelder-Mead",
+ maxit=500,
+ trace=1,
+ verbose=TRUE,
+ seed=384886L,
+ lql.frac = 0.025
+ )
+
+se4 <- m4$se
+
+print(
+ cbind(
+ guess=c(theta.guess[estnames],se=se3,value=m3),
+ truth=c(theta.truth[estnames],se=se1,value=m1),
+ fit.from.guess=c(m4$params[estnames],se=se4,value=m4$value),
+ fit.from.truth=c(m2$params[estnames],se=se2,value=m2$value)
+ )
+ )
Modified: pkg/tests/ou2-nlf.Rout.save
===================================================================
--- pkg/tests/ou2-nlf.Rout.save 2011-05-22 21:01:49 UTC (rev 500)
+++ pkg/tests/ou2-nlf.Rout.save 2011-05-27 19:42:59 UTC (rev 501)
@@ -1,7 +1,8 @@
-R version 2.11.1 (2010-05-31)
+R version 2.12.1 (2010-12-16)
Copyright (C) 2010 The R Foundation for Statistical Computing
ISBN 3-900051-07-0
+Platform: x86_64-unknown-linux-gnu (64-bit)
R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
@@ -16,16 +17,21 @@
Type 'q()' to quit R.
> library(pomp)
+Loading required package: mvtnorm
+Loading required package: subplex
+Loading required package: deSolve
>
> set.seed(594861940L)
>
> data(ou2)
+> estnames=c("alpha.2","alpha.3","tau")
+> theta.truth <- coef(ou2)
+> theta.guess <- theta.truth
+> theta.guess[estnames] <- theta.guess[estnames]*1.5
>
-> p.truth <- coef(ou2)
->
> m1 <- nlf(
+ object=ou2,
-+ start=p.truth,
++ start=theta.truth,
+ lags=c(4,6),
+ nconverge=100,
+ nasymp=2000,
@@ -35,13 +41,14 @@
+ seed=384886L,
+ lql.frac = 0.025
+ )
+> se1 <- rep(NA,length(estnames))
> print(m1)
-[1] -555.547
+[1] -555.0868
>
> m2 <- nlf(
+ object=ou2,
-+ start=p.truth,
-+ est=c("alpha.2","alpha.3","tau"),
++ start=theta.truth,
++ est=estnames,
+ lags=c(4,6),
+ nconverge=100,
+ nasymp=2000,
@@ -53,76 +60,161 @@
+ lql.frac = 0.025
+ )
Nelder-Mead direct search function minimizer
-function value for initial parameters = 555.547033
- Scaled convergence tolerance is 8.2783e-06
+function value for initial parameters = 555.086814
+ Scaled convergence tolerance is 8.27144e-06
Stepsize computed as 0.100000
-BUILD 4 559.966959 555.547033
-HI-REDUCTION 6 556.890993 555.547033
-HI-REDUCTION 8 556.217248 555.181704
-LO-REDUCTION 10 555.617125 555.174895
-REFLECTION 12 555.547033 555.100277
-HI-REDUCTION 14 555.181704 555.100277
-HI-REDUCTION 16 555.174895 555.060649
-HI-REDUCTION 18 555.158167 555.060649
-HI-REDUCTION 20 555.100277 555.060649
-HI-REDUCTION 22 555.073319 555.046470
-HI-REDUCTION 24 555.064516 555.046470
-HI-REDUCTION 26 555.060649 555.046470
-LO-REDUCTION 28 555.049690 555.044469
-REFLECTION 30 555.047655 555.041210
-HI-REDUCTION 32 555.046470 555.041210
-HI-REDUCTION 34 555.044469 555.040900
-HI-REDUCTION 36 555.041431 555.040037
-HI-REDUCTION 38 555.041210 555.039984
-HI-REDUCTION 40 555.040900 555.039847
-LO-REDUCTION 42 555.040037 555.039534
-HI-REDUCTION 44 555.039984 555.039516
-EXTENSION 46 555.039847 555.038378
-HI-REDUCTION 48 555.039534 555.038378
-LO-REDUCTION 50 555.039516 555.038378
-LO-REDUCTION 52 555.039124 555.038378
-REFLECTION 54 555.038819 555.038336
-EXTENSION 56 555.038601 555.037899
-LO-REDUCTION 58 555.038378 555.037899
-LO-REDUCTION 60 555.038336 555.037899
-REFLECTION 62 555.038103 555.037848
-EXTENSION 64 555.038074 555.037593
-EXTENSION 66 555.037899 555.037354
-EXTENSION 68 555.037848 555.037335
-EXTENSION 70 555.037593 555.036537
-LO-REDUCTION 72 555.037354 555.036537
-LO-REDUCTION 74 555.037335 555.036537
-EXTENSION 76 555.037012 555.036094
-LO-REDUCTION 78 555.036938 555.036094
-EXTENSION 80 555.036537 555.035295
-LO-REDUCTION 82 555.036148 555.035295
-LO-REDUCTION 84 555.036094 555.035295
-REFLECTION 86 555.035627 555.035238
-REFLECTION 88 555.035311 555.035071
-LO-REDUCTION 90 555.035295 555.035071
-HI-REDUCTION 92 555.035238 555.035071
-REFLECTION 94 555.035152 555.035045
-LO-REDUCTION 96 555.035150 555.035045
-REFLECTION 98 555.035071 555.035010
-LO-REDUCTION 100 555.035065 555.035010
-LO-REDUCTION 102 555.035045 555.035010
-LO-REDUCTION 104 555.035022 555.035009
-LO-REDUCTION 106 555.035017 555.035009
+BUILD 4 559.254725 555.086814
+HI-REDUCTION 6 556.109178 555.086814
+HI-REDUCTION 8 555.476443 554.544555
+LO-REDUCTION 10 555.105949 554.544555
+LO-REDUCTION 12 555.086814 554.519288
+LO-REDUCTION 14 554.740190 554.519288
+HI-REDUCTION 16 554.618696 554.519288
+LO-REDUCTION 18 554.548158 554.515101
+HI-REDUCTION 20 554.544555 554.508577
+HI-REDUCTION 22 554.519288 554.497436
+REFLECTION 24 554.515101 554.492388
+LO-REDUCTION 26 554.508577 554.492388
+HI-REDUCTION 28 554.497436 554.492388
+REFLECTION 30 554.494893 554.491176
+REFLECTION 32 554.494183 554.489759
+HI-REDUCTION 34 554.492388 554.489759
+REFLECTION 36 554.491176 554.488524
+HI-REDUCTION 38 554.489976 554.488524
+REFLECTION 40 554.489759 554.486675
+HI-REDUCTION 42 554.488704 554.486675
+REFLECTION 44 554.488524 554.486318
+LO-REDUCTION 46 554.487711 554.486318
+EXTENSION 48 554.486945 554.485619
+LO-REDUCTION 50 554.486675 554.485619
+LO-REDUCTION 52 554.486318 554.485607
+HI-REDUCTION 54 554.485830 554.485607
+REFLECTION 56 554.485734 554.485498
+HI-REDUCTION 58 554.485619 554.485498
+HI-REDUCTION 60 554.485607 554.485495
+LO-REDUCTION 62 554.485535 554.485466
+HI-REDUCTION 64 554.485498 554.485466
+HI-REDUCTION 66 554.485495 554.485466
+HI-REDUCTION 68 554.485479 554.485466
Exiting from Nelder Mead minimizer
- 108 function evaluations used
+ 70 function evaluations used
h in NLF = 0.1
-epsilon in NLF = 0.07620834 0.0830427 1.496106
-Fitted param 1 -5.935533 -5.935533 up in 'nlf'
-Fitted param 1 -5.923885 down in 'NLF'
-Fitted param 2 -5.923836 -5.923836 up in 'nlf'
-Fitted param 2 -5.941826 down in 'NLF'
-Fitted param 3 -5.952844 -5.952844 up in 'nlf'
-Fitted param 3 -5.905505 down in 'NLF'
+epsilon in NLF = 0.07715219 0.08600406 1.29893
+Fitted param 1 -5.928762 -5.928762 up in 'nlf'
+Fitted param 1 -5.918791 down in 'NLF'
+Fitted param 2 -5.92007 -5.92007 up in 'nlf'
+Fitted param 2 -5.933666 down in 'NLF'
+Fitted param 3 -5.937785 -5.937785 up in 'nlf'
+Fitted param 3 -5.903043 down in 'NLF'
>
-> se <- m2$se
-> print(cbind(truth=p.truth[names(se)],fit=m2$params[names(se)],se=se),digits=3)
- truth fit se
-alpha.2 -0.5 -0.465 0.0433
-alpha.3 0.3 0.292 0.0508
-tau 1.0 0.917 0.5458
+> se2 <- m2$se
>
+> m3 <- nlf(
++ object=ou2,
++ start=theta.guess,
++ lags=c(4,6),
++ nconverge=100,
++ nasymp=2000,
++ method="Nelder-Mead",
++ maxit=500,
++ trace=1,
++ verbose=TRUE,
++ eval.only=TRUE,
++ seed=384886L,
++ lql.frac = 0.025
++ )
+> se3 <- rep(NA,length(estnames))
+>
+> m4 <- nlf(
++ object=ou2,
++ start=theta.guess,
++ est=estnames,
++ lags=c(4,6),
++ nconverge=100,
++ nasymp=2000,
++ method="Nelder-Mead",
++ maxit=500,
++ trace=1,
++ verbose=TRUE,
++ seed=384886L,
++ lql.frac = 0.025
++ )
+ Nelder-Mead direct search function minimizer
+function value for initial parameters = 10330.933626
+ Scaled convergence tolerance is 0.000153943
+Stepsize computed as 0.150000
+BUILD 4 29742.647014 576.824636
+REFLECTION 6 10330.933626 564.186576
+EXTENSION 8 10330.933626 556.163486
+LO-REDUCTION 10 576.824636 555.749535
+LO-REDUCTION 12 564.186576 555.749535
+LO-REDUCTION 14 557.909279 555.749535
+REFLECTION 16 556.174151 555.518836
+HI-REDUCTION 18 556.163486 555.234932
+REFLECTION 20 555.749535 555.069014
+HI-REDUCTION 22 555.518836 555.069014
+LO-REDUCTION 24 555.234932 555.041059
+REFLECTION 26 555.161158 554.878667
+LO-REDUCTION 28 555.069014 554.878667
+LO-REDUCTION 30 555.041059 554.820011
+HI-REDUCTION 32 554.902564 554.820011
+HI-REDUCTION 34 554.878667 554.820011
+HI-REDUCTION 36 554.863049 554.818271
+REFLECTION 38 554.823167 554.795312
+LO-REDUCTION 40 554.820011 554.795312
+HI-REDUCTION 42 554.818271 554.795007
+EXTENSION 44 554.796015 554.739129
+LO-REDUCTION 46 554.795312 554.739129
+LO-REDUCTION 48 554.795007 554.739129
+EXTENSION 50 554.757998 554.655046
+LO-REDUCTION 52 554.746394 554.655046
+EXTENSION 54 554.739129 554.611468
+EXTENSION 56 554.680717 554.529970
+LO-REDUCTION 58 554.655046 554.529970
+REFLECTION 60 554.611468 554.508046
+REFLECTION 62 554.533844 554.490991
+LO-REDUCTION 64 554.529970 554.490991
+LO-REDUCTION 66 554.508046 554.490991
+HI-REDUCTION 68 554.499927 554.490991
+HI-REDUCTION 70 554.499566 554.490991
+LO-REDUCTION 72 554.495632 554.490991
+LO-REDUCTION 74 554.494415 554.490991
+EXTENSION 76 554.493222 554.488394
+HI-REDUCTION 78 554.492456 554.488394
+EXTENSION 80 554.490991 554.486164
+LO-REDUCTION 82 554.490490 554.486164
+LO-REDUCTION 84 554.488394 554.486164
+REFLECTION 86 554.486677 554.486008
+REFLECTION 88 554.486468 554.485672
+HI-REDUCTION 90 554.486164 554.485662
+HI-REDUCTION 92 554.486008 554.485580
+Exiting from Nelder Mead minimizer
+ 94 function evaluations used
+h in NLF = 0.1
+epsilon in NLF = 0.07703295 0.08568253 1.296531
+Fitted param 1 -5.928954 -5.928954 up in 'nlf'
+Fitted param 1 -5.918596 down in 'NLF'
+Fitted param 2 -5.919574 -5.919574 up in 'nlf'
+Fitted param 2 -5.933857 down in 'NLF'
+Fitted param 3 -5.937754 -5.937754 up in 'nlf'
+Fitted param 3 -5.903045 down in 'NLF'
+>
+> se4 <- m4$se
+>
+> print(
++ cbind(
++ guess=c(theta.guess[estnames],se=se3,value=m3),
++ truth=c(theta.truth[estnames],se=se1,value=m1),
++ fit.from.guess=c(m4$params[estnames],se=se4,value=m4$value),
++ fit.from.truth=c(m2$params[estnames],se=se2,value=m2$value)
++ )
++ )
+ guess truth fit.from.guess fit.from.truth
+alpha.2 -0.75 -0.5000 -0.45749499 -0.45774692
+alpha.3 0.45 0.3000 0.30150315 0.30194588
+tau 1.50 1.0000 1.09790880 1.09526822
+se1 NA NA 0.04374829 0.04397017
+se2 NA NA 0.05423946 0.05463665
+se3 NA NA 0.54596143 0.54625858
+value -10330.93 -555.0868 -554.48554809 -554.48546595
+>
More information about the pomp-commits
mailing list