We know from the previous section that we can restrict our search to prefix codes. For such codes, the entropy \(H(A)\) provides bounds on the expected length. Specifically, we have:

Proof.

For any \(a\in \mathcal{A}\), we define \(l(a) = \lceil - \log_2 p(a) \rceil\). As we have \[ \sum_{a \in \mathcal{A}} 2^{-l(a)} \leq \sum_{a \in \mathcal{A}} 2^{\log_2 p(a)} = \sum_{a \in \mathcal{A}} p(a) = 1 \] there is a prefix code \(c\) that satisfies \(|c(a)| = l(a)\) for any \(a \in \mathcal{A}\). Moreover \[ \mathbb{E} |c(A)| = \sum_{a \in \mathcal{A}} p(a) l(a) < \sum_{a \in \mathcal{A}} p(a) (- \log_2 p(a) + 1) = H(A) + 1. \] ‌\(\blacksquare\)

Probabilities will be used as relative weights: instead we could define the dictionary as

{'a': 5, 'b': 3, 'c': 2}

and the code resulting from the Huffman algorithm would be the same.

The algorithm we implement uses weighted binary trees where:

• terminal nodes are symbol/weight pairs such as

  ('a', 1.0)

• non-terminal nodes are children/weight pairs where children is a list of (terminal or non-terminal) nodes, such as

  ([('c': 0.2), ('a': 0.5)], 0.7)

Handling of these structures is made through a small set of helper functions:

class Node(object):
    "Function helpers to manage nodes as (symbol, weight) pairs"

    @staticmethod
    def symbol(node):
        return node[0]

    @staticmethod
    def weight(node):
        return node[1]

    @staticmethod
    def is_terminal(node):
        return not isinstance(Node.symbol(node), list)

The Huffman algorithm is the following: we create a list of symbol/weight nodes from the initial dictionary and repeat the following steps:

1. pick up two nodes with the least weights, remove them from the list,

2. insert into the list a non-terminal nodes with those two nodes as children and a weight that is the sum of their weights,,

3. stop when there is a single node in the list: this is the root of the binary tree

Note that this algorithm is not entirely deterministic if at any step, the lowest weight correspond to three nodes or more.

Here is the corresponding implementation in make_binary_tree; note that in the other methods of the huffman class, we also create a symbol to code dictionary (self.table) from the binary tree to simplify the coding and decoding of symbols – not all functions are given in this chunk of code ; the symbol_encoder and symbol_decoder create (symbol) coders and decoders from the symbol to code tables and stream_encoder and stream_decoder create stream coder and encoder from those.

class huffman(object):
    def __init__(self, weighted_alphabet):
        self.weighted_alphabet = weighted_alphabet
        self.tree = huffman.make_binary_tree(weighted_alphabet)
        self.table = huffman.make_table(self.tree)
        self.encoder = stream_encoder(symbol_encoder(self.table))
        self.decoder = stream_decoder(symbol_decoder(self.table))

    @staticmethod
    def make_binary_tree(weighted_alphabet):
        nodes = weighted_alphabet.items()
        while len(nodes) > 1:
            nodes.sort(key=Node.weight)
            node1, node2 = nodes.pop(0), nodes.pop(0)
            node = ([node1, node2], Node.weight(node1) + Node.weight(node2))
            nodes.insert(0, node)
        return nodes[0]

    @staticmethod
    def make_table(root, table=None, prefix=None):
        if prefix is None:
            prefix = BitStream()
        if table is None:
            table = {}
        if not Node.is_terminal(root):
            for index in (0, 1):
                new_prefix = prefix.copy()
                new_prefix.write(bool(index))
                new_root = Node.symbol(root)[index]
                huffman.make_table(new_root, table, new_prefix)
        else:
            table[Node.symbol(root)] = prefix
        return table

As usual, the coder and decoder are registered for use with the BitStream instances.

bitstream.register(huffman, reader=lambda h: h.decoder, writer=lambda h: h.encoder)

> \;\; < \;\; + \;\; - \;\; . \;\; , \;\; [ \;\; ]

The language is organised around a single byte pointer p. As an exemple, the instruction “+' corresponds to the C code fragment "(*p)++;"​ and ".’ to “putchar(*p);”​. The classic “Hello World!” program in Brainfuck is the sequence:

