import numpy as np
import matplotlib.pyplot as plt

def self_information(p, base=2):
    """
    计算事件的自信息
    
    参数:
        p: 事件发生的概率 (0 < p <= 1)
        base: 对数的底数，默认为2（比特）
    
    返回:
        自信息值
    """
    if p <= 0 or p > 1:
        raise ValueError("概率必须在(0, 1]范围内")
    return -np.log(p) / np.log(base)

# 计算不同概率事件的自信息
probabilities = [0.01, 0.1, 0.25, 0.5, 0.75, 0.9, 0.99]
info_values = [self_information(p) for p in probabilities]

print("概率 vs 自信息（比特）:")
for p, info in zip(probabilities, info_values):
    print(f"P={p:.2f}: {info:.4f} bits")

# 可视化
plt.figure(figsize=(10, 6))
p_range = np.linspace(0.01, 1, 100)
info_range = [self_information(p) for p in p_range]
plt.plot(p_range, info_range, 'b-', linewidth=2)
plt.xlabel('概率 P', fontsize=12)
plt.ylabel('自信息（比特）', fontsize=12)
plt.title('自信息与概率的关系', fontsize=14)
plt.grid(True, alpha=0.3)
plt.axhline(y=1, color='r', linestyle='--', label='1比特')
plt.legend()
plt.tight_layout()
plt.savefig('self_information.png', dpi=150)
plt.show()

概率 vs 自信息（比特）:
P=0.01: 6.6439 bits
P=0.10: 3.3219 bits
P=0.25: 2.0000 bits
P=0.50: 1.0000 bits
P=0.75: 0.4150 bits
P=0.90: 0.1520 bits
P=0.99: 0.0145 bits

import numpy as np
import matplotlib.pyplot as plt
from scipy.special import entr

def entropy(p, base=2):
    """
    计算离散分布的熵
    
    参数:
        p: 概率分布列表或数组
        base: 对数的底数，默认为2
    
    返回:
        熵值
    """
    p = np.array(p, dtype=float)
    # 过滤掉零概率（避免log(0)）
    p = p[p > 0]
    # 归一化确保概率和为1
    p = p / p.sum()
    # 计算熵
    return -np.sum(p * np.log(p) / np.log(base))

# 测试不同的分布
distributions = {
    "确定性": [1, 0, 0, 0],
    "均匀(4项)": [0.25, 0.25, 0.25, 0.25],
    "偏斜": [0.5, 0.25, 0.125, 0.125],
    "高度偏斜": [0.9, 0.05, 0.03, 0.02],
    "均匀(8项)": [0.125] * 8,
}

print("各种分布的熵:")
print("-" * 40)
for name, dist in distributions.items():
    h = entropy(dist)
    print(f"{name:12s}: H = {h:.4f} bits")

# 研究二元分布的熵随概率的变化
print("\n二元分布的熵:")
print("-" * 40)
binary_probs = [0.01, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.99]
for p in binary_probs:
    h = entropy([p, 1-p])
    print(f"P(X=0)={p:.2f}: H = {h:.4f} bits")

# 可视化二元分布的熵函数
plt.figure(figsize=(10, 6))
p_values = np.linspace(0.001, 0.999, 1000)
entropy_values = [entropy([p, 1-p]) for p in p_values]

plt.plot(p_values, entropy_values, 'b-', linewidth=2.5)
plt.xlabel('P(X=0)', fontsize=13)
plt.ylabel('熵 H(X)（比特）', fontsize=13)
plt.title('二元随机变量的熵函数', fontsize=14)
plt.grid(True, alpha=0.3)
plt.axvline(x=0.5, color='r', linestyle='--', alpha=0.5, label='最大值点')
plt.axhline(y=1, color='g', linestyle='--', alpha=0.5, label='最大熵=1比特')
plt.legend(fontsize=11)
plt.xlim(0, 1)
plt.ylim(0, 1.1)
plt.tight_layout()
plt.savefig('binary_entropy.png', dpi=150)
plt.show()

