pymc.ICAR#

class pymc.ICAR(name, *args, rng=None, dims=None, initval=None, observed=None, total_size=None, transform=UNSET, default_transform=UNSET, **kwargs)[source]#

The intrinsic conditional autoregressive prior.

It is primarily used to model covariance between neighboring areas. It is a special case of the CAR distribution where alpha is set to 1.

The log probability density function is

\[f(\phi| W,\sigma) = -\frac{1}{2\sigma^{2}} \sum_{i\sim j} (\phi_{i} - \phi_{j})^2 - \frac{1}{2}*\frac{\sum_{i}{\phi_{i}}}{0.001N}^{2} - \ln{\sqrt{2\pi}} - \ln{0.001N}\]

The first term represents the spatial covariance component. Each \(\phi_{i}\) is penalized based on the square distance from each of its neighbors. The notation \(i \sim j\) indicates a sum over all the neighbors of \(\phi_{i}\). The last three terms are the Normal log density function where the mean is zero and the standard deviation is N * 0.001 (where N is the length of the vector \(\phi\)). This component imposes a zero-sum constraint by finding the sum of the vector \(\phi\) and penalizing based on its distance from zero.

Parameters:
Wndarray of int

Symmetric adjacency matrix of 1s and 0s indicating adjacency between elements.

sigmascalar, default 1

Standard deviation of the vector of phi’s. Putting a prior on sigma will result in a centered parameterization. In most cases, it is preferable to use a non-centered parameterization by using the default value and multiplying the resulting phi’s by sigma. See the example below.

zero_sum_stdevscalar, default 0.001

Controls how strongly to enforce the zero-sum constraint. The sum of phi is normally distributed with a mean of zero and small standard deviation. This parameter sets the standard deviation of a normal density function with mean zero.

References

Examples

This example illustrates how to switch between centered and non-centered parameterizations.

import numpy as np
import pymc as pm

# 4x4 adjacency matrix
# arranged in a square lattice

W = np.array(
    [
        [0, 1, 0, 1],
        [1, 0, 1, 0],
        [0, 1, 0, 1],
        [1, 0, 1, 0],
    ],
)

# centered parameterization
with pm.Model():
    sigma = pm.Exponential("sigma", 1)
    phi = pm.ICAR("phi", W=W, sigma=sigma)
    mu = phi

# non-centered parameterization
with pm.Model():
    sigma = pm.Exponential("sigma", 1)
    phi = pm.ICAR("phi", W=W)
    mu = sigma * phi

Methods

ICAR.dist(W[, sigma, zero_sum_stdev])

Create a tensor variable corresponding to the cls distribution.