IID Model

Overview

The IID model (Independent and Identically Distributed) is the simplest and most commonly used random effect in hierarchical modeling. It assumes that random effects are:

• Independent - no correlation between different groups
• Identically distributed - all groups share the same variance
• Normally distributed - each effect follows a Gaussian distribution

This makes IID ideal for capturing exchangeable group-level heterogeneity - situations where groups differ from each other, but there's no inherent ordering or spatial relationship between them.

Mathematical Formulation

For \(n\) levels of the index variable, the IID model specifies:

\[u_i \stackrel{\text{iid}}{\sim} \mathcal{N}(0, \tau^{-1}), \quad i = 1, \ldots, n\]

where:

• \(u_i\) is the random effect for level \(i\)
• \(\tau\) is the precision (inverse variance) of the random effects
• \(\sigma = 1/\sqrt{\tau}\) is the standard deviation

Terminology note: throughout this page, level means one value of the index variable (one realization of \(u_i\)). The separate spec key group refers to an additional axis of replication and is unrelated.

The precision matrix for the full vector \(\mathbf{u} = (u_1, \ldots, u_n)^\top\) is simply:

\[\mathbf{Q} = \tau \, \mathbf{I}_n\]

where \(\mathbf{I}_n\) is the \(n \times n\) identity matrix.

When to Use IID

The IID model is appropriate when you have grouped or clustered data where:

Visual Understanding: IID vs Structured Effects

Understanding when to use IID becomes clearer when you see how it differs from structured random effects like RW1. The key difference is in the correlation structure.

IID: Independent Random Effects

Each level \(u_i\) is drawn independently from the same zero-mean Gaussian \(\mathcal{N}(0, \tau^{-1})\). The values are exchangeable: there is no inherent ordering, and the realized effect at one level carries no information about any other.

RW1: Correlated Random Effects

Adjacent positions are correlated - nearby values tend to be similar. This creates smooth trends rather than independent jumps.

Decision Guide

Specification in pyINLA

Allowed Parameters

The IID random-effect spec dict accepts the following keys. The source of truth is model_defs/latent/iid.py's ALLOWED_KEYS set, mirrored 1:1 here.

Any key not listed above raises a clear safety error with the allowed-keys list before the engine runs. Per-observation scaling is done via weights, not scale (the named scale= argument is not a supported IID feature).

Per-Parameter Runnable Examples

Click any parameter to expand a minimal runnable snippet. All snippets assume import numpy as np; import pandas as pd; import math; from pyinla import pyinla and a small dataframe df with columns y, x, group_id (5 levels).

# shared setup for every snippet below
import numpy as np
import pandas as pd
import math
from pyinla import pyinla

rng = np.random.default_rng(42)
n = 60
df = pd.DataFrame({
    'y': rng.standard_normal(n),
    'x': rng.standard_normal(n),
    'group_id': rng.integers(1, 6, size=n),
})

formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid'}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# 'model': 'iid' is the marker that selects this family (required).
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid'}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

formula = {
    'response': 'y',
    'random': [{
        'id': 'group_id', 'model': 'iid',
        'hyper': [{
            'prior': 'pc.prec',
            'param': [1, 0.01],
            'initial': 2,
            'fixed': False,
        }],
    }],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# 'n' must equal the count of unique levels declared by the data or by
# 'values'. The engine refuses to run when they mismatch. To reserve
# room for unobserved levels, declare them with 'values' AND set 'n'
# to that length.
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid',
                 'values': list(range(1, 11)), 'n': 10}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid', 'constr': True}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# start the optimizer at log(tau) = 2  (tau ~ 7.4)
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid', 'initial': 2}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# fix tau = 10  (no estimation of the precision hyperparameter)
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid',
                 'initial': math.log(10), 'fixed': True}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# small jitter for numerical stability
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid', 'diagonal': 1e-4}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

w = rng.uniform(0.5, 1.5, size=n)
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid', 'weights': w}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid',
                 'id_names': ['A', 'B', 'C', 'D', 'E']}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# declare levels 1..5 even if some are absent from this fit
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid',
                 'values': [1, 2, 3, 4, 5]}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

df['rep_id'] = rng.integers(1, 4, size=n)
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid',
                 'nrep': 3, 'replicate': df['rep_id'].values}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