# 验证最大熵原理：均匀分布具有最大熵
print("\n验证最大熵原理（4个取值的情况）:")
print("-" * 40)
n = 4
max_entropy = np.log2(n)
print(f"理论最大熵（均匀分布）: {max_entropy:.4f} bits")

# 生成一些随机分布并计算熵
np.random.seed(42)
for i in range(5):
    random_dist = np.random.dirichlet(np.ones(n))
    h = entropy(random_dist)
    print(f"随机分布 {i+1}: {random_dist.round(3)}, H = {h:.4f} bits")

各种分布的熵:
----------------------------------------
确定性       : H = 0.0000 bits
均匀(4项)    : H = 2.0000 bits
偏斜         : H = 1.7500 bits
高度偏斜     : H = 0.4690 bits
均匀(8项)    : H = 3.0000 bits

二元分布的熵:
----------------------------------------
P(X=0)=0.01: H = 0.0808 bits
P(X=0)=0.10: H = 0.4690 bits
P(X=0)=0.20: H = 0.7219 bits
P(X=0)=0.30: H = 0.8813 bits
P(X=0)=0.40: H = 0.9710 bits
P(X=0)=0.50: H = 1.0000 bits
P(X=0)=0.60: H = 0.9710 bits
P(X=0)=0.70: H = 0.8813 bits
P(X=0)=0.80: H = 0.7219 bits
P(X=0)=0.90: H = 0.4690 bits
P(X=0)=0.99: H = 0.0808 bits

验证最大熵原理（4个取值的情况）:
----------------------------------------
理论最大熵（均匀分布）: 2.0000 bits
随机分布 1: [0.379 0.24  0.213 0.168], H = 1.9717 bits
随机分布 2: [0.118 0.33  0.308 0.244], H = 1.9623 bits
随机分布 3: [0.463 0.158 0.234 0.145], H = 1.8838 bits
随机分布 4: [0.186 0.282 0.354 0.178], H = 1.9768 bits
随机分布 5: [0.276 0.188 0.344 0.192], H = 1.9822 bits

import numpy as np

def joint_entropy(joint_dist):
    """
    计算联合分布的熵
    
    参数:
        joint_dist: 二维数组，表示联合概率分布 P(X,Y)
    
    返回:
        联合熵
    """
    joint_dist = np.array(joint_dist, dtype=float)
    # 展平并过滤零概率
    flat = joint_dist.flatten()
    flat = flat[flat > 0]
    return -np.sum(flat * np.log2(flat))

def conditional_entropy_from_joint(joint_dist, axis=0):
    """
    从联合分布计算条件熵 H(Y|X)
    
    参数:
        joint_dist: 二维联合概率分布
        axis: 条件变量的轴，0表示P(Y|X)，1表示P(X|Y)
    
    返回:
        条件熵
    """
    joint_dist = np.array(joint_dist, dtype=float)
    
    # 计算边缘分布
    if axis == 0:
        marginal_x = joint_dist.sum(axis=1)  # P(X)
    else:
        marginal_x = joint_dist.sum(axis=0)  # P(Y)
    
    # 计算条件熵
    cond_ent = 0
    for i in range(joint_dist.shape[axis]):
        if axis == 0:
            conditional_dist = joint_dist[i, :]
            px = marginal_x[i]
        else:
            conditional_dist = joint_dist[:, i]
            px = marginal_x[i]
        
        if px > 0:
            # 条件分布
            cond_prob = conditional_dist / px
            cond_prob = cond_prob[cond_prob > 0]
            # H(Y|X=i)
            h_y_given_x_i = -np.sum(cond_prob * np.log2(cond_prob))
            cond_ent += px * h_y_given_x_i
    
    return cond_ent

