降维分类方法:以改进质谱流式细胞识别技术

原文:Automatic classification of cellular expression by nonlinear stochastic embedding (ACCENSE)
作者:Karthik Shekhar, Petter Brodin, Mark M.Davis and Arup K.Chakraborty.
来源:Proceedings of the National Academy of Sciences (PNAS), 2014, 111(1): 202-207.
转载:简书不支持公式渲染,便于阅读可参考 原博文

更新

  • 2018.11.22:降维算法部分写得不太明撩,重新整理部分描述。

摘要

质谱流式细胞技术 ( Mass cytometry ) 能够在单细胞水平上识别到近 40 种不同的蛋白质,即提供前所未有的多维信息。由于各式各样的细胞种群数据集的复杂性,要收集有用的生物学知识对计算工具也有新的要求。回顾之前的聚类方法,聚类需要特征 ( 维度 ),而每一种蛋白质 ( 细胞表征由不同蛋白质表示 ) 可当作一种特征,聚类算法就会自动识别不同类型的簇群,即对于不同功能的细胞识别是基于细胞表征相似性来实现区分的。当然,经典方法存在一定局限性,例如单细胞分辨率的损失 ( 特征或维度的减少 );经典方法需要预知簇中的对象数量 ( 本文中指细胞亚群的规模数量 )。

而该论文引入 ACCENSE ( Automatic classification of cellular expression by nonlinear stochastic embedding ) 高维单细胞数据分析工具:

  • 基于密度划分的非线性降维方法,降维方法采用 t-Distributed Stochastic Neighbor Embedding (t-SNE) 算法 ^{[1]}
  • 该算法非常适合于探索性数据分析,同时避免任何手动 阀门(阈值) 的需要,即有别于基于距离的方法 ( 离群点判定 )、基于密度的方法 ( 密度阈值设定 )。
  • 化繁为简,在二维或三维图上展示不同功能的多元细胞群。

再有,本论文将 ACCENSE 应用于 35 参数的质谱流式细胞技术,检测 CD8+ T 细胞的数量 ( 数据来自于特定的无病原和无菌小鼠 ),并将细胞分层到表型亚群中。需要说明的是,对于具体的聚类算法、降维算法中,特定的符号名称会以具体的对象名称替代

正文

背景介绍

  • 免疫系统包含了许多类型的细胞,它们在免疫应答过程中表现出多样化的功能和复杂方式的相互作用,即通过不同蛋白质的表征所定义,故个体细胞的功能与其细胞表型密切相关。这里启示我们,对于不同功能的细胞可通过细胞表型相似性进行聚类区分。

  • 传统流式细胞技术和质谱流式细胞技术的区别:

    • 传统流式细胞技术 ( Flow Cytometry ) ^{[2]} 中,用 荧光基因 标记的抗体染色,其通过单细胞分辨率的光发射信号对靶标蛋白进行量化。且由于有限的光谱和重叠的发射信号,每个细胞限制为 12-16 个参数进行量化。

    • 质谱流式细胞技术 ( Mass Cytometry ) ^{[3]} ,使用 金属螯合探针 可对单个细胞多达 42 个参数的进行量化。

    • 传统流式细胞技术和质谱流式细胞技术相比,主要有两点不同:

      • 标签系统的不同,前者主要使用各种荧光基团作为抗体的标签,后者则使用各种金属元素作为标签;
      • 检测系统的不同,前者使用激光器和光电倍增管,而后者使用 ICP 质谱技术。

聚类算法

  • 质谱流式细胞技术产生的高维数据,以生物学的方式解释是具有挑战性的。然而,很多聚类工具是基于细胞的蛋白表达相似性进行细胞分类的,
  • 例如,SPADE 算法 ^{[4,5]} 使用多元信息定义细胞簇,并在树状结构中显示潜在的表型层次结构。但尚有不足之处:
    • 一是单细胞分辨率的损失;
    • 二是对目标集群数量的需要预知。

降维算法

