2
$\begingroup$

I am stumped as to how to find the one-sided LRT for the Poisson distribution. The majority of the examples I am finding are for two-sided tests, which does not give all that much help. If you have links to any websites, those would be very useful.

Here

$\endgroup$
1
  • $\begingroup$ The procedure remains same, i.e. forming the likelihood ratio statistic $\Lambda(x_1,\ldots,x_n)=\frac{\sup_{\theta\le \theta_0}L(\theta\mid x_1,\ldots,x_n)}{\sup_{\theta}L(\theta\mid x_1,\ldots,x_n)}$ where $L$ is the likelihood and rejecting $H_0$ for small values of $\Lambda$. So first find the restricted and unrestricted maximum likelihood estimators of $\theta$. $\endgroup$ Mar 28 '20 at 7:04
4
$\begingroup$

If anything, the one-sided case is simpler than the two-sided case, since the evidence in favour of the alternative hypothesis is only in one direction with respect to the sufficient statistic. Let's work through the problem to have a look a the test in this case.

Since $X_1,...,X_n \sim \text{IID Pois}(\theta)$ you can easily establish that the log-likelihood function is:

$$\ell_\mathbf{x}(\theta) = T(\mathbf{x}) \log(\theta) -n \theta \quad \quad\ \quad \text{for all } \theta>0,$$

where $T(\mathbf{x}) \equiv \sum x_i$ is the sample total, which is a sufficient statistic for the parameter $\theta$. Moreover, using a well-known result for the Poisson distribution, it can be shown that $T(\mathbf{X}) \sim \text{Pois}(n\theta)$.


Now, suppose you have hypotheses $H_0: \theta \leqslant \theta_0$ and $H_1: \theta > \theta_0$. Then from the Neyman-Pearson lemma we know that the UMP test is the likelihood-ratio test, which in this case is equivalent to a monotone test on the sample mean. Since the alternative hypothesis is for a larger rate parameter, a higher sample mean constitutes evidence in favour of the alternative hypothesis. Hence, the p-value function for the test is:

$$\begin{equation} \begin{aligned} p(\mathbf{x}) \equiv \sup_{H_0} \mathbb{P}(T(\mathbf{X}) \geqslant T(\mathbf{x})|\theta) &= \sup_{\theta \leqslant \theta_0} \mathbb{P}(T \geqslant T(\mathbf{x})|T \sim \text{Pois}(n\theta)) \\[6pt] &= \mathbb{P}(T \geqslant T(\mathbf{x})|T \sim \text{Pois}(n\theta_0)) \\[12pt] &= 1-\mathbb{P}(T < T(\mathbf{x})|T \sim \text{Pois}(n\theta_0)) \\[6pt] &= 1-\sum_{t = 0}^{T(\mathbf{x})-1} \frac{\theta_0^x e^{-\theta_0}}{x!}. \\[6pt] \end{aligned} \end{equation}$$

The p-value function fully determines the test, and so this constitutes a full specification of the test. As can be seen, the p-value is just the upper-tail of the Poisson distribution with parameter $n \theta_0$.


Implementing the test: It is relatively simple to program this test in R using the standard method for programming a hypothesis test. Here we will program a more general version of the test that allows you to stipulate the alternative to be greater or less than the null. (For brevity, we will not include the two-sided test.)

Poisson.rate.test <- function(X, theta0, alternative = "greater") {

#Check validity of inputs
if(!is.numeric(X))      { stop("Error: Data should be numeric"); }
if(length(X) == 0)      { stop("Error: You require at least one observation"); }
if(all(X != as.integer(X)))  { stop("Error: Data should be non-negative integers"); }
if(min(X) < 0)          { stop("Error: Data should be non-negative integers"); }
if(!is.numeric(theta0)) { stop("Error: Parameter theta0 should be numeric"); }
if(length(theta0) != 1) { stop("Error: Parameter theta0 should be a scalar"); }
if(theta0 <= 0)         { stop("Error: Parameter theta0 should be positive"); }
if(!(alternative %in% c("greater", "less"))) 
                        { stop("Error: Alternative must be 'greater' or 'less'"); }

#Set description of test and data
method      <- "Poisson rate test (one-sided)";
data.name   <- paste0(deparse(substitute(X)));

#Set null hypothesis value
null.value  <- theta0;
attr(null.value, "names") <- "rate parameter";

#Calculate test statistics
n           <- length(X);
estimate    <- mean(X);
attr(estimate, "names") <- "estimated rate";
statistic   <- n*estimate;
attr(statistic, "names") <- "sample total";

#Calculate p-value
if (alternative == "less") { 
  p.value   <- ppois(statistic, lambda = n*null.value, 
                     lower.tail = TRUE,  log.p = FALSE); }
if (alternative == "greater") { 
  p.value   <- ppois(statistic, lambda = n*null.value, 
                     lower.tail = FALSE, log.p = FALSE); }
attr(p.value, "names") <- NULL;

#Create htest object
TEST        <- list(method = method, data.name = data.name,
                    null.value = null.value, alternative = alternative,
                    estimate = estimate, statistic = statistic, p.value = p.value);
class(TEST) <- "htest";
TEST; }

We can implement this with some data as an example. I will create a data vector with $n=40$ values that I generated from a Poisson distribution (with a parameter value that I will keep secret). We will test to see if the parameter I used is greater than $\theta_0=15$.

#Create a data vector
DATA <- c(21L, 25L, 13L, 21L,  8L,  9L, 14L, 12L, 14L, 14L, 
          11L, 16L, 19L, 20L, 16L,  9L, 13L, 16L, 17L, 23L, 
          14L, 16L, 16L, 14L, 17L, 14L, 12L,  9L,  8L, 16L, 
          20L, 11L, 10L,  4L, 18L, 15L, 15L, 14L, 16L, 20L);

#Test the alternative hypothesis that theta > 15
TEST <- Poisson.rate.test(DATA, theta0 = 15, alternative = 'greater');
TEST;

        Poisson rate test (one-sided)

data:  DATA
sample total = 590, p-value = 0.6488
alternative hypothesis: true rate parameter is greater than 15
sample estimates:
estimated rate 
         14.75

We can see in this case that there is no evidence that the rate parameter for this data was greater than the hypothesised value of $\theta_0=15$.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.