# 示例1：独立的两个二元变量
print("示例1：独立的二元变量")
print("-" * 40)
# P(X=0)=0.5, P(X=1)=0.5
# P(Y=0)=0.5, P(Y=1)=0.5
# 由于独立，P(X,Y) = P(X)P(Y)
joint_independent = np.array([
    [0.25, 0.25],  # X=0
    [0.25, 0.25],  # X=1
])
h_x = entropy([0.5, 0.5])
h_y = entropy([0.5, 0.5])
h_xy = joint_entropy(joint_independent)
h_y_given_x = conditional_entropy_from_joint(joint_independent, axis=0)

print(f"H(X) = {h_x:.4f} bits")
print(f"H(Y) = {h_y:.4f} bits")
print(f"H(X,Y) = {h_xy:.4f} bits")
print(f"H(Y|X) = {h_y_given_x:.4f} bits")
print(f"验证链式法则: H(X) + H(Y|X) = {h_x + h_y_given_x:.4f}")
print(f"验证独立性: H(X,Y) = H(X) + H(Y)? {np.isclose(h_xy, h_x + h_y)}")

# 示例2：相关的两个二元变量
print("\n示例2：相关的二元变量")
print("-" * 40)
# X和Y高度相关：大多数情况下X=Y
joint_correlated = np.array([
    [0.4, 0.1],  # X=0
    [0.1, 0.4],  # X=1
])
marginal_x = joint_correlated.sum(axis=1)
marginal_y = joint_correlated.sum(axis=0)

h_x = entropy(marginal_x)
h_y = entropy(marginal_y)
h_xy = joint_entropy(joint_correlated)
h_y_given_x = conditional_entropy_from_joint(joint_correlated, axis=0)

print(f"边缘分布 P(X) = {marginal_x}")
print(f"边缘分布 P(Y) = {marginal_y}")
print(f"H(X) = {h_x:.4f} bits")
print(f"H(Y) = {h_y:.4f} bits")
print(f"H(X,Y) = {h_xy:.4f} bits")
print(f"H(Y|X) = {h_y_given_x:.4f} bits")
print(f"验证链式法则: H(X) + H(Y|X) = {h_x + h_y_given_x:.4f}")
print(f"H(X,Y) < H(X) + H(Y)? {h_xy < h_x + h_y} (因为相关)")
print(f"节省的信息量: {h_x + h_y - h_xy:.4f} bits")

示例1：独立的二元变量
----------------------------------------
H(X) = 1.0000 bits
H(Y) = 1.0000 bits
H(X,Y) = 2.0000 bits
H(Y|X) = 1.0000 bits
验证链式法则: H(X) + H(Y|X) = 2.0000
验证独立性: H(X,Y) = H(X) + H(Y)? True

示例2：相关的二元变量
----------------------------------------
边缘分布 P(X) = [0.5 0.5]
边缘分布 P(Y) = [0.5 0.5]
H(X) = 1.0000 bits
H(Y) = 1.0000 bits
H(X,Y) = 1.7219 bits
H(Y|X) = 0.7219 bits
验证链式法则: H(X) + H(Y|X) = 1.7219
H(X,Y) < H(X) + H(Y)? True (因为相关)
节省的信息量: 0.2781 bits

import numpy as np

def mutual_information(joint_dist):
    """
    从联合分布计算互信息 I(X;Y)
    
    参数:
        joint_dist: 二维数组，表示联合概率分布 P(X,Y)
    
    返回:
        互信息值（比特）
    """
    joint_dist = np.array(joint_dist, dtype=float)
    
    # 计算边缘分布
    marginal_x = joint_dist.sum(axis=1)
    marginal_y = joint_dist.sum(axis=0)
    
    # 计算互信息
    mi = 0
    for i in range(joint_dist.shape[0]):
        for j in range(joint_dist.shape[1]):
            pxy = joint_dist[i, j]
            px = marginal_x[i]
            py = marginal_y[j]
            
            if pxy > 0 and px > 0 and py > 0:
                mi += pxy * np.log2(pxy / (px * py))
    
    return mi

# 示例1：独立变量
print("示例1：独立变量")
print("-" * 40)
joint_indep = np.array([[0.25, 0.25], [0.25, 0.25]])
mi = mutual_information(joint_indep)
print(f"I(X;Y) = {mi:.6f} bits (应为0)")

# 示例2：完全相关
print("\n示例2：完全相关 (Y=X)")
print("-" * 40)
joint_perfect = np.array([[0.5, 0], [0, 0.5]])
h_x = entropy([0.5, 0.5])
mi = mutual_information(joint_perfect)
print(f"H(X) = {h_x:.4f} bits")
print(f"I(X;Y) = {mi:.4f} bits")
print(f"I(X;Y) = H(X)? {np.isclose(mi, h_x)}")

# 示例3：部分相关
print("\n示例3：部分相关")
print("-" * 40)
joint_partial = np.array([[0.4, 0.1], [0.1, 0.4]])
marginal_x = joint_partial.sum(axis=1)
marginal_y = joint_partial.sum(axis=0)

h_x = entropy(marginal_x)
h_y = entropy(marginal_y)
h_xy = joint_entropy(joint_partial)
mi = mutual_information(joint_partial)
h_y_given_x = h_y - mi

print(f"H(X) = {h_x:.4f} bits")
print(f"H(Y) = {h_y:.4f} bits")
print(f"H(X,Y) = {h_xy:.4f} bits")
print(f"I(X;Y) = {mi:.4f} bits")
print(f"H(Y|X) = {h_y_given_x:.4f} bits")
print(f"\n验证: I(X;Y) = H(X) + H(Y) - H(X,Y)")
print(f"右边 = {h_x} + {h_y} - {h_xy} = {h_x + h_y - h_xy:.4f}")

# 归一化互信息（用于比较不同规模的变量）
nmi = 2 * mi / (h_x + h_y)
print(f"\n归一化互信息 NMI = {nmi:.4f}")
print("(NMI范围: 0表示独立, 1表示完全相关)")

示例1：独立变量
----------------------------------------
I(X;Y) = 0.000000 bits (应为0)

示例2：完全相关 (Y=X)
----------------------------------------
H(X) = 1.0000 bits
I(X;Y) = 1.0000 bits
I(X;Y) = H(X)? True

示例3：部分相关
----------------------------------------
H(X) = 1.0000 bits
H(Y) = 1.0000 bits
H(X,Y) = 1.7219 bits
I(X;Y) = 0.2781 bits
H(Y|X) = 0.7219 bits

验证: I(X;Y) = H(X) + H(Y) - H(X,Y)
右边 = 1.0 + 1.0 - 1.721928 = 0.2781

归一化互信息 NMI = 0.2781
(NMI范围: 0表示独立, 1表示完全相关)

import numpy as np
import matplotlib.pyplot as plt

def kl_divergence(p, q):
    """
    计算KL散度 D_KL(p || q)
    
    参数:
        p: 真实分布
        q: 近似分布
    
    返回:
        KL散度值（比特）
    """
    p = np.array(p, dtype=float)
    q = np.array(q, dtype=float)
    
    # 归一化
    p = p / p.sum()
    q = q / q.sum()
    
    # 检查支撑集
    if np.any((p > 0) & (q == 0)):
        return np.inf
    
    # 计算KL散度
    mask = p > 0
    return np.sum(p[mask] * np.log2(p[mask] / q[mask]))

# 示例1：相同分布
print("示例1：相同分布")
print("-" * 40)
p = [0.5, 0.3, 0.2]
q = [0.5, 0.3, 0.2]
kl = kl_divergence(p, q)
print(f"D_KL(p || q) = {kl:.6f} bits (应为0)")

# 示例2：不同分布
print("\n示例2：不同分布")
print("-" * 40)
p = [0.5, 0.3, 0.2]
q = [0.2, 0.3, 0.5]
kl_pq = kl_divergence(p, q)
kl_qp = kl_divergence(q, p)
print(f"p = {p}")
print(f"q = {q}")
print(f"D_KL(p || q) = {kl_pq:.4f} bits")
print(f"D_KL(q || p) = {kl_qp:.4f} bits")
print(f"不对称: D_KL(p||q) ≠ D_KL(q||p)? {not np.isclose(kl_pq, kl_qp)}")

