Last night, I encoded some messages in ASCII by hand, including
error-correction coding: like the printer dots Seth Schoen decoded,
matrices of 8 columns of 8 bits each, with an extra column and row for
parity.  (I used even parity, but odd would have been a better
choice.)  Matrix parity lets you correct single-bit errors: a
single-bit error will show up as a parity error in both the row parity
and the column parity.  And it's simple enough and spatial enough that
you can both detect and correct errors by hand.

This got me to thinking about what kind of more powerful
error-correction coding might be feasible to perform by hand, and I
think I came up with one, which at first glance appears to be related
to the Hadamard transform, but in the end seems to be something
simpler.

Existing Hadamard-code work is quite different
----------------------------------------------

For many years, Hadamard matrices have been used for error-correction
coding by the simple expedient of using a row of the matrix for each
codeword; since the rows are perfectly uncorrelated, they have a large
Hamming distance, and the Fast Hadamard Transform of the received
(possibly corrupted) codeword gives you the Hamming distance from all
the N codewords in 2 N lg N additions and subtractions.  (This is
now apparently obsoleted by better codes.)

A problem with this for human use is that it's not just some checksum
symbols you glom onto your data; it makes your data totally
unrecognizable.

I want to discuss a different application of Hadamard matrices for
error-correction codes.  As far as I can tell, this is not related to
"Hadamard-Craigen error correcting codes" (Craigen 2002), and in the
unlikely event that these codes become popular, I hope that nobody
decides to call them "Hadamard-Kragen codes".

The Hadamard transform and its relevance briefly explained
----------------------------------------------------------

The Hadamard transform of a vector is a "holographic representation"
of the vector in the sense that the vector can be recovered from it
(by repeating the Hadamard transform and dividing by the determinant)
and that every value in the Hadamard transform is affected by every
element in the original vector.  This means that any local change in
the original vector results in a global change in the
Hadamard-transformed vector --- in the case of the Hadamard transform,
*every* element changes.

(For simplicity, I'm going to consider only dimension-`2**n` Hadamard
matrices produced by the Sylvester method here, although others
exist.)

The Hadamard transform is particularly interesting among such
holographic representations because the Fast Hadamard Transform
enables computing the HT very efficiently.

The fast Hadamard transform of a vector or characters, illustrated
------------------------------------------------------------------

So suppose you have a vector of characters, say, 'Hadamard'.  You
compute the Fast Hadamard Transform of this vector:

    def hadamard(vec):
        if len(vec) == 1: return vec
        pairs = zip(hadamard(vec[:len(vec)/2]),
                    hadamard(vec[len(vec)/2:]))
        return ([lv + rv for lv, rv in pairs] +
                [lv - rv for lv, rv in pairs])

    >>> hadamard(map(ord, 'Hadamard'))
    [786, 4, -36, -30, -54, -48, -20, -26]

Note that a local change in the input produces a linear global change
in the output:

    >>> hadamard(map(ord, 'Hadamarc'))
    [785, 5, -35, -31, -53, -49, -21, -25]
    >>> hadamard(map(ord, 'Hacamarc'))
    [784, 4, -34, -30, -54, -50, -20, -24]

And you can reconstruct the original input vector from the
Hadamard-transformed version:

    >>> ''.join(chr(c/8) for c in hadamard(hadamard(map(ord, 'Hadamard'))))
    'Hadamard'

And, naturally, local errors in the Hadamard-transformed version
result in linear global errors in the result:

    >>> ''.join(chr(c/8) for c in hadamard([786,4,-36,-30,-54,-48,-20,-26]))
    'Hadamard'
    >>> ''.join(chr(c/8) for c in hadamard([778,4,-36,-38,-54,-48,-20,-26]))
    'Fad_karb'
    >>> ''.join(chr(c/8) for c in hadamard([700,4,-36,-38,-54,-48,-20,-26]))
    '<WZUaWhX'

How to use this for a human-friendly error-correction code
----------------------------------------------------------

