聚类分析方法

第一章 Microarray 介绍

1.1 生物信息处理

基于对生物体“硬件”和“软件”的认识 ,提出暂时地撇开生物的物理属性 ,

着重研究其信息属性 ,从而进入到生物信息处理 (关于生命硬件的信息和软件的信息 ,即生理信息和生命信息 )的一个分支 ,生物信息学。于是 ,为揭开生命之秘、揭示与生命现象相关的复杂系统的运作机制打开一条新的途径。

什么是生物信息处理

生物信息处理的英文是Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ。 1994年初 ,诺贝尔

医学奖获得者美国教授Ｍ·罗德贝尔发表一篇评论 ,题为《生物信息处理 :评估环境卫生的新方法》。他认为生物信息处理是在基因数据库基础上 ,计算机驱动的能快速获得表达基因部分序列的方法。通过ＭＥＤＬＩＮＥ数据库 ,可以查阅到很多与生物信息处理 (Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ)有关的记录，其中Ｊ Ｆ Ａｉｔｏｎ认为生物信息处理是基于计算机的数据库和信息服务;Ｒ Ｐ Ｍｕｒｒａｙ认为生物信息处理包括两方面：第一是大量现存数据的自动化处理 ,第二是新的信息资源的生成；Ｄ Ｂｅｎｔｏｎ在题为《生物信息处理———一个新的多学科工具的原理和潜力》的文章中说 ,生物信息处理的材料是生物学数据 ,其方法来自广泛的各种各样的计算机技术。 其方法来自广泛的各种各样的计算机技术。近年来 ,生物学数据在爆炸式增长 ,新的计算机方法在不断产生。这些方法在结构生物学、遗传学、结构化药品和分子演变学中是研究工作进展的基础。如果生物医学工程要在各个领域都从研究进展中获取最大好处 ,那么生物学数据健全的基础设施的开发与维护是同等重要的。尽管生物信息处理已经作出重要贡献 ,但是在它成熟时就会面临更大的需求在爆炸式增长 ,新的计算机方法在不断产生。这些方法在结构生物学、遗传学、结构化药品和分子演变学中是研究工作进展的基础。如果生物医学工程要在各个领域都从研究进展中获取最大好处 ,那么生物学数据健全的基础设施的开发与维护是同等重要的。尽管生物信息处理已经作出重要贡献 ,但是在它成熟时就会面临更大的需求。

在整体上可以看出 ,生物信息处理的两个基本内容是生物数据库建立和计算

机信息服务 ,也就是生物数据处理的计算机数据库化和程序化。当前这种数据库的内容主要是目录、期刊、遗传基因和细胞三维结构学。服务程序主要用于信息检索和基因序列分析。

所以 ,严格地说 ,当前生物信息处理远未形成独立的学科 ,它同计算机生物

学应用并无重大区别。在1998年第九届世界医药信息学大会上 ,它才作为一个讨论题目被列出来。可以说，生物信息处理技术是一项年轻的研究领域。

1．2 Microarray 技术

1．2．1 Microarray 技术原理

微阵列技术是利用分子杂交的原理，用自动化仪器arrayer把不同的，数以百

计、千计、万计已知部分序列的DNA探针“印”在玻璃片或者尼龙膜上面成阵列。为了比较两份标本中核酸表达的丰度，两份标本中核酸用同位素或者荧光素（红和绿两种）标记，再于微阵列杂交，然后检测杂交信号的强度，通过一定的数据处理系统，把它们转化成两份不同标本中特异基因的丰度，最后对这些数据进行分析。

根据微阵列技术原理，微阵列技术的处理流程如下：

1. 实验设计

2. 样品制备（指mRNA或总RNA样品，包括对照组和实验组）

3. 芯片制备（包括PCR，纯化，点样等步骤）

4. 芯片杂交（将mRNA或总RNA分别进行逆转录生成cDNA，在此步骤中

将对照组和实验组cDNA分别标记CY3和CY5荧光信号）

5. 芯片扫描（采用激光扫描仪，分别用532nm和635nm波长激光扫描芯片，

对于每张芯片，得到CY3和CY5通道两幅图象）

6. 图象处理（采用专门软件，对图象进行分析，提取每个点上的数字信号），

得到原始数据表。

7. 数据校正和筛选（对cy5或cy3信号进行校正，消除实验或扫描等各环节

因素对数据的影响，同时利用筛选规则对数据中的“坏点”，“小点”，“低信号

点”进行筛选，并作标记。）

8. 差异表达基因的确定（采用ratio值对差异基因进行判断，或采用统计方

法如线性回归、主成分分析、调整P值算法等对差异基因进行统计推断）

9. 生物信息学分析（如cluster 算法、差异基因的同源性比对，差异基因的

相关文献检索等）

一个最简单的配置应包括微阵列制作系统 (ａｒｒａｙｅｒ) ,信号收集系统 (ｓｃａｎｎｅｒ) ,计算机和软件 (操作系统和微阵列技术处理的相关软件 )。

1．2．2 Microarray 技术应用领域

Microarray 技术是近几年兴起的新技术，但短短几年中，该技术已经被分子

生物学的很多领域接受，并广泛应用于以下领域：

1、基因表达分析和检测

微阵列技术已经被许多研究小组应用于与基因表达有关的工作中,如对细

菌、动植物和人类的研究。包括：特异性相关的基因、差异表达的基因、基因

功能研究、健康状况的检测、毒理学研究、药物作用机制的研究、定位克隆。

2、功能分析

检测到基因表达差异之后 ,下一步是寻找这些差异的生物学功能。

最近Ｄａｖｉｓ等人[1 7]发明了一种新的方法。主要是应用插入一个独特序

列或标记的突变酵母链。分子标记在特殊的的生长条件下从生存链中扩增 ,

并与高密度微阵列进行杂交。这样不仅可以确定这条链的相对丰度 ,而且可以

在不同时间点反复进行 ,同时还可以精确比较每条缺失链的适应性。

3、基因作图

微阵列技术的应用补充了基因表达研究的方法 ,加强了对疾病易感性和疾

病本质的研究。这种方法无论是在速度上还是在准确性上都远胜于传统方法 ,

它将会改变基因制图的方法。

1．2．3 Microarray 技术发展现状

DNA微阵列技术（DNA microarray technology）是近几年发展起来的应用DNA

微阵列进行基因功能研究的新的生物技术。微阵列自1995年在《Ｓｃｉｅｎｃｅ》上报道后 ,被认为是该年度《Ｓｃｉｅｎｃｅ》上发表的最有影响的文章之一。微阵列是新出现的分子生物学技术 ,是本世纪重要的科学进展 ,它能够高效率、大规模地获取相关生物信息 ,是现代生物技术、微电子技术、机械制造技术、计算机技术的结合。其对科学的深远影响将远胜过ＤＮＡ测序和ＰＣＲ等,使人们更大规模

地获取生物信息 ,使人类基因组计划早日实现。

微阵列技术的迅速发展已经引起了各方面的广泛关注。许多实验室、专业公司

和制药公司都在大力开发与此相关的技术。在制作设备、分析设备、支持软件和探针的构建等方面均投入巨资 ,尤其是一些新兴的从事微阵列相关产业的公司如Ａｆｆｙｍｅｔｒｉｘ ,Ｉｎｃｙｔｅ,Ｓｙｎｔｅｎｉ,Ｃｌｏｎｔｅｃｈ等公司均已研制生产出相关的产品。有供诊断用的芯片如ＨＩＶ ,ｐ53和细胞色素 ｐ450的芯片 ;

有可供研究用的人、大鼠、小鼠不同基因类别的芯片 ;有与不同疾病如肿瘤、

心血管疾病、神经系统疾病相关的芯片也已投入使用。而且很多公司可根据需要定制各种微阵列系统 ,为研究人员提供方便。国内也开展了此项工作 ,清华大学、上海细胞生物所、军事医学科学院放射医学研究所及广州等地正在进行此项研究。微阵列技术的发展为探索生命科学提供了强有力的工具。使一些原本复杂的工作变得简捷。正如ＮＩＨ的主任ＨａｒｏｌｄＶａｒｍｕｓ在旧金山美国细胞生物学年会上指出的 :“应用微阵列技术 ,我们将最终揭示单个细胞的全部基因表达 ,甚至整个机体的基因概况”。同时 ,他还预言 :“微阵列技术将改变我们对生命本质的认识”。

1．3 本次毕业设计的目标

微阵列（Microarray ）技术是一门新兴学科，它是结合了生物学、计算机技

术、电子技术、生物信息学的特征而形成的一门交叉学科。

