[R] Goodness-of-fit for zero-truncated poisson distribution

From: Roger Prins <rogerprins_at_gmail.com>
Date: Wed, 25 Jun 2008 19:53:05 +0200


Hi there,

I am trying to write a function to perform GOF test of data to a zero-truncated Poisson distribution. I am facing 2 problems.

  1. How can I obtain a frequency table for all values within the range of observed values?

For instance if the observations are
obs <- c("A", "A", "A", "A", "B", "C", "C", "D", "E", "E", "F", "G",

         "G", "H", "H", "H", "H")
## counts how many times indiviuals are counted 1, 2, 3, 4 times and so (Fi <- table(table(obs)))

Fi
1 2 4
3 3 2

meaning that 3 individuals have been counted once, 3 counted twice and 2 counted 4 times. Applying table() here does not return the frequency of individuals counted 3 times (which is 0 in this case). How can I achieve such a contingency table where I could end up with something like

Fi
1 2 3 4
3 3 0 2

2) function to estimate lambda with MLE

trunpoismle <- function(xbar, tol = .Machine$double.eps^0.25, nmax = 10000) {

    lambda <- xbar
    i <- 1
    diff <- 1
    while (abs(diff) > tol && i < nmax)
{

        diff <- (lambda - xbar * (1 - exp(-lambda))) /
                (1 - xbar * exp(-lambda))
        lambda <- lambda - diff
        i <- i + 1

    }
    if (i == nmax)
{

        print("Iteration did not converge")     }
    return(lambda)
}

3) Once I get there, I would like to apply a chi-square with (k-2) df, k being the number of classes. Here the difficulty is to reduce the number of classes so that a chi-square test may be applied. For instance:

Fi <- c(8, 5, 0, 4, 0, 1, 2, 0, 0, 1)
muHat <- trunpoismle(mean(Fi))
class <- 1 : 10
Exp <- muHat^class * exp(-muHat) / factorial(class) Exp <- Exp[-1] / (1- Exp[1]) * sum(Fi)
# individuals seen 1, 2, 3, 4, 5, 6, 7, 8, 9, and 10 times # This distribution maybe reduced to:

FiPooled <- c(8, 5, 4, 4)
# individuals seen 1, 2, 3-4 and > 4 times to apply chi-square tests

I would like to collapse the expected values the same way. Does anyone have some tips or point out on a function I could have missed to pool classes of a contingency table (or even to perfom the GOF test itself).

Thanks,

Roger

        [[alternative HTML version deleted]]



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. Received on Wed 25 Jun 2008 - 19:26:11 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Wed 25 Jun 2008 - 19:31:50 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive