C. Hinsley

19 September 2024

The Lempel-Ziv (LZ) algorithm takes a binary string (though it generalizes to arbitrary finite alphabets, not just binary) and passes over it left-to-right, building a dictionary of newly observed substrings, and at the end counts the symbols in the dictionary to obtain a measure of complexity of the binary string. If you normalize the LZ complexity, multiplying by $\log_2(n)/n$, where $n$ is the length of the string, this complexity measure will converge for large $n$ to the entropy rate of a source channel which produced the string.

I was curious if there is a possibility that LZ complexity is sensitive to reversals of the binary string, because it does build a dictionary in a relatively chaotic way — perhaps indicative of a directional preference forward or backward — for sufficiently complex sequences. By sensitivity, I mean such that the discrepancy between the normalized Lempel-Ziv complexity of a binary string and its reversal fails to converge to zero for large $n$.

Below you can find a tabulation of these discrepancies:

Julia code used to calculate the data above:

using Base.Threads

# Lempel-Ziv 1976 complexity (no normalization).
function LZ76_complexity(sequence::Vector{Int})::Float64
    if isempty(sequence)
        return 0.0
    end

    complexity = 0
    i = 1
    n = length(sequence)
    # b = length(unique(sequence))  # Unused variable, can be removed if not needed

    while i <= n
        max_match_length = 0
        
        # Search for the longest match in the window (from the start up to position i-1).
        for j in 1:(i-1)
            match_length = 0
            while (i + match_length <= n) && (sequence[j + match_length] == sequence[i + match_length])
                match_length += 1
                # Prevent overlapping matches.
                if (j + match_length > i - 1)
                    break
                end
            end
            if match_length > max_match_length
                max_match_length = match_length
            end
        end
        
        if max_match_length > 0
            # Found a match; increment complexity and move the pointer.
            complexity += 1
            i += max_match_length + 1
        else
            # No match found; treat the current symbol as a new phrase.
            complexity += 1
            i += 1
        end
    end

    return complexity
end

# Find the binary string of length n (WLOG starting with 0) which has the
# greatest discrepancy in LZ76 complexity from its reverse string.
function maximal_discrepancy(n::Int)::Tuple{Vector{Int}, Int}
    if n < 1
        error("Length of the binary string must be at least 1.")
    end
    
    if n == 1
        return ([0], 0)
    end
    
    if n == 2
        return ([0, 0], 0)
    end

    max_discrepancy = -1.0
    max_string = Int[]

    total_strings = 2^(n-1)  # Since the first bit is fixed to 0.

    # Initialize per-thread storage for maximum discrepancy and corresponding string.
    n_threads = Threads.nthreads()
    per_thread_max_discrepancy = fill(-1.0, n_threads)
    per_thread_max_string = [Int[] for _ in 1:n_threads]

    Threads.@threads for k in 0:(total_strings - 1)
        # Generate the binary string starting with 0.
        bits = zeros(Int, n)
        bits[1] = 0  # First bit is 0.

        for bit_pos in 2:n
            bits[bit_pos] = (k >> (n - bit_pos)) & 1
        end
        
        bits_reversed = reverse(bits)

        # Check for palindromic strings.
        if bits == bits_reversed
            continue  # Skip palindromic strings as discrepancy is zero.
        end

        # If n is even, check for anti-palindromic strings (equal to reverse and inversion).
        if iseven(n) && bits == (1 .- bits_reversed)
            continue  # Skip anti-palindromic strings as discrepancy is zero.
        end

        # Compute complexities.
        c_original = LZ76_complexity(bits)
        c_reversed = LZ76_complexity(bits_reversed)

        # Compute discrepancy.
        discrepancy = abs(c_original - c_reversed)
        
        # Update thread-local maximum discrepancy and corresponding string if needed.
        thread_id = threadid()
        if discrepancy > per_thread_max_discrepancy[thread_id]
            per_thread_max_discrepancy[thread_id] = discrepancy
            per_thread_max_string[thread_id] = copy(bits)
        end
    end

    # After the parallel loop, find the global maximum discrepancy and the corresponding string.
    for thread_id in 1:n_threads
        if per_thread_max_discrepancy[thread_id] > max_discrepancy
            max_discrepancy = per_thread_max_discrepancy[thread_id]
            max_string = copy(per_thread_max_string[thread_id])
        end
    end

    return (max_string, Int(max_discrepancy))
end