4. 时域特征提取（一）：短时能量、短时过零率、短时自相关函数、基频估计（F0）的数学原理与Python实现
各位同学，欢迎来到时域特征提取的第一讲。
很多人一上来就扎进频域，觉得FFT才是王道。我个人习惯是，先把手伸进时域摸一摸。为什么？因为时域特征计算量小、物理意义直观，而且很多场景下——比如端点检测、静音切除——它比频域方法更稳。
这一章我们啃四个基础特征：短时能量、短时过零率、短时自相关函数，以及基频估计（F0）。前两个是轻量级选手，后两个稍微重一点，但都是工程落地里的常客。
核心观点：时域特征不是“低级特征”，而是“高性价比特征”。在资源受限的嵌入式设备上，时域特征往往是第一道防线。
4.1 短时能量（Short-Time Energy）
短时能量，说白了就是衡量一帧语音信号有多“响”。数学上很简单：对一帧内的采样点求平方和。
公式长这样：
E = Σ_{n=0}^{N-1} x²(n)其中 x(n) 是第 n 个采样点，N 是帧长。
我在项目中遇到过一个问题：直接用原始波形算能量，结果被直流偏置搞崩了。录音设备的ADC如果没做好隔直，会有一个固定的直流分量，平方之后能量值虚高。所以我的习惯是——先做去直流，再算能量。
工程小技巧：去直流可以用 x(n) - mean(x)，或者直接加一个高通滤波器（比如截止频率 20Hz 的 IIR）。
Python实现很简单：
import numpy as np
def short_time_energy(frame):
"""
计算一帧语音信号的短时能量
:param frame: 一维numpy数组，一帧语音
:return: 能量值（标量）
"""
# 先去直流
frame = frame - np.mean(frame)
# 计算能量
energy = np.sum(frame ** 2)
return energy嗯，这里要注意：如果你做的是实时流式处理，均值要用滑动窗口算，不能等整帧收齐再算——那会引入一个帧的延迟。
4.2 短时过零率（Zero-Crossing Rate, ZCR）
过零率，顾名思义，就是一帧信号穿过零轴的次数。对于离散信号，我们看相邻两个采样点符号是否相反。
公式：
ZCR = (1/2N) * Σ_{n=1}^{N-1} |sgn(x(n)) - sgn(x(n-1))|其中 sgn 是符号函数，正数返回1，负数返回-1。
你想想看，这个特征有什么用？
清音（比如 /s/、/f/）的波形高频成分多，过零率很高；浊音（比如 /a/、/i/）波形比较平滑，过零率低。所以ZCR常用来做清/浊音判别。
我曾经在做一个噪声环境下的VAD系统时，发现单纯用能量阈值会漏掉弱清音。后来加了一个ZCR辅助判断，准确率直接提了8个点。
def zero_crossing_rate(frame):
"""
计算一帧语音信号的短时过零率
:param frame: 一维numpy数组
:return: 过零率（标量，范围0~1）
"""
# 符号函数
sgn = np.sign(frame)
# 计算相邻符号变化
zcr = 0.5 * np.sum(np.abs(np.diff(sgn)))
return zcr / len(frame)避坑指南：如果信号中有小幅噪声抖动，过零率会虚高。我曾经吃过这个亏——安静环境下ZCR居然飙到0.3。解决办法：加一个死区阈值，比如只统计穿过 ±epsilon 区间的过零。
4.3 短时自相关函数（Short-Time Autocorrelation）
自相关函数，用来衡量信号与自身延迟版本之间的相似度。对于周期信号（比如浊音），自相关会在基频周期处出现峰值。
定义式：
R(k) = Σ_{n=0}^{N-1-k} x(n) * x(n+k)其中 k 是延迟（lag），范围从 0 到 N-1。
我个人的经验是：自相关函数在基频估计中非常经典，但它有两个坑——
- 计算量大：直接算 O(N²)，帧长一长就扛不住。
- 倍频干扰：基频的整数倍延迟处也会出现峰值，容易选错。
工程上常用归一化自相关函数（NCCF）来缓解幅度影响：
def autocorrelation(frame, max_lag=None):
"""
计算一帧语音的短时自相关函数
:param frame: 一维numpy数组
:param max_lag: 最大延迟，默认帧长-1
:return: 自相关序列
"""
n = len(frame)
if max_lag is None:
max_lag = n - 1
# 去直流
frame = frame - np.mean(frame)
# 计算自相关
r = np.correlate(frame, frame, mode='full')
# 只取正延迟部分（从 n-1 到 n-1+max_lag）
r = r[n-1 : n-1+max_lag+1]
return r加速技巧：用 FFT 做自相关，复杂度降到 O(N log N)。原理是 R(k) = IFFT(|FFT(x)|²)。我在嵌入式设备上就是这么干的，帧长512时速度能快10倍。
4.4 基频估计（F0 Estimation）
基频，就是声带振动的频率。成年男性一般在 80-200Hz，女性在 150-350Hz。基频估计是语音分析、语音合成、歌唱评价的核心环节。
基于自相关的基频估计流程：
- 计算一帧语音的短时自相关函数 R(k)。
- 在基频对应的延迟范围内（比如 2ms~20ms，对应 50Hz~500Hz）寻找峰值。
- 峰值位置 k_max 对应的基频 F0 = fs / k_max。
但直接找最大值容易出问题——我遇到过好几次，峰值落在 k=1 处（即延迟1个采样点），那其实是高频噪声，不是基频。
所以工程上要加约束：
- 限制搜索范围：只搜索 k_min 到 k_max 之间的峰值。
- 峰值阈值：只接受超过最大值一定比例（比如 0.3）的峰值。
- 后处理平滑：用中值滤波去掉跳变点。
def estimate_f0(frame, fs, f0_min=50, f0_max=500):
"""
基于自相关的基频估计
:param frame: 一帧语音信号
:param fs: 采样率
:param f0_min: 最小基频 (Hz)
:param f0_max: 最大基频 (Hz)
:return: 基频值 (Hz)，无声时返回0
"""
# 计算自相关
r = autocorrelation(frame)
n = len(frame)
# 计算搜索范围（以采样点为单位）
k_min = int(fs / f0_max) # 最大频率对应最小延迟
k_max = int(fs / f0_min) # 最小频率对应最大延迟
# 边界保护
k_min = max(1, k_min)
k_max = min(n - 1, k_max)
# 在搜索范围内找峰值
search_region = r[k_min:k_max+1]
if len(search_region) == 0:
return 0.0
k_peak = np.argmax(search_region) + k_min
# 峰值阈值检查（防止噪声误判）
peak_val = r[k_peak]
if peak_val < 0.3 * r[0]: # r[0] 是自相关最大值
return 0.0
# 计算基频
f0 = fs / k_peak
return f0我曾经踩过的坑：在低信噪比环境下，自相关峰值会被噪声淹没。后来我改用倒谱法做基频估计，抗噪性好了不少。倒谱法的思路是：对频谱取对数再做FFT，基频会在倒谱的基频周期处出现峰值。这个我们后面章节会细讲。
4.5 四个特征的对比与选型
我把这四个特征放在一起对比一下，方便你根据场景选型：
| 特征 | 计算复杂度 | 主要用途 | 工程落地难点 |
|---|---|---|---|
| 短时能量 | O(N) | 端点检测、静音切除 | 直流偏置、幅度归一化 |
| 短时过零率 | O(N) | 清/浊音判别、VAD辅助 | 噪声抖动、死区阈值设置 |
| 短时自相关 | O(N²) / O(N log N) 用FFT | 基频粗估计、周期检测 | 倍频干扰、计算量 |
| 基频估计(F0) | 依赖自相关/倒谱 | 音高检测、语音合成、歌唱评价 | 低信噪比、倍频/半频错误、后处理平滑 |
我个人在项目中的选型原则是：
- 如果只是做VAD，短时能量 + 过零率就够了，别上自相关。
- 如果需要基频，优先用自相关法，简单可控。
- 如果噪声大（比如车载环境），考虑倒谱法或YIN算法。
一句话总结：时域特征不是“低级”，而是“高效”。用对了地方，它们比频域特征更直接、更省资源。
好了，这一章的内容就到这里。四个特征——短时能量、短时过零率、短时自相关、基频估计——它们的数学原理和Python实现我都给你拆清楚了。代码可以直接拿去用，但记得根据你的实际场景调参。