mathematicskit.information_theory#

mathematicskit.information_theory hub · History · Examples · API reference

Hartley’s measure, Shannon entropy, joint/conditional entropy, mutual information and the Kullback-Leibler divergence via scipy.stats.entropy; the typical set, the Kraft-McMillan inequality, Shannon-Fano, Huffman, Elias gamma, arithmetic and LZ78 codes; channel capacity (binary symmetric, erasure, Gaussian, and Blahut-Arimoto for any discrete memoryless channel), random codes and rate-distortion; and Hamming, convolutional (Viterbi) and Gallager LDPC codes with their decoders.

mathematicskit.information_theory: entropy, source coding, channel capacity, and error-correcting codes.

Hartley’s measure, Shannon entropy, joint/conditional entropy, mutual information, and the Kullback-Leibler divergence via scipy.stats.entropy/scipy.special.entr; the typical set behind Shannon’s source coding theorem, the Kraft-McMillan inequality, Shannon-Fano and Huffman codes, Elias gamma codes, arithmetic coding, and LZ78; channel capacity of the binary symmetric, binary erasure and Gaussian (Shannon-Hartley) channels, the Blahut-Arimoto algorithm, a random-coding illustration of the noisy-channel coding theorem, and rate-distortion functions; Hamming codes, convolutional codes with Viterbi decoding, and Gallager’s LDPC codes with bit-flipping decoding – the coders and decoders hand-rolled, with no scipy/numpy equivalent. Reed-Solomon codes live in mathematicskit.abstract_algebra.

class mathematicskit.information_theory.BitFlipResult(codeword, iterations, converged)[source]#

Bases: object

Container for Gallager’s bit-flipping decoder.

Parameters:
codeword: ndarray#

The decoded word (a codeword when converged).

Type:

ndarray

converged: bool#

Whether every parity check is satisfied.

Type:

bool

iterations: int#

Number of flipping rounds performed.

Type:

int

class mathematicskit.information_theory.ChannelCapacityResult(capacity, input_distribution, iterations, converged, lower_bounds=<factory>, upper_bounds=<factory>)[source]#

Bases: object

Container for a channel-capacity computation (Blahut-Arimoto).

Parameters:
capacity: float#

Capacity in bits per channel use.

Type:

float

converged: bool#

Whether the upper and lower capacity bounds met within the tolerance.

Type:

bool

input_distribution: ndarray#

The capacity-achieving input distribution.

Type:

ndarray

iterations: int#

Number of iterations actually performed.

Type:

int

lower_bounds: ndarray#

Lower bound on the capacity at each iteration, in bits.

Type:

ndarray

upper_bounds: ndarray#

Upper bound on the capacity at each iteration, in bits.

Type:

ndarray

class mathematicskit.information_theory.PrefixCode(codewords=<factory>, probabilities=<factory>)[source]#

Bases: object

A binary prefix code: no codeword is the beginning of another.

Parameters:

  • codewords (dict) – {symbol: codeword}, each codeword a string of "0"/"1".

  • probabilities (dict) – {symbol: probability} for the source the code was built for.

Examples

>>> code = PrefixCode({"a": "0", "b": "10", "c": "11"}, {"a": 0.5, "b": 0.25, "c": 0.25})
>>> code.average_length, code.entropy
(1.5, 1.5)
>>> code.encode("abca")
'010110'
>>> code.decode("010110")
['a', 'b', 'c', 'a']
property average_length: float#

Expected codeword length \(\bar L = \sum_x p(x)\,\ell(x)\) in bits per symbol.

Type:

float

codewords: dict#
decode(bits)[source]#

Split bits into codewords; the prefix property makes this unambiguous.

Parameters:

bits (str)

Return type:

list

Returns:

list – The decoded symbols.

property efficiency: float#

entropy / average_length, at most 1.

Type:

float

encode(message)[source]#

Concatenate the codewords of the symbols in message.

Parameters:

message (iterable) – Symbols of the source alphabet.

Return type:

str

Returns:

str

property entropy: float#

Entropy of the source in bits, the lower bound on average_length.

Type:

float

probabilities: dict#
class mathematicskit.information_theory.TypicalSetResult(n, epsilon, entropy, log2_size, probability)[source]#

Bases: object

Container for the \(\varepsilon\)-typical set of an i.i.d. source.

Parameters:
entropy: float#

Entropy \(H\) of one source symbol, in bits.

Type:

float

epsilon: float#

Tolerance on \(-\tfrac1n \log_2 p(x^n)\) around the entropy.

Type:

float

log2_size: float#

\(\log_2\) of the number of typical sequences.

Type:

float

n: int#

Sequence length.

Type:

int

probability: float#

Total probability of the typical set.

Type:

float

mathematicskit.information_theory.arithmetic_decode(bits, probabilities, length)[source]#

Decode length symbols from an arithmetic codeword.

Parameters:

  • bits (str) – From arithmetic_encode().

  • probabilities (dict or sequence) – The model used to encode.

  • length (int) – Number of symbols in the message.

Return type:

list

Returns:

list

Examples

>>> "".join(arithmetic_decode("01011", {"a": 0.75, "b": 0.25}, 4))
'aaab'
mathematicskit.information_theory.arithmetic_encode(message, probabilities)[source]#

Arithmetic-code message: the shortest binary fraction whose dyadic interval fits in the message’s interval.

The code length is at most \(\lceil -\log_2 p(\text{message})\rceil + 1\) bits, within 2 bits of the ideal for the whole message rather than within 1 bit per symbol. J. Rissanen, “Generalized Kraft Inequality and Arithmetic Coding,” IBM Journal of Research and Development 20(3) (1976), 198-203; R. C. Pasco, Source Coding Algorithms for Fast Data Compression, PhD thesis, Stanford (1976).

Parameters:

  • message (sequence)

  • probabilities (dict or sequence)

Return type:

str

Returns:

str

Examples

>>> arithmetic_encode("aaab", {"a": 0.75, "b": 0.25})  # p = 0.105, 5 bits
'01011'
mathematicskit.information_theory.arithmetic_intervals(message, probabilities)[source]#

The nested intervals of arithmetic coding, one per symbol of message.

Each symbol narrows the current interval \([\ell, h)\) to the sub-interval its probability occupies, so the final width is \(p(\text{message})\). Exact rational arithmetic (fractions.Fraction) avoids the renormalization that practical coders need.

Parameters:

  • message (sequence)

  • probabilities (dict or sequence)

Return type:

list

Returns:

list of tuple of Fraction – [(0, 1), (low_1, high_1), ...], of length len(message) + 1.

Examples

>>> [(float(l), float(h)) for l, h in arithmetic_intervals("ab", {"a": 0.5, "b": 0.5})]
[(0.0, 1.0), (0.0, 0.5), (0.25, 0.5)]
mathematicskit.information_theory.awgn_capacity(snr, bandwidth=1.0)[source]#

Shannon-Hartley capacity \(C = B\log_2(1 + S/N)\) of a band-limited Gaussian channel.

C. E. Shannon, “Communication in the Presence of Noise,” Proceedings of the IRE 37(1) (1949), 10-21, Theorem 2.

Parameters:

  • snr (float or array_like) – Signal-to-noise power ratio \(S/N\) (linear, not dB).

  • bandwidth (float) – \(B\) in Hz; the default 1 gives bits/s/Hz.

Returns:

float or ndarray – Bits per second.

Examples

>>> awgn_capacity(1.0)
1.0
>>> round(awgn_capacity(10 ** (30 / 10), bandwidth=3000))  # a 30 dB telephone line
29902
mathematicskit.information_theory.bec_capacity(epsilon)[source]#

Capacity \(C = 1 - \varepsilon\) of the binary erasure channel.

P. Elias, “Coding for Two Noisy Channels,” Information Theory, Third London Symposium (1955), 61-76.

Parameters:

epsilon (float or array_like) – Erasure probability.

Returns:

float or ndarray

Examples

>>> bec_capacity(0.25)
0.75
mathematicskit.information_theory.binary_entropy(p)[source]#

The binary entropy function \(h(p) = -p\log_2 p - (1-p)\log_2(1-p)\).

Parameters:

p (float or array_like) – Probability in \([0, 1]\).

Returns:

float or ndarray – In bits.

Examples

>>> binary_entropy(0.5)
1.0
>>> round(binary_entropy(0.11), 4)
0.4999
mathematicskit.information_theory.bit_flip_decode(H, received, max_iter=50)[source]#

Gallager’s bit-flipping decoder for a low-density parity-check code.

While some parity check fails, flip every bit that takes part in the largest number of failed checks (Gallager 1963, Section 4.1, decoding algorithm A in its simplest hard-decision form).

Parameters:

  • H (array_like of {0, 1}, shape (m, n)) – Parity-check matrix.

  • received (array_like of {0, 1}, shape (n,))

  • max_iter (int)

Return type:

BitFlipResult

Returns:

BitFlipResult

Examples

>>> H = hamming_parity_check_matrix(3)
>>> r = bit_flip_decode(H, [0, 0, 0, 0, 0, 0, 1])
>>> r.codeword.tolist(), r.converged
([0, 0, 0, 0, 0, 0, 0], True)
mathematicskit.information_theory.blahut_arimoto(transition, tol=1e-10, max_iter=1000)[source]#

Capacity of a discrete memoryless channel by the Blahut-Arimoto algorithm.

Alternates between the output distribution \(q(y) = \sum_x r(x)W(y\mid x)\) and the input update \(r(x) \propto r(x)\exp D\bigl(W(\cdot\mid x)\,\|\,q\bigr)\). At every step \(\log\sum_x r(x)e^{D_x} \le C \le \max_x D_x\), and the iteration stops when the bounds meet within tol. R. Blahut, “Computation of Channel Capacity and Rate-Distortion Functions,” IEEE Transactions on Information Theory 18(4) (1972), 460-473; S. Arimoto, “An Algorithm for Computing the Capacity of Arbitrary Discrete Memoryless Channels,” ibid., 14-20.

Parameters:

  • transition (array_like, shape (m, n)) – transition[x, y] \(= W(y\mid x)\); rows sum to 1.

  • tol (float) – Gap between the capacity bounds, in bits.

  • max_iter (int)

Return type:

ChannelCapacityResult

Returns:

ChannelCapacityResult

Examples

>>> r = blahut_arimoto([[0.9, 0.1], [0.1, 0.9]])
>>> round(r.capacity, 6), r.input_distribution.round(3).tolist()
(0.531004, [0.5, 0.5])
>>> round(blahut_arimoto([[1, 0], [0.5, 0.5]]).capacity, 4)  # the Z channel
0.3219
mathematicskit.information_theory.bsc_capacity(p)[source]#

Capacity \(C = 1 - h(p)\) of the binary symmetric channel with crossover probability \(p\).

Shannon (1948), Section 13: each bit is flipped independently with probability \(p\), and a uniform input achieves capacity.

Parameters:

p (float or array_like)

Returns:

float or ndarray – Bits per channel use.

Examples

>>> bsc_capacity(0.0), bsc_capacity(0.5)
(1.0, 0.0)
>>> round(bsc_capacity(0.11), 3)
0.5
mathematicskit.information_theory.bsc_transmit(bits, p, seed=None)[source]#

Send bits through a binary symmetric channel that flips each bit with probability p.

Parameters:
Return type:

ndarray

Returns:

ndarray of uint8

Examples

>>> bsc_transmit([0, 1, 1, 0], p=0.0).tolist()
[0, 1, 1, 0]
>>> bsc_transmit([0, 1, 1, 0], p=1.0).tolist()
[1, 0, 0, 1]
mathematicskit.information_theory.canonical_code(lengths)[source]#

Binary prefix codewords with the given lengths – the constructive half of Kraft’s inequality.

