DNA序列分类模型

1、 图一基于持怔密刊子分布概率的判别港则对我个人工D岫序审判别从图一中可清楚的看出，在由P(H)、P(R)构成的二维空间内，A类和B类呈现出非常好的聚类倾向。理想的情况是能够找到一条直线将A类和B类的点完全分开，那么自然的，与已知类同侧的点就属于该类。从图上可以看到这条直线位于过原点的角平分线附近，但是A类和B类都分别有一个点位于这条直线附近。所以为了准确起见，我们分别过原点和这两个点做了两条直线作为分界线，将二维空间分成了三个域，待判点位于右下角则属A类；位于左上域则属于B类；如果待判点位于两条分界线之间，则进入了本方案的辨析盲区，方案失效。将待判点描到图上，无须计算就可以将它们归类：人工基因

2、序列的分类结果：A类：34352327293722B类：31213240332426283830不能判别的点：392536这个分类方案简单，计算量小。如果将坐标和界线固定的画成一张坐标纸，就可以象正态概率纸一样方便的对DNA序列进行分类了。但它最主要的缺点就是存在判断盲区，而且，当我们把它应用在182个自然DNA序列的分类时(如图二)，聚类就不那么明显了。所以该模型只能判断序列移位的人工模型，对自然DNA序列判断结果不理想。对自然DNA序列分类的结果：A类：点数143个B类：点数24个不能判别的点数：15个P（R）图二基于特怔密刊子分布慨率的判别淮则判别L皿个自然DN&序模型二：图论最小生成树

3、聚类法从图二中可以看出待判点聚类不明显，为了使以上结果进行精确地分类，为了构造出更好的分类模型，判别率更高，我们以模型一的分类结果（182个点的分类结果）为基础，研究了图论最小生成树的聚类法。基于模型一中对182个点的分类结果，采用贪婪算法先求出这182个点的最小生成树（minimum-weightspanningtree），再求出最小生成树的“颈”按“颈”制定分类的方法和原则。对自然DNA序列分类的结果：A类：点数108个B类：点数74个判别的主要步骤是：1、在最小生成树上，找“主干”，最小生成树上边数最多的没有回路的链，称为“主干”。2、求“子主干”，在最小生成树的主干上，“度”大于2的一

4、点作为起点，肯定能够找到除主干本身外的没有回路的链，边数最多的就称为“子主干”而它的边数就是该点的“深度”。3、找“颈”，在主干上找出“深度”大于2的点，如果这些点之间存在“深度”为零的点，则这些“深度”为零的点之间的边就称为“颈”。4、根据最小生成树，去掉“颈”，就得到基因序列的分类。（如图三所示）通过以上步骤，找出了三个“颈”，分别是（72（3105）、（64（399）、（28（3159），自然分成四类，再将相邻两类合并，得到A、B两类的分类结果，该方法也非常直观，易懂。该方法的缺点是误判概率较难估计，因为不能判别的点数以在“颈”附近相邻点数记，相对来讲，应该很少。PiRi5四三最小生成柯

5、模型三多维向量空间的判别分析模型：PiRi5四三最小生成柯前面两个模型都有自身的局限性，判别分析就是要将与已知类别特征相近的物质归到一类，其关键就是要寻求待分样品和已知类别的相似特征。在这个模型中我们利用多维向量之间的夹角余弦差构造判别函数W(i)，对待判DNA序列片段进行归类。下面我们以二维向量空间为例介绍建模的原理：如图四所示在二维空间有两个类别已知的向量OA(a,a)和向量OB(b,b)，以及待1212判类别的向量OC(ci,c2)，0c与0A、OB之间的夹角分别为a和屋若B，则认为OC更加靠近0A，也就是在性质上OC与0A相近。若其中，向量OA、OB分别代表了A、B两类向量的重心位置，

