《Delta-Sigma数据转换器从入门到精通》笔记之漫谈均值

文章目录

《Delta-Sigma数据转换器从入门到精通》笔记之漫谈均值
前言
一、最小方差无偏估计
二、FIR滤波器
三、最小二乘估计
四、最大似然估计
五、几何
六、协方差
写在最后

前言

“均值”是δ-σ ADC CIC滤波器的重要概念，深刻地理解“均值”是理解δ-σ ADC 以及其他数据处理方法的前提。《Delta-Sigma数据转换器从入门到精通》反复涉及到这一概念，但是作者并未做深入的介绍，本文尝试着从多个视角对此加以解释。

一、最小方差无偏估计

“均值”应该是所有数据处理方法中最常用的也是最朴素的了。试想以下一个情形，一位工程师被要求设计一个测量电池电压（VBAT）系统，这位工程师很快地交出了答卷——连续取样N次，求样本均值，即：
X = [ x [ 1 ] x [ 2 ] x [ 3 ] . . . x [ N ] ] T \mathbf{X}\text{=}{{\left[ \begin{matrix} x[1] & x[2] & x[3] & ... & x[N] \\ \end{matrix} \right]}^{T}} X=[x[1]x[2]x[3]...x[N]]T
x ^ = x ‾ = 1 N ∑ n = 1 N x [ n ] (1) \widehat{x}=\overline{x}=\frac{1}{N}\sum\limits_{\text{n}=1}^{N}{\text{x}[n]}\tag{1} x =x=N1n=1∑Nx[n](1)

问题来了，取样一次不就足够了吗，为什么要取多次？ x ^ \widehat{x} x 为什么不是从 X \mathbf{X} X中随机取一个 x [ n ] x[\text{n}] x[n] ，或者是 x [ n ] x[\text{n}] x[n]的其他线性组合 x ^ = ∑ n = 1 N a n x [ n ] \widehat{x}=\sum\limits_{\text{n}=1}^{N}{{{\text{a}}_{\text{n}}}x[n]} x =n=1∑Nanx[n]?糟糕的工程师会这样解释：这是祖传的，经过实践的，运行的好好的，又没出过问题。一般的工程师会这样解释：测量难免有误差，多测几次才准，至于取均值，也很简单， x [ n ] x[\text{n}] x[n]之间都是平等的，不能厚此薄彼。但是有追求的工程师不满足于这样朴素的回答，我们对此问题进行建模，围绕着模型展开论述。
x [ n ] = θ + ω [ n ] n=1,2, . . . N (2) x[\text{n}]\text{=}\theta \text{+}\omega [\text{n}]\begin{matrix} {} & {} & {} & {} \\ \end{matrix}\text{n=1,2,}...N\tag{2} x[n]=θ+ω[n]n=1,2,...N(2)
其中， θ ( θ ∈ (- ∞ ,+ ∞ )) \theta \begin{matrix} {} \\ \end{matrix}\text{(}\theta \in \text{(-}\infty \text{,+}\infty \text{))} θ(θ∈(-∞,+∞))是要估计（通俗的说，要测量的）参数，误差 ω [ n ] \omega [\text{n}] ω[n]被假设成WGN（高斯白噪声），即 ω [ n ] ∼ N ⁡ (0, σ 2 ) \omega [\text{n}]\sim \operatorname{N}\text{(0,}{{\sigma }^{\text{2}}}\text{)} ω[n]∼N(0,σ2)，这一假设是符合工程师经验的，工程上遇到的许多误差项都可以近似成WGN，有的难以得知误差的CDF，干脆假设成gauss好了。
这一模型，过于理想，因为随着时间的推移 θ \theta θ 总是会减小的。但是，还好在足够短的时间内，模型是够用的。
工程师们常常认为数据是上帝的，但是有的工程师认为，一旦模型建立，论述就要以模型为中心了，这时候，模型就是上帝。我赞同这种观点。
接着式（1）， θ \theta θ 的一个估计 x [ n ] x[\text{n}] x[n]的均值
θ ^ = 1 N ∑ n = 1 N x [ n ] (3) \widehat{\theta }=\frac{1}{N}\sum\limits_{\text{n}=1}^{N}{x[n]}\tag{3} θ =N1n=1∑Nx[n](3)
求 θ ^ \widehat{\theta } θ 期望
E( θ ^ ) = E( 1 N ∑ n = 1 N x [ n ] ) \text{E(}\widehat{\theta })=\text{E(}\frac{1}{N}\sum\limits_{\text{n}=1}^{N}{x[n]}) E(θ )=E(N1n=1∑Nx[n])
E( θ ^ ) = 1 N ∑ n = 1 N E ( x [ n ] ) \text{E(}\widehat{\theta })=\frac{1}{N}\sum\limits_{\text{n}=1}^{N}{\text{E}(x[n]}) E(θ )=N1n=1∑NE(x[n])
E( θ ^ ) = 1 N ∑ n = 1 N θ \text{E(}\widehat{\theta })=\frac{1}{N}\sum\limits_{\text{n}=1}^{N}{\theta } E(θ )=N1n=1∑Nθ
E( θ ^ ) = θ (4) \text{E(}\widehat{\theta })=\theta\tag{4} E(θ )=θ(4)