Assigns codewords in order of increasing length, each the next binary number after the previous one, shifted to its length (the canonical Huffman code of Schwartz and Kallick, 1964).

Parameters:

lengths (sequence of int) – Codeword lengths with Kraft sum at most 1.

Return type:

list

Returns:

list of str – Codewords in the order of lengths.

Examples

>>> canonical_code([2, 1, 3, 3])
['10', '0', '110', '111']
mathematicskit.information_theory.conditional_entropy(joint, base=2)[source]#

Conditional entropy \(H(Y\mid X) = H(X,Y) - H(X)\), with \(X\) indexing the rows.

Parameters:

  • joint (array_like, shape (m, n)) – joint[i, j] = P(X=i, Y=j).

  • base (float)

Return type:

float

Returns:

float

Examples

>>> conditional_entropy([[0.5, 0.0], [0.0, 0.5]])  # Y = X
0.0
mathematicskit.information_theory.convolutional_encode(bits, generators=(7, 5), constraint_length=3)[source]#

Rate-\(1/m\) feed-forward convolutional encoder, terminated with \(K-1\) zeros.

Each input bit shifts into a \(K\)-bit register (newest bit most significant); output \(j\) is the parity of the register masked by generators[j]. The default is the \((7, 5)_8\), \(K=3\) code with free distance 5. P. Elias, “Coding for Noisy Channels,” IRE Convention Record 3(4) (1955), 37-46.

Parameters:

  • bits (array_like of {0, 1})

  • generators (sequence of int) – Tap masks, conventionally written in octal.

  • constraint_length (int) – \(K\).

Return type:

ndarray

Returns:

ndarray of uint8 – Length len(generators) * (len(bits) + K - 1).

Examples

>>> convolutional_encode([1, 0, 1, 1]).tolist()
[1, 1, 1, 0, 0, 0, 0, 1, 0, 1, 1, 1]
mathematicskit.information_theory.elias_gamma_decode(bits)[source]#

Decode a concatenation of Elias gamma codewords.

Parameters:

bits (str)

Return type:

list

Returns:

list of int

Examples

>>> elias_gamma_decode("1" "010" "00101")
[1, 2, 5]
mathematicskit.information_theory.elias_gamma_encode(n)[source]#

Elias’ gamma code of a positive integer: \(\lfloor\log_2 n\rfloor\) zeros, then \(n\) in binary.

A prefix code for all positive integers needing no prior bound on their size, with length \(2\lfloor\log_2 n\rfloor + 1\). P. Elias, “Universal Codeword Sets and Representations of the Integers,” IEEE Transactions on Information Theory 21(2) (1975), 194-203.

Parameters:

n (int) – Positive integer.

Return type:

str

Returns:

str

Examples

>>> [elias_gamma_encode(n) for n in (1, 2, 5, 9)]
['1', '010', '00101', '0001001']
mathematicskit.information_theory.entropy(p, base=2)[source]#

Shannon entropy \(H(X) = -\sum_x p(x)\log p(x)\).

C. E. Shannon, “A Mathematical Theory of Communication,” Bell System Technical Journal 27 (1948), 379-423, Theorem 2. Computed by scipy.stats.entropy().

Parameters:

  • p (array_like) – Probabilities (normalized if they don’t sum to 1). A 2-D array is treated as one joint distribution.

  • base (float) – Logarithm base; 2 gives bits.

Return type:

float

Returns:

float

Examples

>>> entropy([0.5, 0.25, 0.25])
1.5
>>> entropy([1, 0, 0])
0.0
mathematicskit.information_theory.gallager_ldpc_matrix(n, column_weight=3, row_weight=6, seed=None)[source]#

Gallager’s random regular low-density parity-check matrix.

The first band of \(n/w_r\) rows has \(w_r\) consecutive ones each; the other \(w_c - 1\) bands are random column permutations of it. Every column then has weight \(w_c\) and every row weight \(w_r\), for a design rate of \(1 - w_c/w_r\). R. G. Gallager, “Low-Density Parity-Check Codes,” IRE Transactions on Information Theory 8(1) (1962), 21-28; Low-Density Parity-Check Codes (MIT Press, 1963).

