说话人与音频分析

说话人与音频分析识别谁在说话、何时说话，以及音频中存在哪些非语音声音。本文涵盖说话人验证与识别、i-vector、d-vector、x-vector、说话人日志、音频事件分类、音乐信息检索以及语音情感识别。

• 在文件01中，我们建立了信号处理基础：谱图、MFCCs和梅尔滤波器组。在文件02中，我们识别了说了什么内容。现在我们要问：谁说的、何时说的，以及音频中还在发生什么。说话人识别、日志、音频分类和音乐分析有一条共同主线：学习紧凑的嵌入表示，捕捉任务所需的不变性，呼应文件06中的嵌入思想。

• 识别说话人就像在电话中辨认朋友的声音。你不需要理解词语；音色、节奏和嗓音特质对该人是独特的。说话人识别系统学习从原始音频中提取这种"声纹"，忽略说了什么，专注于如何说。

• 说话人识别是两类相关任务的统称：

  • 说话人验证（Speaker Verification, SV）：给定声称的身份和一段音频，判断说话人是否如其声称。这是一个二元决策（接受或拒绝），是语音认证技术的基础（"嘿Siri，这是我的声音吗？"）。
  • 说话人识别（Speaker Identification, SI）：给定一段音频和已知说话人库，确定哪个人产生了该音频。这是一个多分类问题。

• 两项任务共享相同的底层表示：固定维度的说话人嵌入，捕捉说话人身份而不依赖其所说内容。区别仅在于决策阶段：验证比较两个嵌入，识别在候选中寻找最近嵌入。

• 余弦相似度是比较说话人嵌入的标准度量。给定注册嵌入 \(e\) 和测试嵌入 \(t\)：

\[s = \frac{e \cdot t}{\|e\| \, \|t\|}\]

• 阈值 \(\theta\) 决定接受/拒绝：若 \(s > \theta\) 则接受。该阈值在错误接受率（FAR）和错误拒绝率（FRR）之间权衡。等错误率（EER，FAR = FRR时的点）是标准评估指标。EER越低性能越好。最先进系统在标准基准（VoxCeleb）上实现低于1%的EER。

• i-vector（Dehak et al., 2010）是深度学习前主导的说话人嵌入方法。其思想来自因子分析（文件02的矩阵分解和文件04的降维）。通用背景模型（UBM，在多样化说话人上训练的大型GMM）定义超向量空间。每个话语的GMM超向量被投影到低维总变异性空间：

\[M = m + Tw\]

• 其中 \(M\) 是话语的GMM超向量，\(m\) 是UBM均值超向量，\(T\) 是总变异性矩阵（从数据学习），\(w\) 是i-vector，一个低维（通常400-600）表示，同时捕捉说话人和信道变异性。

• 为从i-vector中移除信道变异性，概率线性判别分析（PLDA）将i-vector建模为说话人特定和信道特定潜变量之和。PLDA为验证提供原则性的对数似然比分数：

\[\text{score}(w_1, w_2) = \log \frac{P(w_1, w_2 \mid \text{same speaker})}{P(w_1 \mid \text{speaker}_1) \, P(w_2 \mid \text{speaker}_2)}\]

• d-vector（Variani et al., 2014）是首个神经说话人嵌入。在帧级特征上训练用于说话人分类的DNN，通过对话语中所有帧的最后一层隐藏激活取平均，提取固定维表示。简单但有效，d-vector证明神经网络无需i-vector的复杂统计机制即可学习说话人判别特征。

• x-vector（Snyder et al., 2018）使用时延神经网络（TDNN）架构显著推进了神经说话人嵌入。TDNN是具有特定上下文窗口的1D卷积，与文件03的WaveNet中的空洞卷积相关，但应用于帧级特征而非原始波形样本。

• x-vector架构包含三个阶段：
  • 帧级层：堆叠的TDNN层以逐渐增宽的时间上下文处理MFCC（文件01）。每层看到固定上下文窗口（如第一层 \(\{t-2, t-1, t, t+1, t+2\}\)，后续层更宽）。
  • 统计池化：帧级层之后，计算帧级输出在整个话语上的均值和标准差，生成固定维向量，无论话语长度：

\[ \begin{aligned} \mu &= \frac{1}{T} \sum_{t=1}^{T} h_t \\ \sigma &= \sqrt{\frac{1}{T} \sum_{t=1}^{T} (h_t - \mu)^2} \end{aligned} \]

• 其中 \(h_t\) 是时间 \(t\) 的帧级输出。拼接 \([\mu; \sigma]\) 即为池化表示。

  • 段级层：全连接层处理池化表示。第一个段级层的输出（softmax之前）即为x-vector嵌入。