求 θ ^ \widehat{\theta } θ 方差
var( θ ^ ) = var ⁡ ( 1 N ∑ n = 1 N x [ n ] ) \text{var(}\widehat{\theta })=\operatorname{var}\text{(}\frac{1}{N}\sum\limits_{\text{n}=1}^{N}{x[n]}) var(θ )=var(N1n=1∑Nx[n])
var( θ ^ ) = 1 N 2 ∑ n = 1 N var ( x [ n ] ) \text{var(}\widehat{\theta })=\frac{1}{{{N}^{2}}}\sum\limits_{\text{n}=1}^{N}{\text{var}(x[n]}) var(θ )=N21n=1∑Nvar(x[n])
var( θ ^ ) = 1 N 2 ∑ n = 1 N σ 2 \text{var(}\widehat{\theta })=\frac{1}{{{N}^{2}}}\sum\limits_{\text{n}=1}^{N}{{{\sigma }^{2}}} var(θ )=N21n=1∑Nσ2
var( θ ^ ) = σ 2 N (5) \text{var(}\widehat{\theta })=\frac{{{\sigma }^{2}}}{N}\tag{5} var(θ )=Nσ2(5)
显然的， θ ^ \widehat{\theta } θ 是 θ \theta θ 的一个无偏估计量，并且方差仅为单次取样的1/N，当 N → ∞ \text{N}\to \infty N→∞， var( θ ^ ) → 0 \text{var(}\widehat{\theta })\to 0 var(θ )→0。 θ \theta θ 的无偏估计量当然还有很多，任意的
θ ~ = ∑ n = 1 N a n x [ n ] ( ∑ n = 1 N a n = 1 ) \widetilde{\theta }=\sum\limits_{\text{n}=1}^{N}{{{\text{a}}_{\text{n}}}x[n]}\begin{matrix} {} & {} \\ \end{matrix}(\sum\limits_{\text{n}=1}^{N}{{{\text{a}}_{\text{n}}}\text{=}1}) θ =n=1∑Nanx[n](n=1∑Nan=1)
都是无偏的。那么 θ ^ \widehat{\theta } θ 是最优的吗？如果以最小方差作为判断依据，答案是肯定的，因为它达到了Cramer-Rao下界。
首先求联合PDF
p ( X ; θ ) = ∏ n = 1 N 1 2 π σ 2 e − ( x [ n ] − θ ) 2 2 σ 2 p(\mathbf{X};\theta )=\prod\limits_{n=1}^{N}{\frac{1}{\sqrt{2\pi {{\sigma }^{2}}}}{{e}^{-\frac{{{\left( x[n]-\theta \right)}^{2}}}{2{{\sigma }^{2}}}}}} p(X;θ)=n=1∏N2πσ2 1e−2σ2(x[n]−θ)2
p ( X ; θ ) = 1 ( 2 π σ 2 ) N 2 e − ∑ n= 1 N ( x [ n ] − θ ) 2 2 σ 2 p(\mathbf{X};\theta )=\frac{1}{{{\left( 2\pi {{\sigma }^{2}} \right)}^{\frac{N}{2}}}}{{e}^{-\frac{\sum\limits_{\text{n=}1}^{\text{N}}{{{\left( x[n]-\theta \right)}^{2}}}}{2{{\sigma }^{2}}}}} p(X;θ)=(2πσ2)2N1e−2σ2n=1∑N(x[n]−θ)2
取对数求一阶导
∂ ln ⁡ p ( X ; θ ) ∂ θ = ∂ ∂ θ ( - ln ⁡ [ ( 2 π σ 2 ) N 2 ] − ∑ n= 1 N ( x [ n ] − θ ) 2 2 σ 2 ) \frac{\partial \ln p(\mathbf{X};\theta )}{\partial \theta }=\frac{\partial }{\partial \theta }\left( \text{-}\ln \left[ {{\left( 2\pi {{\sigma }^{2}} \right)}^{\frac{N}{2}}} \right]-\frac{\sum\limits_{\text{n=}1}^{\text{N}}{{{\left( x[n]-\theta \right)}^{2}}}}{2{{\sigma }^{2}}} \right) ∂θ∂lnp(X;θ)=∂θ∂ -ln[(2πσ2)2N]−2σ2n=1∑N(x[n]−θ)2
∂ ln ⁡ p ( X ; θ ) ∂ θ = ∑ n= 1 N ( x [ n ] − θ ) σ 2 \frac{\partial \ln p(\mathbf{X};\theta )}{\partial \theta }=\frac{\sum\limits_{\text{n=}1}^{\text{N}}{\left( x[n]-\theta \right)}}{{{\sigma }^{2}}} ∂θ∂lnp(X;θ)=σ2n=1∑N(x[n]−θ)
∂ ln ⁡ p ( X ; θ ) ∂ θ = N ( x ‾ − θ ) σ 2 \frac{\partial \ln p(\mathbf{X};\theta )}{\partial \theta }=\frac{N\left( \overline{x}-\theta \right)}{{{\sigma }^{2}}} ∂θ∂lnp(X;θ)=σ2N(x−θ)
∂ ln ⁡ p ( X ; θ ) ∂ θ = N x ‾ σ 2 - N θ σ 2 \frac{\partial \ln p(\mathbf{X};\theta )}{\partial \theta }=\frac{N\overline{x}}{{{\sigma }^{2}}}\text{-}\frac{\text{N}\theta }{{{\sigma }^{2}}} ∂θ∂lnp(X;θ)=σ2Nx-σ2Nθ