Parameters:
Return type:

ndarray

Returns:

ndarray of uint8, shape (n * column_weight / row_weight, n)

Examples

>>> H = gallager_ldpc_matrix(12, column_weight=3, row_weight=4, seed=0)
>>> H.shape, set(H.sum(axis=0).tolist()), set(H.sum(axis=1).tolist())
((9, 12), {3}, {4})
mathematicskit.information_theory.hamming_decode(received, r=3)[source]#

Correct up to one error per block by syndrome decoding and return the message bits.

Parameters:

  • received (array_like of {0, 1}) – Length a multiple of \(2^r - 1\).

  • r (int)

Return type:

ndarray

Returns:

ndarray of uint8

Examples

>>> word = hamming_encode([1, 0, 1, 1])
>>> word[4] ^= 1  # one error
>>> hamming_decode(word).tolist()
[1, 0, 1, 1]
mathematicskit.information_theory.hamming_distance(a, b)[source]#

Number of positions where a and b differ.

Parameters:

  • a (array_like) – Equal-length sequences.

  • b (array_like) – Equal-length sequences.

Return type:

int

Returns:

int

Examples

>>> hamming_distance([1, 0, 1, 1], [1, 1, 1, 0])
2
mathematicskit.information_theory.hamming_encode(message, r=3)[source]#

Encode message with the Hamming code, block by block.

Data bits fill the positions that are not powers of two; the parity bit at position \(2^i\) makes every check (row of hamming_parity_check_matrix()) even. R. W. Hamming, “Error Detecting and Error Correcting Codes,” Bell System Technical Journal 29(2) (1950), 147-160.

Parameters:

  • message (array_like of {0, 1}) – Length a multiple of \(k = 2^r - r - 1\).

  • r (int)

Return type:

ndarray

Returns:

ndarray of uint8 – Length len(message) / k * (2**r - 1).

Examples

>>> hamming_encode([1, 0, 1, 1]).tolist()
[0, 1, 1, 0, 0, 1, 1]
mathematicskit.information_theory.hamming_parity_check_matrix(r=3)[source]#

Parity-check matrix of the \((2^r-1,\ 2^r-r-1)\) Hamming code.

Column \(j\) (1-based) is the binary expansion of \(j\), least significant bit in row 0, so the syndrome of a single error spells out its position.

Parameters:

r (int) – Number of parity bits, at least 2.

Return type:

ndarray

Returns:

ndarray of uint8, shape (r, 2**r - 1)

Examples

>>> hamming_parity_check_matrix(3)
array([[1, 0, 1, 0, 1, 0, 1],
       [0, 1, 1, 0, 0, 1, 1],
       [0, 0, 0, 1, 1, 1, 1]], dtype=uint8)
mathematicskit.information_theory.hartley_information(n_messages, length=1, base=2)[source]#

Hartley’s measure \(H_0 = \text{length}\cdot\log(n)\) of a choice among equally likely messages.

R. V. L. Hartley, “Transmission of Information,” Bell System Technical Journal 7(3) (1928), 535-563: a sequence of length selections from an alphabet of n_messages symbols carries \(\text{length}\log n\) units of information.

Parameters:

  • n_messages (int or array_like) – Alphabet size \(n \ge 1\).

  • length (int) – Number of selections.

  • base (float) – Logarithm base: 2 gives bits, 10 gives hartleys.

Return type:

float

Returns:

float or ndarray

Examples

>>> round(hartley_information(26), 4)  # one letter
4.7004
>>> hartley_information(2, length=8)  # one byte
8.0
mathematicskit.information_theory.huffman_code(probabilities)[source]#

Huffman’s optimal binary prefix code.

Repeatedly merge the two least likely nodes into one whose probability is their sum; reading the merges back from the root gives each symbol its codeword. No prefix code has a smaller average length, which lies in \([H, H+1)\). D. A. Huffman, “A Method for the Construction of Minimum-Redundancy Codes,” Proceedings of the IRE 40(9) (1952), 1098-1101.