微阵列技术发展到现在，虽然已经取得了惊人的改变和进步，广泛应用于分子

生物学领域。但是，微阵列技术毕竟是一门新学科、一种新的思维方法，还需要在新的环境和领域下进行试验和完善，特别是在于其它学科、技术的结合方面，还需要研究人员花一定时间来研究和试验。

本次毕业设计要达到以下的目标：

（1） 学习生物信息学的背景知识和微阵列技术的处理流程；

（2） 学习聚类分析中的主要概念和技术方法，并阐述聚类分析在微阵列技术

中的重要地位；

（3） 分析几种常用的聚类方法，将其中存在的一些问题提炼出来加以分析；

（4） 在前几步的基础之上，结合所分析的常用算法设计一种改进算法；

（5） 用C语言实现设计的聚类算法和数据预处理工作。

第二章 聚类分析方法概述

2．1 聚类分析及相关概念

簇（Cluster）是指一个数据对象的集合。聚类就是将物理或抽象对象的集合

分组成为由类似的对象组成的多个簇的过程。由聚类所生成的簇是一组数据对象的集合，这些对象于同一个簇中的对象彼此相似，与其它簇中的对象相异。在许多应用中可以将一个簇中的对象作为一个整体来对待。

聚类是通过对数据对象本身数据的分析，从而将数据对象分成不同的类。聚类

是一种无监督分类法， 没有预先指定的类别。在机器学习领域，聚类是无指导学习（unsupervised clustering）。与分类不同，聚类和无指导学习不依赖预先定义的类和带类标号的训练实例。由于这个原因，聚类是观察式学习，而不是示例式学习。

聚类分析已经广泛的应用在许多领域中，包括模式识别、数据分析、图像处理、以及市场研究。通过聚类，人能够识别密集的和稀疏的区域，因而发现全局的分布模式，以及数据属性之间有趣的相互关系。

“聚类的典型应用是什么？”在商务上，聚类能帮助市场分析人员从客户基本

库中发现不同的客户群，并且用购买模式来刻画不同的客户群的特征。聚类在地球观测数据库中相似地区的确定和汽车保险单持有者的分组上也可以发挥作用。聚类也能用于对Web上的文挡进行分类，以发现信息。在生物学上，聚类能用于推导植物和动物的分类，对基因进行分类，获得数据分布的情况，观察每个簇的特点，集中对特定的某些簇作进一步的分析，获得对种群中固有结构的认识。此外，聚类分析可以作为其他算法（如特征和分类等）的预处理步骤，这些算法再在生成的簇上进行处理。

数据聚类正在蓬勃发展，有贡献的研究领域包括数据挖掘，统计学，机器学习，空间数据库技术，生物学，以及市场营销。由于数据库中收集了大量的数据，聚类分析已经成为生物学种生物信息分析研究领域中一个非常活跃的研究课题。

聚类是一个富有挑战性的研究领域，那么怎样才算是一个好的聚类方法?最重

要的是，一个好的聚类方法要能产生高质量的聚类结果——簇，这些簇要具备以下两个特点：

高的簇内相似性

低的簇间相似性

此外，在不同的领域有一些对聚类更深入的要求，例如在生物信息学中对聚类

算法的更深一步要求如下：

1．能应付脏数据。绝大多数现实世界中的数据库都包含了孤立点、空缺、未

知数据或者错误的数据。一些聚类算法对于这样的数据敏感，可能导致低质量

的聚类结果。

2．对于数据不同的顺序不敏感。一些聚类算法对于输入数据的顺序是敏感的。例如，同一个数据集合，当以不同的顺序提交给同一个算法时，可能生成差别

很大的聚类结果。开发对数据输入顺序不敏感的算法具有重要的意义。

3．模型可解释性和可使用性。用户希望聚类结果是可解释的、可理解的和可

用的。也就是说，聚类可能需要和特定的语义解释和应用相联系。应用目标如

何影响聚类方法的选择也是一个重要的研究课题。

2．2 聚类分析中的数据类型

在本节，我们研究在聚类分析中经常出现的数据类型，以及如何对其进行预处

理。许多基于内存的聚类算法选择如下两种有代表性的数据结构：

·数据矩阵（data matrix，或称为对象与变量结构）：它用p个变量（也称为度量或属性）来表现n个对象，例如用年龄、身高、体重、性别、种族等属性来表现对象“人”。这种数据结构是关系表的形式，或者看成n×p（n个对象×p个变量）

x11

...

xi1

...

xn1

的矩阵。图示如下： ...x1f.........xif.........xnf...x1p.........xip.........xnp·相异度矩阵（dissimilarity matrix，或称为对象－对象结构）：存储n个对象两两之间的近似性，表现形式是一个n×n维的矩阵。图示如下：

0d(2,1)

0d(3,1)d(3,2)0:::d(n,1)d(n,2)......0

在这里d(i ，j)是对象i和对象j 之间相异性的量化表示，通常它是一个非负的数值，当对象i和j越相似或“接近”，其值越接近0；两个对象越不同，其值越大。其中d(i,j)=d(j,i)，而且，d(i,i)=0。关于相异度，我们在下面的内容中进行详细探讨。

聚类分析中的数据类型

区间标度变量

区间标度变量是一个粗略现行标度的连续度量。典型的例子包括重量和高度、

经度和纬度坐标（如聚类房屋），以及大气温度。度量单位的选取对于聚类的结果的很重要的。例如将身高的单位从米变为尺，将体重的单位从公斤变为磅将对聚类的结果产生很大的影响。为了避免出现这种情况，我们必须将数据标准化，将数据的单位化成同样的权重，并将数据中的单位“去掉”。当没有关于数据的先验知识时，这样做是十分有用的。

因为在Microarray分析技术中，涉及到的数据类型是区间标度变量，所以，

下面就对区间标度变量作详细的介绍和讨论，其它的数据类型只做简单的介绍。

“怎样将一个变量的数据标准化？”为了实现度量值的标准化，一种方法是将

原来的度量值转换为无单位的值。给定一个变量f的度量值，可以进行如下的变换：

计算绝对偏差的平均值:

sf(|x1fmf||x2fmf|...|xnfmf|)

其中 mf (xx1f2f ...xnf).

xifmfzif f

计算标准度量值 (z-score)

使用绝对偏差的平均值比使用标准偏差更健壮（robust）。两个对象之间的相

异度（相似度）通常是基于对象间的距离来计算的。

最常用的距离度量方法是欧几里德距离，它的定义如下：

d(i,j)(|xx|2|xx|2...|xx|2)i1j1i2j2ipjp

其中 i = (xi1, xi2, „, xip) 和 j = (xj1, xj2, „, xjp) 是两个p维的数据对象

另一个著名的度量方法是曼哈坦距离，其定义如下：

d(i,j)|xx||xx|...|xx|i1j1i2j2ipjp

其中 i = (xi1, xi2, „, xip) 和 j = (xj1, xj2, „, xjp) 是两个p维的数据对象

上面的距离函数有如下特性：

d(i,j) >= 0：距离是一个非负的数值。d(i,i) = 0：一个对象与自身的距离是零

d(i,j) = d(j,i)：距离函数具有对称性

d(i,j)

明考斯基距离（ Minkowski distance）是欧几里德距离和曼哈坦距离的概化，它的定义如下：

d(i,j)(|xx||xx|...|xx|)i1j1i2j2ipjpq其中 i = (xi1, xi2, „, xip) 和 j = (xj1, xj2, „, xjp) 是两个p维的数据对象, q是一个正整数。 qqq

当q = 1时, d 称为曼哈坦距离（ Manhattan distance）；当q=2时, d 就成为欧几里德距离；当变量的重要性不同时，可以根据每个变量的重要性赋予一个权重。

对于区间标度变量来说，一旦选定一种距离作为衡量相异度（相似度）的基础，对象之间的距离越近，表示两个对象之间越相似，否则，表示两个对象之间越相异。对于类也是这样定义的：两个类之间的距离越近，表示两个类之间越相似，否则，表示两个类之间越相异。但是，类间距离的定义有几种方法，下面给出了四种常用测定两个类之间距离的方法，这儿的mi是类ci的中间值，ni是ci中记录的个数，|p-p’|是两个记录之间的距离。

1．单连通方法（single-link method）：取两个类中最近记录的距离为类的距离，即Minimum distance： Dmin(Ci,Cj)=Min|p-p’|(这儿p在类Ci中p’在Cj中)。此种方法可以生成细长的蛇形类，不适于应用在典型的一堆堆记录集合在一起的情况。