关于 降维算法 在另外一篇博文也有提及,不妨参考学习:利用 t-SNE 降维并可视化数据

  • 同样,降维算法以细胞表征 ( 由不同蛋白质表示 ) 的相似性为依据,把空间组织的细胞群在低维空间上聚类成不同的细胞亚群。以下罗列一些常见的降维算法。

  • PCA 算法:PCA 降维的大致思想就是,挑选特征明显的、显得比较重要的信息保留下来。在本论文中,Newell 等人将主成分分析 ( Principal component analysis,PCA ) 应用于 25 参数的质谱流式细胞技术,检测人的 CD8+ T 细胞的数量,且使用前三种主成分 ( 3D-PCA ) 分离细胞亚群。3D-PCA 以三个汇总变量表示数据,每个汇总变量都是原始维度的 线性组合,并去捕获投影后数据的方差,直至其取值为最大值。然而,PCA 能在数据中所有的可能线性组合中找到最优表达,但也存在限制条件:线性投影可能太严格而不能产生精确的表示 ^{[6]},故作者引入 t-SNE 算法继续展开研究。

  • t-SNE 算法 ^{[7]}:t-Distributed Stochastic Neighbor Embedding,数据降维与可视化的方法,具体的算法细节如下:

    • \{x^{(i)}\} 表示归一化的 N 维蛋白质表达向量编码的细胞表型 i ( i=1, 2, ..., M )。

    • 若在 2D 平面图下,\{y^{(i)}\} 向量是高维向量 \{x^{(i)}\} 对应于低维的映射,它使得具有相似表型的 T 细胞彼此靠近嵌入,表型不相似的则嵌入相对较远的距离。

    • 采用细胞 i 和 j 之间的成对概率 \{p_{i,j}\} 表示 \{x^{(i)}\}\{x^{(j)}\} 之间的相似性。

    • 若在 2D 平面图下,成对概率 \{q_{i,j}\} 表示 \{y^{(i)}\}\{y^{(j)}\} 之间的相似性。

    • 通过最小化 \{p_{i,j}\}\{q_{i,j}\} 的 KL 散度 ( 可理解为代价函数 ),然后找出 "最佳" 嵌入向量 \{y^{(i)}\},即它表示的意义是,高维转低维的表示信息能最大程度被保存下来。

      K-L 散度 ( 详细见附录 1 ),Kullback-Leibler Divergence,又称相对熵,即描述两概率分布 P 和 Q 的差异。KL 散度公式 (1) 如下:

    D_{KL}(\{p_{i,j}\}|\{q_{i,j}\}) = \sum_{i,j} p_{i,j} log \frac{p_{i,j}}{q_{i,j}} \tag{1}

    • \{y^{(i)}\} 可以编码非线性关系,不像 PCA 中被约束为 \{x^{(i)}\} 的线性组合。
    • 最佳嵌入 是通过数值梯度下降法来确定的,即所有数据点的 KL 散度总和减小到最小 ( 详细见附录 2 )。

识别细胞亚群

  • 使用一个高斯核函数,把 t-SNE 的二维细胞散点图加工成 复合图像,如图 1-1 (D) 所示。其中,K_\gamma(y) 通过计算低维空间中所有细胞的位置总和,以表示 t-SNE 二维映射图中细胞的局部密度:

    K_\gamma(y) = (2 \pi \gamma^2)^{-1} \sum_{y' \in Y}exp( -\frac{||y - y'||^{2}}{2\gamma^2}) \tag{2}

  • 在本论文中,K_\gamma(y)局部最大值 表示具有共同表型的 CD8+ T 细胞亚群,且使用了 Matlab 的峰值检测算法识别这些局部最大值。

    当然,也可以在嵌入点上使用 K-Means 聚类算法来识别 T 细胞子集,但其要求事先指定簇的数量。

  • 如何求得 局部最大值,关键是对于公式 (2) 中 \gamma 的参数设定多少有关。即通过比较不同的核-带宽 \gamma 产生的结果,则存在一个 \gamma 值为表型空间中的局部和全局特征提供了准确的粗粒度表示。从图 1-2 中可得,即启示我们可以以数据驱动的方式,求得较合适的 \gamma 值,以近似地识别 CD8+ T 细胞的细胞亚群。

相关图表

  • 如图 1-1 所示,ACCENSE 应用于质谱高维数据。
图 1-1 ACCENSE 应用于质谱高维数据