Parameters:

probabilities (dict or sequence)

Return type:

PrefixCode

Returns:

PrefixCode

Examples

>>> code = huffman_code({"a": 0.4, "b": 0.3, "c": 0.2, "d": 0.1})
>>> sorted((s, len(w)) for s, w in code.codewords.items())
[('a', 1), ('b', 2), ('c', 3), ('d', 3)]
>>> round(code.average_length, 2)
1.9
mathematicskit.information_theory.joint_entropy(joint, base=2)[source]#

Joint entropy \(H(X,Y)\) of a joint distribution joint[i, j] = P(X=i, Y=j).

Parameters:

  • joint (array_like, shape (m, n))

  • base (float)

Return type:

float

Returns:

float

Examples

>>> joint_entropy([[0.25, 0.25], [0.25, 0.25]])
2.0
mathematicskit.information_theory.kl_divergence(p, q, base=2)[source]#

Kullback-Leibler divergence \(D(p\,\|\,q) = \sum_x p(x)\log\frac{p(x)}{q(x)}\).

S. Kullback and R. A. Leibler, “On Information and Sufficiency,” Annals of Mathematical Statistics 22(1) (1951), 79-86. Non-negative (Gibbs’ inequality), zero only when \(p = q\), and not symmetric. Infinite when \(q(x) = 0 < p(x)\). Computed by scipy.stats.entropy().

Parameters:

  • p (array_like) – Distributions on the same alphabet.

  • q (array_like) – Distributions on the same alphabet.

  • base (float)

Return type:

float

Returns:

float

Examples

>>> kl_divergence([0.5, 0.5], [0.5, 0.5])
0.0
>>> round(kl_divergence([0.5, 0.5], [0.9, 0.1]), 4), round(kl_divergence([0.9, 0.1], [0.5, 0.5]), 4)
(0.737, 0.531)
mathematicskit.information_theory.kraft_sum(lengths, radix=2)[source]#

The Kraft sum \(\sum_i D^{-\ell_i}\) of a list of codeword lengths.

Kraft (1949) showed a \(D\)-ary prefix code with lengths \(\ell_i\) exists if and only if the sum is at most 1; McMillan (1956) extended the necessity to every uniquely decodable code. See L. G. Kraft, A Device for Quantizing, Grouping, and Coding Amplitude-Modulated Pulses, MS thesis, MIT (1949), and B. McMillan, “Two Inequalities Implied by Unique Decipherability,” IRE Transactions on Information Theory 2(4) (1956), 115-116.

Parameters:

  • lengths (sequence of int)

  • radix (int) – Code alphabet size \(D\).

Return type:

float

Returns:

float

Examples

>>> kraft_sum([1, 2, 3, 3])
1.0
>>> kraft_sum([1, 1, 2])  # no binary prefix code has these lengths
1.25
mathematicskit.information_theory.lz78_compressed_bits(pairs, alphabet_size)[source]#

Length in bits of an LZ78 parsing: phrase \(i\) costs \(\lceil\log_2 i\rceil + \lceil\log_2 k\rceil\) bits.

Parameters:
Return type:

int

Returns:

int

Examples

>>> lz78_compressed_bits(lz78_encode("ABBABBABBBAABABAA"), alphabet_size=2)
25
mathematicskit.information_theory.lz78_decode(pairs)[source]#

Invert lz78_encode().

Parameters:

pairs (list of tuple)

Return type:

list

Returns:

list

Examples

>>> "".join(lz78_decode(lz78_encode("ABBABBABBBAABABAA")))
'ABBABBABBBAABABAA'
mathematicskit.information_theory.lz78_encode(sequence)[source]#

Lempel-Ziv (LZ78) parsing into (index, symbol) pairs.

Each phrase is the longest previously seen phrase (index into the dictionary, 0 for the empty phrase) extended by one new symbol. For a stationary ergodic source the compressed length per symbol tends to the entropy rate without knowing the source statistics. J. Ziv and A. Lempel, “Compression of Individual Sequences via Variable-Rate Coding,” IEEE Transactions on Information Theory 24(5) (1978), 530-536.