++++++++++[>+++++++>++++++++++>+++>+<<<<-]
>++.>+.+++++++..+++.>++.<<+++++++++++++++.
>.+++.------.--------.>+.>.

Spoon is a variant of Brainfuck created in 1998 by Steven Goodwin. It uses binary encoding of the eight original symbols according to the mapping

The code is a prefix code that was actually designed using Huffman encoding, with probabilities of each command computed from a set of example programs. You may notice that the code appears NOT to be optimal as , has the longest code and has no sibling: it could be shortened to \(001011\). But that’s only because the actual Spoon variant of Brainfuck introduces two extra reserved codes, \(00101110\) to dump the memory (for debug purposes) and \(00101111\) to stop execution. The full binary tree of the code – where the corresponding symbols are * and ! – and the arithmetic representation are displayed in figure .

The Spoon encoding of the “Hello World!” program is a sequence of 245 bits:

11111111110010001011111110101111111111010111010101
10110110110000011010110010100101001010111111100101
00010101110010100101100101001101111111111111111100
10100100010101110010100000000000000000000010100000
000000000000000000000010100101001010010001010

As the commands are exactly \(8\), it is also tempting to define a 3-bit fixed-length encoding of the original symbols and let’s call it Fork ! We adopt the mapping:

The encoding with Fork of the “Hello World!” program is

01001001001001001001001001001011000001001001001001
00100100000100100100100100100100100100100000100100
10000010001001001001011111000010010100000010100010
01001001001001001010010001001001010000001001010000
10010100100100100100100100100100100100100100100101
00000100010010010100011011011011011011100011011011
011011011011011100000010100000100

a stream of 333 binary digits, about 36 % less space-efficient than spoon.

The nodes \(\{3,4,\cdots\}\) and \(2\) have the lowest probability of occurence and therefore shall be grouped. But their probability is the same – \(0.125\) – and we need to make the a decision about which one will be taken as the first node of the group. We decide to go for the infinite set of symbols \(\{3,4,\cdots\}\) first. Somehow it feels the right thing to do, to have the only compound symbol on the left and all the others on the right …

This group has a cumulative probability of \(0.25\), the same as the symbol \(1\). There is only one extra symbol, \(0\) and its probability is higher, therefore no reordering is necessary. We therefore group \(\{2,3,\cdots\}\) (first) and \(1\).

There are only two symbols left, so the algorithm has completed. Let’s look at the results in terms of coding: we ended up with the code \[ 0 \rightarrow 1, \; 1 \rightarrow 01, \; 2 \rightarrow 001, \{3,4,\cdots\} \rightarrow 000 \] From this we can guess what the optimal coding is for the original source, without the reduction attached to the threshold : it is the unary coding of the integer, or more precisely here, the variant that uses \(0\) for the length and \(1\) as an end symbol. \[ 0 \rightarrow 1, \; 1 \rightarrow 01, \; 2 \rightarrow 001, \; 3 \rightarrow 0001, \; 4 \rightarrow 00001, \; 5 \rightarrow 000001, \; \cdots \]

The choice of the variant has no impact on the optimality of the code: unary coding is the optimal solution for the coding of non-negative integers that occur with a probability \(p(i) = 2^{-i-1}\). Here is an implementation of unary coding of integer symbols that uses \(0\) as an end delimiter instead:

def unary_symbol_encoder(stream, symbol):
    return stream.write(symbol * [True] + [False], bool)

def unary_symbol_decoder(stream):
    count = 0
    while stream.read(bool) is True:
        count += 1
    return count

Given those symbol encoder and decoder function, the higher-order function stream_encoder and stream_decoder from the coding module generate a stream encoder and decoder respectively.

unary_encoder = stream_encoder(unary_symbol_encoder)
unary_decoder = stream_decoder(unary_symbol_decoder)

Finally, we define an empty class unary as a type used to register the encoder and decoder in the bitstream module.

class unary(object):
    pass
bitstream.register(unary, reader=unary_decoder, writer=unary_encoder)

After that last step, using unary coding and decoding is as simple as:

