UA MATH567 高维统计III 随机矩阵12 整数环上的区间的应用：DNA序列突变点侦测的统计量及假设检验

UA MATH567 高维统计III 随机矩阵12 整数环上的区间的应用：拐点侦测的统计量及假设检验

整数环上的区间作为随机变量的下标
最大值的概率不等式
应用：DNA序列突变点侦测

整数环上的区间作为随机变量的下标

我们可以在整数上定义区间，虽然只是一个toy definition，但我们可以看到它也是有相关应用的。假设l,r∈Zl,r \in \mathbb{Z}l,r∈Z，称有如下形式的集合为整数环上的区间，
{l+1,l+2,⋯,r}\{l+1,l+2,\cdots,r\}{l+1,l+2,⋯,r}

我们引入两种记号来表示，第一种是(l,r)(l,r)(l,r)，称这种表示方法为区间表示；第二种是(0,⋯,0⏟l个,1,⋯,1⏟(r−l)个,0,⋯)(\underbrace{0,\cdots,0}_{l个},\underbrace{1,\cdots,1}_{(r-l)个},0,\cdots)(l个0,⋯,0,(r−l)个1,⋯,1,0,⋯)，称这种表示方法为向量表示。称r−lr-lr−l是区间的长度，用AnA_nAn表示(0,n)(0,n)(0,n)，我们计数一下它的子区间的数目：

Intervals of length 1: nnn
Intervals of length 2: n−1n-1n−1
Intervals of length kkk: n−(k−1)n-(k-1)n−(k−1)
Intervals of length nnn: 1
empty interval: 1

因此AnA_nAn包含的区间个数为
1+∑k=1n[n−(k−1)]=1+n(n+1)2∝n21+\sum_{k=1}^{n} [n-(k-1)]= 1+\frac{n(n+1)}{2}\propto n^21+k=1∑n[n−(k−1)]=1+2n(n+1)∝n2

我们用In\mathbb{I}_nIn表示AnA_nAn的所有子区间组成的集合。

现在我们引入用区间作为下标的随机变量，假设X1,⋯,Xn∼iidN(0,1)X_1,\cdots,X_n \sim_{iid} N(0,1)X1,⋯,Xn∼iidN(0,1)，∀I∈In\forall I \in \mathbb{I}_n∀I∈In，定义
Z[I]=1∣I∣∑i∈IXiZ[I]=\frac{1}{\sqrt{|I|}}\sum_{i \in I}X_iZ[I]=∣I∣1i∈I∑Xi

其中∣I∣|I|∣I∣表示区间III的长度。下面我们说明几个结论：

结论1 Z[I]∼N(0,1)Z[I] \sim N(0,1)Z[I]∼N(0,1)
For any I∈InI \in \mathbb{I}_nI∈In, let ∣I∣=m+1|I|=m+1∣I∣=m+1, q=inf⁡Iq = \inf Iq=infI, and then the {Xi:i∈I}\{X_i:i \in I\}{Xi:i∈I} can be written as Xq,Xq+1,⋯,Xq+mX_q,X_{q+1},\cdots,X_{q+m}Xq,Xq+1,⋯,Xq+m. So
Z[I]=1m+1∑k=0mXq+kEZ[I]=0,Var(Z[I])=1m+1∑k=0m1=1Z[I]=\frac{1}{\sqrt{m+1}}\sum_{k=0}^{m}X_{q+k} \\ EZ[I]=0,\ Var(Z[I])=\frac{1}{m+1}\sum_{k=0}^{m}1=1Z[I]=m+11k=0∑mXq+kEZ[I]=0, Var(Z[I])=m+11k=0∑m1=1

Since Z[I]Z[I]Z[I] is linear combination of normal random variables, Z[I]∼N(0,1)Z[I] \sim N(0,1)Z[I]∼N(0,1).

