PyGLMs

An implementation of various Generalized Linear Models (GLMs), written in Python.

GLMs are an excellent foundation for predictive modeling and machine learning, so the code is intentionally written for clarity, to help other readers and developers understand the underlying mathematics. Building this package also served as my own personal refresher on the topic.

Overview

The code is packaged as a Python library named turtles-glms (I like turtles).

It is written using numpy for linear algebra operations, scipy for (some) optimization, pandas for displaying tabular results, and matplotlib for plots.

The following model frameworks have been implemented:

1. Multiple Linear Regression (turtles.stats.glms.MLR class)
2. Logistic Regression (turtles.stats.glms.LogReg class, uses GLM parent class)
3. Poisson Regression (turtles.stats.glms.PoissonReg class, uses GLM parent class)

The GLM parent class supports three optimization methods for parameter estimation: Momentum-based Gradient Descent for first-order optimization, Newton's Method for second-order optimization, and Limited-memory Broyden–Fletcher–Goldfarb–Shanno (L-BFGS). The user can specify the desired optimization method during class instantiation.

Momentum-based Gradient Descent and Newton's Method are implemented in Python as part of the code base. L-BFGS is implemented using scipy.optimize; it's a quasi-Newton method that approximates the Hessian (instead of fully computing it, like Newton's Method), so it's quite fast.

Usage

You can pip install the package from PyPI:

pip install turtles-glms

See examples/ in the GitHub repo for example usage of the GLM classes and statistical functions.

Fitting GLMs

You can fit GLMs by instantiating a GLM child class and calling its fit() method.

model = PoissonReg(
    method="newton",
    learning_rate=1.0,
    tolerance=0.00001
)
model.fit(
    X=X, 
    y=y, 
    exposure=exposure
)

A few important notes about fitting turtles GLMs:

1. The fit() method parameters X, y, and (for Poisson) exposure must be numpy arrays. Parameters y and exposure must be of shape (M, 1), where M is the number of rows in the data. The package does not support pandas or polars dataframes at this time. See class / instance method docstrings for exact requirements.
2. Each GLM class has a learning_rate parameter, applicable to Gradient Descent and Newton's optimization methods. The learning rate (or step size) is a hyperparameter that controls the magnitude of parameter updates during the optimization process. If it's too large, the Hessian matrix may become singular, in which case the learning rate should be decreased. This is typically part of the tuning process. (NOTE: learning rate is typically just 1.0 for Newton's method).
3. There are currently no regularization methods implemented in the package. Future versions may include L1, L2, and Elastic Net methods.

Unreleased

Some future updates I'd like to make:

1. Support flattened arrays (in addition to (m, 1))
2. Tweedie GLM
3. Regularization methods