pymc.ICAR#

class pymc.ICAR(name, *args, rng=None, dims=None, initval=None, observed=None, total_size=None, 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

\[\begin{split}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}\end{split}\]

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])

Creates a tensor variable corresponding to the cls distribution.