So suppose you want to be able to correct errors.  You can pick some
subset of the Hadamard transform to transmit along with your raw data.
If there's a single-character error, it will show up in every
transformed value you transmit, but with varying signs; each sign
after the first gives you one bit of information about where the error
is --- *exactly* one bit, I thought, since the rows of the Hadamard
matrix are *perfectly* uncorrelated; so you should be able to localize
the error with lg N bits.  (That turns out to be false, as
demonstrated below; it matters a lot which rows you pick.)  In fact,
since the whole transform is linear, you should be able to just run
the inverse transform on your error signal and get the location of the
error.

Using the inverse transform on arbitrary outputs to localize errors works poorly
--------------------------------------------------------------------------------

Let's try:

    >>> orig = hadamard(map(ord, 'Hadamard'))[:3]
    >>> computed = hadamard(map(ord, 'Hagamard'))[:3]
    >>> [computed_i - orig_i for computed_i, orig_i in zip(computed, orig)]
    [3, 3, -3]
    >>> hadamard([computed_i - orig_i
                  for computed_i, orig_i
                  in zip(computed, orig)] + [0] * 5)
    [3, -3, 9, 3, 3, -3, 9, 3]

Not quite.  We got a nice peak in the transformed error signal at the
location of the actual error, it's true, but we got an equally big one
elsewhere.  It doesn't work with four values either:

    >>> orig = hadamard(map(ord, 'Hadamard'))[:4]
    >>> computed = hadamard(map(ord, 'Hagamard'))[:4]
    >>> hadamard([computed_i - orig_i
                  for computed_i, orig_i
                  in zip(computed, orig)] + [0] * 4)
    [0, 0, 12, 0, 0, 0, 12, 0]

It does work with five:

    >>> computed = hadamard(map(ord, 'Hagamard'))[:5]
    >>> orig = hadamard(map(ord, 'Hadamard'))[:5]
    >>> hadamard([computed_i - orig_i 
                  for computed_i, orig_i
                  in zip(computed, orig)] + [0] * 3)
    [3, 3, 15, 3, -3, -3, 9, -3]

And of course the first value tells us the sign and magnitude of the
actual error:

    >>> computed[0] - orig[0]
    3

What if we pick different values, like the ones from the end?  We
still need five of them to locate the error:

    >>> orig = hadamard(map(ord, 'Hadamard'))[-3:]
    >>> computed = hadamard(map(ord, 'Hagamard'))[-3:]
    >>> computed
    [-45, -23, -29]
    >>> hadamard([0] * 5 + [computed_i - orig_i 
                            for computed_i, orig_i
                            in zip(computed, orig)])
    [-3, -3, 9, -3, 3, 3, -9, 3]
    >>> orig = hadamard(map(ord, 'Hadamard'))[-4:]
    >>> computed = hadamard(map(ord, 'Hagamard'))[-4:]
    >>> hadamard([0] * 4 + [computed_i - orig_i
                            for computed_i, orig_i
                            in zip(computed, orig)])
    [0, 0, 12, 0, 0, 0, -12, 0]
    >>> orig = hadamard(map(ord, 'Hadamard'))[-5:]
    >>> computed = hadamard(map(ord, 'Hagamard'))[-5:]
    >>> hadamard([0] * 3 + [computed_i - orig_i
                            for computed_i, orig_i
                            in zip(computed, orig)])
    [-3, 3, 15, -3, -3, 3, -9, -3]

How about a different error?

    >>> computed = hadamard(map(ord, 'Hadamerd'))[-5:]
    >>> hadamard([0] * 3 + [computed_i - orig_i 
                            for computed_i, orig_i 
                            in zip(computed, orig)])
    [-4, -12, 4, -4, -4, 20, 4, -4]

How about an error in the transformed values?  That produces a
distinctly different pattern:

    >>> computed = [0] + orig[1:]
    >>> hadamard([0] * 3 + [computed_i - orig_i
                            for computed_i, orig_i
                            in zip(computed, orig)])
    [30, -30, -30, 30, 30, -30, -30, 30]

It seems like it ought to be possible to correct multiple errors.

