[R] problem running a function

garciap garciap at usal.es
Sat Mar 19 19:58:31 CET 2011


Dear people,

I'm trying to do some analysis of a data using the models by Royle & Donazio
in their fantastic book, particular the following function:
http://www.mbr-pwrc.usgs.gov/pubanalysis/roylebook/panel4pt1.fn

that applied to my data and in the console is as follows:
> `desman.y` <- structure(c(3L,4L,3L,2L,1L), .Names = c("1", "2", "3", "4",
> "5"))
> 
> RN<-function(y=desman.y,J=5,nsites=39,Nmax=100){
+    ####
+     lik<-function(parms){
+          r<-expit(parms[1])
+          lambda<-exp(parms[2])
+          pvec<-1-(1-r)^(0:Nmax)
+          gN<-dpois(0:Nmax,lambda)
+          gN<-gN/sum(gN)
+          lik<-rep(NA,nsites)
+          for(i in 1:nsites){
+             lik[i]<-sum(dbinom(y[i],J,pvec)*gN)
+          }
+         -1*sum(log(lik))
+        }
+     ####
+ 
+   tmp<-nlm(lik,c(0,0),hessian=TRUE)
+   ests<-tmp$estimate
+   aic<-tmp$minimum*2 + 2*length(ests)
+   se<- sqrt(diag(solve(tmp$hessian)))
+   list(ests=ests,se=se,aic=aic)
+   }

but when I try to see the results (i.e., the list in the function), I've the
follow errors:

> summary(RN)
Error in object[[i]] : object of type 'closure' is not subsettable

or obtains something like:

> print(RN)
function(y=desman.y,J=5,nsites=39,Nmax=100){
   ####
    lik<-function(parms){
         r<-expit(parms[1])
         lambda<-exp(parms[2])
         pvec<-1-(1-r)^(0:Nmax)
         gN<-dpois(0:Nmax,lambda)
         gN<-gN/sum(gN)
         lik<-rep(NA,nsites)
         for(i in 1:nsites){
            lik[i]<-sum(dbinom(y[i],J,pvec)*gN)
         }
        -1*sum(log(lik))
       }
    ####

  tmp<-nlm(lik,c(0,0),hessian=TRUE)
  ests<-tmp$estimate
  aic<-tmp$minimum*2 + 2*length(ests)
  se<- sqrt(diag(solve(tmp$hessian)))
  list(ests=ests,se=se,aic=aic)
  }
> 

this second one seems like the analyses are not done at all.

Any suggestion? Many thanks,

Pablo

--
View this message in context: http://r.789695.n4.nabble.com/problem-running-a-function-tp3390117p3390117.html
Sent from the R help mailing list archive at Nabble.com.



More information about the R-help mailing list