结论2 Cov(Z[I],Z[J])=ρ(I,J),Var(Z[I]−Z[J])=2δ2(I,J)Cov(Z[I],Z[J])=\rho(I,J),Var(Z[I]-Z[J])=2\delta^2(I,J)Cov(Z[I],Z[J])=ρ(I,J),Var(Z[I]−Z[J])=2δ2(I,J)，其中
ρ(I,J)=∣I∩J∣∣I∣∣J∣,δ(I,J)=1−ρ(I,J)\rho(I,J)=\frac{|I \cap J|}{\sqrt{|I||J|}},\delta(I,J)=\sqrt{1-\rho(I,J)}ρ(I,J)=∣I∣∣J∣∣I∩J∣,δ(I,J)=1−ρ(I,J)

Note that
∣I∣Z[I]=∣I∩J∣Z[I∩J]+∣I∖J∣Z[I∖J]∣J∣Z[J]=∣I∩J∣Z[I∩J]+∣J∖I∣Z[J∖I]\sqrt{|I|}Z[I]=\sqrt{|I \cap J|}Z[I \cap J]+\sqrt{|I \setminus J|}Z[I \setminus J] \\ \sqrt{|J|}Z[J]=\sqrt{|I \cap J|}Z[I \cap J]+\sqrt{|J \setminus I|}Z[J \setminus I]∣I∣Z[I]=∣I∩J∣Z[I∩J]+∣I∖J∣Z[I∖J]∣J∣Z[J]=∣I∩J∣Z[I∩J]+∣J∖I∣Z[J∖I]
i) Notice Cov(Z[I],Z[J∖I])=Cov(Z[I∖J],Z[J])=0Cov(Z[I],Z[J \setminus I])=Cov(Z[I \setminus J],Z[J])=0Cov(Z[I],Z[J∖I])=Cov(Z[I∖J],Z[J])=0. Thus,
Cov(Z[I],Z[J])=Cov(∣I∩J∣∣I∣Z[I∩J],∣I∩J∣∣J∣Z[I∩J])=∣I∩J∣∣I∣∣J∣=ρ(I,J)Cov(Z[I],Z[J]) \\=Cov \left( \frac{\sqrt{|I \cap J|}}{\sqrt{|I|}}Z[I \cap J], \frac{\sqrt{|I \cap J|}}{\sqrt{|J|}}Z[I \cap J]\right) \\ = \frac{|I \cap J|}{\sqrt{|I||J|}} = \rho(I,J)Cov(Z[I],Z[J])=Cov(∣I∣∣I∩J∣Z[I∩J],∣J∣∣I∩J∣Z[I∩J])=∣I∣∣J∣∣I∩J∣=ρ(I,J)
ii) Directly calculate,
Var(Z[I]−Z[J])=Var([∣I∩J∣∣I∣−∣I∩J∣∣J∣]Z[I∩J]+∣I∖J∣∣I∣Z[I∖J]−∣J∖I∣∣J∣Z[J∖I])=[∣I∩J∣∣I∣−∣I∩J∣∣J∣]2+∣I∖J∣∣I∣+∣J∖I∣∣J∣=∣I∩J∣+∣I∖J∣∣I∣+∣I∩J∣+∣J∖I∣∣J∣−2ρ(I,J)=2[1−ρ(I,J)]=2δ2(I,J)Var(Z[I]-Z[J]) \\ = Var ( \left[ \frac{\sqrt{|I \cap J|}}{\sqrt{|I|}}- \frac{\sqrt{|I \cap J|}}{\sqrt{|J|}}\right]Z[I \cap J] \\+ \frac{\sqrt{|I \setminus J|}}{\sqrt{|I|}} Z[I \setminus J]- \frac{\sqrt{|J \setminus I|}}{\sqrt{|J|}}Z[J \setminus I]) \\ = \left[ \frac{\sqrt{|I \cap J|}}{\sqrt{|I|}}- \frac{\sqrt{|I \cap J|}}{\sqrt{|J|}}\right]^2+\frac{|I \setminus J|}{|I|}+\frac{|J \setminus I|}{|J|} \\ = \frac{|I \cap J|+|I \setminus J|}{|I|}+\frac{|I \cap J|+|J \setminus I|}{|J|}-2\rho(I,J) \\ = 2[1-\rho(I,J)] =2\delta^2(I,J)Var(Z[I]−Z[J])=Var([∣I∣∣I∩J∣−∣J∣∣I∩J∣]Z[I∩J]+∣I∣∣I∖J∣Z[I∖J]−∣J∣∣J∖I∣Z[J∖I])=[∣I∣∣I∩J∣−∣J∣∣I∩J∣]2+∣I∣∣I∖J∣+∣J∣∣J∖I∣=∣I∣∣I∩J∣+∣I∖J∣+∣J∣∣I∩J∣+∣J∖I∣−2ρ(I,J)=2[1−ρ(I,J)]=2δ2(I,J)

