根据《关于“k-means算法在流式细胞仪中细胞分类的应用”的学习笔记总结》撰写的中期报告

XXXX大学2014届本科毕业设计（论文）中期报告

毕业设计（论文）题目：K-means算法在流式细胞仪中细胞分类的应用

专业（方向）：生物医学工程

学生信息：XXXXXX、XX、生医XXX

指导教师信息：XXXXXX、XX、讲师

报告提交日期：2014年4月18日

 

---

 

下面是在百度网盘中中期报告和实验数据的分享地址，可以直接下载。

http://pan.baidu.com/s/1bnzSer9 

http://pan.baidu.com/s/1hqomDq0

本篇文章中的所有程序、仿真图系亲自在MATLAB2009a得到，方法简单易实践，参考文献也全部给出，而且实验数据资料也自由共享。读者有兴趣可完全照着学习验证！

---

 

一、课题简介

随着信息技术和计算机技术的迅猛发展，人们面临着越来越多的文本、图像、视频以及音频数据，各行各业的信息系统数年来收集了海量数据，而且这些数据还在以指数级增长着。这给我们带来方便的同时，也给我获取有效信息加大了难度。过多无用的信息必然会产生信息距离和有用信息的丢失[1]。为帮助用户从这些大量数据中分析出其间所蕴涵的有价值的知识，数据挖掘（Data Mining，DM）技术应运而生[2]。它能开采出潜在的模式，找出最有价值的信息，指导商业行为，生产控制，或辅助科学研究[3]。数据挖掘技术出现于20世纪80年代后期，在90年代发展迅速，其理论与应用研究已成为国内外的研究热点[4-5]。与此同时，聚类作为数据挖掘的主要方法之一，也越来越引起人们的关注。

迄今为止，聚类还没有一个学术界公认的定义。事实上，聚类是一个无监督的分类，它没有任何先验知识可用。没有任何一种聚类技术可以普遍适用于揭示各种多维数据集所呈现出来的多种多样的结构[6]。聚类是一个将数据集划分为若干组或簇的过程，使得同一类的数据对象之间的相似度较高，而不同类的数据对象之间的相似度较低。聚类问题的关键是把相似的事物聚集在一起，适合探讨样本间的内部关系[3]。1967年，MacQueen[7]首次提出了K均值聚类算法（K-means算法）。该方法虽然不能用于类别属性的数据，但对于数值属性的数据，它能很好地体现聚类在几何和统计学上的意义。在对大型数据集聚类时，K-means算法也比层次聚类算法快得多。迄今为止，很多聚类任务都选择该经典算法或其变种。

流式细胞技术（FCM）是70年代发展起来得一种快速对单细胞定量分析的新技术。它借鉴了荧光显微镜技术，同时利用了荧光染料、激光技术、单抗技术以及计算机技术的发展，将荧光显微镜的激发光源改为激光，使之具有更好的单色性与激发效率，因而大大提高了检测灵敏度，同时将固定的标本台改为流动的单细胞悬液，用计算机进行数据处理，因而大大提高了检测速度与统计精确性，而且能从同一个细胞中可以同时测得多种参数[8]。流式细胞仪是采用流式细胞术的典型仪器，它能迅速的对单个细胞及其群体的化学物质的含量与种类作出分析，对含指定化学物质的细胞进行分离提纯，在一些有关实验室和医院应用广泛[9]。近年来FCM不仅在生物医学与临床检验学中得到了长足发展，而且拓展到生物学的各个领域，微生物学中涉及医学、发酵和环保等的诸多领域，而且它还是细胞学研究中必不可少的工具之一。

在流式细胞仪中使用细胞标志示踪物质为荧光标记物鉴定和分离细胞群和亚群是其非常重要的一个应用，具有深远的意义。例如在免疫学，淋巴细胞是正常机体免疫系统功能最重要的一大细胞群，在免疫应答过程中，未梢血淋巴细胞发育成为功能不同的亚群。各亚群的数量和功能产生异常时，就能导致机体免疫紊乱并产生病理变化[10]。淋巴细胞亚群分析已成为检测细胞免疫和体液免疫功能的重要指标，而FCM可以同时检测一种或几种淋巴细胞膜表面分子（分化群抗原），将不同的淋巴细胞亚群数量的测定来监控病人的免疫状态，指导治疗[11]。对于细胞亚群百分含量的分析，我们一般用单参数直方图、二维散点图、密度图、等高图和三维图来分析结果。而用于显示数据的几种形式的散点图，还可通过设门的方法区分指定的细胞亚群[8-9]。流式细胞术对淋巴细胞亚群的这种鉴定分析已成为医院中最快速、准确和可靠常规检测技术。