• x-vector使用标准的说话人身份交叉熵损失训练。尽管为分类而训练，学习到的中间表示（x-vector）能良好泛化到未见说话人，因为网络学习提取说话人判别特征而非记忆特定说话人。

• ECAPA-TDNN（Desplanques et al., 2020）是当前最先进的基于TDNN的说话人识别架构。相比x-vector引入三项改进：

  • 挤压-激励（SE）：通道注意力（文件08的SENet），根据全局上下文重加权特征通道，使模型能强调与说话人相关的通道。
  • Res2Net风格多尺度特征：在每个TDNN块内，通道被分组并分层处理，在多个时间分辨率上创建特征（类比文件08的多尺度特征提取）。
  • 注意力统计池化：不用等权平均，注意力机制加权每帧对池化统计的贡献。具有更多说话人判别内容的帧（如元音，携带更多说话人信息）获得更高注意力权重：

\[\alpha_t = \frac{\exp(v^T f(h_t))}{\sum_{\tau} \exp(v^T f(h_\tau))}\]

• 其中 \(f\) 是小型神经网络，\(v\) 是学习的注意力向量。注意力加权均值和标准差变为 \(\tilde{\mu} = \sum_t \alpha_t h_t\) 和 \(\tilde{\sigma} = \sqrt{\sum_t \alpha_t (h_t - \tilde{\mu})^2}\)。

• ECAPA-TDNN通常用AAM-Softmax（加性角边际Softmax）训练，该损失在分类损失中添加角边际惩罚，将同一说话人的嵌入在超球面上推得更近、不同说话人的推得更远：

\[L = -\log \frac{e^{s \cos(\theta_{y_i} + m)}}{e^{s \cos(\theta_{y_i} + m)} + \sum_{j \neq y_i} e^{s \cos \theta_j}}\]

• 其中 \(\theta_{y_i}\) 是嵌入与真实类别权重向量的夹角，\(m\) 是边际（通常0.2），\(s\) 是缩放因子（通常30）。该损失来自人脸识别（文件08的ArcFace），对说话人验证非常有效。

• 说话人日志（Speaker Diarisation）回答多说话人录音中"谁在何时说话"。可将其想象为给时间线上色：每种颜色代表不同说话人，系统必须确定每个说话人何时活跃，包括重叠语音。

• 基于聚类的日志是传统流程方法：

  • 分割：使用滑动窗口或说话人变化检测将音频划分为短片段（通常1-2秒）。
  • 嵌入提取：为每个片段提取说话人嵌入（x-vector、ECAPA-TDNN）。
  • 聚类：按说话人分组片段。凝聚层次聚类（AHC）是标准方法：从每个片段为独立簇开始，迭代合并最相似的两个簇，直到满足停止准则（基于距离阈值或目标说话人数）。
  • 重分割：使用基于Viterbi的重对齐细化边界。

• 说话人数量通常事先未知，这使问题比标准聚类更难。使用基于特征值阈值确定 \(k\) 的谱聚类是另一种常见方法。

• 端到端神经日志（EEND）（Fujita et al., 2019）将日志框架化为多标签分类问题。神经网络（通常为基于自注意力的模型，文件07的Transformer）以整段录音为输入，输出每帧每个说话人的二元活跃标签。这直接处理重叠语音，这是基于聚类方法的主要弱点。

• EEND在帧 \(t\) 对 \(S\) 个说话人的输出为：

\[\hat{y}_{t,s} = \sigma(f_s(h_t))\]

• 其中 \(h_t\) 是Transformer在帧 \(t\) 的输出，\(f_s\) 是说话人 \(s\) 的线性投影。训练损失是对说话人和帧求和的二元交叉熵。关键挑战是说话人数量必须固定，或使用可变输出架构处理（EEND-EDA使用带吸引子的编码器-解码器）。

• 日志的排列不变训练（PIT）处理标签歧义问题：由于说话人无固有顺序，损失对所有可能的说话人到输出分配计算并取最小值（这与文件05中源分离使用的相同PIT）。

• 音频分类为整个音频片段分配标签。与转写语音的ASR（文件02）不同，音频分类覆盖更广范围：环境声音（警笛、雨声、狗叫）、音乐流派（摇滚、爵士、古典）和通用音频事件。

• 标准方法遵循文件08的图像分类范式：将音频表示为谱图（2D时频图像），然后应用CNN或Transformer分类器。这种谱图-图像方法利用了计算机视觉数十年的进展。

• 环境声音分类（ESC）使用ESC-50（50类，2000片段）和UrbanSound8K等数据集。典型架构是应用于对数梅尔谱图的CNN（文件06）。数据增强至关重要：时间拉伸、音高偏移、添加背景噪声，以及SpecAugment（文件02的掩码方法应用于谱图）均能提升泛化能力。