# 示例3：KL散度与交叉熵的关系
print("\n示例3：KL散度与交叉熵")
print("-" * 40)
p = [0.1, 0.2, 0.3, 0.4]
q = [0.25, 0.25, 0.25, 0.25]  # 均匀分布

h_p = entropy(p)
cross_entropy = -np.sum(np.array(p) * np.log2(np.array(q)))
kl = kl_divergence(p, q)

print(f"H(p) = {h_p:.4f} bits")
print(f"H(p, q) [交叉熵] = {cross_entropy:.4f} bits")
print(f"D_KL(p || q) = {kl:.4f} bits")
print(f"验证: D_KL(p||q) = H(p,q) - H(p) = {cross_entropy - h_p:.4f}")

# 可视化：KL散度随分布变化
print("\n可视化：二元分布的KL散度")
print("-" * 40)

p_fixed = 0.7  # 固定真实分布
q_values = np.linspace(0.01, 0.99, 100)
kl_values = []

for q in q_values:
    p_dist = [p_fixed, 1 - p_fixed]
    q_dist = [q, 1 - q]
    kl = kl_divergence(p_dist, q_dist)
    kl_values.append(kl if kl != np.inf else 10)  # 截断无穷大

plt.figure(figsize=(10, 6))
plt.plot(q_values, kl_values, 'b-', linewidth=2)
plt.axvline(x=p_fixed, color='r', linestyle='--', label=f'p={p_fixed}')
plt.xlabel('q (近似分布的参数)', fontsize=12)
plt.ylabel('D_KL(p || q)（比特）', fontsize=12)
plt.title(f'二元分布的KL散度 (p固定为{p_fixed})', fontsize=14)
plt.grid(True, alpha=0.3)
plt.legend()
plt.ylim(0, 10)
plt.tight_layout()
plt.savefig('kl_divergence.png', dpi=150)
plt.show()

print("KL散度在q=p时达到最小值0，偏离p时迅速增大")

示例1：相同分布
----------------------------------------
D_KL(p || q) = 0.000000 bits (应为0)

示例2：不同分布
----------------------------------------
p = [0.5, 0.3, 0.2]
q = [0.2, 0.3, 0.5]
D_KL(p || q) = 0.2117 bits
D_KL(q || p) = 0.2117 bits
不对称: D_KL(p||q) ≠ D_KL(q||p)? False

示例3：KL散度与交叉熵
----------------------------------------
H(p) = 1.8464 bits
H(p, q) [交叉熵] = 2.0000 bits
D_KL(p || q) = 0.1536 bits
验证: D_KL(p||q) = H(p,q) - H(p) = 0.1536

import heapq
from collections import Counter

class HuffmanCoder:
    """Huffman编码器"""
    
    def __init__(self, frequencies):
        self.frequencies = frequencies
        self.codebook = {}
        self._build_tree()
    
    def _build_tree(self):
        # 构建优先队列
        heap = [[freq, [char, ""]] for char, freq in self.frequencies.items()]
        heapq.heapify(heap)
        
        while len(heap) > 1:
            lo = heapq.heappop(heap)
            hi = heapq.heappop(heap)
            
            for pair in lo[1:]:
                pair[1] = '0' + pair[1]
            for pair in hi[1:]:
                pair[1] = '1' + pair[1]
            
            heapq.heappush(heap, [lo[0] + hi[0]] + lo[1:] + hi[1:])
        
        # 提取码字
        self.codebook = dict(heapq.heappop(heap)[1:])
    
    def encode(self, text):
        return ''.join(self.codebook[char] for char in text)
    
    def average_code_length(self):
        total_freq = sum(self.frequencies.values())
        avg_len = sum(freq * len(code) for char, code in self.codebook.items() 
                     for ch, freq in self.frequencies.items() if ch == char) / total_freq
        # 更简洁的计算
        avg_len = sum(len(self.codebook[char]) * freq 
                     for char, freq in self.frequencies.items()) / sum(self.frequencies.values())
        return avg_len