现在，数据聚类分析已经成为一个非常活跃的研究课题，但是用聚类的方法对流式细胞仪中检测到的细胞分类数据进行分析还很少有人研究。本课题就是以标记荧光团的标准微球（淋巴细胞）为实验对象，由流式细胞仪得到各细胞的荧光光谱，然后基于此荧光光谱利用MATLAB编程分析，用K-means聚类算法实现不同特性细胞的准确分类。

二、研究方法

本课题研究K-means算法在流式细胞仪中细胞分类的应用，关键是用该算法挖掘出荧光光谱中的有效信息并建立起荧光光谱与细胞分类的正确映射关系。本课题将基于采用FCM对荧光标记的淋巴细胞进行检测而得到的CD3-8-45-4分子的荧光光谱数据，以K-means算法构建聚类模型，并用MATLAB编程仿真实现细胞亚群分类。然后依据仿真结果是否聚类良好，不断改进该算法，直到取得比较满意的结果。

流式细胞术（FCM）是一种快速检测分析单个粒子多物理特性的高技术，通常指细胞通过激光束时在液流中的特性，即粒子的的大小，密度，或是内部结构，以及相对的荧光强度。通过光电系统记录细胞的散射光信号和荧光信号可以得知细胞特性。FCM分析中前向散射光反映颗粒的大小；侧向散射光反映颗粒内部结构复杂程度、表面的光滑程度；荧光反映颗粒被染上荧光部分数量的多少，根据其标记的抗原分子不同，即反映了不同抗原分子的表达情况[12]。下图1为流式细胞仪结构示意图。本课题主要研究荧光标记的淋巴细胞的荧光光谱来实现细胞亚群分类。

在流式细胞仪中，细胞被传送到液流中的激光照射区。任何存在于悬液中的直径为0.2-150微米的粒子或细胞都适用于流式分析。在实际工作中，用实体组织进行流式细胞分析往往是不可能的，分析之前必须对其进行分解。被液滴包绕的粒子称为细胞液柱，当粒子经过激光照射区时，通过激光激发产生散射光。含有荧光的粒子就会表现出其荧光特性。散射光和荧光由光路系统（相应的透镜，滤片和探测器）收集。分光器和滤光片引导散射光和荧光至相应的探测器，把光信号转换为电信号。单个粒子通过其表现出的光散射和荧光属性，通过列表模式（List mode）完成数据采集，并对样本中的细胞亚群进行分析[13]。图2为其工作的基本流程。

 

图1 流式细胞仪结构示意图

荧光物质吸收符合其波长范围的光能量，内部电子受激上升到高能级，然后受激电子迅速衰落回基态，释放过剩能量成为光子，这种能量转换生成的光称为荧光。能够激发荧光物质的波长范围称为激发光谱。荧光物质的发射波长范围叫做发射光谱。光源的谱线愈接近被激发物质的激发光谱的峰值，所产生的荧光信号愈强。而且因为更多的能量消耗在吸收转换而不是荧光转换中，所以发射光波长要高于激发光波长。目前的流式细胞仪大多采用488nm氩离子激光器，原因是488nm的激光器能够激发一种以上的荧光。本课题FCM实验采用的就是488nm的氩离子激光器，检测到三个通道，分别是FL1:53610nm，FL2:59010nm，FL3:62010nm。在FCM对单克隆抗体进行检测时，抗体上标记荧光，然后抗体与被测抗原特异结合，结合后的抗原抗体复合物通过流式细胞仪时就能检测到特异荧光。如果选择不同的单抗及染料，就可同时测定一个细胞上的多个不同特征。每种荧光染料会产生特定波长的荧光和颜色，通过波长选择通透性滤片，可将不同波长的散射光和荧光信号区分开，送入不同的光电倍增管[14]（见图1）。当细胞标记了多色荧光，被激光激发后，得到的荧光信号还可根据需要进行组合分析。本课题即使用该方法对淋巴细胞表面抗原进行检测。

 