• 音频事件检测（声音事件检测，SED）是分类的时间类比：不仅识别存在哪些事件，还要识别其起止时间。AudioSet（Gemmeke et al., 2017）是大规模基准，含527个事件类别和超过200万个来自YouTube的10秒片段，每个片段为弱标注（片段级标签，非帧级）。

• 弱监督SED必须从片段级标签学习帧级预测。标准方法使用生成帧级类别概率的CNN，然后通过注意力池化聚合为片段级预测：

\[\hat{Y}_c = \sigma\left(\sum_t \alpha_{t,c} \cdot f_{t,c}\right)\]

• 其中 \(f_{t,c}\) 是时间 \(t\) 类别 \(c\) 的帧级logit，\(\alpha_{t,c}\) 是注意力权重。片段级预测 \(\hat{Y}_c\) 针对片段级标签训练。

• 声学场景分类（ASC）对整体环境分类："机场"、"公园"、"地铁站"、"办公室"。这是整体性任务：模型必须捕捉一般声学纹理而非特定事件。DCASE挑战系列每年对ASC进行基准测试，获胜系统通常使用多分辨率谱图上CNN集成的方法。

• 音频嵌入是从大规模音频数据学习的通用表示，类比于词嵌入（文件07）或图像特征（文件08），可迁移到下游任务。

• VGGish（Hershey et al., 2017）将图像分类网络VGG（文件08）适配到音频。它处理0.96秒的对数梅尔谱图块，通过预训练于AudioSet的类VGG CNN，为每块生成128维嵌入。VGGish嵌入作为下游任务的通用音频特征，类似于ImageNet预训练CNN提供视觉特征。

• PANNs（Pre-trained Audio Neural Networks, Kong et al., 2020）是一族在完整AudioSet上训练用于音频标注的CNN架构（CNN6、CNN10、CNN14）。最常用的CNN14是14层CNN，使用 \(3 \times 3\) 卷积应用于对数梅尔谱图。PANNs生成2048维嵌入，在多样化音频任务上实现最先进的迁移学习。

• 音频谱图Transformer（AST）（Gong et al., 2021）将Vision Transformer（ViT，文件08）架构直接应用于音频谱图。谱图被划分为 \(16 \times 16\) 块（如同ViT划分图像），每块线性投影为标记嵌入，添加位置嵌入，标准Transformer编码器（文件07）处理该序列。[CLS]标记的输出用于分类。

• AST受益于ImageNet预训练：由于谱图是2D图像，AST从预训练于ImageNet图像的ViT初始化，然后在音频上微调。这种跨模态迁移出人意料地有效，因为两个领域共享低级特征（边缘、纹理），且位置嵌入可插值以适应不同谱图尺寸。

• HTS-AT（Chen et al., 2022）通过分层Swin Transformer架构（文件08的移位窗口注意力）改进AST，在通过多尺度特征提取提升性能的同时降低计算成本。

• BEATs（Chen et al., 2023）使用音频特定的预训练策略：迭代掩码预测配合离散标记器（类似文件02中wav2vec 2.0的方法，但应用于通用音频）。标记器逐步细化，创建语义意义日益增强的离散音频标记。

• 基于嵌入的说话人日志结合说话人嵌入与时间建模。如Pyannote.audio等现代系统使用三阶段流程：(1) 检测说话人轮换和重叠语音的神经分割模型，(2) 对每个检测片段应用嵌入提取（ECAPA-TDNN），(3) 聚类以在整个录音中分配说话人身份。

• 音乐信息检索（MIR）将音频分析应用于音乐。文件01中的谱表示在此特别有用，因为音乐具有丰富的和声结构。

• 节拍跟踪检测音乐的节奏脉冲。标准方法从谱图计算起始强度包络（检测指示音符起始的能量增长），然后使用自相关或时谱图寻找节拍，最后使用动态规划跟踪单个节拍位置，寻找与起始包络最匹配且保持恒定节拍的节拍时间序列。

• 和弦识别识别随时间变化的和声内容。输入通常是色度图（也称音高类剖面）：12维表示，将所有八度折叠在一起，显示12个音高类（C、C#、D、...、B）中每个的能量。CNN或RNN（文件06）将每帧分类为标准和弦标签之一（C大调、A小调、G7等）。

• 色度图通过从STFT（文件01）将每个频率仓映射到其音高类计算：

\[\text{chroma}(p) = \sum_{k : \text{pitch}(k) \bmod 12 = p} |X(k)|^2\]

• 其中 \(p \in \{0, 1, \ldots, 11\}\) 是音高类，\(\text{pitch}(k)\) 将频率仓 \(k\) 映射到其MIDI音符编号。