(A) 质谱细胞计数数据集样本的图示。行对应于不同的细胞,而列对应于测量其表达 (细胞表面抗原和细胞内蛋白) 的不同标记的金属螯合抗体。每一元组对应于指示每个标记的表达水平的质荷比变换值 (反双曲函数)。(C) 来自SPF B6 小鼠的 CD8+ T 细胞的 2D t-SNE 图谱。每个点代表来自训练集的一个细胞 (M = 18304),且数据点是通过对原始数据集进行下采样得到。(D) 通过使用基于内核密度变换 (K_{\gamma}(y)\,{,}\,\gamma = 7),将细胞的局部概率密度嵌入 (C) 的复合图像。并使用标准的峰值检测算法进行识别局部最大值,在二维密度图表示表型亚群的中心。

  • 如图 1-2 所示,展示了峰值随着 \gamma 的增加而变化。
图 1-2 展示了峰值随着 $\gamma$ 的增加而变化

附录

1 t-SNE 中的概率

p_{i,j} 概率

  • 基于蛋白质相似性,设 p_{j|i} (i,j = 1, 2, ..., M) 表示细胞 i 将选择细胞 j 作为其最近邻的概率 ( p_{j|i} 越大,x^{(i)} 和 x^{(j)} 越近 ):

p_{j|i} = \frac{ exp({-||x^{(i)} - x^{(j)}||^2} / { 2\sigma_i^2}) }{ \sum_{k \neq i} exp({-||x^{(i)} - x^{(k)}||^2} / { 2\sigma_i^2}) }, d_{i,j} = ||x^{(i)} - x^{(j)}||_2 \tag{3}

  • 对于概率 p_{j|i} 的几点说明:

    • d_{i,j} 可以使用其他距离范式替代欧式距离范式;
    • 原始的 SNE 算法是不对称的,为简化梯度公式,t-SNE 中让公式 (3) 的条件概率是对称的。即初始化 p_{i|i} = 0 ( 只考虑不同点两两之间的相似度 ),对于任意的 p_{i|j} = p_{j|i},可得:

    p_{i,j} = \frac{ p_{j|i} + p_{i|j} }{2M} = \frac{ exp({-d_{i,j}^2} / { 2\sigma_i^2}) }{ \sum_{k \neq i} exp({-d_{i,k}^2} / { 2\sigma_i^2}) } \tag{4}

  • 不同的点 x_i,方差 \sigma_i 的取值也是不同的。

    • 公式 (3) 中的方差 \sigma_i 是确保对于每一个细胞都有相同的困惑度 ( Complexity )。复杂度可理解为一个点附近的 有效近邻点个数

    • 定义复杂度为 P_i = 2^{H_{j|i}},其近似地解释为细胞 i 的最近邻点的数量。

    • 定义 p_{j|i} 的香农熵 (信息熵) 为 H_{j|i} = - \sum_j p_{j|i} \log_2 p_{j|i},且 H_{j|i} 随着 \sigma_i 的增加而增加。

      在本论文中,t-SNE 图谱的复杂度被设定为 30,即 10-50 范围内的复杂度对最终结果的影响不大 (较好的鲁棒性)。

q_{i,j} 概率

  • 对于低维度下的 \{y_i\},在原始的 SNE 算法 ^{[7]} 中 Hinton 和 Rowers 引用高斯分布函数定义 q_{i,j},但在低维表达中发现了 拥挤问题

    拥挤问题:就是说各个簇聚集在一起,无法区分。譬如,有一高维度数据在降维到 10 维下可以有很好的表达,但是降维到两维后无法得到可信映射。具体情况是,10 维中有数个点之间两两等距离的,在二维下就无法得到可信的映射结果。

    进一步说明,假设一个以数据点 x^i 为中心,半径为 r 的 m 维球(三维空间就是球),其体积是按 r^m 增长的,假设数据点是在 m 维球中均匀分布的,我们来看看其他数据点与 x^i 的距离随维度增大而产生的变化。

  • t-SNE 减轻了拥挤问题,即使用更加偏重长尾分布的方式来将距离转换为概率分布 ^{[8]},故有 q_{i,j}

    q_{i,j} = \frac{ (1 + ||y^{(i)} - y^{(j)}||^2)^{-1} }{ \sum_{k \neq i} (1 + ||y^{(i)} - y^{(k)}||^2)^{-1} }, \Delta_{i,j} = ||y^{(i)} - y^{(j)}||^2 \tag{5}

  • 同样地,对于概率 q_{i,j} 的几点说明:

    • \Delta_{i,j} 可以使用其他距离范式替代欧式距离范式;
    • 原始的 SNE 算法是不对称的,为简化梯度公式,t-SNE 中让公式 (5) 的条件概率是对称的。即初始化 q_{i|i}=0,对于任意的 q_{i|j} = q_{j|i}

