(I think this was published previously on kragen-tol.)
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.
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 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.
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'
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.
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.
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]
...
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!
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.