• 源分离基础（文件05详述）将音乐录音分离为独立乐器（人声、鼓、贝斯、其他）。这对混音、卡拉OK和音乐转录等MIR应用至关重要。如Demucs（文件05）等模型在标准MUSDB18基准上实现出色的分离质量。

• 音乐标注为歌曲分配标签（流派、情绪、乐器、年代）。本质上是应用于音乐的音频分类，使用相同的谱图上CNN方法。Million Song Dataset和MagnaTagATune是标准基准。

• 音频指纹从短片段识别特定录音，即使存在噪声、混响或压缩伪影。经典系统是Shazam，它哈希谱图中的星座点（显著峰值）。神经方法学习鲁棒嵌入，对声学退化不变，同时在不同录音间保持判别力，呼应文件06和文件08中的不变特征学习。

编程任务（使用CoLab或Notebook）

• 任务1：使用统计池化的说话人嵌入提取。构建简单的x-vector风格模型，通过TDNN层和统计池化处理帧级特征以生成说话人嵌入。

import jax
import jax.numpy as jnp
import jax.random as jr
import matplotlib.pyplot as plt

# 为多个说话人生成模拟帧级MFCC特征
def generate_speaker_data(key, n_speakers=5, utterances_per_speaker=20,
                          n_frames=100, n_features=40):
    """生成具有说话人依赖模式的合成说话人数据."""
    keys = jr.split(key, 3)
    all_features = []
    all_labels = []

    # 每个说话人有特征性频谱模式
    speaker_patterns = jr.normal(keys[0], (n_speakers, n_features)) * 0.5

    for spk in range(n_speakers):
        for utt in range(utterances_per_speaker):
            k = jr.fold_in(keys[1], spk * utterances_per_speaker + utt)
            noise = jr.normal(k, (n_frames, n_features)) * 0.3
            features = speaker_patterns[spk][None, :] + noise
            all_features.append(features)
            all_labels.append(spk)

    perm = jr.permutation(keys[2], len(all_features))
    features = jnp.stack(all_features)[perm]
    labels = jnp.array(all_labels)[perm]
    return features, labels

key = jr.PRNGKey(42)
features, labels = generate_speaker_data(key)
n_speakers = 5
n_features = 40

# x-vector风格模型
def init_xvector(key, n_features=40, hidden=128, embed_dim=64, n_speakers=5):
    keys = jr.split(key, 8)
    params = {
        # TDNN层1: 上下文 [-2, 2]
        'tdnn1_w': jr.normal(keys[0], (5, n_features, hidden)) * jnp.sqrt(2.0 / (5 * n_features)),
        'tdnn1_b': jnp.zeros(hidden),
        # TDNN层2: 上下文 [-2, 2]
        'tdnn2_w': jr.normal(keys[1], (5, hidden, hidden)) * jnp.sqrt(2.0 / (5 * hidden)),
        'tdnn2_b': jnp.zeros(hidden),
        # TDNN层3: 上下文 [-3, 3]
        'tdnn3_w': jr.normal(keys[2], (7, hidden, hidden)) * jnp.sqrt(2.0 / (7 * hidden)),
        'tdnn3_b': jnp.zeros(hidden),
        # 段级层（池化后: 2*hidden -> embed_dim）
        'seg1_w': jr.normal(keys[3], (2 * hidden, embed_dim)) * jnp.sqrt(2.0 / (2 * hidden)),
        'seg1_b': jnp.zeros(embed_dim),
        # 分类头
        'cls_w': jr.normal(keys[4], (embed_dim, n_speakers)) * jnp.sqrt(2.0 / embed_dim),
        'cls_b': jnp.zeros(n_speakers),
    }
    return params

def xvector_forward(params, x, return_embedding=False):
    """x: (batch, frames, features) -> logits 或嵌入."""
    # TDNN层（1D卷积）
    h = jax.lax.conv_general_dilated(
        x.transpose(0, 2, 1), params['tdnn1_w'].transpose(2, 1, 0),
        window_strides=(1,), padding='SAME'
    ).transpose(0, 2, 1) + params['tdnn1_b']
    h = jax.nn.relu(h)

    h = jax.lax.conv_general_dilated(
        h.transpose(0, 2, 1), params['tdnn2_w'].transpose(2, 1, 0),
        window_strides=(1,), padding='SAME'
    ).transpose(0, 2, 1) + params['tdnn2_b']
    h = jax.nn.relu(h)

    h = jax.lax.conv_general_dilated(
        h.transpose(0, 2, 1), params['tdnn3_w'].transpose(2, 1, 0),
        window_strides=(1,), padding='SAME'
    ).transpose(0, 2, 1) + params['tdnn3_b']
    h = jax.nn.relu(h)

    # 统计池化: 时间上的均值和标准差
    mu = jnp.mean(h, axis=1)
    sigma = jnp.std(h, axis=1)
    pooled = jnp.concatenate([mu, sigma], axis=-1)

    # 段级层 -> 嵌入
    embedding = jax.nn.relu(pooled @ params['seg1_w'] + params['seg1_b'])

    if return_embedding:
        return embedding

    # 分类
    logits = embedding @ params['cls_w'] + params['cls_b']
    return logits