6、则可以断定OC属于A类；反之，则属于B类。当a0图4二维空间向量夹角示意图当aP时，cosacosP图4二维空间向量夹角示意图根据平面几何的知识：ac+acbc+bccosa=,1122cosp=1122aa2+a2、c2+c2bb2+b2,cc2+c21212V1212其中a,b,c(i=1,2)分别为A、B、C各点的平面坐标值。iii令W=cosa-cosp则oC的判别式可写为:a2+a2工c2+c2当W0时，ap，OC属于A类;当Wp，OC属于B类。若有i个待归类的向量则W可以表示为W(i)，它是向量坐标的函数，称之判别函数。对于DNA序列片段的分类问题，虽然没有二维向量分类这么简单，但

7、如果我们能用向量表示每个片段的结构，并且找出A、B两类片段的重心向量，就可以用上面同样的方法对DNA片段进行归类了。我们知道，DNA的排序是千变万化的，但是反映遗传信息的是由4种碱基所组成的64种密码子的有效排列，据此可以用一个64维向量表示一个DNA片段，向量中的每个元素对应一个密码子在这个片段中出现的频率，第i个片段的向量表示为：P;p(k)k=1,2,364iiPi(k)表示第k种密码子在第i序列片段中出现的频率；那么所有属于A类和属于B类的DNA片段向量的集合就将整个64维空间分为两个互不重叠的域，我们姑且称之为A域和B域。若已知这两个域的重心向量，则可以分别求出待判别向量和A域重心向

8、量间的夹角a以及它和B域重心向量间的夹角P，夹角越小表明它离所对应的重心向量越近。如果ap则表明这个向量在空间上更加接近A域，也就是这个向量所代表的DNA片段在特征上更类似于A类DNA片段，则可以认为这个片段是属于A类的。反之，则是属于B类的。值得注意的是，在这里“夹角”已经失去了它在二维和三维空间里的具体含义，而仅仅是表示两个向量间关系的“关系系数”，因此用夹角余弦所表达的意义和用夹角表达的意义是完全一致的，夹角越小，夹角余弦值越大，二者的关系也就越密切。所以我们用关系函数W(j)=cosa-cosp来判断待判向量与哪一类序列的关系更加密切。现在我们就以给出的A、B两类片段来估计这两个域的中