求二阶导
∂ 2 ln ⁡ p ( X ; θ ) ∂ θ 2 = ∂ ∂ θ ( N x ‾ σ 2 - N θ σ 2 ) \frac{{{\partial }^{2}}\ln p(\mathbf{X};\theta )}{\partial {{\theta }^{2}}}=\frac{\partial }{\partial \theta }\left( \frac{N\overline{x}}{{{\sigma }^{2}}}\text{-}\frac{\text{N}\theta }{{{\sigma }^{2}}} \right) ∂θ2∂2lnp(X;θ)=∂θ∂(σ2Nx-σ2Nθ)
∂ 2 ln ⁡ p ( X ; θ ) ∂ θ 2 = - N σ 2 \frac{{{\partial }^{2}}\ln p(\mathbf{X};\theta )}{\partial {{\theta }^{2}}}=\text{-}\frac{\text{N}}{{{\sigma }^{2}}} ∂θ2∂2lnp(X;θ)=-σ2N
即
1 -E ⁣ ⁣ [ ⁣ ⁣ ∂ 2 ln ⁡ p ( X ; θ ) ∂ θ 2 ⁣ ⁣ ] ⁣ ⁣ = σ 2 N \frac{1}{\text{-E }\!\![\!\!\text{ }\frac{{{\partial }^{2}}\ln p(\mathbf{X};\theta )}{\partial {{\theta }^{2}}}\text{ }\!\!]\!\!\text{ }}=\frac{{{\sigma }^{2}}}{\text{N}} -E [ ∂θ2∂2lnp(X;θ) ] 1=Nσ2
据此，均值作为该模型参数 θ \theta θ的最佳估计得到了充分的证明，对于一名工程师而言，可以放心地使用了。

二、FIR滤波器

方差减小，意味着直观上，输出更为平滑了，暗示着我们均值是一个低通滤波器。审视式（2）给出的模型，可以将噪声 ω [ n ] \omega [\text{n}] ω[n]看成交流分量，而 θ \theta θ 是直流分量，显然，我们需要一个低通滤波器提取 θ \theta θ。仔细审视式（1），可以发现它同时就是FIR滤波器，其Z变换为
H(Z)= ∑ n = 0 N − 1 1 N Z − n \text{H(Z)=}\sum\limits_{n=0}^{N-1}{\frac{1}{\text{N}}{{Z}^{-n}}} H(Z)=n=0∑N−1N1Z−n
H(Z)= 1 N ∑ n = 0 N − 1 Z − n (6) \text{H(Z)=}\frac{1}{\text{N}}\sum\limits_{n=0}^{N-1}{{{Z}^{-n}}}\tag{6} H(Z)=N1n=0∑N−1Z−n(6)
当N确定时，利用MATLAB可以画出幅频/相频特性曲线。