图2 流式细胞仪基本流程

人体内的淋巴细胞并不是功能单一的群体，但其在光学显微镜下的形态基本上是一样的，都是以致密的核及少量细胞质组成的小而圆的细胞为特征，因此要对其进一步的分类观察就不能采用形态学的方法，而是要对其表面标志的检测。包括对T细胞、B细胞和NK细胞及其有关亚群的检测，实质上是对其相应的表面标志进行检测，据此可以建立起相应的细胞计数方法，以判断机体的免疫水平[15]。常用于鉴定和检测淋巴细胞表面标志的是分化抗原簇（cluster of defferentiation，CD）。CD抗原的鉴定和检测依赖与其相应的单克隆抗体，国际制定了细胞分化抗原的CD统一编号，统一把识别同一抗原的两种以上的单抗划分为同一抗体组。T淋巴细胞表达CD3；B细胞表达CD19和CD20；NK细胞表达CD56和CD16，不表达CD3。利用各种单克隆抗体与淋巴细胞表面抗原结合，再配多色荧光染料，即可以把淋巴细胞区分为各亚群。本课题即检测淋巴细胞表面CD3-8-45-4四中分子的特征以划分细胞亚群。

K-means算法（K均值聚类算法）是一种动态聚类方法，由James MacQueen于1976年提出，但算法本身思想就已经由Stuart Lloyd于1957年给出，是一种最广泛使用的基于划分的聚类算法。它把n个对象划分为k个簇(cluster)，使簇内相似度较高，簇间相似度较低，抑或使簇间相异度较高，簇内相异度较低[16]。它与处理混合正态分布的最大期望算法（Expectation-maximization algorithm，EM）很相似，因为他们都试图找到数据中自然聚类的中心。但前者倾向于寻找空间内容易区分的簇，而后者则允许该簇拥有各种形态。

K均值聚类属于无监督学习，样本集中没有给定类别，只有特征。而监督学习中，训练样本集却已经给出了类别标签，标签信息可以协助机器学习同类样本之间存在的共性，在预测时只需判定给定样本与哪个类别的训练样本最相似即可，典型的有回归分析和分类[17]。K-means算法只能通过学习寻找样本集数据所体现的规律性，即计算出样本数据间的相似度或相异度，将样本集划分为k个簇。

给定一系列N维向量样本集和聚类参数k，运用K-means算法将这些向量聚为k类。其基本步骤[18]是

（1）任意选择k个点作为初始聚类质心；

（2）将每个点指派到最近的质心，形成k个簇；

（3）重新计算每个簇的质心；

（4）重复（2）（3）直到质心不再发生变化。

算法很简单，但很多复杂算法的实际效果也许还不如它，而且它的局部性较好，容易并行化，对大规模数据集很有意义。该算法的时间复杂度为O(tkmn)，空间复杂度为O((m+k)n)，其中，t为迭代次数，k为簇数，m为点数，n为维数。MATLAB工具箱中的kmeans函数就能直接实现这些步骤[19]。本课题拟直接调用该函数。

在k-means算法中，误差的平法和（Sum of the Squared Error，SSE）经常作为度量聚类质量的准则函数[20]。对于任何一个样本点，其准则函数为

 

其中，dist函数是样本空间中样本点x与其对应的簇质心Ci的距离，Ci取第i个簇所有样本点的均值。准则函数计算的是所有样本点与其对应的簇质心的距离平方和，使准则函数最小的簇划分则为最优的聚类。

在k-means算法中，相异度（dissimilarity）或相似度（similarity）常被用来描述“最近”这个概念，采用的方法就是计算样本间的“距离”(Distance)。采用什么样的方法计算距离是很讲究，甚至关系到分类的正确与否。下面是常用的是三种距离的定义公式[21-23]。

①　欧几里得距离（Euclidean Distance）

 

②　曼哈顿距离（Manhattan Distance）

 

③　余弦相似度（Cosine Similarity）

 

其中，X，Y是样本集中任意两个元素项，各自具有n个可度量特征属性。前两种计算样本间相异度，后一种计算相似度。当采用曼哈顿距离时，质心的计算改求均值为求中位数。