df['week'] = rng.integers(1, 5, size=n)
formula = {
    'response': 'y',
    'random': [{
        'id': 'group_id', 'model': 'iid',
        'group': df['week'].values,
        'ngroup': 4,
        'control.group': {'model': 'ar1'},
    }],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# enforce sum of all 5 random effects equals zero
formula = {
    'response': 'y',
    'random': [{
        'id': 'group_id', 'model': 'iid',
        'extraconstr': {'A': np.ones((1, 5)),
                         'e': np.zeros(1)},
    }],
}
result = pyinla(formula=formula, family='gaussian', data=df)

df['group2_id'] = rng.integers(1, 6, size=n)
formula = {
    'response': 'y',
    'random': [
        {'id': 'group_id',  'model': 'iid'},
        {'id': 'group2_id', 'model': 'iid', 'copy': 'group_id'},
    ],
}
result = pyinla(formula=formula, family='gaussian', data=df)

# enable VB correction on this component  (or pass a list of node indices)
formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid', 'vb.correct': True}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

formula = {
    'response': 'y',
    'random': [{'id': 'group_id', 'model': 'iid', 'debug': True}],
}
result = pyinla(formula=formula, family='gaussian', data=df)

Level Identifier Handling

The id column should contain integer values identifying which level of the random effect each observation belongs to. pyINLA handles indices flexibly:

Key points:

• Indices are mapped to consecutive positions internally via sort(unique(values))
• np.nan in the index column means "this row has no contribution from this random effect" (useful for mixed designs)
• String/factor columns must be converted to integers first (e.g. pd.factorize)

# These all work:
df['idx'] = [1, 2, 3, 1, 2, 3]       # Standard
df['idx'] = [10, 20, 30, 10, 20, 30]   # Non-consecutive (mapped to 0,1,2)
df['idx'] = [1, np.nan, 2, np.nan, 1, 2]  # np.nan drops the RE contribution for that row

# Convert string groups to integers
df['group_idx'] = pd.factorize(df['group_name'])[0] + 1

NaN in the weights vector

When weights is supplied, individual entries may also be np.nan / NA. pyINLA accepts these without error: the engine treats those rows as having no weighted contribution from this random effect. There is no auto-imputation or row-dropping side effect.

This behaviour is covered by the nan_in_weights case in the iid sweep.

Hyperparameters

The IID model has one hyperparameter: the precision \(\tau\) of the random effects. In pyINLA, we work with the log-transformed precision \(\theta = \log(\tau)\) for numerical stability.

Default Configuration

Customizing Hyperparameters

Override defaults using the hyper key in your random effect specification:

'random': [
    {
        'id': 'group_id',
        'model': 'iid',
        'hyper': [{
            'prior': 'loggamma',
            'param': [1, 0.01],      # More informative prior
            'initial': 2,             # Different starting point
            'fixed': False
        }]
    }
]

Available Priors

Log-Gamma Prior (Default)

The default loggamma prior places a gamma distribution on the precision \(\tau\):

\[\tau \sim \text{Gamma}(a, b)\]

With default param=[1, 0.00005], this is a very weakly informative prior that allows \(\tau\) to range over several orders of magnitude.

# Default loggamma prior
'hyper': [{'prior': 'loggamma', 'param': [1, 0.00005]}]

# More informative: expect precision around 10-100
'hyper': [{'prior': 'loggamma', 'param': [1, 0.01]}]

PC Prior (Penalised Complexity)

The pc.prec prior is a principled choice that penalizes deviation from a base model (no random effect). It's parameterized in terms of the standard deviation \(\sigma = 1/\sqrt{\tau}\):

\[P(\sigma > u) = \alpha\]

With param=[u, alpha], you specify that the probability of \(\sigma\) exceeding \(u\) is \(\alpha\).

# PC prior: P(sigma > 1) = 0.01
# "I believe standard deviation is unlikely to exceed 1"
'hyper': [{'prior': 'pc.prec', 'param': [1, 0.01]}]

# PC prior: P(sigma > 0.5) = 0.05
# More restrictive: effects should be small
'hyper': [{'prior': 'pc.prec', 'param': [0.5, 0.05]}]

Other Supported Priors

Beyond loggamma and pc.prec, the precision hyperparameter (internally \(\theta = \log \tau\)) accepts several more priors. Pick by what you want to place the prior on: the log-precision \(\theta\), the precision \(\tau\), or the standard deviation \(\sigma = 1/\sqrt{\tau}\).

# Gaussian on log-precision, centred at 4 with precision 1
'hyper': [{'prior': 'normal', 'param': [4, 1]}]

# Improper flat prior on log-precision
'hyper': [{'prior': 'flat'}]

# Custom log-density as an R expression (here: standard normal on theta)
'hyper': [{'prior': 'expression: -theta * theta / 2'}]

# Tabulated log-density on a grid: theta values | log pi(theta) values
'hyper': [{'prior': 'table: -3 -1 1 3 | -4.5 -1.1 -1.1 -4.5'}]

Unknown prior names are rejected up front with a suggestion of near matches. The full registry is in pyinla/api/priors.py.

Fixed Precision

To fix the precision at a known value (no estimation), set fixed=True and specify the initial value on the log scale:

import math

# Fix precision at tau = 10 (variance = 0.1)
'hyper': [{
    'initial': math.log(10),  # log(tau) = 2.303
    'fixed': True
}]

Interpreting Results

Hyperparameter Summary

After fitting, access the precision estimate via result.summary_hyperpar:

print(result.summary_hyperpar)
#                      mean        sd  0.025quant  0.5quant  0.975quant     mode
# Precision for z  6.995828  1.913305    3.986632  6.748838   11.421926  6.294978

To convert precision to standard deviation:

tau_mean = result.summary_hyperpar['mean'].values[0]
sigma_mean = 1 / np.sqrt(tau_mean)
print(f"Estimated SD of random effects: {sigma_mean:.3f}")

Random Effect Estimates

The posterior estimates for each group's random effect are available in result.summary_random:

# Random effects for each group
print(result.summary_random['group_id'])
#    ID      mean        sd  0.025quant  0.5quant  0.975quant      mode
# 0   1  0.234567  0.089123    0.059876  0.234567    0.409258  0.234567
# 1   2 -0.123456  0.091234   -0.302345 -0.123456    0.055433 -0.123456
# ...