最大值的概率不等式

定义
Z∗[In]=max⁡{Z[I]:I∈In}Z^*[\mathbb{I}_n]=\max\{Z[I]:I \in \mathbb{I}_n\}Z∗[In]=max{Z[I]:I∈In}

我们的目标是导出下面的概率不等式：
P(Z∗[In]>2.01log⁡n)→0,asn→∞P(Z^*[\mathbb{I}_n]>\sqrt{2.01 \log n}) \to 0,\ as\ n \to \inftyP(Z∗[In]>2.01logn)→0, as n→∞

这个不等式看上去让人摸不着头脑，但我们介绍完这个不等式的推导后会介绍它到底有什么用。

推导这个概率不等式的思路与推导随机矩阵的概率不等式的思路类似，我们用ϵ\epsilonϵ-net作为工具，将和式放大为ϵ\epsilonϵ-net的绝对值乘以概率的最值即可。构造ϵ\epsilonϵ-net需要度量空间，因此我们首先要在In\mathbb{I}_nIn上定义度量，然后再推导不等式。

结论3 (In,δ)(\mathbb{I}_n,\delta)(In,δ)是一个度量空间。

证明
i) δ(I,J)=0⇔ρ(I,J)=1⇔I=J\delta(I,J)=0 \Leftrightarrow \rho(I,J)=1 \Leftrightarrow I = Jδ(I,J)=0⇔ρ(I,J)=1⇔I=J，后半句的证明如下：
⇒\Rightarrow⇒: If ρ(I,J)=1\rho(I,J)=1ρ(I,J)=1, assume I≠JI \ne JI=J. We can decompose JJJ as
J=(I∩J)⊔(J∖I)J = (I \cap J) \sqcup (J \setminus I)J=(I∩J)⊔(J∖I)where ∣J∖I∣>0|J \setminus I|>0∣J∖I∣>0. Besides, I∩J⊂II \cap J \subset II∩J⊂I, ∣I∩J∣≤I|I \cap J| \le I∣I∩J∣≤I. So
ρ(I,J)=∣I∩J∣∣I∣∣J∣≤∣I∩J∣∣I∩J∣(∣I∩J∣+∣J∖I∣)<1\rho(I,J)=\frac{|I \cap J|}{\sqrt{|I||J|}}\le \frac{|I \cap J|}{\sqrt{|I \cap J|(|I \cap J|+|J \setminus I|)}}<1ρ(I,J)=∣I∣∣J∣∣I∩J∣≤∣I∩J∣(∣I∩J∣+∣J∖I∣)∣I∩J∣<1

This is contradict to ρ(I,J)=1\rho(I,J)=1ρ(I,J)=1. So I=JI = JI=J.

⇐\Leftarrow⇐: If I=JI = JI=J, ∣I∩J∣=∣I∣|I \cap J| = |I|∣I∩J∣=∣I∣, ∣I∣∣J∣=∣I∣2=∣I∣\sqrt{|I||J|} = \sqrt{|I|^2} = |I|∣I∣∣J∣=∣I∣2=∣I∣, thus
ρ(I,J)=1\rho(I,J)=1ρ(I,J)=1