# 英文字母频率（近似）
english_freq = {
    'e': 12.7, 't': 9.1, 'a': 8.2, 'o': 7.5, 'i': 7.0,
    'n': 6.7, 's': 6.3, 'h': 6.1, 'r': 6.0, 'd': 4.3,
    'l': 4.0, 'c': 2.8, 'u': 2.8, 'm': 2.4, 'w': 2.4,
    'f': 2.2, 'g': 2.0, 'y': 2.0, 'p': 1.9, 'b': 1.5,
    'v': 1.0, 'k': 0.8, 'j': 0.2, 'x': 0.2, 'q': 0.1, 'z': 0.1
}

# 归一化为概率
total = sum(english_freq.values())
probs = [f / total for f in english_freq.values()]

# 计算熵
h = entropy(probs)
print(f"英文字母的熵: {h:.4f} bits/symbol")

# Huffman编码
coder = HuffmanCoder(english_freq)
avg_len = coder.average_code_length()
print(f"Huffman编码的平均码长: {avg_len:.4f} bits/symbol")
print(f"冗余度: {avg_len - h:.4f} bits/symbol")
print(f"编码效率: {h / avg_len * 100:.2f}%")

# 显示部分码字
print("\n部分Huffman码字:")
for char in 'etaoin':
    print(f"  '{char}': {coder.codebook[char]}")

from sklearn.datasets import load_iris
from sklearn.feature_selection import mutual_info_classif

# 加载数据
iris = load_iris()
X, y = iris.data, iris.target

# 计算每个特征与标签的互信息
mi_scores = mutual_info_classif(X, y, discrete_features=False, random_state=42)

print("Iris数据集特征的互信息得分:")
for name, score in zip(iris.feature_names, mi_scores):
    print(f"  {name:20s}: {score:.4f}")

import numpy as np
from itertools import product

class InformationTheoryCalculator:
    """信息论度量计算器"""
    
    @staticmethod
    def entropy(p, base=2):
        """计算熵"""
        p = np.array(p, dtype=float)
        p = p[p > 0]
        p = p / p.sum()
        return -np.sum(p * np.log(p) / np.log(base))
    
    @staticmethod
    def joint_entropy(joint_dist):
        """计算联合熵"""
        joint_dist = np.array(joint_dist, dtype=float)
        flat = joint_dist.flatten()
        flat = flat[flat > 0]
        return -np.sum(flat * np.log2(flat))
    
    @staticmethod
    def conditional_entropy(joint_dist, axis=0):
        """计算条件熵 H(Y|X)，axis=0表示条件在行上"""
        joint_dist = np.array(joint_dist, dtype=float)
        
        if axis == 0:
            marginal = joint_dist.sum(axis=1)
        else:
            marginal = joint_dist.sum(axis=0)
        
        cond_ent = 0
        for i in range(joint_dist.shape[axis]):
            if axis == 0:
                cond_dist = joint_dist[i, :]
            else:
                cond_dist = joint_dist[:, i]
            
            px = marginal[i]
            if px > 0:
                cond_prob = cond_dist / px
                cond_prob = cond_prob[cond_prob > 0]
                h = -np.sum(cond_prob * np.log2(cond_prob))
                cond_ent += px * h
        
        return cond_ent
    
    @staticmethod
    def mutual_information(joint_dist):
        """计算互信息"""
        joint_dist = np.array(joint_dist, dtype=float)
        marginal_x = joint_dist.sum(axis=1)
        marginal_y = joint_dist.sum(axis=0)
        
        mi = 0
        for i in range(joint_dist.shape[0]):
            for j in range(joint_dist.shape[1]):
                pxy = joint_dist[i, j]
                px = marginal_x[i]
                py = marginal_y[j]
                
                if pxy > 0 and px > 0 and py > 0:
                    mi += pxy * np.log2(pxy / (px * py))
        
        return mi
    
    @staticmethod
    def kl_divergence(p, q):
        """计算KL散度"""
        p = np.array(p, dtype=float)
        q = np.array(q, dtype=float)
        p = p / p.sum()
        q = q / q.sum()
        
        if np.any((p > 0) & (q == 0)):
            return np.inf
        
        mask = p > 0
        return np.sum(p[mask] * np.log2(p[mask] / q[mask]))
    
    @staticmethod
    def cross_entropy(p, q):
        """计算交叉熵"""
        p = np.array(p, dtype=float)
        q = np.array(q, dtype=float)
        p = p / p.sum()
        q = q / q.sum()
        
        mask = p > 0
        return -np.sum(p[mask] * np.log2(q[mask]))


