LDPC-Style Decoding

This example builds a small low-density parity-check style decoding problem. The transmitted codeword is observed through a binary symmetric channel. Unary factors encode the channel likelihoods, and parity-check factors enforce even-parity constraints.

The graph contains cycles, so iterative sum-product belief propagation is used to estimate the posterior probability of each bit.

Variable Nodes

Each code bit is a binary discrete variable with state order [0, 1]:

using FactorGraph

states = [0, 1]
bitIds = [:x1, :x2, :x3, :x4, :x5, :x6]

variables = [
    DiscreteVariable(bitId, 2; label = string(bitId), states = states)
    for bitId in bitIds
]

Channel Factors

Assume a binary symmetric channel with bit-flip probability p. The example transmits a valid nonzero codeword and flips one received bit. If the received bit is 0, the likelihood over states [0, 1] is [1 - p, p]. If the received bit is 1, the likelihood is [p, 1 - p].

p = 0.08
transmitted = [1, 0, 1, 1, 1, 0]
received = [1, 0, 1, 0, 1, 0]

channelLikelihood(y, p) = y == 0 ? [1.0 - p, p] : [p, 1.0 - p]

channelFactors = [
    DiscreteFactor(
        bitIds[index], channelLikelihood(received[index], p);
        label = "channel_$(bitIds[index])", initialize = true
    )
    for index in eachindex(bitIds)
]

Parity-Check Factors

The parity-check matrix defines which bits participate in each even-parity constraint:

H = [
    1 1 0 1 0 0
    0 1 1 0 1 0
    1 0 1 0 0 1
]

function evenParityTable(degree)
    table = zeros(Float64, ntuple(_ -> 2, degree))

    for assignment in Iterators.product(ntuple(_ -> 0:1, degree)...)
        table[map(bit -> bit + 1, assignment)...] = iseven(sum(assignment)) ? 1.0 : 0.0
    end

    return table
end

parityFactors = [
    begin
        checkVariables = bitIds[findall(H[checkIndex, :] .== 1)]
        DiscreteFactor(
            checkVariables..., evenParityTable(length(checkVariables));
            label = "check_$checkIndex"
        )
    end
    for checkIndex in axes(H, 1)
]

For each parity factor, table axes follow the order of checkVariables.

Factor Graph Construction

Build the factor graph:

graph = factorGraph(variables, vcat(channelFactors, parityFactors))

The graph can be rendered as an SVG factor graph figure:

saveGraphFigure("../ldpcd.svg", graph; label = (showEdgeIds = true, tooltipDetail = :full))

Running Belief Propagation

Run loopy sum-product belief propagation on the graph:

inference = sumproduct(graph)

gbp!(graph, inference; iterations = 25, tolerance = 1e-8, schedule = :flooding)

Posterior Marginals

The posterior marginal of each bit gives the decoder's soft information:

posterior = [
    marginal(graph, inference, bitId)
    for bitId in bitIds
]

6-element Vector{Vector{Float64}}:
 [0.1228176070393919, 0.8771823929606081]
 [0.8771823929606081, 0.1228176070393919]
 [0.06337378100400985, 0.9366262189959902]
 [0.24533915400373224, 0.7546608459962677]
 [0.07106698431820907, 0.928933015681791]
 [0.928933015681791, 0.07106698431820907]

Hard Decision

A hard decision chooses the most likely state for each bit:

decoded = [
    states[argmax(marginal(graph, inference, bitId))]
    for bitId in bitIds
]

6-element Vector{Int64}:
 1
 0
 1
 1
 1
 0

Parity-Check Validation

The decoded word satisfies the parity checks:

mod.(H * decoded, 2)

3-element Vector{Int64}:
 0
 0
 0