从含噪时间序列估计离散化步长及高斯噪声参数

问题 #

我们有一款用于测量物理过程的传感器，通过PWM传输数据。该传感器会对数据进行离散化，理想情况下读取PWM信号时应仅能看到离散步长的数值——比如信号在200到400之间变化时，只会出现200 225 250 ... 375 400这类数值，对应离散化步长为25。但实际的PWM读取器存在高斯噪声，得到的数值是198 224 261 275这类含噪数据。

现在需要正确估计离散化步长及噪声参数。我曾考虑枚举所有可能的步长，计算每个步长下数值取余的结果，认为当余数与高斯分布拟合度最佳时，对应的步长即为目标步长，拟合度采用基于余数的高斯分布最小二乘和衡量。但这种朴素方法存在明显问题：step与step / 2可能产生相同的拟合效果，虽可通过缩小步长猜测范围解决，但不确定这种基于余数的方法是否为最优方案。

另外，因怀疑步长与噪声参数受环境条件影响，我们不希望混合不同运行的数据，导致可用数据量有限。典型数据点直方图如下：

value   count
682     2
775     8
776     14
807     7
838     3
868     1
869     2
900     1
931     3
993     2
1024    2
1055    2
1117    3
1179    2
1210    1
1241    1
1272    2
1303    1
1334    1
1365    2
1397    1
1427    1
1428    2
1458    1
1521    2
1551    1
1552    1
1583    1
1614    1
1645    1
1676    1
1707    3
1738    1
1769    1
1800    2
1831    1
1862    3
1893    1
1924    2
1955    1
1986    2
2047    1
2048    2
2079    1
2111    1
2141    1
2142    1
2173    1
2204    2
2234    1
2235    2
2266    1
2297    1
2325    1
2359    2
2390    3
2483    3
2514    2
2545    1
2575    1
2607    2
2638    3
2669    1
2731    4
2762    2
2794    2
2823    1
2825    1
2854    1
2858    1
2889    1
2918    2
2980    3
3011    2
3042    2
3071    1
3073    1
3104    2
3134    2
3135    1
3197    3
3228    1
3259    1
3287    1
3289    1
3290    2
3321    2
3351    1
3352    2
3383    1
3414    1
3445    4
3506    1
3508    3
3539    1
3570    3
3601    2
3630    1
3632    1
3661    1
3663    2
3694    2
3723    1
3725    2
3754    2
3756    2
3785    16
3787    7
3818    45
3820    1
3822    1
3825    2
3849    3
3880    2
3911    1
3942    3
3971    1
3973    2
4002    1
4004    2
4035    5
4095    1
4097    3
4128    2
4152    1
4157    2
4191    2
4220    1
4222    1
4251    1
4253    3
4276    1
4313    1
4315    4
4344    1
4346    3
4375    3
4377    4

解决方案 #

1. 改进版余数拟合方法 #

你的初始思路方向正确，可以通过以下方式解决step和step/2混淆的问题：

• 约束步长范围：先从直方图的峰值间距提取候选步长的合理区间。比如观察数据中的峰值（如775/776、807、838），相邻间距在31~32左右，先把步长候选框定在物理合理的范围（比如结合传感器手册给出的可能离散步长，或峰值间距的众数区间）。
• 引入额外拟合指标：除了余数的高斯拟合误差，计算余数分布的峰度——当步长正确时，余数是围绕0（或步长一半，取决于离散化对齐方式）的高斯分布，峰度接近3；若步长为真实值的一半，余数分布会呈现双峰，峰度会明显偏离3。
• 加权最小二乘：对计数多的数值点赋予更高拟合权重，这类点更接近真实离散值，噪声影响更小。

2. 基于峰值间距的统计方法 #

理想离散数据含噪后，直方图会呈现多个对应真实离散值的峰值：

• 提取峰值：先对直方图做平滑处理（比如移动平均），再找出所有局部峰值的位置。
• 计算峰值间距：统计所有相邻峰值的距离，取间距的众数作为候选步长。从你给出的数据看，3785与3818间距33、3818与3849间距31，这类接近32的间距占比最高，大概率是真实步长。
• 验证步长：用候选步长将所有峰值对齐到最近的步长倍数，检查是否大部分峰值能匹配离散点；同时计算每个峰值周围数据的高斯噪声参数（均值为离散点，标准差为噪声标准差）。

3. 贝叶斯估计方法（适配小数据量） #

当数据量有限时，贝叶斯方法能更高效利用先验信息：

• 定义模型：假设每个观测值x_i = k_i * step + ε_i，其中k_i为整数，ε_i ~ N(0, σ²)。
• 设置先验：给step设置基于传感器规格的先验分布（比如在合理步长范围内的均匀分布），给σ设置半正态先验（噪声标准差非负）。
• MCMC采样：用马尔可夫链蒙特卡洛（MCMC）方法采样得到step和σ的后验分布，取后验均值或众数作为估计值。这种方法能自然处理小数据量的不确定性，同时避免步长倍数混淆（真实步长的后验概率会远高于其约数）。

4. 噪声参数估计 #

确定步长step后，按以下步骤估计噪声参数：

• 对每个观测值，计算其到最近step倍数的距离，得到噪声样本ε_i = x_i - round(x_i / step) * step。
• 计算这些噪声样本的均值（应接近0，因高斯噪声对称）和标准差，即为噪声参数；若均值偏离0，需额外估计离散化的偏移量（比如真实离散值为base + k*step）。

内容的提问来源于stack exchange，提问作者Igor Bukanov