工程师们如果对你说，均值可以对数据进行滤波，输出更为平滑，请不要怀疑，他们是在用通俗的语言在向你解释FIR滤波器。
和FIR相对应的，还有一个叫IIR的滤波器，其差分方程如下式。
y ( n ) + ∑ p = 1 P a p y ( n − p ) = ∑ q = 0 Q b q x ( n − q ) y(n)+\sum\limits_{p=1}^{P}{{{a}_{p}}y(n-p)}=\sum\limits_{q=0}^{Q}{{{b}_{q}}x(n-q)} y(n)+p=1∑Papy(n−p)=q=0∑Qbqx(n−q)
y ( n ) = ∑ p = 1 P - a p y ( n − p ) + ∑ q = 0 Q b q x ( n − q ) (7) y(n)=\sum\limits_{p=1}^{P}{\text{-}{{a}_{p}}y(n-p)}\text{+}\sum\limits_{q=0}^{Q}{{{b}_{q}}x(n-q)}\tag{7} y(n)=p=1∑P-apy(n−p)+q=0∑Qbqx(n−q)(7)
Z变换可得
H(Z) = ∑ q = 0 Q b q Z - q 1 - ∑ p = 1 P a p Z - p \text{H(Z)}=\frac{\sum\limits_{q=0}^{Q}{{{b}_{q}}{{\text{Z}}^{\text{-}q}}}}{1\text{-}\sum\limits_{p=1}^{P}{{{a}_{p}}{{\text{Z}}^{\text{-}p}}}} H(Z)=1-p=1∑PapZ-pq=0∑QbqZ-q

对于式（7），我常常喜欢把它写成
y ( n ) = ∑ p = 1 P a p y ( n − p ) + ∑ q = 0 Q b q x ( n − q ) y(n)=\sum\limits_{p=1}^{P}{{{a}_{p}}y(n-p)}\text{+}\sum\limits_{q=0}^{Q}{{{b}_{q}}x(n-q)} y(n)=p=1∑Papy(n−p)+q=0∑Qbqx(n−q)
前者是差分方程的理解，后者意味着输出是历史输入和历史输出的加权组合，ARM公司提供的库也是这样设计的。当然，这只是个人理解，更多的场合，还是使用式（7）更为方便。相比于FIR，IIR意味着什么呢？多了分母，意味着反馈？很多工程师就是这么理解的。我们重新审视式（6），对其乘以（1-Z-1）再除以（1-Z-1）
H(Z)= 1 N ( ∑ n = 0 N − 1 Z − n ) ( 1 − Z − 1 ) 1 − Z − 1 \text{H(Z)=}\frac{1}{\text{N}}\frac{\left( \sum\limits_{n=0}^{N-1}{{{Z}^{-n}}} \right)\left( 1-{{Z}^{-1}} \right)}{1-{{Z}^{-1}}} H(Z)=N11−Z−1(n=0∑N−1Z−n)(1−Z−1)
H(Z)= 1 N ( 1 + Z − 1 + Z − 1 + . . . + Z − ( N − 1 ) ) ( 1 − Z − 1 ) 1 − Z − 1 \text{H(Z)=}\frac{1}{\text{N}}\frac{\left( 1+{{Z}^{-1}}+{{Z}^{-1}}+...+{{Z}^{-(N-1)}} \right)\left( 1-{{Z}^{-1}} \right)}{1-{{Z}^{-1}}} H(Z)=N11−Z−1(1+Z−1+Z−1+...+Z−(N−1))(1−Z−1)
H(Z)= 1 N 1 − Z − N 1 − Z − 1 (8) \text{H(Z)=}\frac{1}{\text{N}}\frac{1-{{Z}^{-N}}}{1-{{Z}^{-1}}}\tag{8} H(Z)=N11−Z−11−Z−N(8)
FIR便写成了IIR的形式。用MATLAB绘制一下N=10的幅频/相频特性曲线

完全没有变化，只是matlab也误以为它是IIR

由此工程师们应该明白，判断一个滤波器到底是FIR还是IIR，其Z变换是否有分母并不是一个判据。
式（6）是均值滤波（具体的讲，是滑动均值滤波），和它等价的式（8）可以理解成输入累加（分母）减去N项前的输入累加（分子），即数据先经过一个累加器，再经过一个减去N阶差分的减法器，也即式（8）进一步变换成下式
H(Z)= 1 N × 1 1 − Z − 1 × ( 1 − Z − N ) \text{H(Z)=}\frac{1}{\text{N}}\times \frac{1}{1-{{Z}^{-1}}}\times \left( 1-{{Z}^{-N}} \right) H(Z)=N1×1−Z−11×(1−Z−N)
H(Z)= 1 N × H 1 ( Z ) × H 2 ( Z ) \text{H(Z)=}\frac{1}{\text{N}}\times {{\text{H}}_{\text{1}}}\left( Z \right)\times {{\text{H}}_{\text{2}}}\left( Z \right) H(Z)=N1×H1(Z)×H2(Z)
当然，也可以反过来理解，先经过减法器后经过累加器，即
H(Z)= 1 N × H 2 ( Z ) × H 1 ( Z ) \text{H(Z)=}\frac{1}{\text{N}}\times {{\text{H}}_{\text{2}}}\left( Z \right)\times {{\text{H}}_{\text{1}}}\left( Z \right) H(Z)=N1×H2(Z)×H1(Z)
理解了这一点，δ-σ ADC 中的CIC滤波器就可以理解一半了。再次审视式（6），它有一项“ 1 N \frac{1}{\text{N}} N1”，我们在阅读ADI提供的PDM麦克风驱动源码的时候发现，并没有这一项，这不稀奇，只要把公式乘以N，等效成输出放大N倍即可，这样做的好处是没有了除法，不会出现浮点运算带来的误差累积，也方便后续的各种定点滤波。