ii) δ(I,J)=1−ρ(I,J)=1−ρ(J,I)=δ(J,I)\delta(I,J)=\sqrt{1-\rho(I,J)}=\sqrt{1-\rho(J,I)}=\delta(J,I)δ(I,J)=1−ρ(I,J)=1−ρ(J,I)=δ(J,I)

iii) ∀I,J,K∈In\forall I,J,K \in \mathbb{I}_n∀I,J,K∈In，要说明
δ(I,K)≤δ(I,J)+δ(J,K)\delta(I,K) \le \delta(I,J)+\delta(J,K)δ(I,K)≤δ(I,J)+δ(J,K)

参考结论2，用以区间为下标的随机变量之差的方差表示δ\deltaδ，则不等式等价于
12Var(Z[I]−Z[K])≤12Var(Z[I]−Z[J])+12Var(Z[J]−Z[K])\sqrt{\frac{1}{2}Var(Z[I]-Z[K])} \\ \le \sqrt{\frac{1}{2}Var(Z[I]-Z[J])}+\sqrt{\frac{1}{2}Var(Z[J]-Z[K])}21Var(Z[I]−Z[K])≤21Var(Z[I]−Z[J])+21Var(Z[J]−Z[K])

引入X=Z[I]−Z[J],Y=Z[J]−Z[K]X=Z[I]-Z[J],Y=Z[J]-Z[K]X=Z[I]−Z[J],Y=Z[J]−Z[K]，
EX=EY=0EX=EY=0EX=EY=0

则上式等价于
E[(X+Y)2]≤E[X2]+E[Y2]\sqrt{E[(X+Y)^2]} \le \sqrt{E[X^2]}+\sqrt{E[Y^2]}E[(X+Y)2]≤E[X2]+E[Y2]

也就是等价于L2L^2L2中的范数满足的三角不等式：
∥X+Y∥L2≤∥X∥L2+∥Y∥L2\left\| X+Y \right\|_{L^2} \le \left\| X \right\|_{L^2}+\left\| Y \right\|_{L^2}∥X+Y∥L2≤∥X∥L2+∥Y∥L2

所以
δ(I,K)≤δ(I,J)+δ(J,K)\delta(I,K) \le \delta(I,J)+\delta(J,K)δ(I,K)≤δ(I,J)+δ(J,K)

引理极值理论的概率不等式
假设X1,⋯,Xn∼iidN(0,1)X_1,\cdots,X_n \sim_{iid} N(0,1)X1,⋯,Xn∼iidN(0,1)，则
P(X(n)>2log⁡n)≤14πlog⁡nP(X_{(n)}>\sqrt{2 \log n}) \le \frac{1}{\sqrt{4\pi \log n}}P(X(n)>2logn)≤4πlogn1

这是极值理论的一个重要结论，极值理论是研究一列随机变量的最大值/最小值的分布的理论，最大值X(n)X_{(n)}X(n)可能收敛到一些特殊的分布，比如Gumbel分布等，这个不等式说明这种收敛的收敛速率为O(1/log⁡n)O(1/\sqrt{\log n})O(1/logn)，这是一个非常小的收敛速率，也就是说极值的收敛是非常慢的。另外一类重要的收敛是均值的收敛，研究均值收敛的理论是中心极限定理，经典的中心极限定理给出的均值的收敛速率是O(1/n)O(1/n)O(1/n)。

现在我们可以在(In,δ)(\mathbb{I}_n,\delta)(In,δ)上定义一个ϵ\epsilonϵ-net，记为N\mathcal{N}N，它的cardinality满足
∣N∣≤4nϵ4|\mathcal{N}|\le \frac{4n}{\epsilon^4}∣N∣≤ϵ44n