>>> stream = BitStream()
>>> stream.write([0,1,2,3], unary)
>>> print stream
0101101110
>>> stream.read(unary, 4)
[0, 1, 2, 3]

• the unary coding of \(\lfloor i/l \rfloor\), -the Huffman coding of \(i \mbox{ mod } l\).

The probability distribution needed to perform the second part of the encoding is: \[ P(i \mbox{ mod } l) = \frac{(1-\theta)\theta^i}{1-\theta^l}. \] Note that this distribution is quite flat with respect to the original one. For that reason, the length of the optimal coding are almost constant for the Huffman part. Precisely, the length of the coding of \(0 \leq i < l\) is \(\lfloor \log_2 l \rfloor\) (if \(i< 2^{\lfloor \log_2 l + 1\rfloor} - l\)) or \(\lfloor \log_2 l \rfloor +1\) (otherwise). In particular, if \(l\) is a power of two, we have \(2^{\lfloor \log_2 l + 1\rfloor} - l = l\) is therefore every \(0 \leq i < l\) is coded with the same length of \(\log_2 l\).

Rice Coding

The parameter n in the rice constructor is the Golomb parameter. The optional signed option may be used to enable the coding of negative integers.

class rice(object):
    def __init__(self, n, signed=False):
        self.n = n
        self.signed = signed

    @staticmethod
    def select_parameter(mean):
        golden_ratio = 0.5 * (1.0 + numpy.sqrt(5))
        theta = mean / (mean + 1.0)
        log_ratio = log(golden_ratio - 1.0) / log(theta))
        return int(maximum(0, 1 + floor(log2(log_ratio))))

The following encoder and decoder implementations demonstrate how we deal with signed values: by prefixing the code stream with a bit sign (0 for \(+\), 1 for \(-\)) before we encode the absolute value of the integer. More complex schemes – that intertwin negative and positive values – are possible and may be useful to deal with two-sided geometric distribution that are not centered around \(0\) (see for example ). We then encode the fixed-length part of the code and follow with the unary code.

In the following code, options is meant to be an instance of the rice class, stream is an instance of BitStream and symbol an integer.

def rice_symbol_encoder(options):
    def encoder(stream, symbol):
        if options.signed:
            stream.write(symbol < 0)
        symbol = abs(symbol)
        remain, fixed = divmod(symbol, 2 ** options.n)
        fixed_bits = []
        for _ in range(options.n):
            fixed_bits.insert(0, bool(fixed % 2))
            fixed = fixed >> 1
        stream.write(fixed_bits)
        stream.write(remain, unary)
    return encoder


def rice_symbol_decoder(options):
    def decoder(stream):
        if options.signed and stream.read(bool):
            sign = -1
        else:
            sign = 1
        fixed_number = 0
        for _ in range(options.n):
            fixed_number = (fixed_number << 1) + int(stream.read(bool))
        remain_number = 2 ** options.n * stream.read(unary)
        return sign * (fixed_number + remain_number)
    return decoder

Finally, the rice class, its encoder and decoder are registered for integration with the BitStream instances.

rice_encoder = lambda r: stream_encoder(rice_symbol_encoder(r))
rice_decoder = lambda r: stream_decoder(rice_symbol_decoder(r))
bitstream.register(rice, reader=rice_decoder, writer=rice_encoder)

The following Python session demonstrates the basic usage:

>>> data = [0, 8, 0, 8, 16, 0, 32, 0, 16, 8, 0, 8]
>>> rice.select_parameter(mean(data))
3
>>> stream = BitStream()
>>> stream.write(data, rice(3))
>>> stream
000000010000000010000110000000011110000000011000010000000010
>>> stream.read(rice(3), 12)
[0, 8, 0, 8, 16, 0, 32, 0, 16, 8, 0, 8]

On this particular data, the selection of the Golomb parameter was effective if we consider the following test of possible parameter values between \(0\) (unary coding) and \(6\).

>>> for i in range(7):
...     stream = BitStream(data, rice(i))
...     print "rice n={0}: {1} bits".format(i, len(stream)) 
... 
rice n=0: 108 bits
rice n=1: 72 bits
rice n=2: 60 bits
rice n=3: 60 bits
rice n=4: 64 bits
rice n=5: 73 bits
rice n=6: 84 bits