三、最小二乘估计

我们再次审视式（2）所给出的模型，概率论学得好的工程师会发现这是一个平稳随机过程，LS（最小二乘）和ML（最大似然）的方法便在我们脑海中产生。本节先介绍LS。
构造目标函数
J( θ )= ∑ n= 1 N ( x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ - θ ) 2 \text{J(}\theta \text{)=}\sum\limits_{\text{n=}1}^{\text{N}}{{{\left( \text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ -}\theta \right)}^{2}}} J(θ)=n=1∑N(x [ n ] -θ)2
J( θ )= ∑ n= 1 N ( x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ 2 -2x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ θ + θ 2 ) \text{J(}\theta \text{)=}\sum\limits_{\text{n=}1}^{\text{N}}{\left( \text{x }\!\![\!\!\text{ n}{{\text{ }\!\!]\!\!\text{ }}^{\text{2}}}\text{-2x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }\theta +{{\theta }^{2}} \right)} J(θ)=n=1∑N(x [ n ] 2-2x [ n ] θ+θ2)
可见是一个二阶函数，可以求极值。令
∂ J( θ ) ∂ θ = 0 \frac{\partial \text{J(}\theta \text{)}}{\partial \theta }\text{=}0 ∂θ∂J(θ)=0
即
∂ θ ∑ n= 1 N ( x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ 2 -2x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ θ + θ 2 ) = 0 \frac{\partial }{\theta }\sum\limits_{\text{n=}1}^{\text{N}}{\left( \text{x }\!\![\!\!\text{ n}{{\text{ }\!\!]\!\!\text{ }}^{\text{2}}}\text{-2x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }\theta +{{\theta }^{2}} \right)}\text{=}0 θ∂n=1∑N(x [ n ] 2-2x [ n ] θ+θ2)=0
∑ n= 1 N ( -2x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ + 2 θ ) = 0 \sum\limits_{\text{n=}1}^{\text{N}}{\left( \text{-2x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }+2\theta \right)}\text{=}0 n=1∑N(-2x [ n ] +2θ)=0
N θ = ∑ n= 1 N x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ \text{N}\theta \text{=}\sum\limits_{\text{n=}1}^{\text{N}}{\text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }} Nθ=n=1∑Nx [ n ]
θ = 1 N ∑ n= 1 N x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ \theta \text{=}\frac{\text{1}}{\text{N}}\sum\limits_{\text{n=}1}^{\text{N}}{\text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }} θ=N1n=1∑Nx [ n ]

四、最大似然估计

ML的方法和上文求Cramer-Rao下界差不多，许多推导步骤因此跳过。由于 ω [ n ] ∼ N ⁡ (0, σ 2 ) \omega [\text{n}]\sim \operatorname{N}\text{(0,}{{\sigma }^{\text{2}}}\text{)} ω[n]∼N(0,σ2)，因此 x [ n ] x[n] x[n]的联合PDF为
p ( X ; θ ) = ∏ n = 1 N 1 2 π σ 2 e − ( x [ n ] − θ ) 2 2 σ 2 p(\mathbf{X};\theta )=\prod\limits_{n=1}^{N}{\frac{1}{\sqrt{2\pi {{\sigma }^{2}}}}{{e}^{-\frac{{{\left( x[n]-\theta \right)}^{2}}}{2{{\sigma }^{2}}}}}} p(X;θ)=n=1∏N2πσ2 1e−2σ2(x[n]−θ)2
p ( X ; θ ) = 1 ( 2 π σ 2 ) N 2 e − ∑ n= 1 N ( x [ n ] − θ ) 2 2 σ 2 p(\mathbf{X};\theta )=\frac{1}{{{\left( 2\pi {{\sigma }^{2}} \right)}^{\frac{N}{2}}}}{{e}^{-\frac{\sum\limits_{\text{n=}1}^{\text{N}}{{{\left( x[n]-\theta \right)}^{2}}}}{2{{\sigma }^{2}}}}} p(X;θ)=(2πσ2)2N1e−2σ2n=1∑N(x[n]−θ)2
取对数并且令其一阶导为0
∂ ln ⁡ p ( X ; θ ) ∂ θ = 0 \frac{\partial \ln p(\mathbf{X};\theta )}{\partial \theta }=0 ∂θ∂lnp(X;θ)=0
N x ‾ σ 2 - N θ σ 2 = 0 \frac{N\overline{x}}{{{\sigma }^{2}}}\text{-}\frac{\text{N}\theta }{{{\sigma }^{2}}}\text{=}0 σ2Nx-σ2Nθ=0
θ = x ‾ = 1 N ∑ n= 1 N x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ \theta \text{=}\overline{\text{x}}\text{=}\frac{\text{1}}{\text{N}}\sum\limits_{\text{n=}1}^{\text{N}}{\text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }} θ=x=N1n=1∑Nx [ n ]