2．完全连通方法（complete-link method）：取两个类中最远记录的距离为类的距离，即Maximum distance：Dmax(Ci,Cj)=Max|p-p’|(这儿p在类Ci中p’在Cj中)。同第1中方法相反，此种方法易生成很小的记录都聚集在一起的类。

3．平均连通方法（group-average link）：计算两个类中所有记录对的距离平均值，即Average distance： Davg(Ci,Cj)=(1/ni*nj)∑p∈Ci∑p’ ∈Cj。效果介于1、2种算法之间。

4．Ward方法（Ward’s method）：计算两个类中所有记录的距离的和，即Sum distance： Dsum(Ci,Cj)=∑pЄCi∑p’ ЄCj。易于用在生成类层次的情况，对例外的数据（outliers）很敏感，很难应用于生成蛇形类。

5．类平均值连通方法：取两个类中间值的距离为类的距离，即Mean distance： Dmean(Ci,Cj)=|mi-mj|（mi是类ci的中间值）。此种方法是一种结合上面各种方法的措施，生成的聚类效果各项指标比较平均。

二元变量

一个二元变量只有两个状态：0或1，0表示该变量为空，1表示该变量存在。例如，给出一个描述病人的变量smoker，1表示病人抽烟，而0表示病人不抽烟。二元变量要采用特定的方法来计算其相异度。

标称变量

标称变量是二元变量的推广，它可以具有多于两个的状态，比如 变量

map_color可以有 red, yellow, blue, green四种状态。

序数型变量

一个序数型变量可以是离散的也可以是连续的。离散的序数型变量类似于标称变量，除了它的M个状态是以有意义的序列排序的，比如职称；连续的序数型变量类似于区间标度变量，但是它没有单位，值的相对顺序是必要的，而其实际大小并不重要。

比例标度型变量

比例标度型变量是总取正的度量值，有一个非线性的标度，近似的遵循指数标度，比如 AeBt or Ae-Bt。

混合类型的变量

在许多真实的数据库中。对象是被混合类型的变量描述的。一般来说，一个数据库可能包含上面列出的全部六种变量类型。

上面介绍的非区间标度变量是不能像处理区间标度变量来对待的，否则，会误导聚类结果，所以要采用特定的方法来计算其相异度。

2．3 主要聚类分析方法分类

目前在文献中存在大量的聚类算法。算法的选择取决于数据的类型、聚类的目的和应用。如果聚类分析被用作描述或探查的工具，可以对同样的数据尝试多种算法，已发现数据可能揭示的结果。

大体上，主要的聚类算法可以划分为如下几类：

（1）划分方法（partitioning method）：给定一个n个对象或元组的数据库，一个划分方法构建数据的k个划分，每个划分表示一个聚簇，并且k

个对象；（ii）每个对象必须属于且只属于一个组。注意在某些模糊划分技术中的二个要求可以放宽。

给定要构件的划分的数目k，划分方法首先创建一个初始划分。然后采用一种迭代的重定位技术，尝试通过对象在划分间移动来改进划分。一个好的划分的一般准则是：在同一个类中的对象之间进可能“接近”或相关，而不同类中的对象之间尽可能“远离”或“不同”。还有许多其它划分质量的评判准则。 为了达到全局最优，基于划分的聚类会要求穷举所有可能的划分。实际上，绝大多数应用采用了以下两个比较流行的启发式方法：（i）k－平均算法，在该算法中，每个簇用该簇中对象的平均值来表示。（ii）k－中心点算法，在该算法中，每个簇用接近聚类中心的一个对象来表示。这些启发式聚类方法对中小规模的数据库中发现球状簇很适用。为了对大规模的数据集进行聚类，以及处理复杂形状的聚类，基于划分的方法需要进一步的扩展。

（2）层次的方法（hierarchical method）：层次的方法对给定数据对象集合进行层次的分解。根据层次的分解如何形成，层次的方法可以分为凝聚的和分裂的。凝聚的方法，也称为自低向上的方法，一开始将每个对象作为单独的一个组，然后相继的合并相近的对象或组，直到所有的组合并为一个（层次的最上层），或者达到一个终止条件。分裂的方法，也称为自顶向下的方法，一开始将所有的对象置于一个簇中。在迭代的每一步中，一个簇被分裂为更小的簇，直到最终每个对象在单独的一个簇中，或者达到一个终止条件。

层次的方法的缺陷在于，一旦一个步骤（合并或分裂）完成，它就不能被撤消。这个严格规定是有用的，由于不用担心组合数目的不同选择，计算代价会较小。但是，该技术的一个主要问题是它不能更正错误的决定。

（3）基于密度的方法（density-based method）：绝大多数划分方法基于对象之间的距离进行聚类。这样的方法只能发现球状的簇，而在发现任意形状的簇上遇到了困难。随之提出了基于密度的另一类聚类方法，其主要思想是：只要临近区域的密度（对象或数据点的数目）超过某个阀值，就继续聚类。也就是说，对给定类中的每个数据点，在一个给定范围的区域中必须至少包含某个数目的点。这样的方法可以用来过滤“噪声”孤立点数据，发现任意形状的簇。 DBSCAN是一个由代表性的基于密度的方法，它根据一个密度阀值来控制簇的增长。

（4）基于网络的方法（grid-based method）：基于网络的方法把对象空间尽量化为有限数目的单元，形成一个网络结构。所有的聚类操作都在这个网络结

构（即量化的空间）上进行。这种方法的主要优点是它的处理速度很快，其处理时间独立于数据对象的数目，只与量化空间中每一维单元数目有关。

（5）基于模型的方法（model-based method）：基于模型的方法为每个簇假定了一个模型，寻找数据对给定模型的最佳拟合。一个基于模型的算法可能通过构建反映数据点空间分布的密度函数来定位聚类。它也基于标准的统计数字自动决定聚类的数目，考虑“噪声”数据或孤立点，从而产生健壮的聚类方法。 一些聚类算法集成了多种聚类方法的思想，所以有时将某个给定的算法划分为属于某类聚类方法是很困难的。此外，某些应用可能有特定的聚类标准，要求综合多个聚类技术。

2．4 聚类分析在Microarray分析中的位置及应用

微阵列（Microarray）技术在基因表达分析中的工作流程简单概括起来，分为以下几个过程：基因数据的收集、数据预处理、数据聚类分析、预测。聚类分析在整个过程中是承上启下的作用，聚类分析的对象（基因数据）来源于前期收集的数据，通过聚类分析，可以发现已经收集数据的正误，同时，也能根据数据本身的特征对数据进行分类，最后的聚类成果将为后面的预测工作提供材料。可以说，聚类分析在整个Microarray分析工程中非常重要，聚类分析算法的选择将影响到数据预处理过程，也将影响到最后的预测结果。所以，选择合适的聚类分析算法将有益于及时发现源数据中的“脏数据（有误的数据）”并修正源数据，对源数据进行高质量的聚类，产生高质量的聚类结果，为后面的预则过程提供高质量的材料，从而改善预测结果，这样，就提高了整个微阵列（Microarray）分析技术的工作质量！

第三章 聚类分析方法的分析

3．1 评价聚类分析方法好坏的准则

在微阵列（Microarray）分析中，保证最后聚类的质量是最重要的，以下几个方面是Microarray分析中对聚类算法选择的要求：

1．高的簇内相似性。

2．低的簇间相似性。

聚类定义的内涵就是：将一个对象的集合分割成几个类，每个类内的对象之间是相似的，但与其他类的对象是不相似的。所以说，这两点是一个聚类算法最根本也是最重要的要求。但是，一般来说一个数据库没有一种最好的分类方法。这是因为，聚类算法中的“相似度”是人为主观定义的，“相似性”也是根据算法中“相似度”的定义来判断的。因此，好的聚类算法要在类中对象的相似程度和类的数目之间找到一个最佳的结合点。

3．能应付脏数据。绝大多数现实世界中的数据库都包含了孤立点、空缺、未知数据或者错误的数据。一些聚类算法对于这样的数据敏感，可能导致低质量的聚类结果。

4．对于数据不同的顺序不敏感。一些聚类算法对于输入数据的顺序是敏感的。例如，同一个数据集合，当以不同的顺序提交给同一个算法时，可能生成差别很大的聚类结果。开发对数据输入顺序不敏感的算法具有重要的意义。

5．模型可解释性和可使用性。用户希望聚类结果是可解释的、可理解的和可用的。也就是说，聚类可能需要和特定的语义解释和应用相联系。应用目标如何影响聚类方法的选择也是一个重要的研究课题。

