POSTERIOR_LOGSUMEXP

This function computes the logarithm of summed exponentials in a numerically stable way, which is essential for posterior normalization in log space and aggregation of log-weights.

The core operation is:

\log\left(\sum_i e^{a_i}\right)

and with optional weights b_i:

\log\left(\sum_i b_i e^{a_i}\right)

Working in log space helps avoid overflow and underflow when posterior log densities span large ranges.

Excel Usage

=POSTERIOR_LOGSUMEXP(a, axis, b, keepdims, return_sign)
  • a (list[list], required): 2D array of log-space values.
  • axis (int, optional, default: null): Axis along which to aggregate (dimensionless index).
  • b (list[list], optional, default: null): Optional 2D array of linear-space scaling weights.
  • keepdims (bool, optional, default: false): Keep reduced dimensions in output shape.
  • return_sign (bool, optional, default: false): Return sign information for weighted sums that may be negative.

Returns (list[list]): Log-sum-exp result as a 2D array, or paired log-value/sign rows when return_sign is true.

Example 1: Stable aggregate for a single row of log values

Inputs:

a
0 1 2

Excel formula:

=POSTERIOR_LOGSUMEXP({0,1,2})

Expected output:

2.40761

Example 2: Weighted aggregate using posterior scaling factors

Inputs:

a b
0 1 2 1 2 1

Excel formula:

=POSTERIOR_LOGSUMEXP({0,1,2}, {1,2,1})

Expected output:

2.62652

Example 3: Axis reduction with retained dimensions

Inputs:

a axis keepdims
0 1 1 true
2 3

Excel formula:

=POSTERIOR_LOGSUMEXP({0,1;2,3}, 1, TRUE)

Expected output:

Result
1.31326
3.31326
Example 4: Signed output for mixed-sign weighted terms

Inputs:

a b return_sign
1 2 1 -1 true

Excel formula:

=POSTERIOR_LOGSUMEXP({1,2}, {1,-1}, TRUE)

Expected output:

Result
1.54132
-1

Python Code

import numpy as np
from scipy.special import logsumexp as scipy_logsumexp

def posterior_logsumexp(a, axis=None, b=None, keepdims=False, return_sign=False):
    """
    Compute stable log-sum-exp aggregates for posterior normalization and evidence calculations.

    See: https://docs.scipy.org/doc/scipy/reference/generated/scipy.special.logsumexp.html

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

    Args:
        a (list[list]): 2D array of log-space values.
        axis (int, optional): Axis along which to aggregate (dimensionless index). Default is None.
        b (list[list], optional): Optional 2D array of linear-space scaling weights. Default is None.
        keepdims (bool, optional): Keep reduced dimensions in output shape. Default is False.
        return_sign (bool, optional): Return sign information for weighted sums that may be negative. Default is False.

    Returns:
        list[list]: Log-sum-exp result as a 2D array, or paired log-value/sign rows when return_sign is true.
    """
    try:
        def to2d(v):
            return [[v]] if not isinstance(v, list) else v

        def to_excel_2d(value):
            arr = np.asarray(value)
            if arr.ndim == 0:
                return [[float(arr)]]
            if arr.ndim == 1:
                return [[float(val) for val in arr.tolist()]]
            if arr.ndim == 2:
                return [[float(val) for val in row] for row in arr.tolist()]
            flat = arr.ravel()
            return [[float(val) for val in flat.tolist()]]

        a = to2d(a)
        if not isinstance(a, list) or not all(isinstance(row, list) for row in a):
            return "Error: a must be a 2D list"

        try:
            a_arr = np.array(a, dtype=float)
        except Exception:
            return "Error: a must contain numeric values"

        b_arr = None
        if b is not None:
            b = to2d(b)
            if not isinstance(b, list) or not all(isinstance(row, list) for row in b):
                return "Error: b must be a 2D list"
            try:
                b_arr = np.array(b, dtype=float)
            except Exception:
                return "Error: b must contain numeric values"

        axis_arg = None if axis is None else int(axis)

        if return_sign:
            log_val, sign_val = scipy_logsumexp(
                a_arr,
                axis=axis_arg,
                b=b_arr,
                keepdims=bool(keepdims),
                return_sign=True
            )
            log_row = to_excel_2d(log_val)
            sign_row = to_excel_2d(sign_val)
            width = max(len(log_row[0]), len(sign_row[0])) if log_row and sign_row else 1
            padded_log = (log_row[0] + [""] * (width - len(log_row[0]))) if log_row else [""] * width
            padded_sign = (sign_row[0] + [""] * (width - len(sign_row[0]))) if sign_row else [""] * width
            return [padded_log, padded_sign]

        result = scipy_logsumexp(
            a_arr,
            axis=axis_arg,
            b=b_arr,
            keepdims=bool(keepdims),
            return_sign=False
        )
        return to_excel_2d(result)
    except Exception as e:
        return f"Error: {str(e)}"

Online Calculator

2D array of log-space values.
Axis along which to aggregate (dimensionless index).
Optional 2D array of linear-space scaling weights.
Keep reduced dimensions in output shape.
Return sign information for weighted sums that may be negative.