How about a longer text?  Correcting one error in eight characters by
appending five error-check symbols is not at all impressive.  Can we
correct an error in a 64-character text with 11 Hadamard results?
(Not the way I've been doing it, but see later sections for
successfully doing it with 7.)

    >>> text='correct an error in a 64-character text with 11 Hadamard results'
    >>> len(text)
    64
    >>> orig = hadamard(map(ord, text))[:11]
    >>> errtext = text[:20] + 'x' + text[21:]
    >>> errtext
    'correct an error in x 64-character text with 11 Hadamard results'
    >>> computed = hadamard(map(ord, errtext))[:11]
    >>> orig
    [5806, -74, 100, -304, 68, -576, -166, -78, -170, 310, 212]
    >>> computed
    [5829, -51, 123, -281, 45, -599, -189, -101, -147, 333, 235]
    >>> hadamard([computed_i - orig_i
                  for computed_i, orig_i
                  in zip(computed, orig)] + [0] * 53)
    [69, 23, 23, -23, 253, 23, 23, -23, -69, -23, -23, 23, 115, -23,
    -23, 23, 69, 23, 23, -23, 253, 23, 23, -23, -69, -23, -23, 23,
    115, -23, -23, 23, 69, 23, 23, -23, 253, 23, 23, -23, -69, -23,
    -23, 23, 115, -23, -23, 23, 69, 23, 23, -23, 253, 23, 23, -23,
    -69, -23, -23, 23, 115, -23, -23, 23]

No, it identifies four equally likely locations for the error.  13?

    >>> orig = hadamard(map(ord, text))[:13]
    >>> computed = hadamard(map(ord, errtext))[:13]
    >>> hadamard([computed_i - orig_i
                  for computed_i, orig_i
                  in zip(computed, orig)] + [0] * 51)
    [69, -23, -23, -23, 299, 23, 23, 23, -69, 23, 23, 23, 69, -23,
    -23, -23, 69, -23, -23, -23, 299, 23, 23, 23, -69, 23, 23, 23, 69,
    -23, -23, -23, 69, -23, -23, -23, 299, 23, 23, 23, -69, 23, 23,
    23, 69, -23, -23, -23, 69, -23, -23, -23, 299, 23, 23, 23, -69,
    23, 23, 23, 69, -23, -23, -23]

No, mysteriously the new Hadamard results didn't give us any new
information about the location of the error; we still have four
equally likely locations.  It's as if the relevant Hadamard rows are
somehow correlated with the basis comprised of the rows we have
already, though not any of those rows individually.

It turns out that better choices are available:

    >>> orig_full = hadamard(map(ord, text))
    >>> computed_full = hadamard(map(ord, errtext))
    >>> orig = orig_full[:11] + [0] * 51 + orig_full[-2:]
    >>> computed = computed_full[:11] + [0] * 51 + computed_full[-2:]
    >>> hadamard([computed_i - orig_i 
                  for computed_i, orig_i
                  in zip(computed, orig)])
    [115, 23, -23, -23, 207, 23, 69, -23, -115, -23, 23, 23, 161, -23,
    -69, 23, 23, 23, 69, -23, 299, 23, -23, -23, -23, -23, -69, 23,
    69, -23, 23, 23, 23, 23, 69, -23, 299, 23, -23, -23, -23, -23,
    -69, 23, 69, -23, 23, 23, 115, 23, -23, -23, 207, 23, 69, -23,
    -115, -23, 23, 23, 161, -23, -69, 23]

That reduces the choices to two.  What if we choose by a different
approach?

    >>> orig = orig_full[::5]
    >>> len(orig)
    13
    >>> computed = computed_full[::5]
    >>> orig_expanded = [0] * 64
    >>> orig_expanded[::5] = orig
    >>> orig_expanded
    [5806, 0, 0, 0, 0, -576, 0, 0, 0, 0, 212, 0, 0, 0, 0, 86, 0, 0, 0,
    0, 438, 0, 0, 0, 0, -108, 0, 0, 0, 0, -220, 0, 0, 0, 0, 4, 0, 0,
    0, 0, -358, 0, 0, 0, 0, 96, 0, 0, 0, 0, -102, 0, 0, 0, 0, 120, 0,
    0, 0, 0, 406, 0, 0, 0]
    >>> computed_expanded = [0] * 64
    >>> computed_expanded[::5] = computed
    >>> hadamard([computed_i - orig_i 
                  for computed_i, orig_i
                  in zip(computed_expanded, orig_expanded)])
    [69, 161, -23, 161, 23, -69, 23, 23, 23, -161, 23, 115, -23, -23,
    69, -23, -23, 69, -23, -23, 299, 23, 23, 23, 23, 23, -69, 23, -23,
    -23, -23, 69, -23, 161, -23, -115, 23, 23, -69, 23, 23, 115, 115,
    115, 69, -23, 69, -115, -23, -23, 69, -23, 23, 23, 23, -69, -69,
    23, -69, 115, -23, 69, 253, 69]