3．2 常见聚类分析方法的分析

3．2．1 划分方法

划分方法: 将一个包含n个数据对象的数据库组织成k个划分（k

典型的划分方法有：k-平均和k-中心点

1．k-平均方法（基于质心的技术）

k-平均方法以k为参数，把n个对象分为n个簇，以使簇内具有较高的相似度，而簇间的相似度较低。相似度的计算根据一个簇中对象的平均值（被看做簇的重心）来进行。

给定k，算法的处理流程如下:

算法：k-平均。划分的k-平均算法基于簇中的平均值。

输入：簇的数目k和包含个n对象的数据库

输出：k个簇，效果最好的聚类

方法：

1． 随机的选择k个对象作为初始的簇中心；

2．repeat

3．计算每个簇的平均值，并用该平均值代表相应的簇；

4. 将每个对象根据其与各个簇中心的距离，重新分配到与它最近的簇中；

5．更新簇的平均值，即计算每个簇中对象的平均值；

6. until直到不再有新的分配发生

该算法流程的图示如下：

对k-平均方法的分析如下：

优点

1、当结果簇是密集的，而簇于簇之间区别明显时，它的效果较好；

2、对于处理大数据集，该算法是相对高效的: 算法复杂度O(tkn), 其中n 是数据

对象的个数, k 是簇的个数, t是迭代的次数，通常k, t

3、算法通常终止于局部最优解。

缺点

1、只有当平均值有意义的情况下才能使用，对于类别字段不适用；

2、必须事先给定要生成的簇的个数k；

3、对“噪声”和异常数据敏感，少量的该类数据能够对平均值产生极大的影响；

4、不能发现非凸面形状的数据；

5、对初始簇的随机选择算法对聚类结果也会产生随机的影响。

由于上述的几个缺点，简单的用k-平均方法来处理Microarray分析中的聚类问题是不满足要求的，需要改进其中的步骤、流程，或者是和其它的方法结合起来用。

2．k-中心点方法（基于有代表性的对象的技术）

从上面对k-平均方法的分析我们可以看出：一个有极大值的对象可能相当程度上扭曲数据的分布，k-平均算法对于脏数据是敏感的。那么，“如何修改这个算法来消除这种敏感性？”不采用簇中对象的平均值作为参照点，可以选用簇中位置最中心的对象，即中心点。这样的划分方法仍然是基于最小化所有对象与其参照点之间的相异度之和的原则来执行的。这是k-中心点方法的基础。

k-中心点聚类算法的基本策略如下所示： 算法：k-中心点。对基于中心点或中心对象的划分的典型k-中心点算法。 输入：簇的数目k和包含个n对象的数据库。

输出：k个簇，使得所有对象与其最近中心点的相异度总和最小。

方法：

1． 随机的选择k个对象作为初始的簇中心；

2．repeat

3．指派每个剩余的对象给离它最近的中心点所代表的簇；

4. 随机的选择一个非中心点对象Orandom；

5．计算用Orandom代替Oj的总代价S；

6．If S

可以说，k-中心点方法是k-平均方法的改进，这种改进是为了防止有极大值脏数据。但是，k-中心点方法的执行代价比k-平均方法高。此外，以下的两个缺点没有得到解决：

1、必须事先给定要生成的簇的个数k；

2、对初始簇的随机选择算法对聚类结果也会产生随机的影响。

3．2．2 层次方法

层次聚类，按照聚类的层次结构可以把所有的记录层次聚类可以分为两种：凝聚的方式和分割的方式，取决于聚类层次结构的形成是自顶向下的还是自底向上的。

凝聚的方式：这是一种至底向上的方法，将每一条记录看作一个类，然后根据一些规则将他们聚合成越来越大的类，直到满足一些预先设定的条件。大多数的层次聚类方法属于这一类。

分割的方式：这种自顶向下的方法是一个与凝聚的方式相反的过程，将整个数据库作为一个大的类，然后按照一些规则将这个类分成小的类，直到满足一些预定的条件，例如类的数目到了预定值，最近的两个类之间的最小距离大于设定值。

下图给出了对于集合{a,b,c,d,e}层次聚类两种方式的聚类过程。从这个图我们可以看出，凝聚的方式是将每一个记录看作一个类，再按照一定的规则逐步将这些类合并。举个例子，如果类C1和类C2之间的距离小于预定的最小距离，那么他们就会被合并为一个类，这儿两个类的距离是由两个类中距离最近的一对记录来确定的。

分割的方式则是先将所有的记录作为一个大的类，然后再根据一些规则将它进行分割，例如最近的两个记录之间的距离。

无论凝聚的方式还是分割方式，用户都可以根据自己的要求来设定所得类的个数，也可以不需要指定聚类的个数，但用户可以指定希望得到的簇的数目作为一个结束条件。

层次方法的主要缺点：

1、没有良好的伸缩性：时间复杂度至少是 O(n2)；

2、一旦一个合并或分裂被执行，就不能修复，在选择合并或分裂的时候特别值得

注意。前面介绍的划分方法和层次方法都是用距离来衡量相似度的，下面给出了四种常用测定两个类之间距离的方法，这儿的mi是类Ci的中间值，ni是Ci中记录的个数，|p-p’|是两个记录之间的距离。

Minimum distance： Dmin(Ci,Cj)=Min|p-p’|;(这儿p在类Ci中p’在Cj中) Maximum distance：Dmax(Ci,Cj)=Max|p-p’|;(这儿p在类Ci中p’在Cj中) Mean distance： Dmean(Ci,Cj)=|mi-mj|;

Average distance： Davg(Ci,Cj)=(1/ni*nj)∑p∈Ci∑p’ ∈ Cj;

层次聚集虽然比较简单，但是在选择凝聚或者分割点的时候经常会遇到一些困难，这个是非常关键的，因为一旦记录被凝聚或者分割以后，下一步的工作是建立在新形成的类的基础之上的。因此，如果其中任何一步没有做好的话，就会影响最终聚集的结果。这个方法并不是太好，因为要牵涉到很大数量的类和记录。

一个比较有前途的能够提高聚集质量的方向是将层次聚集和其它的聚集结合起来进行。

3．2．3 基于密度的方法

上面介绍的两种聚类方法都是基于距离来判断对象之间、类之间相似度的，这两种聚类方法有一个共同点，就是只能发现球状的聚类结果，为了发现任意形状的聚类结果，提出了基于密度的聚类方法。这类方法将簇看做是数据空间中被低密度区域分割开的高密度对象区域。下面就介绍一种基于密度的聚类方法。

DBSCAN(基于高密度连接区域的密度聚类方法)

DBSCAN（Density-Based Spatial Clustering of Applications with Noise）是一个基于密度的聚类算法。该算法将于有足够高密度的区域划分为簇，并可以在带有“噪声”的空间数据库中发现任意形状的聚类。它定义簇为密度相连的点的最大集合。

基于密度的聚类的基本想法涉及一些新的定义。我们先给出这些定义，然后用一个例子加以说明。

