Decoding the Reed–Muller $RM(m, m)$ code, and a matrix related to Pascal's triangle

Question

The Reed–Muller $RM(m, m)$ code sends the message $p(X_0, \ldots , X_{m - 1}) = \sum_{S \subset \{0, \ldots , m - 1\}} \alpha_S \cdot X_S$ to its set of evaluations $\left\{ p(x_0, \ldots , x_{m - 1}) \right\}_{(x_0, \ldots , x_{m - 1}) \in \{0, 1\}^m}$. This code is perhaps not very useful in practice, because both its message length and block size are $2^m$, its distance is 1, and its encoding function is a bijection. On the other hand, (efficiently) decoding (i.e., inverting) this map seems rather interesting. Decoding a codeword $\left\{ \beta_{x_0, \ldots , x_{m - 1}} \right\}_{(x_0, \ldots , x_{m - 1}) \in \{0, 1\}^m}$ amounts to finding its algebraic normal form. Though various methods apparently exist, I want to propose one:

It's fairly easy to show that decoding $\left\{ \beta_{x_0, \ldots , x_{m - 1}} \right\}_{(x_0, \ldots , x_{m - 1}) \in \{0, 1\}^m}$ amounts to multiplying it by the inverse of a certain $2^m \times 2^m$ matrix related to Pascal's triangle. This matrix appears in OEIS A047999; it's also discussed by Massey in connection with RM codes (see the matrix $\mathbf{M}_m$ in the section "two useful matrices"). It also appears in a number of places throughout the literature (see e.g. Preparata, Saluja and Ong, and an extremely interesting treatment by Callan which shows that the Thue–Morse sequence shows up in its inverse!). It also apparently shows up in the LU decomposition of the Walsh matrix; see this image. As an example, $\mathbf{M}_3$ is reproduced here:

$$\mathbf{M}_3 = \begin{bmatrix} 1 & & & & & & & \\ 1 & 1 & & & & & & \\ 1 & & 1 & & & & & \\ 1 & 1 & 1 & 1 & & & \\ 1 & & & & 1 & & & \\ 1 & 1 & & & 1 & 1 & & \\ 1 & & 1 & & 1 & & 1 & \\ 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 \end{bmatrix} $$

You can see that its 1 entries take the shape of Sierpinski's triangle.

Naïvely multiplying $\left\{ \beta_{x_0, \ldots , x_{m - 1}} \right\}_{(x_0, \ldots , x_{m - 1}) \in \{0, 1\}^m}$ by the inverse of this matrix would take $\Theta((2^m)^2)$ time. It seems like there should be an FFT-like recursive procedure which instead takes $O(2^m \cdot m)$ time. In fact, I'm sure I could work it out. But someone else must have done this already! Any reference would be appreciated. Thanks.

Answer 1

Let $C$ be the coefficients of the algebraic normal form: then we have the linear matrix equation [you may need to check if your coordinate order is the reverse of mine] $$ C=fA_m $$ where $A_m$ is defined via a tensor power mod 2: $$ A_m=\left[ \begin{array}{cc} 1 & 0 \\ -1 & 1 \end{array}\right]^{\otimes m},\quad A_0=[1]. $$ Then the recursive algorithm you are looking for starts with the truth table $f=f(a),$ where $0\leq a\leq 2^m-1$ namely lets $$ f_{0,a}=f(a),\quad 0\leq a\leq 2^m-1. $$ Then,

for $k=0$ to $m-1$ do

$\quad$ for $b=0$ to $2^{m-k-1}-1$ do

$\quad\quad f_{(k+1),b \cdot 2^{k + 1}:(b+1) \cdot 2^{k + 1}}=\left[f_{k,2b \cdot 2^k:(2b + 1) \cdot 2^k}~~-f_{k,2b \cdot 2^k:(2b + 1) \cdot 2^k}+f_{k,(2b + 1) \cdot 2^k:2(b + 1) \cdot 2^k} \right]$

$\quad$ end for

end for

will yield $C=f_{m,0:2^m}.$

Recall the fast Hadamard transform over reals which scans pairs of entries from the left and generates the new vector as sums and differences.

This iteration works modulo 2 and scans pairs of entries from the left and generates the new vector as the left element of the pair followed by the sum of the elements in the pair, see below. Its complexity is clearly $m2^m.$

References: This was folklore for a long time, unsure where it first appeared in an academic publication. The terms ANF transform or Mobius transform are used. One relevant reference which is available online is here.