N\mathcal{N}N上的元素满足：∀I∈In\forall I \in \mathbb{I}_n∀I∈In，∃J∈N\exists J \in \mathcal{N}∃J∈N，δ(I,J)<ϵ\delta(I,J)<\epsilonδ(I,J)<ϵ，根据结论2，用方差表示就是
Var(Z[I]−Z[J])≤2ϵ2Var(Z[I]-Z[J]) \le 2\epsilon^2Var(Z[I]−Z[J])≤2ϵ2

引入Z∗[N]=max⁡{Z[I]:I∈N}Z^*[\mathcal{N}]=\max\{Z[I]:I \in \mathcal{N}\}Z∗[N]=max{Z[I]:I∈N}，根据上面的引理
P(Z∗[N]>2log⁡∣N∣)≤14πlog⁡∣N∣P(Z∗[N]>2log⁡4nϵ4)≤14πlog⁡4nϵ4P(Z^*[\mathcal{N}]>\sqrt{2 \log |\mathcal{N}|}) \le \frac{1}{\sqrt{4\pi \log |\mathcal{N}|}} \\ P(Z^*[\mathcal{N}]>\sqrt{2 \log \frac{4n}{\epsilon^4}}) \le \frac{1}{\sqrt{4\pi \log \frac{4n}{\epsilon^4}}} P(Z∗[N]>2log∣N∣)≤4πlog∣N∣1P(Z∗[N]>2logϵ44n)≤4πlogϵ44n1

这是关于Z∗[N]Z^*[\mathcal{N}]Z∗[N]，我们来说明它可以逼近Z∗[In]Z^*[\mathbb{I}_n]Z∗[In]，
Z∗[In]=max⁡I∈InZ[I]=max⁡I∈In(Z[J]+(Z[I]−Z[J]))≤Z∗[N]+max⁡I∈In,J∈N(Z[I]−Z[J])Z^*[\mathbb{I}_n] = \max_{I \in \mathbb{I}_n} Z[I] = \max_{I \in \mathbb{I}_n} (Z[J]+(Z[I]-Z[J])) \\ \le Z^{*}[\mathcal{N}]+\max_{I \in \mathbb{I}_n,J \in \mathcal{N}}(Z[I]-Z[J])Z∗[In]=I∈InmaxZ[I]=I∈Inmax(Z[J]+(Z[I]−Z[J]))≤Z∗[N]+I∈In,J∈Nmax(Z[I]−Z[J])

max⁡I∈In,J∈N(Z[I]−Z[J])\max_{I \in \mathbb{I}_n,J \in \mathcal{N}}(Z[I]-Z[J])maxI∈In,J∈N(Z[I]−Z[J])可以理解成∣In∣∣N∣|\mathbb{I}_n||\mathcal{N}|∣In∣∣N∣个标准正态的最大值，根据引理，
P(max⁡I∈In,J∈N(Z[I]−Z[J])>2log⁡∣In∣∣N∣)≤14πlog⁡∣In∣∣N∣P(max⁡I∈In,J∈N(Z[I]−Z[J])>2ϵ2log⁡n4)≤14πlog⁡n4P(\max_{I \in \mathbb{I}_n,J \in \mathcal{N}}(Z[I]-Z[J])>\sqrt{2 \log |\mathbb{I}_n||\mathcal{N}|}) \le \frac{1}{\sqrt{4\pi \log |\mathbb{I}_n||\mathcal{N}|}} \\ P(\max_{I \in \mathbb{I}_n,J \in \mathcal{N}}(Z[I]-Z[J])>\sqrt{2}\epsilon\sqrt{2 \log n^4}) \le \frac{1}{\sqrt{4\pi \log n^4}} P(I∈In,J∈Nmax(Z[I]−Z[J])>2log∣In∣∣N∣)≤4πlog∣In∣∣N∣1P(I∈In,J∈Nmax(Z[I]−Z[J])>2ϵ2logn4)≤4πlogn41