定义

 邻域N(q): 给定对象半径 内的区域，即{q  D | distance(p,q)

对象p 从对象q 出发是直接密度可达 ：给定义个对象集合D，如果p是在q的- 邻域内，而q是一个核心对象，我们说对象p从对象q出发是直接密度可达的，用数学表达法就是：pN(q) 且|N(q)|  MinPts

对象p从对象q关于  和MinPts密度可达：存在对象链p1，p2，„，pn，p1=q，pn=p，pi∈D，pi+1是从pi关于  和MinPts直接密度可达的（非对称）

对象p和q关于  和MinPts密度相连：存在对象o ∈D，使得对象p和q 从o关于  和MinPts密度可达

在上图中，给定  圆的半径，MinPts=3。基于上述的定义：

·在被标记的点中，P、M、O、R及S都是核心对象，因为它们在 -领域内至少包含了三个点。

·Q是从M直接密度可达的，M是从P值直接密度可达的。

·因为Q是从M值直接密度可达的。M是从P直接密度可达的，所以Q是从P密度可达的。

·O、R和S都是密度相连的。

DBSCAN基本思想如下：

簇：基于密度可达性的最大的密度相连对象的集合

噪音：不在任何簇中的对象

边界对象：不是核心对象，但在簇中，即至少从一个核心对象直接可达 如下图所示：

DBSCAN算法流程如下： 算法：DBSCAN。基于高密度连接区域的密度聚类方法。

输入：阀值 、 MinPts和包含个n对象的数据库。

输出：基于密度的聚类结果。

方法：

1）任意选择没有加簇标签的点 p；

2）找到从p关于 和 MinPts 密度可达的所有点 ；

3）如果|N(q)|MinPts ，则p是核心对象,形成一个新的簇，给簇内所有的对象

点加簇标签；

4）如果p 是边界点, 则处理数据库的下一点；

5）重复上述过程，直到所有的点处理完毕。

示意图如下：

不足之处：

1、只能发现密度相仿的簇；

2、对用户定义的参数（  和 MinPts ）敏感；

3、计算复杂度至少为O(n)，采用R-树等空间索引技术，计算复杂度为o(nlogn)。

3．2．4 小节

本章介绍了Microarray分析中对聚类算法选择的要求，并对聚类分析中常用的几种算法（划分方法、层次方法、基于密度的方法）进行了详细的分析。从上面的分析我们可以看出，各种聚类方法各有特色，并且也都有各自的优点和缺点。划分方法和层次方法都是基于距离来判断相似度的；而基于密度的方法是基于密度这一个概念来判断相似度的，最终而言，相似度是根据用户输入的参数  和 MinPts来判断的。正是因为判断相似度方法的不同，所以它们聚类的思想大相径庭，得出的结果也是很不一样的。

除了这几种聚类方法，还有基于网络的方法、基于模型的方法等聚类方法，但是，这些聚类方法与本课题无关，就不再具体研究和分析了。对上述几种聚类方法的分析是为了：①研究和分析算法的思想，并总结其优缺点；②研究对算法性能的改进方法；③将算法的思想应用于具体的事物，解决具体问题。在下一章我们将就一个具体的问题进行分析、讨论，最后设计一种适合解决具体问题的算法，并实现。

2

第四章 Microarray 分析中聚类分析算法的设计

4．1 问题的提出

在微阵列（Microarray）分析中，概括来说，可以分为以下几个阶段：获取数据、分析数据、预测。要求通过对常用聚类方法的研究和分析，选择或者设计一个算法，对获取的数据进行聚类。要求能读取各种不同的数据格式，对于聚类算法，要求能达到以下几点：1、在簇内相似性和簇间相异性之间协调，尽量达到高簇内相似性和簇间相异性；2、对脏数据不敏感；3、对数据顺序不敏感；4、选择或者设计的算法模型可解释和可使用；5、聚类的结果质量高。

获取的数据文件格式如下：

1、该数据文件存储成文本文件的格式，在这种数据文件中，数据是文本格式的，除第一行外，每一行代表一条基因纪录。第一行存储的是每一条基因纪录的属性个数和各个属性的名称。行于行之间是新行符（Newline）分隔开的，一条基因记录第一个数据存储的是该基因纪录的名称，后面的就是相应属性的值。图示如下：

2、该数据文件存储成文本文件的格式，在这种数据文件中，数据是文本格式的，每一行代表一条基因纪录。一条基因记录第一个数据存储的是该基因纪录的名称，后面的就是相应属性的值，行于行之间是用新行符分隔开的。相当于数据格式1去掉了第一行之后的文件格式。图示如下：

4．2 问题分析和算法设计

根据用户提出的需求，将任务分为两大部分。1、实现数据的预处理（包括数据读取、输出和规范化）；2、应用合适的聚类算法，将读取出的数据进行聚类分析。

4．2．1 数据预处理

数据预处理的任务是实现数据读取、规范化和将数据分解成为矢量数据。

1、 数据读取。是指从数据文件中将数据读出来，并存放在一个数组中。对三种不

同的数据文件格式分别编写不同的函数：

Get_data1(char *filename,int rows，double data[MAX_ROWS][MAX_COLUMNS])、Get_data2(char *filename,int rows，double data[MAX_ROWS][MAX_COLUMNS])、其中参数filename表示数据文件名，参数rows表示要读取数据的行数，数组data是用来返回读取的数据，数组data中，每一行代表一条基因记录。如果程序中出现读取文件的错误，将会错误退出。

2、数据规范化。是指对从数据文件中读取的数据按照一定的算法进行规范化处理，以减少数据的运算量。 规范化的算法如下： （1） 读入数据；

（2） 求每一行记录数据中的极小值元素；

（3） 该行中的每一个元素减去极小值，得到一新的行纪录。 以后的聚类及分析，均用此规范后的数据 函数格式：

Standard(double data[MAX_ROWS][MAX_COLUMNS],int rows,int columns) 其中参数data是数据数组，也用来最后返回规范化的数据；参数rows和columns分别表示数据数组的行和列。

3、分解成为矢量数据。是指将数据文件中的数据分解，将每一行的记录数据（也称为矢量数据）分离开来，使得一个文件中只有一个矢量数据，这种文件的扩展名为“.dis”，此外，还要产生一个矢量数据的索引文件，文件的扩展名为“*.ml”。

为什么要对数据文件进行分解成为矢量数据呢？主要是为了以后对聚类结果的分析提供材料。

4．2．2 聚类算法设计

通过上面的数据预处理，我们可以得到一个数据数组data[MAX_ROWS][MAX_COLUMNS]，该数组每一行代表一条数据记录，也就是参加聚类的对象，下面的聚类算法就是对该数组进行聚类分析。

根据用户的要求，对于聚类算法，要求能达到以下几点：1、在簇内相似性和簇间相异性之间协调，尽量达到高簇内相似性和簇间相异性；2、对脏数据不敏感；3、对数据顺序不敏感；4、选择或者设计的算法模型可解释和可使用；5、聚类的结果质量高。

考虑前面分析过的聚类算法：

1、 如果采用k-平均方法，则有以下问题需要解决：

(1) 必须事先给定要生成的簇的个数k；

(2) 对“噪声”和异常数据敏感，少量的该类数据能够对平均值产生极大的影

响；

(3) 对初始簇的随机选择算法对聚类结果也会产生随机的影响。

由于上述的几个缺点，简单的用k-平均方法来处理Microarray分析中的聚类问题是不满足要求的，需要改进其中的不足之处，或者是和其它的方法结合起来用。

2、如果采用k-中心点方法，虽然，k-中心点方法是k-平均方法的改进，这种改进防止有极大值脏数据的干扰。但是以下的两个缺点没有得到解决：

(1) 必须事先给定要生成的簇的个数k；

(2) 对初始簇的随机选择算法对聚类结果也会产生随机的影响。

所以，简单的使用k-中心点方法来处理Microarray分析中的聚类问题也是不满足要求的，需要改进其中的不足之处，或者是和其它的方法结合起来用。

3、如果采用层次聚类方法，可以用凝聚的层次划分方法，但是在选择凝聚的时候经常会遇到一些困难，在什么情况下才可以凝聚呢？这个是非常关键的，一个比较有希望的能够提高聚集质量的方向是将层次聚集和其它的聚集结合起来进行。解决“在什么情况下才可以凝聚？”这个问题。

4、如果采用基于密度的聚类方法，只能发现密度相仿的簇，而且需要用户输入参数 和 MinPts，算法对用户定义的参数（  和 MinPts ）很敏感，所以，如果采用基于密度的聚类方法，需要对输入参数这一个问题改进。

根据上面的分析，简单的采用四种聚类方法的任意一种都是不满足用于需求的，需要对上述的算法进行改进或者是和其它算法结合使用。有以下几个问题有待于解决：

（1）尽量不要用户输入对算法影响很大的参数； （2）避免对“噪声”和异常数据敏感；

（3）避免初始簇的随机选择算法对聚类结果也会产生随机的影响；

（4）如果采用层次聚类方法，要通过合适的选择来决定凝聚的时候。 根据用户的需求，结合层次聚类方法和划分方法的聚类思想，设计如下的算法，可以成功的解决以上的四个问题。算法的聚类设计思想如下：

Microarray分析聚类算法设计

[初始条件]

设有S个向量（也就是数组data[][]中的行数），维数为t（也就是数组data[][]中的列数），构成一个组G， ︱G︱=S，按照下面的算法聚成若干类，限定类数不超过W。

[算法流程] 一、初始化

（1） 读入数据组Ｇ至数组 data[][]，同时纪录数组的行数Rｏｗｓ和

列数Cｏｌｕｍｎｓ

（2） 计算G的平均向量Ｓ０ （3） 计算每个向量Ｓｋ与Ｓ０的距离

xk(SklS0l)

l1

1

2

t

2

k=1,2,„„S

令x=max xk x=min xk k=1,2,„„S

二、对剩余组的向量，依次操作以下步骤，直到剩余组为空

假设已经有G1,G2,„„,Gn个类，

（一）在剩余组中（选择剩余组中参加下一步聚类的数据向量）

1、计算剩余组的平均向量S0 2、计算每个向量Sk与S0的距离

(4) 将{ x}和{ x}作为初始的两个类，下面对剩下的向量进行聚类

1

2

xk

j

2

(SS)kl0ll1

t

令 x=min xk

（二）计算每个类的平均值X1,X2,„„,Xn（均为向量）

计算类之间的距离

dij

(x

l1n

t

il

xjl)

2

i

计算类之间的平均距离

D1(dij)/(n(n1)/2) i

jii1

j

j

（三）判断x可能属于某个类q

1、计算x和各个类平均值之间的距离

xm

j

j

j

2(xx)lml

j

l1

t

,m=1,„„,n

2、xq=min xm m=1,„„,n 设有较小值的那个类为第q类

j

j

（四）如果n>=W（即已经达到了最大的类数），则x 属于第q类，x 从剩余

组中删除。若剩余组为空，则聚类结束退出；否则，对剩余组中的其它元素，转（一）。

如果n

计算类之间的距离

dij

(x

l1n1

t

il

xjl)

2

i

D2(dij)/((n1)n/2)

jii1

i

如果D2>D1，则x 属于新类，即Gn+1={x }，x 从剩余组中删除，若剩余组为空，则聚类结束退出；否则，对剩余组中的其它元素，转（一）。 如果D2

j

j

j

j

j

j

该算法是对层次聚类算法的改进。该算法是用距离来衡量相似度的，间隔距离越小表示相似度越高，否则，相似度就低。该算法需要一个参数W，表示限制类数不超过W，这个参数有两个用处：1、对于非专业用户，参数W可以设置为一个很大的数，这样对算法本身没有实际上的影响；2、对于专业用户，参数W可以设置为一个专业化参数，用来限制最大类数。初始化的时候是选择离中心最远和最近的两个向量形成类，不会有随机因素的影响。而且，在剩余向量组中选择下一步参加聚类的向量时，应用了一定算法，挑选和剩余组中心距离最近的那个向量，这样，避免了对“噪声”数据和异常数据的敏感。在最后决定剩余组的向量属于哪个类的时候，除了之前考虑了类内部的相似性，而且还考虑了类之间的相异性（即类间平均距离要增大）。通过以上的考虑和设计，前面提到的几个问题就得到了较好解决。

第五章 Microarray分析中聚类算法的实现

5．1 Microarray 分析中的数据预处理

在本次的算法分析和研究中，使用的编程语言是C语言。本节介绍C语言中对Microarray 分析中的数据文件格式操作的相应知识和处理方法，并应用C语言实现数据预处理。

5．1．1 Microarray 分析中的数据预处理技术

Microarray 分析中的数据文件格式在前面已经介绍了，两种数据文件格式都是存储成文本文件的格式，在这两种数据文件中，数据是文本格式的，除了第一行有区别之外，其它的每一行代表一条基因纪录。行于行之间是新行符（Newline）分隔开的，一条基因记录第一个数据存储的是该基因纪录的名称，后面的就是相应属性的值。对于这种文本文件格式，需要用到C语言中相应的处理技术，下面就对其处理技术做简单介绍。

C文件系统是为适应各种设备而设计的，设备包括终端、磁盘和磁带驱动器等。各种设备的差别很大，但缓冲文件系统把每种设备都转换成称为流（stream）的逻辑设备。所有流的性质完全类似。因为流基本上与设备无关，所以能写入磁盘文件的同一函数也能写入另一类型设备，如控制台终端等。

文件指针（file pointer）是贯穿缓冲I/O系统的主线。文件指针是指向定义文件操作信息的指针，信息中包括文件的名字、状态和当前读写位置。本质上，文件指针标识一个特定磁盘文件，被相连的流用来指导缓冲文件函数的操作。

文件指针是FILE型指针变量（类型FILE在stdio.h中定义）。为了读写文件，程序需要使用文件指针。为得到一个文件指针变量，应书写的语句结构如下所示：

FILE *fp; //这样就创建了一个文件指针；

下面介绍一些主要的文件操作函数，我编写源程序时用到了其中的函数。 1、fopen()

函数fopen()打开一个流，将该流和一个文件关联，然后返回有关的文件指针。最常见的文件是磁盘文件。函数fopen()的原型是：

FILE *fopen(const char *filename,const char *mode);

其中，filename指向构成文件名的串，串中可以有路径说明。Mode指向的串决定文件的打开方式，其合法值见下表所示：

函数fopen()返回一个文件指针，空文件指针表示fopen()失败。程序中该更

2、fclose()关闭文件

函数fclose()关闭fopen()打开的流。fclose()把遗留在缓冲区的数据写入文fclose()的原型是：

int fclose(FILE * fp);

其中，fp是fopen()返回的文件指针，fclose()返回零值表示关闭成功，返回其它值表示关闭错误。函数fclose()判断并报告错误内容。

件，实施操作系统级的关闭操作。

改函数fopen()返回的文件指针。

3、fgetc()、fputc()读、写字符

函数fputc()向预先用fopen()“写打开”的文件写一个字符，其原型是： int fputc(int ch,FILE * fp);

其中，fp是fopen()返回的文件指针，ch是被写入文件的字符。文件指针告诉fputc()哪个文件要写入。尽管ch定义为int型，但质写入低字节。

fputc()成功时返回写入的字符，否则返回eof。

函数fgetc()向预先用fopen()“读打开”的文件读一个字符，其原型是： int fgetc(FILE * fp);

其中，fp是fopen()返回的文件指针。fgetc()返回整数，但字符包含在低字节中。除非发生错误，否则高字节为零。

到达文件尾的时候，fgetc()返回eof。

4、fgets()、fputs()读、写字符串

C语言还支持函数fgets()、fputs()，用来向文件读写字符串而不是单字符。两个函数的原型是：

Int fputs(const char * str,FILE * fp); Char *fgets(char * str,int length,FILE * fp);

函数fputs()向指定的流写串（由str指定）。出错时，fputs()返回eof。 函数fgets()从指定的流读串，读到新行字符（newline）或已经读入length-1个字符时结束。读到新行时，把新行符也作为读入串的一部分(与gets()的处理不同)。结果串用null结尾。成功时，函数返回指针str的值，失败时返回空指针(null)。

5、feof()

6、rewind() 函数feof()判断读入文件数据时是否到达文件尾。函数原型是： int feof(FILE * fp); 到达文件尾时，feof（）返回真值，否则返回零。

函数rewind()重置文件的位置指示，将其指到函数变元所指定的文件开始处，即该函数具有“重绕”文件的功能。函数rewind()的原型是：

7、ferror()

函数ferror()判断文件操作时是否出错，其原型是： int ferror(FILE * fp);

其中，fp是有效文件指针。在最近一次文件操作中出错时，函数返回真值；否则，返回假值。

8、remove()

5．1．2 Microarray 分析中的数据预处理实现

根据上面对Microarray 分析中的数据预处理的分析，将数据预处理分为以下三个部分：1、从数据文件中读取数据；2、对读取的数据进行规范化处理；3、将规范化的数据分解成为矢量数据，为以后的工作提供材料。

1、从数据文件中读取数据的函数。

（1）读取第一种数据文件格式的函数Get_data1(char *filename,double

data1[MAX_ROWS][MAX_COLUMNS])；

用途：从第一种数据文件中读取数据至数组data1中，同时纪录数据的实际行数ROWS和列数COLUMNS。

（2）读取第二种数据文件格式的函数Get_data2(char *filename,double

data1[MAX_ROWS][MAX_COLUMNS])； 函数remove()删除指定文件，其原型是： int remove(const char *filename); 成功时，remove()返回零，否则返回非零。 void rewind(FILE * fp); 其中，fp是有效文件指针。

用途：从第二种数据文件中读取数据至数组data1中，同时纪录数据的实际行数ROWS和列数COLUMNS。

2、对读取的数据进行规范化处理的函数。

Standard(double data1[MAX_ROWS][MAX_COLUMNS]，int rows)

说明：对读取的数组data1按照前面的规范化算法进行规范化操作，参数rows是数组的列数（即向量维数）。输出结果数组data1是规范化之后的数据。

3、将规范化的数据分解成为矢量数据。

具体的说明如下：

（1）对第一种数据格式分解的程序名称是：“nordata1.c”,说明：

nordata1.c 是对第一类的数据进行操作的，命令行格式为：

命令 文件 行/列 单个/连续输出 一系列的记录号

如：

nordata normal.dat r s 1 6 4 9 100 200 „„

输出第1、6、4、9、100、200、„„行数据 输出第1--6、81--90、100--200、„„行数据，如果参数个数为奇数，侧最nordata normal.dat r a 1 6 81 90 100 200 „„ 后一个参数单个输出

nordata normal.dat r a 1 6 81 90 100 200 205„„

输出到文件时，产生no_r1.dis,no_r6.dis,„„

关于*.dis文件的命名规则：

关于*.ml文件的命名规则：

文件名的前三个字符（不足补下划线‘_’）+‘_’+第几次使用命令+‘.ml’ 文件名前两个字符（不足补下划线‘_’）+‘_’+行/列+数据序号+‘.dis’ nor_1.ml,里面的内容为产生的*.dis文件列表 normal.tol里面是*.ml文件所对应的命令行，供以后参考用 输出第1--6、81--90、100--200、205„„行数据

关于*.tol文件的命名规则：

如果一次生成*1.ml，接着是*2.ml„„*100.ml„如果中途删除*1.ml，则下次生成*.ml，其它文件形式类似。

（2）对第二种数据格式分解的程序名称是：“nordata2.c”,说明：

nordata2.c 是对第二类的数据进行操作的，命令行格式为：

命令 文件 行/列 单个/连续输出 一系列的记录号

如：

nordata2 tumor.dat r s 1 6 4 9 100 200 „„

输出第1、6、4、9、100、200、„„行数据 输出第1--6、81--90、100--200、„„行数据，如果参数个数为奇数，侧最nordata2 tumor.dat r a 1 6 81 90 100 200 „„ 后一个参数单个输出

nordata2 tumor.dat r a 1 6 81 90 100 200 205„„

输出到文件时，产生tu_r1.dis,tu_r6.dis,„„

其它的命名说明和方法与程序nordata1.c的处理方法相同

tum_1.ml,里面的内容为产生的*.dis文件列表 tumor.tol里面是*.ml文件所对应的命令行，供以后参考用 输出第1--6、81--90、100--200、205„„行数据 文件名+‘.tol’

5．2 Microarray 分析中聚类算法的实现

由于用户要求为了以后的工作需要，在编码的时候尽量不要用C++的扩展功能，所以，我的程序编码基本上都是通俗易懂的C语句。程序中用到了一些接口程序，先对接口函数做一定的说明，以下是接口函数：

*********************************************************************** Create_Class(double *data1,int r，int Columns)

说明： 新建一个序号为 r 的类，新类中包含一个矢量数据 data1[]，参数Columns是矢量数据data1的维数，该函数没有返回值；

*********************************************************************** double Distance_Class(int p,int q，int Columns)

说明：求出序号为p和q的类之间的距离，参数Columns是矢量数据data1的维数，并将距离值返回；

*********************************************************************** Ave_Class(int r，double *average1, int Columns)

说明：求出序号为ｒ的类的均值，参数Columns是矢量数据data1的维数，并返回均值矢量至矢量average1中；

*********************************************************************** Add_Class(double *data1,int r, int Columns)

说明：将矢量数据data1[]加入到序号为 r 的类中，参数Columns是矢量数据data1的维数，该函数没有返回值；

*********************************************************************** Delete_Class(int r)

说明：删除序号为ｒ的类中最后加入的矢量元素，该函数没有返回值。

***********************************************************************

此外，为了提高程序的可读性，还编写了一些需要的对矢量数据的处理函数，说明如下：

1、 数void exchange(double *data1,double *data2,int C)

用途：交换C维矢量数据data1 和 data2的函数，该函数无返回值。

/* the function is used to exchange the vectors(data1 and data2) */ void exchange(double *data1,double *data2,int C)

{

int i;

double temp;

for (i=0;i

{

temp=data1[i];

}

data1[i]=data2[i]; data2[i]=temp; }

2、函数void average(double data1[MAX_ROWS][MAX_COLUMNS],int R,int C,double *data2)

用途：获取数组平均矢量的函数，其中参数二维数组data1源数据数组，参数R和C表示源数组的要求平均矢量的实际行数和列数，参数矢量data2用来返回平均矢量。

/* the function is used to get the average vector of the array data MAX_ROWS---the max rows of the array

MAX_COLUMNS---the max columns of the array

*/

void average(double data1[MAX_ROWS][MAX_COLUMNS],int R,int C,double *data2)

{

int i,j;

for (j=0;j

{

}

3、函数double distance(double *data1,double *data2,int C)

用途：计算两个矢量欧几里德距离的函数，其中参数矢量data1 和 data2表示要求距离的两个矢量，参数C表示矢量数据的维数，返回值是double型数据。

/* the function is used to compute the distance between the vectors(data1 data2[j]=0; for(i=0;i

and data 2) */

double distance(double *data1,double *data2,int C)

{

int i;

double d;

d=0;

for (i=0;i

d=d+(data1[i]-data2[i])*(data1[i]-data2[i]);

d=sqrt(d);

return d;

}

在实现算法的时候，将该算法编成一个函数，使用时通过参数来调用该函数，MAX_ROWS 和 MAX_COLUMNS分别表示源矩阵数据的最大行数和最大列数。参数二维数组data表示进行聚类的源数组数据，参数Rows 和 Columns分别表示源矩阵的实际行数和实际列数，参数W表示的是聚类的最大类数目。在实现该聚类算法时，要用到一些接口函数和矢量处理函数，上面已经对聚类函数中用到的接口函数和矢量处理函数做了介绍，下面是该函数的主要语句：

/* the function is used to cluster the vector */

void cluster(double data[MAX_ROWS][MAX_COLUMNS],int Rows,int Columns,int W)

{

int i,j,k,q;

int max_no,min_no;

double max,min;

double ave[MAX_COLUMNS];

double c_ave[MAX_ROWS][MAX_COLUMNS],d[MAX_ROWS][MAX_COLUMNS]; double x[MAX_ROWS],D1,D2;

/* the explain of variables :

q---record the no of the selected class max_no,min_no---record the no of the max and min

num---the amounts of the rest vectors

c_num---the amounts of classes

max,min---record the average,max,min

ave---record the average of vectors

c_ave[i]---record the average of classes in no i

d[i][j]---record the distance between class no i and class no j x[]---record distances

D1,D2---record the average of the distances in classes

*/

/* initialization */

average(data,Rows,Columns,ave); //get the average vector

for(i=0;i

{

x[i]=distance(data[i],ave,Columns); //get the distance

{

max=x[i];max_no=i;

min=x[i];min_no=i;

}

else

{

if (x[i]>max)

{ max_no=i; max=x[i]; }

if (x[i]

{ min_no=i; min=x[i]; }

}

}

exchange(data[max_no],data[Rows-1],Columns);

exchange(data[min_no],data[Rows-2],Columns);

create_class(data[Rows-1],0,Columns); if (i==0)

create_class(data[Rows-2],1,Columns);

num=Rows-2;

c_num=2;

/* cluster the rest vectors,handle a vector one time */ while (num>0)

{

/*get the vector to handle*/

average(data,num,Columns,ave); //get the average vector for (i=0;i

{

}

exchange(data[min_no],data[num-1],Columns); //exchange vectors

/* get the average of the classes */

for(i=0;i

{

ave_class(i,c_ave[i],Columns);

//compute the distance to the classes,and get the min x[i]=distance(data[i],ave,Columns); //get the distance if (i==0) { } else if (x[i]>max) min=x[i]; { min_no=i; } min=x[i]; min_no=i;

x[i]=distance(data[num-1],c_ave[i],Columns);

if (i==0)

{ min=x[i];min_no=i; }

else

}

q=min_no; //q--- no of the selected class

if (c_num>=W)

{

add_class(data[num-1],q,Columns);

num--;

}

else

{

/* compute the the old average distance among classes */ D1=0;

for (i=0;i

{

for(j=i+1;j

D1=D1+d[i][j];

}

}

D1=D1/(c_num*(c_num-1)/2);

add_class(data[num-1],q,Columns); //need to modify

/* compute the the new distance among classes after add vector*/ D2=0; if (x[i]

for (i=0;i

{

for(j=i+1;j

D2=D2+d[i][j];

}

}

D2=D2/(c_num*(c_num-1)/2);

if (D2>D1)

{

delete_class(q);

create_class(data[num-1],c_num,Columns);

c_num++;

}

num--;

}

} // the last circle

} //the function

第六章 算法结果分析和展望

6．1 算法结果分析

在该算法中，参数W是唯一有可能影响聚类结果的用户输入参数，在算法设计的时候，是这样考虑的：W是对聚类过程的一个限制参数，要求聚类结果类的数目不能超过W。这个参数也是用来区分专业用户和非专业用户的，非专业用户不知道聚类的可能结果，只是依靠数据的本身属性和关联来进行聚类；专业用户可以在聚类之前估计或是猜测聚类的数目最多的数值，这样将参数W输入，就可以约束聚类结果不至于太散。

下面是用户输入参数W和不输入参数的聚类结果比较：

1、不输入参数W

2、输入参数W=3

对算法进行聚类结果示例如下：

该算法聚类结果整体上提高了聚类的质量，在处理“脏数据”方面，大部分时候效果是不错的，如下图所示：

通过对本算法运行结果的分析，发现：该算法通过对参数W的灵活使用，较好的解决了用户输入参数对算法聚类结果产生影响的问题；选择将要参加聚类的矢量的步骤，避免了初始簇的随机选择算法对聚类结果产生随机影响的问题；选择将要参加聚类的矢量的步骤，同时也解决了对数据顺序的敏感性； 在决定凝聚的时候，即考虑了簇内部的相似性，也考虑了簇之间的相异性。这些都是比较好的改进之处。该算法聚类结果整体上提高了聚类的质量，在处理“脏数据”方面，大部分时候效果是不错的，但是，有时候对“脏数据”还是敏感的，所以，虽然该算法通过结合传统聚类算法的特点，改进不足之处，从整体上提高了聚类算法的质量，为以后的工作提供了较高质量的聚类结果。

但是，该算法毕竟是层次聚类方法和k-平均划分方法的结合应用，在表示类中

心的时候用的是“质心”来表示一个类，这样如果一个类中有极大值的对象元素，那么这个元素对“质心”的影响也是很大的，尽管在进行聚类之前对参加聚类的元素进行了算法挑选，但是，仍然不能排除有“脏数据”的可能性。这样会在一定程度上影响聚类结果。

为了处理这种情况，可以考虑以下改进方案：

（一）改进类中心的表示方法。

选择基于质心和基于代表对象方法之间的中间策略。不用单个质心或单

个对象来代表一个类，而是选择数据空间中固定数目的具有代表性的点。一个类的代表点通过如下方式产生：首先选择类中分散的对象，然后根据一个特定的分数或收缩因子向类中心“收缩”或移动它们。如下图所示：

上图大圆圈表示一个类，里面的点表示该类中的对象元素，很明显，其

中有两个点是有极大值的点或者是“脏数据”。这样如果用平均值来代表该类，将会产生较大的误差；如果用中心点代表该类，又会显得太过于个体化。为了避免这两种极端的现象，我们可以使用一种折衷的方法，就是上面那种选择基于质心和基于代表对象方法之间的中间策略。这样，我们小圆圈中的对象元素平均值来代表该类，这样就比用平均值代表和中心点代表的方法效果要好。

那么怎样来实现这种策略呢？我们可以按下面的步骤进行：

（1） 计算类的直径；

（2） 选择一个特定的分数或收缩因子向类中心“收缩”或移动它们（按

照一个特定的分数或收缩因子构造上图中的小圆圈；

（3） 计算小圆圈中对象点的平均值，用来代表该类。

特定的分数或收缩因子可以选择统计意义的数，例如：1/3、1/2、2/3等。

（二）在聚类过程中剔除孤立点。如下图所示：

在聚类的过程中剔除上图的几个孤立点。可以采用这样的方法实现：在

聚类过程中探测类的增长速度，如果一个类增长速度太慢，就去掉它。

6．2 展望

聚类分析在Microarray分析的整个过程中起到了承上启下的作用，聚类分析的对象（基因数据）来源于前期收集的数据，通过聚类分析，可以发现已经收集数据的正误，同时，也能根据数据本身的特征对数据进行分类，最后的聚类结果将为后面的预测工作提供材料。可以说，聚类分析在整个Microarray分析工程中非常重要，聚类分析算法的选择将影响到数据预处理过程，也将影响到最后的预测结果。聚类分析方法在Microarray分析中的广泛和深入应用将会进一步的改善Microarray分析的质量和效果。

我们深信，随着Microarray技术及其应用范围的不断扩大，以基因类型和分子描绘的电子计算机预示和论断的临床数据以及生物信息学分析（包含前面提到的聚类分析），将为基因工程药物的发展、新药筛选和基因治疗以及逐步实现个性化治疗开拓广泛的前景。

结束语

微阵列（Microarray ）技术是一门新兴学科，它是结合了多门学科的特征而形成的一门交叉学科，国际上有许多研究人员正在该领域做研究工作，在国内也正处于探讨和发展阶段。目前，聚类分析这个领域的内容在国内一般是研究生教育开设的课程，将聚类分析应用于具体的科学领域是更深的研究方向。对于只有一个本科生水平的我来说，领域背景和研究内容都是全新的，需要在短短的毕业设计规定期间内进行学习和探索。

经过三个多月的学习和努力，我按期完成了任务书中规定的工作：

（1）学习微阵列（Microarray）技术的处理流程和相关领域知识；

（2）分析和研究常用的聚类算法，指出其中的不足之处，并提供改进意见；

（3）在前两步工作的基础上，对微阵列（Microarray）技术中收集的原数据

文件格式进行分析和研究，设计相适应的聚类算法；

（4）用C语言编程实现设计的算法，并实现对数据的预处理；

（5）对算法的运行结果进行采集和分析，并探索进一步的改进策略。

本次毕业设计给我提供了解先进科学领域的机会，也给我提供了将大学所学习的理论知识、专业知识同科研活动相结合的机会。通过本次毕业设计，我学到了很多新知识，对聚类分析在微阵列技术分析的应用做了一定探索，并且历经了一项科研开发的全过程（分析、设计、实现、测试和结果分析）。我亲身感受到了科研开发工作的艰辛和对投入精力的需要；同时，也感受到在科研开发过程中与指导教师之间交流和沟通的重要性。

聚类分析是一个应用很广泛的领域，理论程度深，需要一定的数学基础，以往在这方面的研究多是对常用的聚类方法进行分析和改进，很少有提出新的聚类思路和方法，本次毕业设计对聚类方法做了一定的研究和探索，并形成了自己的想法和思路。值得庆幸的是，我考取了北京大学计算机系硕士研究生，聚类分析领域和我的专业方向联系非常紧密，在以后的研究生学习和研究期间，我将对该领域做进一步的学习和探索。

三个多月的毕业设计即将结束，在这里，我要特别感谢我的指导教师方伟武老师。方老师是中科院应用数学研究所的研究员，方老师的研究工作和学术交流工作特别忙，但是方老师经常在百忙之中抽出时间对我进行关心和指导，正是在方老师的指引下，我对领域知识和研究内容有了详细的了解和较深的认识，也锻炼了自己

的独立工作能力和自学能力，方老师对科学一丝不苟的科学作风给我留下了深刻的印象。

同时，我要感谢胡矩老师。胡老师在我的毕业论文选题、毕业设计过程中的注意事项和毕业论文的成文思路方面都给予了很大的帮助。

此外，我要感谢那些曾经在毕业设计过程中给于我关心和帮助的人，谢谢你们！ 就要毕业了，就要离开这令人回味无穷的母校了。我要感谢母校——北京信息工程学院，是您给了我丰富多彩的大学生活，是您给了我步入社会的必备知识，是您为我插上了飞向理想的翅膀！

参考文献

（1）Peter Gwynne and Guy Page.《MICROARRAY ANALYSIS:the next revolution

in molecular biology》.Internet

（2）田爱景.《关于生物信息处理于生命信息科学》.湖北大学学报.1998年.第20

卷第4期

（3）杨劲松，陈诗书.《微阵列（microarrays）技术及其应用》.生命科学.2001

年4月.第13卷第2期

（4）陈执中.《DNA微阵列技术的研究应用进展》.药物生物技术.2001,8(6):357-360

（5）罗田.《微阵列计划(Microarray Project)》.生物工程进展.2001,Vol.21,No.3

（6）（英）P.史尼斯（美）R.索卡尔著，赵铁桥译.《数值分类学----数值分类的

原理和应用》.科学出版社.1984

（7）（加）Jiawei Han、Micheline Kamber著，范明 孟小峰等译.《数据挖掘概

念与技术》.机械工业出版社.2001

（8）方开泰 等著.《聚类分析》.北京 地质出版社.1982

（9）谭浩强 著.《C程序设计（第二版）》.清华大学出版社.2000

（10）许士良 主编.《C常用算法程序集（第二版）》.清华大学出版社.2001

（11）（美）赫伯特 希尔特著，王子恢等译.《C语言大全（第四版）》.电子工业

出版社.2001