GLM_POISSON

Overview

The GLM_POISSON function fits a Generalized Linear Model (GLM) with a Poisson distribution family, designed specifically for modeling count data. Poisson regression is appropriate when the response variable represents counts of events—such as the number of customer arrivals, defect occurrences, or species observations—where values are non-negative integers.

This implementation uses the statsmodels library’s GLM framework. For complete documentation, see the statsmodels GLM reference and the GLM class API. The source code is available on GitHub.

The Poisson distribution assumes that the mean and variance of the response variable are equal. The model relates the expected count \mu to predictor variables through a link function. The canonical link for Poisson regression is the log link:

\log(\mu_i) = \beta_0 + \beta_1 x_{i1} + \beta_2 x_{i2} + \cdots + \beta_p x_{ip}

This function supports three link functions: log (default), identity, and sqrt. The log link ensures predicted counts are always positive and provides easily interpretable coefficients—each unit increase in a predictor multiplies the expected count by e^{\beta}, known as the incidence rate ratio (IRR).

Parameters are estimated via maximum likelihood using iteratively reweighted least squares (IRLS). The function returns coefficient estimates, standard errors, z-statistics, p-values, and confidence intervals. Model fit is assessed using deviance, Pearson chi-squared, AIC, BIC, and log-likelihood statistics. The incidence rate ratio for each coefficient is also computed as \text{IRR} = e^{\beta}.

A key assumption of Poisson regression is equidispersion—that the variance equals the mean. When observed variance exceeds the mean (overdispersion), consider using negative binomial regression or quasi-Poisson methods instead. The Pearson chi-squared statistic divided by residual degrees of freedom provides a rough overdispersion diagnostic; values substantially greater than 1 suggest overdispersion.

This example function is provided as-is without any representation of accuracy.

Excel Usage

=GLM_POISSON(y, x, glm_poisson_link, fit_intercept, alpha)

y (list[list], required): Column vector of count data (non-negative values). Each element represents an observation of the dependent variable.
x (list[list], required): Matrix of predictor variables. Each column represents a different predictor variable; each row represents an observation.
glm_poisson_link (str, optional, default: “log”): The link function to use for the Poisson GLM model.
fit_intercept (bool, optional, default: true): If True, adds an intercept term to the model.
alpha (float, optional, default: 0.05): Significance level for confidence intervals. Must be between 0 and 1.

Returns (list[list]): 2D list with GLM results and statistics, or error string.

Examples

Example 1: Poisson log link with single predictor

Inputs:

y	x
2	1
3	2
5	3
7	4
11	5

Excel formula:

=GLM_POISSON({2;3;5;7;11}, {1;2;3;4;5})

Expected output:

parameter	coefficient	std_error	z_statistic	p_value	ci_lower	ci_upper	incidence_rate_ratio
intercept	0.2836	0.5951	0.4766	0.6337	-0.8828	1.45	1.328
x1	0.4224	0.1491	2.833	0.004604	0.1302	0.7146	1.526
deviance	0.0252
pearson_chi2	0.02546
aic	21.17
bic	20.39
log_likelihood	-8.585

Example 2: Poisson identity link with single predictor

Inputs:

y	x	glm_poisson_link
5	1	identity
8	2
10	3
15	4

Excel formula:

=GLM_POISSON({5;8;10;15}, {1;2;3;4}, "identity")

Expected output:

parameter	coefficient	std_error	z_statistic	p_value	ci_lower	ci_upper	incidence_rate_ratio
intercept	1.775	3.055	0.581	0.5612	-4.212	7.762	5.9
x1	3.09	1.316	2.349	0.01885	0.5113	5.669	21.98
deviance	0.158
pearson_chi2	0.1558
aic	20.29
bic	19.06
log_likelihood	-8.145

Example 3: Poisson sqrt link with single predictor

Inputs:

y	x	glm_poisson_link
3	1	sqrt
6	2
9	3
12	4

Excel formula:

=GLM_POISSON({3;6;9;12}, {1;2;3;4}, "sqrt")

Expected output:

parameter	coefficient	std_error	z_statistic	p_value	ci_lower	ci_upper	incidence_rate_ratio
intercept	1.226	0.6124	2.003	0.04521	0.02616	2.427	3.409
x1	0.5744	0.2236	2.569	0.01021	0.1361	1.013	1.776
deviance	0.06545
pearson_chi2	0.06541
aic	19.1
bic	17.88
log_likelihood	-7.552

Example 4: Poisson log link without intercept

Inputs:

y	x	glm_poisson_link	fit_intercept	alpha
4	2	log	false	0.1
6	3
8	4
10	5

Excel formula:

=GLM_POISSON({4;6;8;10}, {2;3;4;5}, "log", FALSE, 0.1)

Expected output:

parameter	coefficient	std_error	z_statistic	p_value	ci_lower	ci_upper	incidence_rate_ratio
x1	0.4971	0.04623	10.75	5.703e-27	0.4211	0.5732	1.644
deviance	1.453
pearson_chi2	1.57
aic	18.47
bic	17.86
log_likelihood	-8.236

Python Code

import math
import statsmodels.api as sm
from statsmodels.genmod.families import Poisson as sm_Poisson
from statsmodels.genmod.families.links import Log as sm_Log, Identity as sm_Identity, Sqrt as sm_Sqrt