欧几里得距离是最容易理解的一种距离计算方法，源自于欧氏空间中两点间的距离公式。曼哈顿距离则是由计算在曼哈顿从一个十字路口开车到另外一个十字路口的驾驶距离而来的，它考虑了分布的形状，因为来源如此，也被称为城市街区距离(City Block distance)[24]。而余弦相似度是用两向量的余弦来度量的一种距离计算方法。本课题中则拟采用曼哈顿距离计算其相异度。

在kmeans算法中取k个初始聚类质心进行划分时，参数k是事先给定的，初始质心是随机确定，一旦初始值选择得不好，可能无法得到较好的聚类结果。本课题拟直接调用subclust函数，设置参数，采用减法聚类法[25]确定聚类参数k的值，再以该值用kmeans函数聚类划分。减法聚类法是将每个数据点都作为可能的聚类中心，并根据各个数据点周围的数据点密度来计算该点作为聚类中心的可能性。被选为聚类中心的数据点周围具有最高的数据点密度，同时该数据点附近的数据点被排除作为聚类中心的可能性。在选出第一个聚类中心后，继续采用这样的方式对剩余的可能作为聚类中心的数据点选择下一个聚类中心。

三、阶段性成果（仿真实验与结果）

本课题实验拟通过使用减法聚类发获得可能的聚类参数k，然后应用k-means算法将四维数据样本集划分为k类来获得分类结果。由于这些数据是四维的，很难通过可视化辨析分类结果的合理程度，本实验拟使用kmeans函数的分类指数调用silhouette函数绘出轮廓图，来演示类中的每个点与相邻类的中点的接近程度以分析分类是否合理。获得较好的分类结果后，随机5次重复运行，通过参数replicates选出SSE全局最小的最少迭代次数，确保最优初始聚类中心。

1．减法聚类

装载数据，直接调用subclust函数。其调用格式为

[C,S]=subclust(X,radii,Xbounds,options)

设置参数radii为0.5，Xbounds为默认值，options为[1.25 0.5 0.15 0]。其中，函数的返回值C为聚类中心向量，它的每一行代表一个聚类中心的位置。向量S包含了聚类中心在每一维坐标上的影响范围，所有聚类中心在同一方向上具有相同的影响范围。

程序为

[C,S]=subclust(A,0.5,[],[1.25 0.5 0.15 0])

运行后

C =

   75.0000   67.0000   50.0000  197.0000

   91.0000   71.0000  479.0000  272.0000

  563.0000   78.0000  541.0000  434.0000

S =

  128.8702  132.9361  136.6484  123.0366

显然，减法聚类法将样本集分成3类。

2．K-means聚类

装载数据，设定分类参数k为3，采用曼哈顿距离计算方式直接调用kmeans函数聚类，并用可选参数display打印每次迭代信息，使用silhouette函数绘制轮廓图。

程序

idx3=kmeans(A,3,'dist','city','display','iter');

[silh3,h]=silhouette(A,idx3,'city');

迭代信息

 iter  phase      num          sum

     1      1    15608 4.64974e+006

     2      1     1304 4.55743e+006

     3      1      673 4.53688e+006

     4      1      415 4.52699e+006

     5      1      356 4.51966e+006

     6      1      436 4.51107e+006

     7      1      509 4.49984e+006

     8      1      684 4.47984e+006

     9      1     1691 4.34833e+006

    10      1     3142 3.62312e+006

    11      1     1104 3.43982e+006

    12      1       96 3.43766e+006

    13      1       13 3.43762e+006

    14      2        0 3.43762e+006

14 iterations, total sum of distances = 3.43762e+006

kmeans函数重新指定类间的点并重新计算类的质心时，距离的总和会随迭代次数的增加逐渐减小。本次聚类迭代过程中，算法设置在第二阶段没有任何改变，表示经过13次迭代后，第一阶段已经达到了最小。

轮廓图如下

 

图3 三分类的轮廓图

silouette值的取值范围是[-1,1]，样本点的silouette值越大，表明该点的分类越合理。当silouette值为负值时，表明该点的分类不合理，还有比目前更好的方案。从本次聚类的轮廓图中可以看出，第2类中的大部分点的silouette值比较大（大于0.6），表示该类在某种程度上区分于相邻的类。第1、3类中虽然仍包含了少数silouette值较低的点，但大部分点的silouette值大于0.4，区分度也还可以。不过，第1类中有几个点具有负值，表示这些点的分类不是很好。