2 数值梯度下降法

  • 在 [7] 中的概述过程,获得优化的梯度公式,如下所示:

    \frac{ \partial D_{KL}(\{p_{i,j}\} | \{q_{i,j}\}) }{ \partial_{y_t}^{(i)} } = 4 \sum_j \frac{ (p_{i,j} - q_{i,j}) }{ (1 + ||y_t^{(i)} - y_t^{(j)}||^2) } (y_t^{(i)} - y_t^{(j)}) \tag{6}

  • 通过梯度下降法迭代计算局部最大值:

    y_{t+1}^{(i)} = y_{t}^{(i)} + \eta(t) \frac{ \partial D_{KL}(\{p_{i,j}\} | \{q_{i,j}\}) }{ \partial_{y_t}^{(i)} } + \alpha(t)(y_{t}^{(i)} - y_{t-1}^{(i)}) \tag{7}

    • y_t^{(i)} 表示迭代 t 次的解,\eta(t) 表示学习速率,\alpha(t) 表示迭代 t 次的动量。
    • 学习速率初始值为 \eta(t) = 100\,^{[9]},且动能量 \alpha(t) 设定为:

    \alpha(t) = \begin{cases} 0.8, & t < 300 \\ 0.5, & t \geq 300 \end{cases}

不足

  • t-SNE 主要用于可视化,很难用于其他目的。譬如测试集合降维,因为他没有显式的预估部分,不能在测试集合直接降维。
  • 关于核-带宽 \gamma 参数设定问题:文中展示了 \gamma 参数的大小与识别细胞亚群能力的数量关系。然而,数据驱动方式虽能实现自动聚类,但缺乏对于 \gamma 参数设定范围该如何控制的说明。

参考

  • [1] Maaten L, Hinton G. Visualizing data using t-SNE [J]. Journal of machine learning research, 2008, 9(Nov): 2579-2605.
  • [2] Cantor H, Simpson E, Sato V L, et al. And functional studies of peripheral t-cells binding different amounts of fluorescent anti-thy 1.2 (theta) Antibody using a fluorescence--activated cell sorter (FACS) [J]. 1975.
  • [3] Bendall S C, Nolan G P, Roederer M, et al. A deep profiler's guide to cytometry [J]. Trends in immunology, 2012, 33(7): 323-332.
  • [4] Qiu P, Simonds E F, Bendall S C, et al. Extracting a cellular hierarchy from high-dimensional cytometry data with SPADE [J]. Nature biotechnology, 2011, 29(10): 886.
  • [5] Bendall S C, Simonds E F, Qiu P, et al. Single-cell mass cytometry of differential immune and drug responses across a human hematopoietic continuum [J]. Science, 2011, 332(6030): 687-696.
  • [6] Van Der Maaten L, Postma E, Van den Herik J. Dimensionality reduction: a comparative [J]. J Mach Learn Res, 2009, 10: 66-71.
  • [7] Maaten L, Hinton G. Visualizing data using t-SNE [J]. Journal of machine learning research, 2008, 9(Nov): 2579-2605.
  • [8] Chrispher. t-SNE 完整笔记 [OL]. www.datakit.cn. 2017.
  • [9] Jacobs R A. Increased rates of convergence through learning rate adaptation[J]. Neural networks, 1988, 1(4): 295-307.
最后编辑于
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 148,827评论 1 317
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 63,511评论 1 266
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 99,318评论 0 218
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 42,108评论 0 189
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 50,112评论 1 266
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 39,387评论 1 185
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 30,905评论 2 283
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 29,657评论 0 177
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 33,144评论 0 223
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 29,744评论 2 225
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 31,100评论 1 236
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 27,565评论 2 222
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 32,041评论 3 216
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 25,769评论 0 9
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 26,301评论 0 178
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 34,173评论 2 239
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 34,313评论 2 242

推荐阅读更多精彩内容