# 测试
calc = InformationTheoryCalculator()

# 创建一个联合分布
joint = np.array([
    [0.4, 0.1, 0.0],
    [0.1, 0.3, 0.1]
])

print("联合分布:")
print(joint)
print(f"\nH(X,Y) = {calc.joint_entropy(joint):.4f} bits")
print(f"H(Y|X) = {calc.conditional_entropy(joint, axis=0):.4f} bits")
print(f"I(X;Y) = {calc.mutual_information(joint):.4f} bits")

from collections import Counter
import re

class TextEntropyAnalyzer:
    """文本熵分析器"""
    
    def __init__(self, text):
        self.text = text.lower()
        self.chars = list(self.text)
        self.words = re.findall(r'\b[a-z]+\b', self.text)
    
    def char_entropy(self):
        """字符级别的熵"""
        counter = Counter(self.chars)
        total = len(self.chars)
        probs = [count / total for count in counter.values()]
        return InformationTheoryCalculator.entropy(probs)
    
    def word_entropy(self):
        """单词级别的熵"""
        counter = Counter(self.words)
        total = len(self.words)
        probs = [count / total for count in counter.values()]
        return InformationTheoryCalculator.entropy(probs)
    
    def bigram_entropy(self):
        """二元字符组的熵"""
        bigrams = [self.text[i:i+2] for i in range(len(self.text)-1)]
        counter = Counter(bigrams)
        total = len(bigrams)
        probs = [count / total for count in counter.values()]
        return InformationTheoryCalculator.entropy(probs)
    
    def redundancy(self):
        """计算冗余度"""
        h_char = self.char_entropy()
        h_bigram = self.bigram_entropy()
        # 理论最大熵（假设所有字符等概率）
        unique_chars = len(set(self.chars))
        h_max = np.log2(unique_chars) if unique_chars > 1 else 0
        
        return {
            'char_entropy': h_char,
            'bigram_entropy': h_bigram,
            'max_entropy': h_max,
            'redundancy_ratio': 1 - h_char / h_max if h_max > 0 else 0
        }


# 测试
sample_text = """
Information theory is a mathematical framework for quantifying information.
It was founded by Claude Shannon in 1948. The central concept is entropy,
which measures the uncertainty or randomness of a random variable.
"""

analyzer = TextEntropyAnalyzer(sample_text)
stats = analyzer.redundancy()

print("文本熵分析结果:")
print(f"字符熵: {stats['char_entropy']:.4f} bits")
print(f"二元组熵: {stats['bigram_entropy']:.4f} bits")
print(f"最大熵: {stats['max_entropy']:.4f} bits")
print(f"冗余度: {stats['redundancy_ratio']*100:.2f}%")

import numpy as np
import matplotlib.pyplot as plt