def cross_entropy_loss(params, features, labels):
    logits = xvector_forward(params, features)
    one_hot = jax.nn.one_hot(labels, n_speakers)
    log_probs = jax.nn.log_softmax(logits)
    return -jnp.mean(jnp.sum(one_hot * log_probs, axis=-1))

grad_fn = jax.jit(jax.value_and_grad(cross_entropy_loss))

# 训练
params = init_xvector(jr.PRNGKey(0))
lr = 1e-3
losses = []

for epoch in range(300):
    loss_val, grads = grad_fn(params, features, labels)
    params = jax.tree.map(lambda p, g: p - lr * g, params, grads)
    losses.append(float(loss_val))

# 提取嵌入并用t-SNE风格2D投影可视化（使用PCA）
embeddings = xvector_forward(params, features, return_embedding=True)

# 简单PCA到2D
emb_centered = embeddings - jnp.mean(embeddings, axis=0)
_, _, Vt = jnp.linalg.svd(emb_centered, full_matrices=False)
proj_2d = emb_centered @ Vt[:2].T

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

axes[0].plot(losses, color='#3498db', linewidth=1.5)
axes[0].set_xlabel(' epoch')
axes[0].set_ylabel('交叉熵损失')
axes[0].set_title('说话人分类训练')
axes[0].set_yscale('log')

colors = ['#3498db', '#e74c3c', '#27ae60', '#f39c12', '#9b59b6']
for spk in range(n_speakers):
    mask = labels == spk
    axes[1].scatter(proj_2d[mask, 0], proj_2d[mask, 1], c=colors[spk],
                    label=f'说话人 {spk}', alpha=0.7, s=30)
axes[1].set_xlabel('主成分 1')
axes[1].set_ylabel('主成分 2')
axes[1].set_title('说话人嵌入（PCA投影）')
axes[1].legend()

plt.tight_layout()
plt.show()

# 验证演示: 余弦相似度
emb_norm = embeddings / jnp.linalg.norm(embeddings, axis=-1, keepdims=True)
sim_matrix = emb_norm @ emb_norm.T
print(f"嵌入形状: {embeddings.shape}")
print(f"平均同说话人相似度: {jnp.mean(sim_matrix[labels[:, None] == labels[None, :]]):.4f}")
print(f"平均异说话人相似度: {jnp.mean(sim_matrix[labels[:, None] != labels[None, :]]):.4f}")

• 任务2：使用余弦相似度评分的说话人验证。给定预计算的说话人嵌入，实现验证系统以计算EER（等错误率）并绘制DET曲线。

import jax
import jax.numpy as jnp
import jax.random as jr
import matplotlib.pyplot as plt

def generate_verification_pairs(key, n_speakers=20, dim=64, n_pairs=2000):
    """生成说话人嵌入和验证试验对."""
    keys = jr.split(key, 5)

    # 说话人质心，带一定方差
    centroids = jr.normal(keys[0], (n_speakers, dim))
    centroids = centroids / jnp.linalg.norm(centroids, axis=-1, keepdims=True)

    # 生成注册和测试嵌入，带说话人内方差
    enroll_embs = []
    test_embs = []
    trial_labels = []  # 1 = 同说话人（目标），0 = 不同（冒名）

    for i in range(n_pairs):
        k1, k2, k3 = jr.split(jr.fold_in(keys[1], i), 3)
        is_target = jr.bernoulli(k1).astype(int)

        spk1 = jr.randint(k2, (), 0, n_speakers)
        emb1 = centroids[spk1] + jr.normal(jr.fold_in(k3, 0), (dim,)) * 0.15

        if is_target:
            spk2 = spk1
        else:
            spk2 = (spk1 + jr.randint(jr.fold_in(k3, 1), (), 1, n_speakers)) % n_speakers

        emb2 = centroids[spk2] + jr.normal(jr.fold_in(k3, 2), (dim,)) * 0.15

        enroll_embs.append(emb1)
        test_embs.append(emb2)
        trial_labels.append(int(is_target))

    return (jnp.stack(enroll_embs), jnp.stack(test_embs),
            jnp.array(trial_labels))