五、几何

我不是一个聪明的工程师，常常难以接受纯粹的公式推导，总是试图探究其背后的物理意义，特别是其几何意义。也正是这个原因，我更喜欢“均值”、“方差”的另一种充满直观的称呼：“一阶原点矩”、“二阶中心距”，也喜欢将方差看成正方形，将协方差看成矩形。也只有借助几何的概念，我才能理解微分方程的通解，特解，齐次解、理解傅里叶变换、理解Hilbert空间、理解线性代数、理解…。确实是笨了点，这或许注定了我只是一名普通的工程师，为了一口饭四处奔波。
在这里，我也从几何的视角解释均值。考察观测向量
X = [ x [ 1 ] x [ 2 ] x [ 3 ] . . . x [ N ] ] T \mathbf{X}\text{=}{{\left[ \begin{matrix} x[1] & x[2] & x[3] & ... & x[N] \\ \end{matrix} \right]}^{T}} X=[x[1]x[2]x[3]...x[N]]T
和目标向量
θ = [ θ θ θ . . . θ ] T \mathbf{\theta }\text{=}{{\left[ \begin{matrix} \theta & \theta & \theta & ... & \theta \\ \end{matrix} \right]}^{T}} θ=[θθθ...θ]T
θ = θ [ 1 1 1 . . . 1 ] T \mathbf{\theta }\text{=}\theta {{\left[ \begin{matrix} 1 & 1 & 1 & ... & 1 \\ \end{matrix} \right]}^{T}} θ=θ[111...1]T
在Euclidean空间 R N {{\text{R}}^{\text{N}}} RN当中，向量 θ \mathbf{\theta } θ张成了1维子空间 span ( θ ) ∈ R N \text{span}(\mathbf{\theta })\in {{\text{R}}^{\text{N}}} span(θ)∈RN，我们的目标是在该空间当中寻找一个向量 θ ^ \mathbf{\hat{\theta }} θ^，该向量最接近 X \mathbf{X} X。所谓最接近，是指误差向量范数 ∥ e ∥ = ∥ X - θ ^ ∥ \left\| \mathbf{e} \right\|\text{=}\left\| \mathbf{X}\text{-}\mathbf{\hat{\theta }} \right\| ∥e∥= X-θ^ 最小，显然的，只有span( θ \mathbf{\theta } θ)和 e \mathbf{e} e正交时， ∥ e ∥ \left\| \mathbf{e} \right\| ∥e∥最小，这时候 θ ^ \mathbf{\hat{\theta }} θ^就是 X \mathbf{X} X在span( θ \mathbf{\theta } θ)上的投影，一个内积公式便可以解决了。
首先求向量内积。
< θ , X > = θ T X <\mathbf{\theta },\mathbf{X}>={{\mathbf{\theta }}^{T}}\mathbf{X} <θ,X>=θTX
< θ , X > = θ ( x ⁣ ⁣ [ ⁣ ⁣ 1 ⁣ ⁣ ] ⁣ ⁣ +x ⁣ ⁣ [ ⁣ ⁣ 2 ⁣ ⁣ ] ⁣ ⁣ + . . . +x ⁣ ⁣ [ ⁣ ⁣ N ⁣ ⁣ ] ⁣ ⁣ ) <\mathbf{\theta },\mathbf{X}>=\theta \left( \text{x }\!\![\!\!\text{ 1 }\!\!]\!\!\text{ +x }\!\![\!\!\text{ 2 }\!\!]\!\!\text{ +}...\text{+x }\!\![\!\!\text{ N }\!\!]\!\!\text{ } \right) <θ,X>=θ(x [ 1 ] +x [ 2 ] +...+x [ N ] )
< θ , X > = θ ∑ n = 1 N x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ <\mathbf{\theta },\mathbf{X}>=\theta \sum\limits_{n=1}^{N}{\text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }} <θ,X>=θn=1∑Nx [ n ]