That has a max, 299, at the right spot, with "only" 13
Hadamard-transformed values to give us six bits of error-location
information.

So, although this can be made to work, it doesn't perform as well as I
had hoped, and I didn't really understand why.  I did find out how to
stop it; see later sections.

Using the assumption that there's just one error doesn't help
-------------------------------------------------------------

But that approach is, in essence, attempting to find an estimate of an
*arbitrary* difference between the two vectors.  But we're assuming
that the difference here comes from a transmission error corrupting a
single letter.

What if we just use the signs of the differences, and try to figure
out which single input position could have caused that pattern of
signs?  That doesn't work either with poor choices of rows:

    sgn = lambda x: 1 if x > 0 else 0 if x == 0 else -1
    matrix = lambda n: [hadamard([i==j for i in range(n)]) for j in range(n)]
    def find(prefix, items):
        for ii, item in enumerate(items):
            if item[:len(prefix)] == prefix: yield ii

    for i in range(64):
         print i, list(find([sgn(computed_i - orig_i) 
                             for computed_i, orig_i 
                             in zip(computed_full[:i], orig_full[:i])], 
                            matrix(64)))

    0 [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17,
       18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33,
       34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49,
       50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63]
    ...
    2 [0, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24, 26, 28, 30, 32,
       34, 36, 38, 40, 42, 44, 46, 48, 50, 52, 54, 56, 58, 60, 62]
    3 [0, 4, 8, 12, 16, 20, 24, 28, 32, 36, 40, 44, 48, 52, 56, 60]
    4 [0, 4, 8, 12, 16, 20, 24, 28, 32, 36, 40, 44, 48, 52, 56, 60]
    5 [4, 12, 20, 28, 36, 44, 52, 60]
    ...
    9 [4, 20, 36, 52]
    ...
    17 [20, 52]
    ...
    33 [20]
    ...

Picking particularly helpful rows of the matrix *does* help
-----------------------------------------------------------

Somehow, almost all of the additional rows we get are adding no new
information.  Rows 0, 1, 2, 4, 8, 16, and 32 apparently added a bit
each.  (I say row 0 because we need the information about which
direction the change perturbed the input in order to interpret the
others.)  What if we use just those rows?

    select = lambda x: [x[i] for i in [0, 1, 2, 4, 8, 16, 32]]
    >>> len(set(tuple(select(row)) for row in matrix(64)))
    64

This looks promising!  Compare to:

    >>> len(set(tuple(row[:5]) for row in matrix(64)))
    8

And it works, see?

    >>> list(find([sgn(computed_i - orig_i)
                   for computed_i, orig_i 
                   in zip(select(computed_full), select(orig_full))],
                  [select(row) for row in matrix(64)]))
    [20]

So there we have seven Hadamard values that successfully localized the
error using only their signs.  (If the change were negative, we'd've
had to try `-sgn` as well as `sgn`; you could presumably avoid that by
using the pattern of sign *changes* instead of signs).  What if we use
the reverse-transform approach with these values instead of the
arbitrary ones we were using before?

    def expand(vals):
        rv = [0] * 64
        for val, pos in zip(vals, [0, 1, 2, 4, 8, 16, 32]): rv[pos] = val
        return rv

    hadamard(expand(computed_i - orig_i
                    for computed_i, orig_i 
                    in zip(select(computed_full), select(orig_full))))

    [69, 23, 23, -23, 115, 69, 69, 23, 23, -23, -23, -69, 69, 23, 23,
    -23, 115, 69, 69, 23, 161, 115, 115, 69, 69, 23, 23, -23, 115, 69,
    69, 23, 23, -23, -23, -69, 69, 23, 23, -23, -23, -69, -69, -115,
    23, -23, -23, -69, 69, 23, 23, -23, 115, 69, 69, 23, 23, -23, -23,
    -69, 69, 23, 23, -23]