key = jr.PRNGKey(42)
enroll, test, labels = generate_verification_pairs(key)

# 计算余弦相似度分数
enroll_norm = enroll / jnp.linalg.norm(enroll, axis=-1, keepdims=True)
test_norm = test / jnp.linalg.norm(test, axis=-1, keepdims=True)
scores = jnp.sum(enroll_norm * test_norm, axis=-1)

# 计算不同阈值下的FAR和FRR
thresholds = jnp.linspace(-1.0, 1.0, 500)

target_scores = scores[labels == 1]
impostor_scores = scores[labels == 0]

fars = []
frrs = []
for thresh in thresholds:
    far = jnp.mean(impostor_scores >= thresh)  # 错误接受
    frr = jnp.mean(target_scores < thresh)     # 错误拒绝
    fars.append(float(far))
    frrs.append(float(frr))

fars = jnp.array(fars)
frrs = jnp.array(frrs)

# 寻找EER: FAR ≈ FRR处
eer_idx = jnp.argmin(jnp.abs(fars - frrs))
eer = float((fars[eer_idx] + frrs[eer_idx]) / 2)
eer_threshold = float(thresholds[eer_idx])

print(f"等错误率（EER）: {eer:.4f} ({eer*100:.2f}%)")
print(f"EER阈值: {eer_threshold:.4f}")

fig, axes = plt.subplots(1, 3, figsize=(18, 5))

# 分数分布
bins = jnp.linspace(-0.5, 1.0, 60)
axes[0].hist(target_scores, bins=bins, alpha=0.6, color='#27ae60',
             label='目标（同说话人）', density=True)
axes[0].hist(impostor_scores, bins=bins, alpha=0.6, color='#e74c3c',
             label='冒名（不同说话人）', density=True)
axes[0].axvline(eer_threshold, color='#f39c12', linestyle='--', linewidth=2,
                label=f'EER阈值 = {eer_threshold:.3f}')
axes[0].set_xlabel('余弦相似度分数')
axes[0].set_ylabel('密度')
axes[0].set_title('分数分布')
axes[0].legend()

# FAR vs FRR
axes[1].plot(thresholds, fars, color='#e74c3c', linewidth=2, label='FAR')
axes[1].plot(thresholds, frrs, color='#3498db', linewidth=2, label='FRR')
axes[1].axvline(eer_threshold, color='#f39c12', linestyle='--', linewidth=1.5)
axes[1].scatter([eer_threshold], [eer], color='#f39c12', s=100, zorder=5,
                label=f'EER = {eer:.4f}')
axes[1].set_xlabel('阈值')
axes[1].set_ylabel('错误率')
axes[1].set_title('FAR和FRR随阈值变化')
axes[1].legend()

# DET曲线（FAR vs FRR）
axes[2].plot(fars, frrs, color='#9b59b6', linewidth=2)
axes[2].plot([0, 1], [0, 1], 'k--', alpha=0.3)
axes[2].scatter([eer], [eer], color='#f39c12', s=100, zorder=5,
                label=f'EER = {eer:.4f}')
axes[2].set_xlabel('错误接受率')
axes[2].set_ylabel('错误拒绝率')
axes[2].set_title('DET曲线')
axes[2].set_xlim([0, 0.5])
axes[2].set_ylim([0, 0.5])
axes[2].legend()
axes[2].set_aspect('equal')

plt.tight_layout()
plt.show()

• 任务3：音频谱图块嵌入（AST风格）。实现音频谱图Transformer的块提取和嵌入层，可视化谱图如何被标记化。

import jax
import jax.numpy as jnp
import jax.random as jr
import matplotlib.pyplot as plt

# 生成合成谱图（谐波结构 + 噪声）
def generate_spectrogram(key, n_time=128, n_freq=128):
    """创建具有谐波模式的合成谱图."""
    k1, k2 = jr.split(key)
    spec = jr.normal(k1, (n_time, n_freq)) * 0.1

    # 添加谐波带（模拟语音共振峰）
    for f0 in [15, 30, 45, 70]:
        width = 3
        envelope = jnp.exp(-0.5 * ((jnp.arange(n_freq) - f0) / width) ** 2)
        time_mod = 0.5 + 0.5 * jnp.sin(2 * jnp.pi * jnp.arange(n_time) / 40)
        spec += jnp.outer(time_mod, envelope)

    return jnp.clip(spec, 0, None)

key = jr.PRNGKey(42)
spectrogram = generate_spectrogram(key)
n_time, n_freq = spectrogram.shape

# 块提取参数
patch_h = 16  # 时间
patch_w = 16  # 频率
stride_h = 16
stride_w = 16
embed_dim = 192  # ViT-Small维度

