[Pomp-commits] r840 - in branches/mif2: . R inst/examples src tests

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Tue Mar 26 12:50:29 CET 2013


Author: kingaa
Date: 2013-03-26 12:50:29 +0100 (Tue, 26 Mar 2013)
New Revision: 840

Modified:
   branches/mif2/DESCRIPTION
   branches/mif2/R/authors.R
   branches/mif2/inst/examples/bbs.R
   branches/mif2/inst/examples/blowflies.R
   branches/mif2/inst/examples/euler.sir.R
   branches/mif2/inst/examples/gillespie.sir.R
   branches/mif2/inst/examples/gompertz.R
   branches/mif2/inst/examples/ou2.R
   branches/mif2/inst/examples/ricker.R
   branches/mif2/inst/examples/rw2.R
   branches/mif2/inst/examples/verhulst.R
   branches/mif2/src/SSA_wrapper.c
   branches/mif2/src/blowfly.c
   branches/mif2/src/dprocess.c
   branches/mif2/src/partrans.c
   branches/mif2/tests/filtfail.R
   branches/mif2/tests/filtfail.Rout.save
   branches/mif2/tests/ou2-mif2.R
   branches/mif2/tests/ou2-mif2.Rout.save
Log:
- some updates from main branch of pomp


Modified: branches/mif2/DESCRIPTION
===================================================================
--- branches/mif2/DESCRIPTION	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/DESCRIPTION	2013-03-26 11:50:29 UTC (rev 840)
@@ -2,9 +2,20 @@
 Type: Package
 Title: Statistical inference for partially observed Markov processes
 Version: 0.44-1
-Date: 2013-02-26
-Author: Aaron A. King, Edward L. Ionides, Carles Breto, Steve Ellner, Bruce Kendall, Helen Wearing, Matthew J. Ferrari, Michael Lavine, Daniel C. Reuman
+Date: 2013-03-26
+Author: Aaron A. King, Edward L. Ionides, Carles Breto, Steve Ellner, Bruce Kendall, Dao Nguyen, Helen Wearing, Matthew J. Ferrari, Michael Lavine, Daniel C. Reuman
 Maintainer: Aaron A. King <kingaa at umich.edu>
+Author at R: c(person(given=c("Aaron","A."),family="King",role=c("aut","cre"),email="kingaa at umich.edu"),
+	  person(given=c("Edward","L."),family="Ionides",role=c("aut")),
+	  person(given=c("Carles"),family="Breto",role=c("aut")),
+	  person(given=c("Stephen","P."),family="Ellner",role=c("ctb")),
+	  person(given=c("Matthew","J."),family="Ferrari",role=c("ctb")),
+	  person(given=c("Bruce","E."),family="Kendall",role=c("ctb")),
+	  person(given=c("Michael"),family="Lavine",role=c("ctb")),
+	  person(given=c("Dao"),family="Nguyen",role=c("ctb")),
+	  person(given=c("Daniel","C."),family="Reuman",role=c("ctb")),
+	  person(given=c("Helen"),family="Wearing",role=c("ctb")),
+	  person(given=c("Simon","N."),family="Wood",role=c("ctb")))
 URL: http://pomp.r-forge.r-project.org
 Description: Inference methods for partially-observed Markov processes
 Depends: R(>= 2.14.1), stats, methods, graphics, mvtnorm, subplex, deSolve

