DNA序列.docx
《DNA序列.docx》由会员分享，可在线阅读，更多相关《DNA序列.docx（40页珍藏版）》请在冰豆网上搜索。
DNA序列
DNA序列分类
摘要
2000年6月，人类基因组计划中DNA全序列草图完成，DNA全序列结构的研究成为生物信息学的一个重要课题，而DNA序列分类是研究DNA全序列结构的基础。
本文采用Fisher判别对题目给出的DNA序列进行分类。
问题一是在20个已知类别的人工制造DNA序列的条件下，对20个未标明类别的人工DNA序列进行分类。
我们首先根据题中所给的20个已知类别的人工制造的序列（其中序列标号1—10为A类，11-20为B类），可算出A、C、G、T在序列中出现的频率，然后建立三种分类模型，分别为K-均值聚类，系统聚类和Fisher判别模型，并分析比较模型分类的正确性和稳定性，用于确定最优的分类模型，即Fisher判别模型。
最后用该方法把20个未标明类别的人工序列进行分类。
利用SPSS和EXCEL软件得出以下结果：
22、23、25、27、29、34、35、36、37为A类，20、21、24、26、28、30、31、32、37、38、39、40为B类。
对于问题二是对题目给出的自然序列，利用问题一中的分类方法进行分类，它是问题一的推广。
给出的自然序列的长度发生了变化，与问题一类似，现在本文用问题一中选出的正确率比较高的分类方法即Fisher判别法对182个自然DNA序列进行分类即可。
运用SPSS和VC++，得到下面的结果：
3、20、45、70、101、136、4、21、47、71、104、139、5、25、49、73、105、141、6、27、52、77、106、142、8、31、53、79、109、145、9、32、55、82、112、147、10、33、58、89、113、148、13、359、90、115、149、14、36、60、91、117、154、15、38、61、93、118、155、16、39、62、97、120、158、17、41、64、98、124、171、18、42、67、99、132、172、19、44、69、100、134、176为A类，其余为B类。
关键词：
频率，K-均值聚类模型，系统聚类模型，Fisher判别模型，正确率
一、问题重述
人类基因组计划中DNA全序列草图是由4个字符A，T，C，G按一定顺序排成的长约30亿的序列，其中没有“断句”也没有标点符号。
虽然人类对它知之甚少，但也发现了其中的一些规律性和结构。
例如，在全序列中有一些是用于编码蛋白质的序列片段，即由这4个字符组成的64种不同的3字符串，其中大多数用于编码构成蛋白质的20种氨基酸。
又例如，在不用于编码蛋白质的序列片段中，A和T的含量特别多些，于是以某些碱基特别丰富作为特征去研究DNA序列的结构也取得了一些结果。
此外，利用统计的方法还发现序列的某些片段之间具有相关性，等等。
这些发现让人们相信，DNA序列中存在着局部的和全局性的结构，充分发掘序列的结构对理解DNA全序列是十分有意义的。
目前在这项研究中最普通的思想是省略序列的某些细节，突出特征，然后将其表示成适当的数学对象。
作为研究DNA序列的结构的尝试，提出以下对序列集合进行分类的问题：
1）请从20个已知类别的人工制造的序列（其中序列标号1—10为A类，11-20为B类）中提取特征，构造分类方法，并用这些已知类别的序列，衡量你的方法是否足够好。
然后用你认为满意的方法，对另外20个未标明类别的人工序列（标号21—40）进行分类，把结果用序号（按从小到大的顺序）标明它们的类别（无法分类的不写入）
2）同样方法对182个自然DNA序列（它们都较长）进行分类，像1）一样地给出分类结果。
二、问题分析
对于问题一，属于分类问题，为了解决这一问题，本文首先从已知类别的人工制造序列来提取特征（本文是根据A、C、G、T在序列中出现的频率来提取特征的），然后学习K-均值聚类，系统聚类和Fisher判别这三种模型，并用SPSS、EXCEL和VC++软件实现分类，最后根据已知类别的人工制造序列来检验这三种分类方法的正确性，通过比较三种模型的正确率与稳定率的值，从中选出正确率与稳定率较高的一种分类模型，即Fisher判别模型，并用该方法把20个未标明类别的人工序列进行分类即可。
对于问题二，序列的长度发生了变化，与问题一类似，本文用问题一中选出的正确率与稳定率较高的分类模型对182个自然DNA序列进行分类即可。
三、模型假设
3.1模型假设
1.较长的182个自然序列与已知类别的20个样本序列具有共同的特征。
2.认为所讨论的序列都是从DNA序列中任意截取的一部分.
3.所研究的这些DNA序列都是稳定的,至于在极少数情况下发生的DNA变异不予考虑
四、与符号说明
:
表示题目给出的DNA的序列号
：
分别表示第
个DNA序列中的A、C、G、T碱基的出现频率
:
表示A碱基在序列中出现的频率
:
表示C碱基在序列中出现的频率
:
表示G碱基在序列中出现的频率
:
表示T碱基在序列中出现的频率
：
表示A类中A碱基的初始聚类中心
：
表示A类中C碱基的初始聚类中心
：
表示A类中G碱基的初始聚类中心
：
表示A类中T碱基的初始聚类中心
：
表示B类中A碱基的初始聚类中心
：
表示B类中C碱基的初始聚类中心
：
表示B类中G碱基的初始聚类中心
：
表示B类中T碱基的初始聚类中心
：
表示Fisher判别法中的线性判别函数
：
表示欧式距离判别法中的线性判别系数
、
、
、
：
表示Fisher判别函数的系数
：
表示Fisher判别法的判别临界值
：
表示已分类的20个序列间的距离
：
表示A,B两个聚类间的最短距离
：
表示A类中的样本xu和B类中的样本xv之间的距离
：
表示A类总体的样本离差阵
：
表示B类总体的样本离差阵
五、模型的建立与求解
5.1问题一的模型建立于求解
观察并分析数据预处理结果，归纳总结出A类和B类的有效特征，将其表示成适当的数学对象，并选择适当的分类方法，建立普遍意义下聚类数学模型，再用得到的模型对其余未知类别的20个人工制造的DNA序列进行分类，并运用SPSS和EXCEL进行求解。
由题意，建立的数学模型应该保证分类结果具有以下特点：
（1）类别间差异尽量大；
（2）类别内差异尽量小；
（3）样品能够尽可能的落入A、B范围，且只能落入其中的一个。
5.1.1问题一的模型准备—提取特征
DNA序列分类是一个复杂的统计分析问题，数据量大，影响因素多，无法直接从20条已知类别的人工制造的DNA序列中提取出所有的有效特征，因此有必要对这20条DNA序列进行预处理。
因为不同段的DNA序列中，碱基总数不同且差异大，所以不可根据碱基总数判别，但是易发现不同碱基出现的频率不同，且易分类，故本文是根据A、C、G、T在序列中出现的频率来提取特征的。
根据题中所给的20个已知类别的人工制造的序列（其中序列标号1—10为A类，11-20为B类），可算出A、C、G、T在序列中出现的频率，利用word和excel软件，可以得出下面的结果，如表4.1.1：
表4.1.1.20个样本中A、C、G、T碱基出现的次数和频率
类别
序列号
A
C
G
T
AF
CF
GF
TF
A类
1
33
19
44
15
0.297297
0.171171
0.396396
0.135135
2
30
18
46
17
0.27027
0.162162
0.414414
0.153153
3
30
24
50
7
0.27027
0.216216
0.45045
0.063063
4
47
12
20
32
0.423423
0.108108
0.18018
0.288288
5
26
26
47
12
0.234234
0.234234
0.423423
0.108108
6
39
14
44
14
0.351351
0.126126
0.396396
0.126126
7
39
11
40
21
0.351351
0.099099
0.36036
0.189189
8
31
18
41
21
0.279279
0.162162
0.369369
0.189189
9
23
23
48
17
0.207207
0.207207
0.432432
0.153153
10
20
30
45
15
0.181818
0.272727
0.409091
0.136364
B类
11
39
5
11
55
0.354545
0.045455
0.1
0.5
12
36
3
16
55
0.327273
0.027273
0.145455
0.5
13
28
11
14
57
0.254545
0.1
0.127273
0.518182
14
33
9
13
55
0.3
0.081818
0.118182
0.5
15
32
0
7
71
0.290909
0
0.063636
0.645455
16
40
9
10
51
0.363636
0.081818
0.090909
0.463636
17
39
27
15
29
0.354545
0.245455
0.136364
0.263636
18
32
13
10
55
0.290909
0.118182
0.090909
0.5
19
24
16
8
62
0.218182
0.145455
0.072727
0.563636
20
22
19
7
62
0.2
0.172727
0.063636
0.563636
5.1.2判别函数
5.1.2.1K-均值聚类模型
K均值聚类属于划分聚类算法，十分简洁和效率。
现已知20个样本中A、C、G、T碱基出现的频率的数据集合和确定的聚类数目k=2，根据距离函数不断迭代将数据分入A、B俩个聚类中。
首先从所给样本数据对象中随机选取8个对象作为初始聚类中心点，然后对于所剩下的其它对象，则根据它们与所选8个中心点的相似度（距离）分别分配给与其最相似的聚类，然后在重新计算所获聚类的聚类中心（该聚类中所有对象的均值），不断重复这一过程直到标准测度函数开始收敛为止，其基本流程如下：
1.从20个样本中选择8个对象作为初始聚类中心，即
，
，
，
，
，
，
，
。
2.根据每个聚类对象的均值（中心对象），计算每个对象与这些中心对象的距离，距离公式如下：
（1）
并根据最小距离对相应对象进行划分。
3.将所有样本分配到相应的聚类后，重新计算聚类中心，得到最终聚类中心，如下图5.1.2所示。
表5.1.2最终聚类中心
聚类
A
B
a
32
32
g
40
11
c
20
9
t
18
58
表5.1.3迭代历史记录a
迭代
聚类中心内的更改
A
B
1
16.930
16.396
2
.000
.000
a.由于聚类中心内没有改动或改动较小而达到收敛。
任何中心的最大绝对坐标更改为.000。
当前迭代为2。
初始中心间的最小距离为74.993。
4.循环执行第2和3步，直到每个聚类不再发生变化或者标准测度函数开始收敛为止。
运用SPSS软件分析得出结果，如下表5.1.4所示：
表5.1.420个样本DNA序列分类结果
DNA序列
判别类型
原类型
距离
1
A
A
5.274
2
A
A
6.941
3
A
A
15.673
4
A
A
29.487
5
A
A
12.756
6
A
A
10.703
7
A
A
11.623
8
A
A
3.977
9
A
A
12.756
10
A
A
16.930
11
B
B
9.039
12
B
B
9.873
13
B
B
5.389
14
B
B
4.101
15
B
B
16.396
16
B
B
10.900
17
A
B
28.833
18
B
B
4.776
19
B
B
11.212
20
B
B
14.679
表5.1.5每个聚类中的案例数
聚类
A
11.000
B
9.000
有效
20.000
缺失
.000
上述判别结果表明，原B类总体中只有第17个序列判归类别为A，与原类别不同，其余序列与原类别相同，可以猜测第17个序列有可能是属于原类别时的错分序列。
总的正确率达95%。
5.1.2.2系统聚类模型
系统聚类首先根据一批数据或指标找出能度量这些数据或指标之间相似程度的统计量；然后以统计量作为划分类型的依据，把一些相似程度大的变量首先聚合为一类，而把另一些相似程度较小的变量聚合为另一类，直到所有的变量都聚合完毕，最后根据各类之间的亲疏关系，逐步画成一张完整的分类系统图。
其相似程度由距离或者相似系数定义。
进行类别合并的准则是使得类间差异最大，而类内差异最小。
其基本流程如下：
第一步：
设初始DNA序列样本共有n=20个，每个样本自成一类，即建立20类，
。
计算各类之间的距离（初始时即为各样本间的距离），得到一个N*N维的距离矩阵D（0）。
这里，标号（0）表示聚类开始运算前的状态。
第二步：
假设前一步聚类运算中已求得距离矩阵D（n），n为逐次聚类合并的次数，则求D（n）中的最小元素。
如果它是Gi（n）和Gj（n）两类之间的距离，则将Gi（n）和Gj（n）两类合并为一类
，由此建立新的分类：
。
第三步：
计算合并后新类别之间的距离，得D（n+1）。
计算
与其它没有发生合并的
之间的距离，可采用多种不同的距离计算准则进行计算，公式如下：
（2）
第四步：
返回第二步，重复计算及合并，直到得到满意的分类结果。
过程中，我们选择最短距离法：
设A和B是两个聚类，则两类间的最短距离定义为：
（3）
其中，du,v表示A类中的样本xu和B类中的样本xv之间的距离，DA,B表示A类中的所有样本和B类中的所有样本之间的最小距离。
具体流程图如下图5.1.2.1所示：
图5.1.2.1迭代次数为n的系统聚类算法流程图
否
运用SPSS软件分析得出结果，如下表5.1.4所示：
聚类表
阶
群集组合
系数
首次出现阶群集
下一阶
群集1
群集2
群集1
群集2
1
19
20
.001
0
0
16
2
1
2
.001
0
0
7
3
14
18
.002
0
0
5
4
11
16
.003
0
0
9
5
13
14
.003
0
3
11
6
5
9
.004
0
0
8
7
1
8
.004
2
0
12
8
5
10
.005
6
0
13
9
11
12
.006
4
0
11
10
6
7
.006
0
0
12
11
11
13
.010
9
5
16
12
1
6
.010
7
10
14
13
3
5
.012
0
8
14
14
1
3
.025
12
13
19
15
4
17
.026
0
0
18
16
11
19
.028
11
1
17
17
11
15
.035
16
0
18
18
4
11
.101
15
17
19
19
1
4
.243
14
18
0
图5.1.2.1系统分析冰柱图
由上图划线分成两类，可以分析得出1、2、3、5、6、7、8、9、10分为一类，4、11、12、13、14、15、16、17、18、19、20分为另一类，具体如下表5.1.2.1所示。
表5.1.2.120个样本DNA序列分类结果
DNA序列
判别类型
原类型
1
A
A
2
A
A
3
A
A
4
B
A
5
A
A
6
A
A
7
A
A
8
A
A
9
A
A
10
A
A
11
B
B
12
B
B
13
B
B
14
B
B
15
B
B
16
B
B
17
B
B
18
B
B
19
B
B
20
B
B
上述判别结果表明，原A类总体中只有第4个序列判归类别为B，与原类别不同，其余序列与原类别相同，可以猜测第4个序列有可能是属于原类别时的错分序列。
总的正确率达95%。
5.1.2.3Fisher判别模型
利用前面计算的结果，可按下列步骤进行计算：
1）确定
、
、
、
及判别函数
由
可求得Fisher判别函数的系数
、
、
、
，利用SPSS软件，可得到判别函数如下：
（4）
2）确定判别临界值y0
由
计算判别临界值
其中，
（5）
（6）
3）确定判别准则
因
，故判别准则为：
4）对已知类别的序列判别分类
对20个已知类别的人工制造的序列用线性判别函数进行判别归类，
为A类，
为B类。
具体过程结果如下表5.1.2.2所示。
表5.1.2.2按照案例顺序的统计量
案例数目
最高组
第二最高组
判别式得分
P（D>d|G=g）
实际组
预测组
p
df
P（G=g|D=d）
到质心的平方Mahalanobis距离
组
P（G=g|D=d）
到质心的平方Mahalanobis距离
函数1
初始
1
1
1
.743
1
1.000
.108
2
.000
32.827
3.029
2
1
1
.656
1
1.000
.198
2
.000
34.180
3.146
3
1
1
.229
1
1.000
1.446
2
.000
43.613
3.903
4
1
2**
.011
1
.716
6.400
1
.284
8.246
-.171
5
1
1
.669
1
1.000
.182
2
.000
33.971
3.128
6
1
1
.486
1
1.000
.484
2
.000
37.179
3.397
7
1
1
.980
1
1.000
.001
2
.000
28.899
2.675
8
1
1
.733
1
1.000
.116
2
.000
25.609
2.360
9
1
1
.727
1
1.000
.122
2
.000
33.062
3.049
10
1
1
.834
1
1.000
.044
2
.000
26.958
2.491
11
2
2
.707
1
1.000
.141
1
.000
25.260
-2.325
12
2
2
.317
1
1.000
1.001
1
.000
19.370
-1.700
13
2
2
.861
1
1.000
.031
1
.000
27.311
-2.525
14
2
2
.735
1
1.000
.114
1
.000
25.638
-2.363
15
2
2
.378
1
1.000
.776
1
.000
39.465
-3.581
16
2
2
.747
1
1.000
.104
1
.000
25.792
-2.378
17
2
2
.182
1
.999
1.779
1
.001
16.547
-1.367
18
2
2
.838
1
1.000
.042
1
.000
31.426
-2.905
19
2
2
.276
1
1.000
1.188
1
.000
42.137
-3.791
20
2
2
.171
1
1.000
1.878
1
.000
45.858
-4.071
交叉验证a
1
1
1
.983
3
1.000
.165
2
.000
31.249
2
1
1
.885
3
1.000
.649
2
.000
33.124
3
1
1
.562
3
1.000
2.050
2
.000
45.544
4
1
2**
.006
3
1.000
12.392
1
.000
75.334
5
1
1
.736
3
1.000
1.271
2
.000
33.566
6
1
1
.530
3
1.000
2.212
2
.000
38.199
7
1
1
.560
3
1.000
2.059
2
.000
28.933
8
1
1
.973
3
1.000
.229
2
.000
24.420
9
1
1
.541
3
1.000
2.157
2
.000
33.494
10
1
1
.187
3
1.000
4.799
2
.000
28.857
11
2
2
.655
3
1.000
1.619
1
.000
24.967
12
2
2
.230
3
1.000
4.314
1
.000
20.603
13
2
2
.850
3
1.000
.797
1
.000
26.389
14
2
2
.951
3
1.000
.345
1
.000
24.520
15
2
2
.118
3
1.000
5.876
1
.000
45.825
16
2
2
.661
3
1.000
1.592
1
.000
25.465
17
2
1**
.000
3
1.000
28.595
2
.000
46.523
18
2
2
.982
3
1.000
.170
1
.000
29.869
19
2
2
.279
3
1.000
3.839
1
.000
46.474
20
2
2
.058
3
1.000
7.475
1
.000
57.584
升级会员 