由于 ∥ θ ∥ ≠ 1 \left\| \mathbf{\theta } \right\|\ne 1 ∥θ∥=1，因此需要做单位化
< θ , θ > = θ T θ =N θ 2 <\mathbf{\theta },\mathbf{\theta }>={{\mathbf{\theta }}^{T}}\mathbf{\theta }\text{=N}{{\theta }^{2}} <θ,θ>=θTθ=Nθ2
由此可得
θ ^ = θ < θ , X > < θ , θ > \mathbf{\hat{\theta }}\text{=}\mathbf{\theta }\frac{<\mathbf{\theta },\mathbf{X}>}{<\mathbf{\theta },\mathbf{\theta }>} θ^=θ<θ,θ><θ,X>
θ ^ = θ [ 1 1 1 . . . 1 ] T θ ∑ n = 1 N x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ N θ 2 \mathbf{\hat{\theta }}\text{=}\theta {{\left[ \begin{matrix} 1 & 1 & 1 & ... & 1 \\ \end{matrix} \right]}^{T}}\frac{\theta \sum\limits_{n=1}^{N}{\text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }}}{\text{N}{{\theta }^{2}}} θ^=θ[111...1]TNθ2θn=1∑Nx [ n ]
θ ^ = 1 N ∑ n = 1 N x ⁣ ⁣ [ ⁣ ⁣ n ⁣ ⁣ ] ⁣ ⁣ [ 1 1 1 . . . 1 ] T \mathbf{\hat{\theta }}\text{=}\frac{1}{\text{N}}\sum\limits_{n=1}^{N}{\text{x }\!\![\!\!\text{ n }\!\!]\!\!\text{ }}{{\left[ \begin{matrix} 1 & 1 & 1 & ... & 1 \\ \end{matrix} \right]}^{T}} θ^=N1n=1∑Nx [ n ] [111...1]T
嘿，真形象、真直观、真简单、真有趣！

六、协方差