Modified: branches/mif2/R/authors.R
===================================================================
--- branches/mif2/R/authors.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/R/authors.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -1,9 +1,9 @@
 list(
      aak=person(given=c("Aaron","A."),family="King",role=c("aut","cre"),email="kingaa at umich.edu"),
-     eli=person(given=c("Edward","L."),family="Ionides",role=c("aut"),email="ionides at umich.edu"),
-     cb=person(given=c("Carles"),family="Breto",role=c("aut")),
-     spe=person(given=c("Stephen","P."),family="Ellner",role=c("aut")),
-     bek=person(given=c("Bruce","E."),family="Kendall",role=c("aut")),
+     eli=person(given=c("Edward","L."),family="Ionides",role=c("ctb")),
+     cb=person(given=c("Carles"),family="Breto",role=c("ctb")),
+     spe=person(given=c("Stephen","P."),family="Ellner",role=c("ctb")),
+     bek=person(given=c("Bruce","E."),family="Kendall",role=c("ctb")),
      mf=person(given=c("Matthew","J."),family="Ferrari",role=c("ctb")),
      ml=person(given=c("Michael"),family="Lavine",role=c("ctb")),
      dcr=person(given=c("Daniel","C."),family="Reuman",role=c("ctb")),

Modified: branches/mif2/inst/examples/bbs.R
===================================================================
--- branches/mif2/inst/examples/bbs.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/bbs.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -1,6 +1,6 @@
 require(pomp)
 
-tc <- textConnection("
+flu <- read.csv2(text="
 day;reports
 1;3
 2;8
@@ -18,9 +18,6 @@
 14;5
 ")
 
-flu <- read.csv2(file=tc)
-close(tc)
-
 po <- pomp(
            data=flu,
            times="day",

Modified: branches/mif2/inst/examples/blowflies.R
===================================================================
--- branches/mif2/inst/examples/blowflies.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/blowflies.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -868,14 +868,13 @@
 718;8103;4
 720;6803;4
 '
-
 raw.data <- subset(
-                   read.csv2(textConnection(blowfly.data),comment.char="#"),
-                   set==4
+                   read.csv2(text=blowfly.data,comment.char="#"),
+                   set==4,
+                   select=-set
                    )
-
 pomp(
-     data=subset(raw.data[c("day","y")],day>14&day<400),
+     data=subset(raw.data,day>14&day<400),
      times="day",
      t0=14,
      rprocess=discrete.time.sim(

Modified: branches/mif2/inst/examples/euler.sir.R
===================================================================
--- branches/mif2/inst/examples/euler.sir.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/euler.sir.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -213,7 +213,7 @@
 
 
 po <- pomp(
-           data=read.csv2(textConnection(dat)),
+           data=read.csv2(text=dat),
            times="time",
            t0=0,
            params=c(

Modified: branches/mif2/inst/examples/gillespie.sir.R
===================================================================
--- branches/mif2/inst/examples/gillespie.sir.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/gillespie.sir.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -525,7 +525,7 @@
 '
 
 pomp(
-     data=read.csv2(textConnection(dat)),
+     data=read.csv2(text=dat),
      times="time",
      t0=0,
      params=c(

Modified: branches/mif2/inst/examples/gompertz.R
===================================================================
--- branches/mif2/inst/examples/gompertz.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/gompertz.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -105,7 +105,7 @@
 '
 
 po <- pomp(
-           data=read.csv2(textConnection(dat)),
+           data=read.csv2(text=dat),
            times="time",
            t0=0,
            params=c(K=1,r=0.1,sigma=0.1,tau=0.1,X.0=1),

Modified: branches/mif2/inst/examples/ou2.R
===================================================================
--- branches/mif2/inst/examples/ou2.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/ou2.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -104,7 +104,7 @@
 '
 
 pomp( 
-     data=read.csv2(textConnection(dat)),
+     data=read.csv2(text=dat),
      times="time",
      t0=0,
      rprocess=discrete.time.sim("ou2_step",PACKAGE="pomp"),

Modified: branches/mif2/inst/examples/ricker.R
===================================================================
--- branches/mif2/inst/examples/ricker.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/ricker.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -55,7 +55,7 @@
 '
 
 pomp(
-     data=read.csv2(textConnection(dat)),
+     data=read.csv2(text=dat),
      times="time",
      t0=0,
      params=c(r=exp(3.8),sigma=0.3,phi=10,N.0=7,e.0=0), # originally used to generate the data

Modified: branches/mif2/inst/examples/rw2.R
===================================================================
--- branches/mif2/inst/examples/rw2.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/rw2.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -104,7 +104,7 @@
 '
 
 pomp(
-     data=read.csv2(textConnection(dat)),
+     data=read.csv2(text=dat),
      times="time",
      t0=0,
      params=c(x1.0=0,x2.0=0,s1=1,s2=3,tau=1), # parameters at which data were generated

Modified: branches/mif2/inst/examples/verhulst.R
===================================================================
--- branches/mif2/inst/examples/verhulst.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/inst/examples/verhulst.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -1004,7 +1004,7 @@
 '
 
 pomp(
-     data=read.csv2(textConnection(dat)),
+     data=read.csv2(text=dat),
      times="time",
      t0=0,
      params=c(n.0=10000,K=10000,r=0.9,sigma=0.4,tau=0.1),

Modified: branches/mif2/src/SSA_wrapper.c
===================================================================
--- branches/mif2/src/SSA_wrapper.c	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/src/SSA_wrapper.c	2013-03-26 11:50:29 UTC (rev 840)
@@ -3,6 +3,10 @@
 #include "pomp_internal.h"
 #include <R_ext/Constants.h>
 
+typedef double (*_pomp_rxnrate) (const int *j, const double *t, const double *x, const double *p,
+			   int *stateindex, int *parindex, int *covarindex, 
+			   int *ncovar, double *covar);
+
 void F77_SUB(rndstart)(void) { GetRNGstate(); }
 void F77_SUB(rndend)(void) { PutRNGstate(); }
 double F77_SUB(unifrnd)(void) { return unif_rand(); }
@@ -11,7 +15,7 @@
 double F77_SUB(gammarnd)(double shape, double scale) { return rgamma(shape,scale); }
 void F77_SUB(multinomrnd)(int N, double *p, int ncat, int *ix) { rmultinom(N,p,ncat,ix); }
 
-void F77_NAME(driverssa)(void *fprob, int *nvar, int *nevent, int *npar, int *nreps, int *ntimes, 
+void F77_NAME(driverssa)(_pomp_rxnrate *fprob, int *nvar, int *nevent, int *npar, int *nreps, int *ntimes, 
 			 int *kflag, double *xstart, double *times, double *params, double *xout,
 			 double *e, double *v, double *d, int *nzero, int *izero, int *istate, 
 			 int *ipar, int *ncov, int *icov, int *lcov, int *mcov, double *tcov, double *cov);

Modified: branches/mif2/src/blowfly.c
===================================================================
--- branches/mif2/src/blowfly.c	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/src/blowfly.c	2013-03-26 11:50:29 UTC (rev 840)
@@ -5,7 +5,7 @@
 #include "pomp.h"
 
 #define P       (p[parindex[0]]) // growth rate
-#define NZERO   (p[parindex[1]]) // density-dependence parameter
+#define N0      (p[parindex[1]]) // density-dependence parameter
 #define DELTA   (p[parindex[2]]) // survival parameter
 #define SIGMAP  (p[parindex[3]]) // recruitment noise SD
 #define SIGMAD  (p[parindex[4]]) // survivorship noise SD
@@ -29,7 +29,7 @@
   double eps = rgammawn(SIGMAD,dt)/dt;
   int k;
 
-  R = rpois(P*N[tau]*exp(-N[tau]/NZERO)*dt*e);
+  R = rpois(P*N[tau]*exp(-N[tau]/N0)*dt*e);
   S = rbinom(N[0],exp(-DELTA*dt*eps));
   E = e;
   EPS = eps;
@@ -60,7 +60,7 @@
 #undef EPS
 
 #undef P
-#undef NZERO
+#undef N0
 #undef DELTA
 #undef SIGMAP
 #undef SIGMAD

Modified: branches/mif2/src/dprocess.c
===================================================================
--- branches/mif2/src/dprocess.c	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/src/dprocess.c	2013-03-26 11:50:29 UTC (rev 840)
@@ -12,7 +12,7 @@
   int nprotect = 0;
   int *xdim, npars, nvars, nreps, nrepsx, ntimes;
   SEXP X, fn, fcall, rho;
-  SEXP dimP, dimF;
+  SEXP dimF;
 
   PROTECT(gnsi = duplicate(gnsi)); nprotect++;
 

Modified: branches/mif2/src/partrans.c
===================================================================
--- branches/mif2/src/partrans.c	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/src/partrans.c	2013-03-26 11:50:29 UTC (rev 840)
@@ -12,7 +12,7 @@
   int nprotect = 0;
   SEXP fn, fcall, rho, ans, nm;
   SEXP pdim, pvec;
-  SEXP tparams;
+  SEXP tparams = R_NilValue;
   int mode = -1;
   char direc;
   int qmat;

Modified: branches/mif2/tests/filtfail.R
===================================================================
--- branches/mif2/tests/filtfail.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/tests/filtfail.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -5,7 +5,8 @@
 ### the following example tests to make sure that states are updated properly
 ### upon filtering failures
 
-"time,admissions,discharges,patients,cases
+records <- read.csv(text="
+time,admissions,discharges,patients,cases
 0,4,2,8,
 1,0,1,10,2
 2,2,0,9,1
@@ -27,12 +28,8 @@
 18,4,0,7,1
 19,0,0,11,0
 20,1,4,11,
-" -> csvtext
+")
 
-tc <- textConnection(csvtext)
-records <- read.csv(tc)
-close(tc)
-
 po <- pomp(
            data=subset(records[c("time","cases")],!is.na(cases)),
            times="time",

Modified: branches/mif2/tests/filtfail.Rout.save
===================================================================
--- branches/mif2/tests/filtfail.Rout.save	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/tests/filtfail.Rout.save	2013-03-26 11:50:29 UTC (rev 840)
@@ -1,6 +1,6 @@
 
-R version 2.15.2 (2012-10-26) -- "Trick or Treat"
-Copyright (C) 2012 The R Foundation for Statistical Computing
+R version 2.15.3 (2013-03-01) -- "Security Blanket"
+Copyright (C) 2013 The R Foundation for Statistical Computing
 ISBN 3-900051-07-0
 Platform: x86_64-unknown-linux-gnu (64-bit)
 
@@ -26,7 +26,8 @@
 > ### the following example tests to make sure that states are updated properly
 > ### upon filtering failures
 > 
-> "time,admissions,discharges,patients,cases
+> records <- read.csv(text="
++ time,admissions,discharges,patients,cases
 + 0,4,2,8,
 + 1,0,1,10,2
 + 2,2,0,9,1
@@ -48,12 +49,8 @@
 + 18,4,0,7,1
 + 19,0,0,11,0
 + 20,1,4,11,
-+ " -> csvtext
++ ")
 > 
-> tc <- textConnection(csvtext)
-> records <- read.csv(tc)
-> close(tc)
-> 
 > po <- pomp(
 +            data=subset(records[c("time","cases")],!is.na(cases)),
 +            times="time",
@@ -121,4 +118,4 @@
 > 
 > proc.time()
    user  system elapsed 
-  0.456   0.052   0.526 
+  0.464   0.044   0.524 

Modified: branches/mif2/tests/ou2-mif2.R
===================================================================
--- branches/mif2/tests/ou2-mif2.R	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/tests/ou2-mif2.R	2013-03-26 11:50:29 UTC (rev 840)
@@ -76,9 +76,26 @@
              )  
 mif2b <- continue(mif2b,Nmif=50)
 
+mif2c <- mif(ou2,Nmif=50,start=guess1,
+             pars=c('alpha.2','alpha.3'),ivps=c('x1.0','x2.0'),
+             rw.sd=c(
+               x1.0=0.5,x2.0=.5,
+               alpha.2=0.1,alpha.3=0.1),
+             transform=F,
+             Np=1000,
+             var.factor=1,
+             cooling.type="hyperbolic",
+             cooling.fraction=0.05,
+             max.fail=100,
+             method="mif2"
+             )  
+mif2c <- continue(mif2c,Nmif=50)
+
 compare.mif(list(mif1b,mif2b))
 
 compare.mif(list(mif1a,mif1b))
 compare.mif(list(mif2a,mif2b))
 
+compare.mif(list(mif1b,mif2c))
+
 dev.off()

Modified: branches/mif2/tests/ou2-mif2.Rout.save
===================================================================
--- branches/mif2/tests/ou2-mif2.Rout.save	2013-03-21 14:35:10 UTC (rev 839)
+++ branches/mif2/tests/ou2-mif2.Rout.save	2013-03-26 11:50:29 UTC (rev 840)
@@ -1,6 +1,6 @@
 
-R version 2.15.2 (2012-10-26) -- "Trick or Treat"
-Copyright (C) 2012 The R Foundation for Statistical Computing
+R version 2.15.3 (2013-03-01) -- "Security Blanket"
+Copyright (C) 2013 The R Foundation for Statistical Computing
 ISBN 3-900051-07-0
 Platform: x86_64-unknown-linux-gnu (64-bit)
 
@@ -102,15 +102,32 @@
 See '?mif' for instructions on specifying the cooling schedule. 
 > mif2b <- continue(mif2b,Nmif=50)
 > 
+> mif2c <- mif(ou2,Nmif=50,start=guess1,
++              pars=c('alpha.2','alpha.3'),ivps=c('x1.0','x2.0'),
++              rw.sd=c(
++                x1.0=0.5,x2.0=.5,
++                alpha.2=0.1,alpha.3=0.1),
++              transform=F,
++              Np=1000,
++              var.factor=1,
++              cooling.type="hyperbolic",
++              cooling.fraction=0.05,
++              max.fail=100,
++              method="mif2"
++              )  
+> mif2c <- continue(mif2c,Nmif=50)
+> 
 > compare.mif(list(mif1b,mif2b))
 > 
 > compare.mif(list(mif1a,mif1b))
 > compare.mif(list(mif2a,mif2b))
 > 
+> compare.mif(list(mif1b,mif2c))
+> 
 > dev.off()
 null device 
           1 
 > 
 > proc.time()
    user  system elapsed 
- 42.138   0.064  42.548 
+ 49.343   0.076  49.753 



More information about the pomp-commits mailing list