And that does indeed have a clear peak of 161 at position 20, which is
where the error is.

What if we have a different error?  That seems to work too:

    >>> errtext = text[:35] + '!' + text[36:]
    >>> errtext
    'correct an error in a 64-character !ext with 11 Hadamard results'
    >>> computed_full = hadamard(map(ord, errtext))
    >>> select(computed_full)
    [5723, 9, 183, -15, -253, -55, 5]
    >>> select(orig_full)
    [5806, -74, 100, 68, -170, 28, -78]

    ...     hadamard(expand(computed_i - orig_i
    ...                     for computed_i, orig_i 
    ...                     in zip(select(computed_full), select(orig_full))))

    [-83, -249, -249, -415, 83, -83, -83, -249, 83, -83, -83, -249,
    249, 83, 83, -83, 83, -83, -83, -249, 249, 83, 83, -83, 249, 83,
    83, -83, 415, 249, 249, 83, -249, -415, -415, -581, -83, -249,
    -249, -415, -83, -249, -249, -415, 83, -83, -83, -249, -83, -249,
    -249, -415, 83, -83, -83, -249, 83, -83, -83, -249, 249, 83, 83,
    -83]

And indeed there's a negative peak at position 35, where all seven
differences of 83 ('t' - '!') line up with the same sign, producing
-581.  And using the approach of just looking at the signs?  That
works too, but since it's a negative change, we have to negate the
signs:

    >>> list(find([sgn(computed_i - orig_i)
    ...                    for computed_i, orig_i 
    ...                    in zip(select(computed_full), select(orig_full))],
    ...                   [select(row) for row in matrix(64)]))
    []
    >>> list(find([-sgn(computed_i - orig_i)
    ...                    for computed_i, orig_i 
    ...                    in zip(select(computed_full), select(orig_full))],
    ...                   [select(row) for row in matrix(64)]))
    [35]

This is very far from a proof that the method works, but it's
suggestive.  Here are the rows of `select(matrix(64))`:

    [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, 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, 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, -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, 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, -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, -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, 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, 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, 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, -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, -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, 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, -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, -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, 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, -1, -1, -1, -1, -1, -1, -1, -1, -1,
    -1, -1, -1, -1, -1, -1, -1]

These appear to be precisely the most regular rows, the ones with
sequencies 0, 63, 31, 15, 7, 3, and 1.  You certainly don't need a
Hadamard matrix to compute those!  It's closely related to the matrix
parity code I started with --- but instead of XOR, you're using
addition, and instead of two dimensions, you're in 6, and instead of
separately recording the totals modulo whatever of two 5-dimensional
hyperplanes, you're recording only the difference between them.

So maybe this is already known under some other name!

Using finite fields instead of the integers?
--------------------------------------------

I did all the above in the ring Z, the integers, which has the problem
that you need a potentially unbounded number of bits to represent its
members.

All of this also works over finite fields such as Z/7, the integers
mod 7, as well as the integers themselves.  In fact, it isn't even
necessary to have a full field, which is why it's possible to do it in
Z at all; you only need to be able to divide by the power of two that
is your vector length.  Unfortunately, this eliminates the most
desirable group for computation, Z/256 or in general `Z/2**n`, since
it will lose the upper bits of each character.

One option would be to use Z/257, the integers mod 257.  Then you can
represent the Hadamard transform results in a byte each, plus a
usually-empty list of which transform results are equal to 256; or in
9 bits each.

Another somewhat more computationally expensive option, which might
avoid that additional unnecessary message expansion, is to represent
the entire message in some other finite field: start by interpreting
it as a large number in base 256, then transform it into a sequence of
digits in e.g. Z/7 or Z/257, then do the Hadamard computation in that
field.

Using a field with fewer members makes the message longer, which
improves the ratio lg N/N, which I strongly suspect is proportional to
the number of Hadamard results you need to include to localize and
correct a single error.

How do you find the peaks in the transformed error signal in Z/p,
where big numbers can wrap around to be small ones again, and division
can easily produce big numbers?  I have no idea.  But the
sign-sequence-lookup approach should work, even though "sign" itself
is not a well-defined concept in Z/p.
-- 
To unsubscribe: http://lists.canonical.org/mailman/listinfo/kragen-tol

Reply via email to