class InformationBottleneck:
    """
    信息瓶颈方法演示
    
    信息瓶颈原理：在压缩输入X的同时，保留关于目标Y的最大信息
    优化目标：min I(X;T) - β*I(T;Y)
    其中T是中间表示，β是权衡参数
    """
    
    def __init__(self, beta=0.5):
        self.beta = beta
    
    def simulate_channel(self, p_x, p_y_given_x, compression_rate):
        """
        模拟信息瓶颈通道
        
        参数:
            p_x: 输入分布
            p_y_given_x: 条件分布 P(Y|X)
            compression_rate: 压缩率（保留的信息比例）
        
        返回:
            I(X;T), I(T;Y)
        """
        p_x = np.array(p_x)
        p_y_given_x = np.array(p_y_given_x)
        
        n_x = len(p_x)
        n_y = p_y_given_x.shape[1]
        
        # 简化的瓶颈模型：T是X的有损版本
        # I(X;T) ≈ compression_rate * H(X)
        h_x = InformationTheoryCalculator.entropy(p_x)
        i_xt = compression_rate * h_x
        
        # I(T;Y) 的上界是 min(I(X;Y), I(X;T))
        # 计算 I(X;Y)
        p_y = p_y_given_x.T @ p_x
        joint_xy = np.zeros((n_x, n_y))
        for i in range(n_x):
            for j in range(n_y):
                joint_xy[i, j] = p_x[i] * p_y_given_x[i, j]
        
        i_xy = InformationTheoryCalculator.mutual_information(joint_xy)
        i_ty = min(i_xy, i_xt)  # 数据处理不等式
        
        return i_xt, i_ty
    
    def plot_tradeoff_curve(self, p_x, p_y_given_x):
        """绘制信息权衡曲线"""
        compression_rates = np.linspace(0, 1, 100)
        i_xt_values = []
        i_ty_values = []
        
        for rate in compression_rates:
            i_xt, i_ty = self.simulate_channel(p_x, p_y_given_x, rate)
            i_xt_values.append(i_xt)
            i_ty_values.append(i_ty)
        
        plt.figure(figsize=(10, 8))
        plt.plot(i_xt_values, i_ty_values, 'b-', linewidth=2.5)
        plt.xlabel('I(X;T) - 输入与表示的互信息', fontsize=12)
        plt.ylabel('I(T;Y) - 表示与目标的互信息', fontsize=12)
        plt.title(f'信息瓶颈权衡曲线 (β={self.beta})', fontsize=14)
        plt.grid(True, alpha=0.3)
        
        # 标记几个关键点
        for rate in [0.25, 0.5, 0.75]:
            i_xt, i_ty = self.simulate_channel(p_x, p_y_given_x, rate)
            plt.plot(i_xt, i_ty, 'ro', markersize=8)
            plt.annotate(f'rate={rate}', (i_xt, i_ty), 
                        xytext=(10, 10), textcoords='offset points', fontsize=9)
        
        plt.tight_layout()
        plt.savefig('information_bottleneck.png', dpi=150)
        plt.show()


# 测试
p_x = [0.3, 0.3, 0.2, 0.2]
p_y_given_x = np.array([
    [0.7, 0.3],
    [0.6, 0.4],
    [0.2, 0.8],
    [0.1, 0.9]
])

ib = InformationBottleneck(beta=0.5)
ib.plot_tradeoff_curve(p_x, p_y_given_x)

print("信息瓶颈分析完成")
print("蓝色曲线显示了压缩率与信息保留之间的权衡")

import cv2
import numpy as np

def block_entropy(image, block_size=8):
    h, w = image.shape
    entropies = np.zeros((h // block_size, w // block_size))
    
    for i in range(0, h, block_size):
        for j in range(0, w, block_size):
            block = image[i:i+block_size, j:j+block_size]
            # 计算直方图
            hist = cv2.calcHist([block], [0], None, [256], [0, 256])
            hist = hist.flatten() / block.size
            # 计算熵
            hist = hist[hist > 0]
            entropies[i//block_size, j//block_size] = -np.sum(hist * np.log2(hist))
    
    return entropies

# 使用示例
image = cv2.imread('test.jpg', cv2.IMREAD_GRAYSCALE)
ent_map = block_entropy(image)

plt.imshow(ent_map, cmap='hot')
plt.colorbar(label='Entropy (bits)')
plt.title('Local Entropy Map')
plt.show()