随机5次重复运行，通过参数replicates选出全局最小SSE，确保最优初始质心。

程序为

[idx3,cent3,sumdist]=kmeans(A,3,'dist','city','display','final','replicates',5);

运行后

11 iterations, total sum of distances = 4.25664e+006

19 iterations, total sum of distances = 4.25991e+006

4 iterations, total sum of distances = 3.43762e+006

12 iterations, total sum of distances = 4.27462e+006

11 iterations, total sum of distances = 3.43762e+006

显然保证SSE全局最小时，最少的迭代次数为4次。这时的初始聚类中心及其与类内各点的距离和分别是

cent3 =

    98    78   470   277

   552    97   553   394

    78    78    54   205

sumdist =

      894441

      957646

     1585528

可以看到，这里的最优初始聚类中心与减法聚类得到的初始聚类中心是不一样的，但差别很小。

四、下一步的工作及研究思路

聚类完成后，大部分样本点的分类还是比较合理，但仍有极少数点显然没有划分到最合适的类。猜测很可能是FCM检测数据的过程中产生了噪声和异常值，导致样本集中出现了孤立点（outlier），不能产生很好的分类。本课题还需要用样本数据预处理，去除噪声，并借以检测出孤立点，删除它们，消除对结果产生的不良影响。对此，可以采用基于统计学的方法。

 

参考文献

1．王立伟.数据挖掘研究现状综述.图书与情报，2008(5)：41~46

2．王慧中，彭安群.数据挖掘研究现状及发展趋势.工矿自动化，2011(2)：29~32

3．王桂芹，黄道.数据挖掘技术综述.电脑应用技术，2007(69)：9~14

4．李菁菁，邵培基，黄亦潇.数据挖掘在中国的现状和发展研究.管理工程学报，2004(3)：10~15

5．Tipawan Silwattananusarn，KulthidaTuamsuk.Data Mining and Its Applications for Knowledge Management：A Literature Review from 2007 to 2012.International Journal of Data Mining & Knowledge Management Process，2012(9)：13~24

6．孙吉贵，刘杰，赵连宇.聚类算法研究.软件学报，2008(1)：48~61

7．J. MacQueen.Some methods for classification and analysis of multivariate observations.Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability，1967(1)：281~297

8．宋平根，李素文.流式细胞术的原理和应用.北京：北京师范大学出版社，1992

9．M. R. Melamed，T. Lindmo，M． L． Mendelsohn．Flow Cytometry and Sorting．New York：John Wiley and Sons，1990

10．曹雪涛.医学免疫学.北京：人民卫生出版社，2013

11．张瑞华，王进，兰文君.流式细胞术及其在生物医学领域中的应用.山东轻工业学院学报，2011(8)：20~24

12．张艺.流式细胞仪构成与工作原理.医疗设备信息，2005(8)：25~26

13．J. Paul Robinson.Flow Cytometry.Encyclopedia of Biomaterials and Biomedical Engineering，2004：630~640

14．瑞菲尔·努纳兹.流式细胞术原理与科研应用简明手册.北京：化学工业出版社，2005

15．王兰兰，吴建民.临床免疫学与检验.北京：人民卫生出版社，2007

16．http://en.wikipedia.org/wiki/K-means_clustering

17．邵峰晶，于忠清.数据挖掘原理与算法.北京：中国水利水电出版社，2003

18．韩家炜，坎布尔.数据挖掘：概念与技术.北京：机械工业出版社，2001

19．张德丰等.MATLAB概率与数理统计分析.北京：机械工业出版社，2010

20．Pang-Ning Tan，Michael Steinbach，Vipin Kumar.数据挖掘导论.北京：人民邮电出版社，2011

21．李雄飞，李军.数据挖掘与知识发现.北京：高等教育出版社，2003

22．陈燕.数据挖掘技术与应用.北京：清华大学出版社，2011

23．杨杰，姚莉秀.数据挖掘技术及其应用.上海：上海交通大学出版社，2011

24．http://zh.wikipedia.org/wiki/曼哈頓距離

25．周品.MATLAB概率与数理统计.北京：清华大学出版社，2012