n_patches_h = n_time // stride_h
n_patches_w = n_freq // stride_w
n_patches = n_patches_h * n_patches_w

print(f"谱图: {n_time} x {n_freq}")
print(f"块大小: {patch_h} x {patch_w}")
print(f"块数量: {n_patches_h} x {n_patches_w} = {n_patches}")

# 提取块
def extract_patches(spec, patch_h, patch_w, stride_h, stride_w):
    """从谱图提取非重叠块."""
    patches = []
    positions = []
    for i in range(0, spec.shape[0] - patch_h + 1, stride_h):
        for j in range(0, spec.shape[1] - patch_w + 1, stride_w):
            patch = spec[i:i+patch_h, j:j+patch_w]
            patches.append(patch.flatten())
            positions.append((i, j))
    return jnp.stack(patches), positions

patches, positions = extract_patches(spectrogram, patch_h, patch_w, stride_h, stride_w)
print(f"块形状: {patches.shape}")  # (n_patches, patch_h * patch_w)

# 线性投影（块嵌入）
patch_dim = patch_h * patch_w
k1, k2 = jr.split(jr.PRNGKey(0))
W_embed = jr.normal(k1, (patch_dim, embed_dim)) * jnp.sqrt(2.0 / patch_dim)
b_embed = jnp.zeros(embed_dim)

# 可学习位置嵌入
pos_embed = jr.normal(k2, (n_patches + 1, embed_dim)) * 0.02  # +1用于CLS

# CLS标记
cls_token = jnp.zeros((1, embed_dim))

# 前向传播
patch_tokens = patches @ W_embed + b_embed  # (n_patches, embed_dim)
tokens = jnp.concatenate([cls_token, patch_tokens], axis=0)  # (n_patches+1, embed_dim)
tokens = tokens + pos_embed  # 添加位置嵌入

print(f"标记序列形状: {tokens.shape}")
print(f"每个标记维度: {embed_dim}")

# 可视化
fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# 带块网格的原始谱图
axes[0, 0].imshow(spectrogram.T, aspect='auto', origin='lower', cmap='magma')
for i in range(0, n_time + 1, stride_h):
    axes[0, 0].axvline(i - 0.5, color='white', linewidth=0.5, alpha=0.5)
for j in range(0, n_freq + 1, stride_w):
    axes[0, 0].axhline(j - 0.5, color='white', linewidth=0.5, alpha=0.5)
axes[0, 0].set_title(f'谱图，{patch_h}x{patch_w}块网格')
axes[0, 0].set_xlabel('时间帧')
axes[0, 0].set_ylabel('频率仓')

# 可视化单个块
n_show = min(16, n_patches)
patch_grid = patches[:n_show].reshape(n_show, patch_h, patch_w)
combined = jnp.concatenate([patch_grid[i] for i in range(min(8, n_show))], axis=1)
axes[0, 1].imshow(combined.T, aspect='auto', origin='lower', cmap='magma')
axes[0, 1].set_title(f'前{min(8, n_show)}个块（拼接）')
axes[0, 1].set_xlabel('块索引（水平）')
axes[0, 1].set_ylabel('块内频率')

# 标记嵌入相似度矩阵
token_norms = tokens / jnp.linalg.norm(tokens, axis=-1, keepdims=True)
sim = token_norms @ token_norms.T
im = axes[1, 0].imshow(sim, cmap='RdBu_r', vmin=-1, vmax=1)
axes[1, 0].set_title('标记相似度矩阵（余弦）')
axes[1, 0].set_xlabel('标记索引')
axes[1, 0].set_ylabel('标记索引')
plt.colorbar(im, ax=axes[1, 0], fraction=0.046)

# 位置嵌入相似度
pos_norms = pos_embed / jnp.linalg.norm(pos_embed, axis=-1, keepdims=True)
pos_sim = pos_norms @ pos_norms.T
im2 = axes[1, 1].imshow(pos_sim, cmap='RdBu_r', vmin=-1, vmax=1)
axes[1, 1].set_title('位置嵌入相似度')
axes[1, 1].set_xlabel('位置索引')
axes[1, 1].set_ylabel('位置索引')
plt.colorbar(im2, ax=axes[1, 1], fraction=0.046)

plt.tight_layout()
plt.show()

• 任务4：用于和弦分析的简单色度图计算。从合成谐波信号计算并可视化色度图，演示音乐信息检索中使用的音高类折叠。

import jax
import jax.numpy as jnp
import matplotlib.pyplot as plt

# 生成合成音乐信号: C大调和弦 -> G大调和弦
sr = 16000
duration = 2.0
t = jnp.linspace(0, duration, int(sr * duration))

