使用Dirichlet分布进行随机初始化

应用场景其实很多，比如说HMM模型中的3个参数矩阵在baum-welch期望迭代中，需要设置1个初始化的起点，假如说我有1个强先验的假设，比如说我手头上有3个隐状态，我统计训练数据频数，或者是先验想要将每一个隐状态在初始状态分布（概率）设置为[0.8, 0.1, 0.1]，即我希望第1个隐状态有80%左右的期望是状态1（或者0，在0-index下），然后状态2、3基本上不出现；

但是我们是随机初始化，所以如果限制死了参数又会和“随机”相互矛盾，我就希望：

大致围绕[0.8, 0.1, 0.1]参数可以任意次数随机生成1个概率，这个概率基本上有很大概率每次都是第1个状态给出很大的p值，或者说我随机生成100次，其中有绝大多数生成机会中（也就是频数意义中），状态1的取值（也就是p）是1个比2、3状态远远大的概率值。

那就得使用概率的概率，生成概率分布的概率，也就是狄利克雷分布。

关键在浓度参数的设置，即拟合概率的概率，python实现方式主要有2种：

1，numpy模块：

np.random.dirichlet

比如说：

import numpy as np
# 定义α参数
for i in range(10):
    alpha = np.array([0.8, 0.1, 0.1])

    # 生成一个概率行向量
    probability_vector = np.random.dirichlet(alpha)

    print(probability_vector)

import numpy as np
import matplotlib.pyplot as plt

def plot_dirichlet_probabilities(alpha, n_samples=10, figsize=(15, 8)):
    """
    绘制狄利克雷分布生成的概率向量柱状图
    
    Args:
        alpha (array-like): 狄利克雷分布的α参数
        n_samples (int): 生成的样本数量，默认为10
        figsize (tuple): 图形大小，默认为(15, 8)
    """
    # 生成概率向量
    alpha = np.array(alpha)
    probability_vectors = []
    
    for i in range(n_samples):
        prob_vector = np.random.dirichlet(alpha)
        probability_vectors.append(prob_vector)
        print(f"Sample {i+1}: {prob_vector}")
    
    # 转换为numpy数组便于处理
    prob_matrix = np.array(probability_vectors)
    
    # 设置图形
    fig, ax = plt.subplots(figsize=figsize)
    
    # 设置柱状图参数
    x = np.arange(n_samples)
    width = 0.25
    state_labels = ['State 0 (+)', 'State 1 (-)', 'State 2 (Neutral)']
    colors = ['red', 'blue', 'gray']
    
    # 绘制每个状态的柱状图
    for i in range(3):
        bars = ax.bar(x + i * width, prob_matrix[:, i], width, 
                     label=state_labels[i], color=colors[i], alpha=0.7)
        
        # 在柱子上添加数值标签
        for j, bar in enumerate(bars):
            height = bar.get_height()
            ax.text(bar.get_x() + bar.get_width()/2., height + 0.01,
                   f'{height:.3f}', ha='center', va='bottom', fontsize=8)
    
    # 设置图形属性
    ax.set_xlabel('Sample Index')
    ax.set_ylabel('Probability')
    ax.set_title(f'Dirichlet Distribution Probability Vectors (α={alpha})')
    ax.set_xticks(x + width)
    ax.set_xticklabels([f'Sample {i+1}' for i in range(n_samples)])
    ax.legend()
    ax.set_ylim(0, 1.1)
    
    # 添加网格
    ax.grid(True, alpha=0.3)
    
    plt.xticks(rotation=45 if n_samples > 20 else 0)
    plt.tight_layout()
    plt.show()
    
    return prob_matrix

# 使用示例
if __name__ == "__main__":
    # 定义α参数
    alpha = [0.8, 0.1, 0.1]
    
    # 生成并绘制10个样本
    print("=== 10个样本 ===")
    prob_matrix_10 = plot_dirichlet_probabilities(alpha, n_samples=10)
    
    # 生成并绘制100个样本
    print("\n=== 100个样本 ===")
    prob_matrix_100 = plot_dirichlet_probabilities(alpha, n_samples=100, figsize=(20, 10))

可以很直观看出来，基本上（频数意义上），生成的序列中都是第1个状态的值最大

如果是100个样本的话，

规律其实就比较明显了，当然堆叠比并行要好可视化点：

import numpy as np
import matplotlib.pyplot as plt

def plot_dirichlet_stacked_bars(alpha, n_samples=10, figsize=(15, 8)):
    """
    绘制狄利克雷分布生成的概率向量堆叠柱状图
    
    Args:
        alpha (array-like): 狄利克雷分布的α参数
        n_samples (int): 生成的样本数量，默认为10
        figsize (tuple): 图形大小，默认为(15, 8)
    """
    # 生成概率向量
    alpha = np.array(alpha)
    probability_vectors = []
    
    for i in range(n_samples):
        prob_vector = np.random.dirichlet(alpha)
        probability_vectors.append(prob_vector)
        print(f"Sample {i+1}: {prob_vector}")
    
    # 转换为numpy数组便于处理
    prob_matrix = np.array(probability_vectors)
    
    # 设置图形
    fig, ax = plt.subplots(figsize=figsize)
    
    # 设置参数
    x = np.arange(n_samples)
    state_labels = ['State 0 (+)', 'State 1 (-)', 'State 2 (Neutral)']
    colors = ['red', 'blue', 'gray']
    
    # 绘制堆叠柱状图
    bottom = np.zeros(n_samples)
    
    for i in range(3):
        bars = ax.bar(x, prob_matrix[:, i], bottom=bottom, 
                     label=state_labels[i], color=colors[i], alpha=0.8)
        
        # 在每个堆叠部分添加数值标签
        for j, bar in enumerate(bars):
            height = bar.get_height()
            if height > 0.05:  # 只有当高度足够大时才显示标签
                ax.text(bar.get_x() + bar.get_width()/2., 
                       bottom[j] + height/2,
                       f'{height:.3f}', 
                       ha='center', va='center', fontsize=8, 
                       color='white', fontweight='bold')
        
        # 更新底部位置
        bottom += prob_matrix[:, i]
    
    # 设置图形属性
    ax.set_xlabel('Sample Index')
    ax.set_ylabel('Probability')
    ax.set_title(f'Dirichlet Distribution Stacked Probability Vectors (α={alpha})')
    ax.set_xticks(x)
    ax.set_xticklabels([f'Sample {i+1}' for i in range(n_samples)])
    ax.legend()
    ax.set_ylim(0, 1.05)
    
    # 添加网格
    ax.grid(True, alpha=0.3, axis='y')
    
    # 添加水平线表示总和为1
    ax.axhline(y=1, color='black', linestyle='--', alpha=0.5, linewidth=1)
    
    plt.xticks(rotation=45 if n_samples > 15 else 0)
    plt.tight_layout()
    plt.show()
    
    return prob_matrix

# 使用示例，替换原来的调用
if __name__ == "__main__":
    # 定义α参数
    alpha = [0.8, 0.1, 0.1]
    
    # 生成并绘制10个样本
    print("=== 10个样本 ===")
    prob_matrix_10 = plot_dirichlet_stacked_bars(alpha, n_samples=10)
    
    # 生成并绘制100个样本
    print("\n=== 100个样本 ===")
    prob_matrix_100 = plot_dirichlet_stacked_bars(alpha, n_samples=100, figsize=(20, 10))

这个就很清晰了，虽然有个别蓝色的生成不符合我们的单次期望，但是从频数意义上讲，我们大多数随机化的设置中，都生成了期望的向量序列：

2，scipy模块

from scipy.stats import dirichlet

dirichlet.rvs