综上，
P(Z∗[In]>2log⁡4nϵ4+2ϵ2log⁡n4)<14πlog⁡n4+14πlog⁡4nϵ4→0,asn→∞P(Z^*[\mathbb{I}_n]>\sqrt{2 \log \frac{4n}{\epsilon^4}}+\sqrt{2}\epsilon\sqrt{2 \log n^4})\\<\frac{1}{\sqrt{4\pi \log n^4}} + \frac{1}{\sqrt{4\pi \log \frac{4n}{\epsilon^4}}} \to 0,\ as\ n \to \inftyP(Z∗[In]>2logϵ44n+2ϵ2logn4)<4πlogn41+4πlogϵ44n1→0, as n→∞

最后考虑一下
2log⁡4nϵ4+2ϵ2log⁡n4\sqrt{2 \log \frac{4n}{\epsilon^4}}+\sqrt{2}\epsilon\sqrt{2 \log n^4}2logϵ44n+2ϵ2logn4

我们希望得到一个最大的下界，于是选择一个ϵ\epsilonϵ使得这个值最小，比如取ϵ=0.001/2\epsilon=0.001/\sqrt{2}ϵ=0.001/2，就可以得到待证的不等式了。

应用：DNA序列突变点侦测

先从维基百科扒来一个标准密码子表，简单介绍一下DNA序列突变点的含义，如果一段DNA序列为ATGATCACTATGTAG，则它翻译的结果是甲硫氨酸-异亮氨酸-苏氨酸-甲硫氨酸，如果这个序列突变为ATGATCACTATCTAG，它翻译的结果就变成了甲硫氨酸-异亮氨酸-苏氨酸-异亮氨酸，这两种结果中硫元素的含量具有显著差异，根据这个性质，我们可以设计试验寻找突变点。要检验某段DNA中是否有使ATG变成ATC的突变点，翻译得到一个很长的肽链，用XiX_iXi表示一小段肽链的硫元素含量，Xi∼N(μi2,1),i=1,⋯,nX_i \sim N(\mu_i^2,1),i=1,\cdots,nXi∼N(μi2,1),i=1,⋯,n，这里方差被标准化为1，我们可以用假设检验建模：
H0:μ1=⋯=μn=0Ha:μ1=⋯=μl=μr+1=⋯=μn=0,μl+1=⋯=μr≠0H_0:\mu_1=\cdots=\mu_n=0 \\ H_a:\mu_1 = \cdots = \mu_l = \mu_{r+1} = \cdots = \mu_n = 0,\mu_{l+1}=\cdots=\mu_r \ne 0H0:μ1=⋯=μn=0Ha:μ1=⋯=μl=μr+1=⋯=μn=0,μl+1=⋯=μr=0

定义
X[l,r]=1l−r∑i=l+1rXiX[l,r]=\frac{1}{\sqrt{l-r}}\sum_{i=l+1}^rX_iX[l,r]=l−r1i=l+1∑rXi

在原假设下，X[l,r]∼N(0,1)X[l,r] \sim N(0,1)X[l,r]∼N(0,1)，要检验这段DNA中是否包含使ATG变成ATC的突变点，引入统计量
T=max⁡l<r∣X[l,r]∣T = \max_{l<r}|X[l,r]|T=l<rmax∣X[l,r]∣

这个统计量的作用非常直观，如果这段DNA不含突变，那么TTT的值就会比较小，否则TTT是比较大的，为了让这个方法具有实用价值，我们需要确定多大才算大。根据上面的不等式：
P(T=Z∗[In]>2.01log⁡n)→0,asn→∞P(T=Z^*[\mathbb{I}_n]>\sqrt{2.01 \log n}) \to 0,\ as\ n \to \inftyP(T=Z∗[In]>2.01logn)→0, as n→∞

这里的nnn表示的是肽链的数目，这个概率不等式说明，用一个比较严格的标准，当TTT比2.01log⁡n\sqrt{2.01 \log n}2.01logn时，我们可以显著地拒绝原假设，并得出这段DNA存在使ATG变成ATC的突变点。