9、心向量。对已给出的两类序列中的各个密码子的频率求平均,由这些平均频率所组成的向量就是这两个域的中心向量估计。两个域的重心向量勺和Pb的估计为:p(1),10p(2),.,10p(64)mmmm=1m=1mm=1p(1)napp(1)napn=11n=11p(64)n=11第j(j=21,22,23，40)个待判序列的向量和A域重心向量的夹角余弦:圮P(k)P(k)jAcosacosa圮P2(k)P2(k)jAk=1k=1第j(j=21,22,23，40)个待判序列的向量和B域重心向量的夹角余弦:圮p(k)P(k)jBcosp=1兄p2(k3jk=兄p2(k3jk=1k=1第j个待判序列的判别函

10、数为：圮p(k)P(k)PB2(k)W(j)=jAk=1PA2(k)k=1k=1圮P(k)P(k)jBk=1圮p2(k)P2(k)jBk=1k=1(j=21,22，40)时，aB,向量Prp(k)与PB的夹角较小，第j个DNA片段属于B类。四、结果分析模型一的求解及结果分析：问题(1)的分类结果：A类：22232729343537B类：21242628303132333840不能判断的点：25363917分辨率X=x100%=85%，由于没有发生错判的情况，所以误判率为零。4问题（2）的分类结果：（如图二）TOC o 1-5 h zA类：1234568910111314151617181920

11、212223252627293132333536373839404142444546474849515253545556575859606162646566676869707172737476777879818283848788899091939495969798991001011041051061081091111121131151171181201241251261291301321331341351361371381391401411421451461471481491501521531541551571581601641651671681691711721741751761771781

12、79180181B类：7124375858692102107110114116119121122123128151159161162163170182不能判断的点：24283034506380103127131143144156166173168入=x100%=92.3%182这个方案简单、灵活，使用方便，计算量小。但由于存在判断盲区，而使它分辨率降低。所以这种方案适合于对基因片段的粗略的估计。模型二的求解及结果分析：（1）自然基因中的分类结果为:TOC o 1-5 h zA类：1235689111314151617181920212527313233353738394041424445464

13、749515255565758596061626364666768697071727377787981828789909193949596100101104105106108109111112113115117118120124132133134135136137139141145147148150152153154155157158164165167171172176177178B类：471012222324262829303436434850535465747576808384858688929798991021031071101141161191211221231251261271281

14、29130131138140142143144146149151156159160161162163166168169170173174175179180181182仔细分析模型一和模型二的结果，发现在模型二被判别为B类的点的集合中都包括了模型一中被判别为B类的所有点。说明模型二是在模型一的基础上更进一步的精确。这种方法简单易懂，非常直观。模型三的求解及结果分析：在这个模型中，我们将64维空间分成了两个分别以户和户为重心的、互斥的域，BB从而对基因片段进行分类。利用Matlab软件，我们编制计算机程序20个人工合成基因片段以及182个自然DNA序列片段分别进行了判类计算。同时为了评估这种分类方

15、法的优劣，我们还对已知类别的20个基因片段进行了重新判类分析，结果如下：判别函数的值：W（i）=0.65260.63170.7481-0.09230.71980.61450.51940.53310.63130.6246-0.6162-0.5975-0.6681-0.6274-0.7024-0.5414-0.0072-0.6082-0.6681-0.6309根据判类准则：当W（i）0时，第i个序列片段属于A类;当W（i）0时，第i个序列片段属于B类已知人工序列片段的重新分类结果如下：A类：1,2,3,5,6,7,8,9,10；B类：4,11,12,13,14,15,16,17,18,19,20；

16、假设原来给定的分类是准确的话，显然新方案在第四段处发生了误判。由判别函数值我们也看到W(4)和W(17)的绝对值相对其他值而言要小的多，因此，虽然第十七片段并没有发生误判，但我们仍然认为是不可靠的。同时我们也必须看到，由于样本非常有限，每类给出的样本分别只有10个，从这么少的样本中抽取的共同特征是有限的。因此，这种分类方法的误判率=x100%=10%。和模型二相比较，这个误判率是可以接受的。而且我们相信，随着样本数的增加，重心向量的估计值会更靠近真实值，从而使误判率进一步降低，因此，多维空间的判别分析不失是一种对基因片段分类好的方法。采用与上述相同的计算程序我们对另外二十个不知类别的人工基因片

17、段进行了分类计算。为了更直观的反映判类函数在分类计算中的作用，我们在坐标纸上描出了每一片段的判别函数值，如图五所示，凡是处于零线以上的均属于A类，以下的则属于B类，零线两侧非常近的点代表的则是那些可能发生误判的片段。由图可知第21到第40人工片段2D22243口32配3日334口AHDWA片段序列号图5人工DNA判别函数的点阵图分类结果：A类：22,23,25,27,29,34,35,36,37,39；B类：21,24,26,28,30,31,32,33,38,40。对182种自然DNA序列片段的分类结果如下：A类共120个片段：12345689111314151617181920212527

18、293132333536373839414244454647495152535556575859606162636466676869707172737778798182838788899091939495969799100101104105106108109111112113115117118120124125126132133134135136137139141145147148149150152153154155157158160164165167169171172176177178B类共62个片段:7101222232426283034404348505465747576808485869298102103107110114116119121122123127128129130131138140142143144146151156159161162163166168170173174175179180181182鉴于篇幅有限，其判别函数值列于附录上，这里仅将数据点阵绘出，如图六所示。同样位于零线以下的点代表属于B类的片段，而位于零线以上的点则代表了属于A类的点。由图上看出这些点中绝大多离零线较远，只有少数几个点位于零线附近，是可能的误判点。我们取0.02作为判别可靠性的阈值，即当|