Parameters:

sequence (sequence)

Return type:

list

Returns:

list of tuple – (index, symbol); the last pair’s symbol is None when the sequence ends inside an already-known phrase.

Examples

>>> lz78_encode("ABBABBABBBAABABAA")
[(0, 'A'), (0, 'B'), (2, 'A'), (2, 'B'), (1, 'B'), (4, 'A'), (5, 'A'), (3, 'A')]
>>> lz78_encode("ABA")  # ends inside the known phrase "A"
[(0, 'A'), (0, 'B'), (1, None)]
mathematicskit.information_theory.minimum_ebn0(spectral_efficiency)[source]#

Smallest energy per bit \(E_b/N_0 = (2^\eta - 1)/\eta\) that supports \(\eta\) bits/s/Hz.

As \(\eta \to 0\) this tends to \(\ln 2\), the Shannon limit of \(-1.59\) dB below which no code communicates reliably.

Parameters:

spectral_efficiency (float or array_like) – \(\eta = C/B > 0\).

Returns:

float or ndarray – Linear ratio (10 * log10 for dB).

Examples

>>> minimum_ebn0(1.0)
1.0
>>> round(float(10 * np.log10(minimum_ebn0(1e-9))), 2)
-1.59
mathematicskit.information_theory.mutual_information(joint, base=2)[source]#

Mutual information \(I(X;Y) = H(X) + H(Y) - H(X,Y)\).

Shannon (1948), Section 20: the rate of transmission of a noisy channel is \(H(X) - H(X\mid Y)\). It is zero exactly when \(X\) and \(Y\) are independent.

Parameters:

  • joint (array_like, shape (m, n)) – joint[i, j] = P(X=i, Y=j).

  • base (float)

Return type:

float

Returns:

float

Examples

>>> mutual_information([[0.5, 0.0], [0.0, 0.5]])
1.0
>>> mutual_information(np.outer([0.3, 0.7], [0.6, 0.4]))  # independent
0.0
mathematicskit.information_theory.random_code_error_rate(n, rate, p, trials=1000, seed=None)[source]#

Block error rate of random codes with maximum-likelihood decoding on a binary symmetric channel.

Shannon’s proof of the noisy-channel coding theorem (1948, Theorem 11) draws \(2^{nR}\) codewords of length \(n\) uniformly at random. Averaged over that ensemble, the error probability of minimum-Hamming-distance decoding tends to 0 as \(n\) grows when \(R < C = 1 - h(p)\), and to 1 when \(R > C\). Each trial draws a fresh codebook, sends its first codeword, and breaks ties in distance at random.

Parameters:

  • n (int) – Block length.

  • rate (float) – Code rate \(R\); the codebook has \(\mathrm{round}(2^{nR})\) codewords.

  • p (float) – Crossover probability.

  • trials (int)

  • seed (int or numpy.random.Generator, optional)

Return type:

float

Returns:

float – Estimated block error probability.

Examples

>>> random_code_error_rate(16, 0.25, p=0.0, trials=200, seed=0) < 0.05
True
mathematicskit.information_theory.rate_distortion_binary(D, p=0.5)[source]#

Rate-distortion function \(R(D) = h(p) - h(D)\) of a Bernoulli(\(p\)) source under Hamming distortion.

Zero for \(D \ge \min(p, 1-p)\). C. E. Shannon, “Coding Theorems for a Discrete Source with a Fidelity Criterion,” IRE National Convention Record 7(4) (1959), 142-163.

Parameters:

  • D (float or array_like) – Allowed fraction of bit errors.

  • p (float)

Returns:

float or ndarray – Bits per source symbol.

Examples

>>> rate_distortion_binary(0.0), rate_distortion_binary(0.5)
(1.0, 0.0)
>>> round(rate_distortion_binary(0.11), 3)
0.5
mathematicskit.information_theory.rate_distortion_gaussian(D, variance=1.0)[source]#

