R help - Feb 2011 - question mle again

A few day ago, I was looking for an answer to my question but didn't
get one. Anybody who can help now?

Hello,


I tried to use mle to fit a distribution(zero-inflated negbin for
count data). My call is very simple:

mle(ll)

ll() takes the three parameters, I'd like to be estimated (size, mu
and prob). But within the ll() function I have to judge if the current
parameter-set gives a nice fit or not. So I have to apply them to
observation data. But how does the method know about my observed data?
The mle()-examples define this data outside of this method and it
works. For a simple example, it was fine but when it comes to a loop
(tapply) providing different sets of observation data, it doesn't work
anymore. I'm confused - is there any way to do better?

Here is a little example which show my problem:

# R-code ---------------------------------

lambda.data <- runif(10,0.5,10)

ll <- function(lambda = 1) {
       cat("x in ll()",x,"\n")
       y.fit <- dpois(x, lambda)

       sum( (y - y.fit)^2 )

       }

lapply(1:10, FUN = function(x){

       raw.data <- rpois(100,lambda.data[x])

       freqTab <- count(raw.data)
       x <- freqTab$x
       y <- freqTab$freq / sum(freqTab$freq)
       cat("x in lapply", x,"\n")
       fit <- mle(ll)

       coef(fit)
       })

Can anybody help?

Antje

Hi,

On Mon, Feb 7, 2011 at 8:15 AM, Antje Niederlein
<niederlein-rstat at yahoo.de> wrote:
> A few day ago, I was looking for an answer to my question but didn't
> get one. Anybody who can help now?
>
> Hello,
>
>
> I tried to use mle to fit a distribution(zero-inflated negbin for
> count data). My call is very simple:
>
> mle(ll)
>
> ll() takes the three parameters, I'd like to be estimated (size, mu
> and prob). But within the ll() function I have to judge if the current
> parameter-set gives a nice fit or not. So I have to apply them to
> observation data. But how does the method know about my observed data?
> The mle()-examples define this data outside of this method and it
> works. For a simple example, it was fine but when it comes to a loop
> (tapply) providing different sets of observation data, it doesn't work
> anymore. I'm confused - is there any way to do better?
When a function cannot find a variable inside its own environment, it
will look to its parent environment.  If you define a function in the
global environment, the global environment is its parent environment.
However, if you define a function in the global environment, but then
proceed to use lapply() with another function, the actual variable
ll() needs to access is neither passed to II (so it is not in its
environment) nor is it in the global environment (II's parent
environment).  It is in the function in lapply's environment, which is
inaccessible to II.  I have made some small changes to your code that
gets around this, but I am still not convinced this is really doing
what you want, but that is a whole other question/problem.

Also, for future reference, you are more likely to get a response/help
if you mention the required packages.  I made educated guesses, that
you are using mle() from "stats4" and count() from "plyr" (I
realize
you may not even be aware that those functions came from non-default
loading packages).

HTH,

Josh

Here are my edits to your code:

foo <- function(x) {
  ## load required packages (guessing here)
  require(stats4)
  require(plyr)

  ## define ll function _inside_ foo
  ## this is important if you want it to have access
  ## to arguments in foo
  ll <- function(lambda = 1) {
    cat("x in ll()", x, "\n")
    y.fit <- dpois(x, lambda)
    sum( (y - y.fit)^2 )
  }

  ## Your calculations
  ## (though I'm not convinced this is what you really want)
  raw.data <- rpois(100, lambda.data[x])
  freqTab <- count(raw.data)
  x <- freqTab$x
  y <- freqTab$freq / sum(freqTab$freq)

  cat("x in lapply", x, "\n")
  fit <- mle(ll)
  coef(fit)
}

## Data
lambda.data <- runif(10, 0.5, 10)
## Run it through lapply for x = 1:10
lapply(1:10, FUN = foo)

>
> Here is a little example which show my problem:
>
> # R-code ---------------------------------
>
> lambda.data <- runif(10,0.5,10)
>
> ll <- function(lambda = 1) {
> ? ? ? cat("x in ll()",x,"\n")
> ? ? ? y.fit <- dpois(x, lambda)
>
> ? ? ? sum( (y - y.fit)^2 )
>
> ? ? ? }
>
> lapply(1:10, FUN = function(x){
>
> ? ? ? raw.data <- rpois(100,lambda.data[x])
>
> ? ? ? freqTab <- count(raw.data)
> ? ? ? x <- freqTab$x
> ? ? ? y <- freqTab$freq / sum(freqTab$freq)
> ? ? ? cat("x in lapply", x,"\n")
> ? ? ? fit <- mle(ll)
>
> ? ? ? coef(fit)
> ? ? ? })
>
> Can anybody help?
>
> Antje
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>


-- 
Joshua Wiley
Ph.D. Student, Health Psychology
University of California, Los Angeles
http://www.joshuawiley.com/

Antje Niederlein <niederlein-rstat <at> yahoo.de> writes:
> 
> A few day ago, I was looking for an answer to my question but didn't
> get one. Anybody who can help now?
> 
> Hello,
> 
> I tried to use mle to fit a distribution(zero-inflated negbin for
> count data). My call is very simple:
> 
  I will point out that this is one of the reasons I wrote mle2
(in the bbmle package), which differs from mle in taking an explicit
'data' argument.  I *think* the following does what you want (although
I admit I haven't looked at the output closely):

library(plyr)
library(bbmle)
lambda.data <- runif(10,0.5,10)

ll <- function(lambda = 1) {
  cat("x in ll()",x,"\n")
  y.fit <- dpois(x, lambda)
  sum( (y - y.fit)^2 )
}

lapply(1:10, FUN = function(x){

  raw.data <- rpois(100,lambda.data[x])
  
  freqTab <- count(raw.data)
  x <- freqTab$x
  y <- freqTab$freq / sum(freqTab$freq)
  cat("x in lapply", x,"\n")
  fit <- mle2(ll,data=data.frame(x,y))

  coef(fit)
})