我们再把视线放回δ-σ ADC的滑动均值滤波中。正如上文所述，滤波器输出变得平滑了，因为输出序列的每一项方差都为 σ 2 N \frac{{{\sigma }^{2}}}{N} Nσ2，每一项都更可能地接近 θ \theta θ，我们自然地想到，输出序列项之间相互更为接近了。这当然是朴素的想法，作为一名工程师，应该学会用数学的语言去阐述。
对于任意的s和t， ω [ n ] \omega [\text{n}] ω[n]序列的自协方差
γ ω ( s , t ) = cov ⁡ ( ω s , ω t ) {{\gamma }_{\omega }}(s,t)=\operatorname{cov}({{\omega }_{s}},{{\omega }_{t}}) γω(s,t)=cov(ωs,ωt)
γ ω ( s , t ) = { σ 2 s=t 0 s ≠ t {{\gamma }_{\omega }}(s,t)=\left\{ \begin{matrix} {{\sigma }^{2}}\text{ s=t} \\ 0\text{ s}\ne \text{t} \\ \end{matrix} \right. γω(s,t)={σ2 s=t0 s=t
滤波后的序列
γ ω ‾ ( s , t ) = cov ⁡ ( ω s ‾ , ω t ‾ ) {{\gamma }_{\overline{\omega }}}(s,t)=\operatorname{cov}({{\omega }_{\overline{s}}},{{\omega }_{\overline{\text{t}}}}) γω(s,t)=cov(ωs,ωt)
γ ω ‾ ( s , t ) = cov ⁡ ( 1 N ∑ n= 0 N-1 ω s-n , 1 N ∑ n= 0 N-1 ω t-n ) {{\gamma }_{\overline{\omega }}}(s,t)=\operatorname{cov}(\frac{1}{\text{N}}\sum\limits_{\text{n=}0}^{\text{N-1}}{{{\omega }_{\text{s-n}}}},\frac{1}{\text{N}}\sum\limits_{\text{n=}0}^{\text{N-1}}{{{\omega }_{\text{t-n}}}}) γω(s,t)=cov(N1n=0∑N-1ωs-n,N1n=0∑N-1ωt-n)
γ ω ‾ ( s , t ) = 1 N 2 cov ⁡ ( ∑ n= 0 N-1 ω s-n , ∑ n= 0 N-1 ω t-n ) {{\gamma }_{\overline{\omega }}}(s,t)=\frac{1}{{{\text{N}}^{\text{2}}}}\operatorname{cov}(\sum\limits_{\text{n=}0}^{\text{N-1}}{{{\omega }_{\text{s-n}}}},\sum\limits_{\text{n=}0}^{\text{N-1}}{{{\omega }_{\text{t-n}}}}) γω(s,t)=N21cov(n=0∑N-1ωs-n,n=0∑N-1ωt-n)
γ ω ‾ ( s , t ) = 1 N 2 E ( ∑ n= 0 N-1 ω s-n × ∑ n= 0 N-1 ω t-n ) {{\gamma }_{\overline{\omega }}}(s,t)=\frac{1}{{{\text{N}}^{\text{2}}}}E\left( \sum\limits_{\text{n=}0}^{\text{N-1}}{{{\omega }_{\text{s-n}}}}\times \sum\limits_{\text{n=}0}^{\text{N-1}}{{{\omega }_{\text{t-n}}}} \right) γω(s,t)=N21E(n=0∑N-1ωs-n×n=0∑N-1ωt-n)
γ ω ‾ ( s , t ) = 1 N 2 E [ ( ω s + ω s-1 + . . . + ω s-(N-1) ) ( ω t + ω t-1 + . . . + ω t-(N-1) ) ] {{\gamma }_{\overline{\omega }}}(s,t)=\frac{1}{{{\text{N}}^{\text{2}}}}E\left[ \left( {{\omega }_{\text{s}}}\text{+}{{\omega }_{\text{s-1}}}\text{+}...\text{+}{{\omega }_{\text{s-(N-1)}}} \right)\left( {{\omega }_{\text{t}}}\text{+}{{\omega }_{\text{t-1}}}\text{+}...\text{+}{{\omega }_{\text{t-(N-1)}}} \right) \right] γω(s,t)=N21E[(ωs+ωs-1+...+ωs-(N-1))(ωt+ωt-1+...+ωt-(N-1))]
γ ω ‾ ( s , t ) = { 1 N 2 ∑ n= 0 N-1- ∣ s-t ∣ σ 2 ∣ s − t ∣ ≤ N − 1 0 ∣ s − t ∣ > N − 1 {{\gamma }_{\overline{\omega }}}(s,t)=\left\{ \begin{matrix} \frac{1}{{{\text{N}}^{\text{2}}}}\sum\limits_{\text{n=}0}^{\text{N-1-}\left| \text{s-t} \right|}{{{\sigma }^{2}}}\begin{matrix} {} & {} \\ \end{matrix}\left| s-t \right|\le N-1 \\ 0\begin{matrix} {} & {} & {} & {} \\ \end{matrix}\left| s-t \right|>N-1 \\ \end{matrix} \right. γω(s,t)=⎩ ⎨ ⎧N21n=0∑N-1-∣s-t∣σ2∣s−t∣≤N−10∣s−t∣>N−1
由式(9)可知，对于 ∣ s − t ∣ ≤ N − 1 \left| s-t \right|\le N-1 ∣s−t∣≤N−1， ω ‾ [ n ] \overline{\omega }[\text{n}] ω[n]序列的自协方差 γ ω ‾ ( s , t ) {{\gamma }_{\overline{\omega }}}(s,t) γω(s,t)＞0，而且 ∣ s − t ∣ \left| s-t \right| ∣s−t∣越小， γ ω ‾ ( s , t ) {{\gamma }_{\overline{\omega }}}(s,t) γω(s,t)越大，越相关，自然地， ω ‾ [ n ] \overline{\omega }[\text{n}] ω[n]相邻项之间越接近，直观的来看，输出越平滑。

写在最后

夜深了，该结束了，而“均值”的故事还远没有结束，比如上文提及的 θ ∈ (- ∞ ,+ ∞ ) \theta \in \text{(-}\infty \text{,+}\infty \text{)} θ∈(-∞,+∞)，如果我们先验地认为 θ ∈ (a,b) \theta \in \text{(a,b)} θ∈(a,b)，会怎么样呢？手机上的锂离子电池电压，它的值就可以认为是 θ ∈ [ 0 ,4.2 ⁣ ⁣ ] ⁣ ⁣ \theta \in [0\text{,4}\text{.2 }\!\!]\!\!\text{ } θ∈[0,4.2 ] ，这种条件下数据该如何处理？，工程师们常常使用的中值滤波、去除离群值后再取均值，它的统计参数又该发生什么变化？
这些和δ-σ ADC的关系就不太大了，等《Delta-Sigma数据转换器从入门到精通》系列笔记完成后找机会再聊！
晚安，迟哥。