def glm_poisson(y, x, glm_poisson_link='log', fit_intercept=True, alpha=0.05):
    """
    Fits a Generalized Linear Model with Poisson family for count data.

    See: https://www.statsmodels.org/stable/generated/statsmodels.genmod.generalized_linear_model.GLM.html

    This example function is provided as-is without any representation of accuracy.

    Args:
        y (list[list]): Column vector of count data (non-negative values). Each element represents an observation of the dependent variable.
        x (list[list]): Matrix of predictor variables. Each column represents a different predictor variable; each row represents an observation.
        glm_poisson_link (str, optional): The link function to use for the Poisson GLM model. Valid options: log, identity, sqrt. Default is 'log'.
        fit_intercept (bool, optional): If True, adds an intercept term to the model. Default is True.
        alpha (float, optional): Significance level for confidence intervals. Must be between 0 and 1. Default is 0.05.

    Returns:
        list[list]: 2D list with GLM results and statistics, or error string.
    """
    def to2d(value):
        return [[value]] if not isinstance(value, list) else value

    try:
      # Normalize inputs
      y = to2d(y)
      x = to2d(x)

      # Validate y is a column vector
      if not isinstance(y, list) or len(y) == 0:
        return "Error: Invalid input: y must be a non-empty 2D list."
      for row in y:
        if not isinstance(row, list) or len(row) != 1:
          return "Error: Invalid input: y must be a column vector (2D list with single column)."

      # Validate x is a 2D matrix
      if not isinstance(x, list) or len(x) == 0:
        return "Error: Invalid input: x must be a non-empty 2D list."
      num_rows = len(x)
      num_cols = len(x[0]) if isinstance(x[0], list) else 0
      if num_cols == 0:
        return "Error: Invalid input: x must be a 2D list with at least one column."
      for row in x:
        if not isinstance(row, list) or len(row) != num_cols:
          return "Error: Invalid input: x must have consistent number of columns."

      # Check dimensions match
      if len(y) != num_rows:
        return "Error: Invalid input: y and x must have the same number of rows."

      # Validate link function
      valid_links = ['log', 'identity', 'sqrt']
      if not isinstance(glm_poisson_link, str):
        return "Error: Invalid input: glm_poisson_link must be a string."
      if glm_poisson_link not in valid_links:
        return f"Error: Invalid input: glm_poisson_link must be one of {valid_links}."

      # Validate fit_intercept
      if not isinstance(fit_intercept, bool):
        return "Error: Invalid input: fit_intercept must be a boolean."

      # Validate alpha
      if not isinstance(alpha, (int, float)):
        return "Error: Invalid input: alpha must be a number."
      if math.isnan(alpha) or math.isinf(alpha):
        return "Error: Invalid input: alpha must be finite."
      if alpha <= 0 or alpha >= 1:
        return "Error: Invalid input: alpha must be between 0 and 1."

      # Convert y and x to numeric values
      try:
        y_values = [float(row[0]) for row in y]
      except (ValueError, TypeError):
        return "Error: Invalid input: y must contain numeric values."

      try:
        x_values = [[float(val) for val in row] for row in x]
      except (ValueError, TypeError):
        return "Error: Invalid input: x must contain numeric values."

      # Check for non-finite values
      for val in y_values:
        if math.isnan(val) or math.isinf(val):
          return "Error: Invalid input: y must contain finite values."
      for row in x_values:
        for val in row:
          if math.isnan(val) or math.isinf(val):
            return "Error: Invalid input: x must contain finite values."

      # Check that y values are non-negative (count data)
      for val in y_values:
        if val < 0:
          return "Error: Invalid input: y must contain non-negative values (count data)."

      # Add intercept if requested
      if fit_intercept:
        x_values = [[1.0] + row for row in x_values]

      # Select link function
      if glm_poisson_link == 'log':
        link = sm_Log()
      elif glm_poisson_link == 'identity':
        link = sm_Identity()
      elif glm_poisson_link == 'sqrt':
        link = sm_Sqrt()

      # Create Poisson family
      family = sm_Poisson(link=link)

      # Fit GLM
      try:
        model = sm.GLM(y_values, x_values, family=family)
        results = model.fit()
      except Exception as e:
        return f"Error: statsmodels.genmod.generalized_linear_model.GLM error: {e}"

      # Extract results
      params = results.params
      std_errors = results.bse
      z_stats = results.tvalues
      p_values = results.pvalues
      conf_int = results.conf_int(alpha=alpha)

      # Build results table
      output = [['parameter', 'coefficient', 'std_error', 'z_statistic', 'p_value', 'ci_lower', 'ci_upper', 'incidence_rate_ratio']]

      # Add parameter results
      for i in range(len(params)):
        if fit_intercept and i == 0:
          param_name = 'intercept'
        else:
          predictor_idx = i if not fit_intercept else i - 1
          param_name = f'x{predictor_idx + 1}'

        coef = float(params[i])
        std_err = float(std_errors[i])
        z_stat = float(z_stats[i])
        p_val = float(p_values[i])
        ci_low = float(conf_int[i, 0])
        ci_high = float(conf_int[i, 1])
        irr = math.exp(coef)

        output.append([param_name, coef, std_err, z_stat, p_val, ci_low, ci_high, irr])

      # Add model statistics
      output.append(['deviance', float(results.deviance), '', '', '', '', '', ''])
      output.append(['pearson_chi2', float(results.pearson_chi2), '', '', '', '', '', ''])
      output.append(['aic', float(results.aic), '', '', '', '', '', ''])
      output.append(['bic', float(results.bic_llf), '', '', '', '', '', ''])
      output.append(['log_likelihood', float(results.llf), '', '', '', '', '', ''])

      return output
    except Exception as e:
      return f"Error: {str(e)}"

Overview

Excel Usage

Examples

Python Code

Online Calculator