Rate-distortion function \(R(D) = \max\bigl(0, \tfrac12\log_2(\sigma^2/D)\bigr)\) of a Gaussian source under squared error.

Shannon (1959). Each extra bit per sample divides the achievable mean squared error by 4 (6.02 dB).

Parameters:

  • D (float or array_like) – Allowed mean squared error, positive.

  • variance (float) – Source variance \(\sigma^2\).

Returns:

float or ndarray

Examples

>>> rate_distortion_gaussian(0.25)
1.0
>>> rate_distortion_gaussian(2.0)
0.0
mathematicskit.information_theory.shannon_fano_code(probabilities, method='fano')[source]#

A Shannon-Fano prefix code.

method="fano" (R. M. Fano, The Transmission of Information, MIT RLE Technical Report 65, 1949): sort by probability and split the list recursively where the two halves’ totals are most nearly equal, prefixing 0 and 1. method="shannon" (Shannon 1948, Section 9): give the \(i\)-th most likely symbol the first \(\lceil -\log_2 p_i\rceil\) bits of the binary expansion of \(\sum_{j<i} p_j\). Both are within one bit of the entropy, but neither is always optimal (compare huffman_code()).

Parameters:

  • probabilities (dict or sequence)

  • method (str)

Return type:

PrefixCode

Returns:

PrefixCode

Examples

>>> shannon_fano_code({"a": 0.4, "b": 0.3, "c": 0.2, "d": 0.1}).codewords
{'a': '0', 'b': '10', 'c': '110', 'd': '111'}
>>> shannon_fano_code({"a": 0.4, "b": 0.3, "c": 0.2, "d": 0.1}, method="shannon").codewords
{'a': '00', 'b': '01', 'c': '101', 'd': '1110'}
mathematicskit.information_theory.typical_set(p, n, epsilon)[source]#

Size and probability of the \(\varepsilon\)-typical set of an i.i.d. source.

A sequence \(x^n\) is typical when \(\bigl|-\tfrac1n\log_2 p(x^n) - H\bigr| \le \varepsilon\). The asymptotic equipartition property (Shannon 1948, Theorem 3) says the typical set has probability tending to 1 and about \(2^{nH}\) members, so \(nH\) bits suffice to index the likely sequences – the source coding theorem. Computed exactly by summing over sequence types (symbol counts) with multinomial coefficients, so the cost grows as \(\binom{n+k-1}{k-1}\), not \(k^n\).

Parameters:

  • p (array_like) – Symbol probabilities of the source (no zeros).

  • n (int) – Sequence length.

  • epsilon (float) – Tolerance in bits per symbol.

Return type:

TypicalSetResult

Returns:

TypicalSetResult

Examples

>>> r = typical_set([0.9, 0.1], n=100, epsilon=0.1)
>>> round(r.probability, 3), round(r.log2_size / 100, 3), round(r.entropy, 3)
(0.759, 0.529, 0.469)
mathematicskit.information_theory.viterbi_decode(received, generators=(7, 5), constraint_length=3)[source]#

Maximum-likelihood decoding of a terminated convolutional code by the Viterbi algorithm.

Dynamic programming over the \(2^{K-1}\)-state trellis: keep, for every state, the input path with the smallest Hamming distance to the received bits, then trace back from the all-zero final state. A. J. Viterbi, “Error Bounds for Convolutional Codes and an Asymptotically Optimum Decoding Algorithm,” IEEE Transactions on Information Theory 13(2) (1967), 260-269.

Parameters:

  • received (array_like of {0, 1}) – Hard-decision channel output of convolutional_encode().

  • generators (sequence of int)

  • constraint_length (int)

Return type:

ndarray

Returns:

ndarray of uint8 – The decoded input bits (termination bits removed).

Examples

>>> word = convolutional_encode([1, 0, 1, 1])
>>> word[[2, 7]] ^= 1  # two errors, far apart
>>> viterbi_decode(word).tolist()
[1, 0, 1, 1]