# C大调（C4=261.6, E4=329.6, G4=392.0）前半段
# G大调（G3=196.0, B3=246.9, D4=293.7）后半段
half = len(t) // 2

c_major = (0.5 * jnp.sin(2 * jnp.pi * 261.63 * t[:half]) +
           0.4 * jnp.sin(2 * jnp.pi * 329.63 * t[:half]) +
           0.3 * jnp.sin(2 * jnp.pi * 392.00 * t[:half]))

g_major = (0.5 * jnp.sin(2 * jnp.pi * 196.00 * t[:half]) +
           0.4 * jnp.sin(2 * jnp.pi * 246.94 * t[:half]) +
           0.3 * jnp.sin(2 * jnp.pi * 293.66 * t[:half]))

signal = jnp.concatenate([c_major, g_major])

# 计算STFT
n_fft = 4096  # 高音高精度
hop_length = 512
window = jnp.hanning(n_fft)

def stft(signal, n_fft, hop_length, window):
    n_frames = 1 + (len(signal) - n_fft) // hop_length
    frames = jnp.stack([
        signal[i * hop_length : i * hop_length + n_fft] * window
        for i in range(n_frames)
    ])
    return jnp.fft.rfft(frames, n=n_fft)

S = stft(signal, n_fft, hop_length, window)
power_spec = jnp.abs(S) ** 2
freqs = jnp.fft.rfftfreq(n_fft, 1.0 / sr)

# 通过将频率仓映射到音高类计算色度图
# 从频率到MIDI音符编号: 69 + 12 * log2(f / 440)
note_names = ['C', 'C#', 'D', 'D#', 'E', 'F', 'F#', 'G', 'G#', 'A', 'A#', 'B']

def freq_to_chroma(freq):
    """将频率映射到音高类（0-11）。freq <= 0时返回-1."""
    midi = 69 + 12 * jnp.log2(jnp.clip(freq, 1e-10, None) / 440.0)
    return jnp.round(midi).astype(int) % 12

# 构建色度图: 为每个音高类求和功率谱能量
chromagram = jnp.zeros((power_spec.shape[0], 12))
valid_freqs = freqs[1:]  # 跳过直流
valid_power = power_spec[:, 1:]

for p in range(12):
    # 找到属于该音高类的频率仓
    chroma_bins = freq_to_chroma(valid_freqs)
    mask = (chroma_bins == p).astype(jnp.float32)
    chromagram = chromagram.at[:, p].set(
        jnp.sum(valid_power * mask[None, :], axis=1)
    )

# 归一化每帧
chromagram = chromagram / (jnp.max(chromagram, axis=1, keepdims=True) + 1e-8)

# 可视化
fig, axes = plt.subplots(3, 1, figsize=(14, 10))

# 波形
axes[0].plot(t[:3000], signal[:3000], color='#3498db', linewidth=0.5,
             label='C大调')
axes[0].plot(t[half:half+3000], signal[half:half+3000], color='#e74c3c',
             linewidth=0.5, label='G大调')
axes[0].set_title('波形: C大调 → G大调')
axes[0].set_ylabel('振幅')
axes[0].set_xlabel('时间 (s)')
axes[0].legend()

# 谱图（对数尺度）
time_axis = jnp.arange(power_spec.shape[0]) * hop_length / sr
axes[1].imshow(jnp.log1p(power_spec[:, :500].T), aspect='auto', origin='lower',
               cmap='magma', extent=[0, time_axis[-1], 0, freqs[500]])
axes[1].set_title('功率谱图')
axes[1].set_ylabel('频率 (Hz)')
axes[1].set_xlabel('时间 (s)')

# 色度图
im = axes[2].imshow(chromagram.T, aspect='auto', origin='lower', cmap='YlOrRd',
                     extent=[0, time_axis[-1], -0.5, 11.5])
axes[2].set_yticks(range(12))
axes[2].set_yticklabels(note_names)
axes[2].set_title('色度图（随时间变化的音高类能量）')
axes[2].set_ylabel('音高类')
axes[2].set_xlabel('时间 (s)')
plt.colorbar(im, ax=axes[2], fraction=0.046, label='归一化能量')

# 标记预期的活跃音高类
mid_frame = chromagram.shape[0] // 2
print(f"C大调区域 - 预期: C, E, G")
print(f"  色度值: {dict(zip(note_names, [f'{v:.2f}' for v in chromagram[mid_frame//2]]))}")
print(f"G大调区域 - 预期: G, B, D")
print(f"  色度值: {dict(zip(note_names, [f'{v:.2f}' for v in chromagram[mid_frame + mid_frame//2]]))}")

plt.tight_layout()
plt.show()