基于数值对生物序列进行聚类

Question

我正在尝试根据 Atchley 因子（代表每个氨基酸的 5 个数字）将几个固定长度 (13) 的氨基酸序列聚类成 K 个簇。

例如，我有一个字符串输入向量，如下所示：

  key <- HDMD::AAMetric.Atchley


sequences <- sapply(1:10000, function(x) paste(sapply(1:13, function (X) sample(rownames(key), 1)), collapse = ""))

然而，我的实际序列列表超过 10^5（指定需要计算效率）。

然后我通过以下方式将这些序列转换为数字向量：

  key <- HDMD::AAMetric.Atchley

  m1 <- key[strsplit(paste(sequences, collapse = ""), "")[[1]], ]
  p = 13
  output <-
    do.call(cbind, lapply(1:p, function(i)
      m1[seq(i, nrow(m1), by = p), ]))

我想output（现在是 65 维向量）以一种有效的方式。

我最初使用的是 Mini-batch kmeans，但我发现重复时结果非常不一致。我需要一致的聚类方法。

我也担心维数灾难，考虑到在 65 维时，欧氏距离不起作用。

我看到的很多高维聚类算法都假设数据中存在离群值和噪声，但由于这些是生物序列转换为数值，所以没有噪声或离群值。

除此之外，特征选择将不起作用，因为每个氨基酸的每个属性和每个氨基酸在生物学背景下都是相关的。

您建议如何对这些向量进行聚类？

Answer 1

我认为自组织地图在这里可以提供帮助 - 至少实施速度非常快，所以您很快就会知道它是否有帮助：

使用操作中的数据以及：

rownames(output) <- 1:nrow(output)
colnames(output) <- make.names(colnames(output), unique = TRUE)

library(SOMbrero)

你提前定义簇数

fit <- trainSOM(x.data=output , dimension = c(5, 5), nb.save = 10, maxit = 2000, 
               scaling="none", radius.type = "gaussian")

nb.save 用作进一步探索训练在迭代过程中如何发展的中间步骤：

plot(fit, what ="energy")

似乎需要更多迭代

检查集群的频率：

table(my.som$clustering)
  1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25 
428 417 439 393 505 458 382 406 271 299 390 303 336 358 365 372 332 268 437 464 541 381 569 419 467

根据新数据预测聚类：

predict(my.som, output[1:20,])
#output
 1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 
19 12 11  8  9  1 11 13 14  5 18  2 22 21 23 22  4 14 24 12

检查哪些变量对聚类很重要：

summary(fit)
#part of output
Summary

      Class :  somRes 

      Self-Organizing Map object...
         online learning, type: numeric 
         5 x 5 grid with square topology
         neighbourhood type: gaussian 
         distance type: euclidean 

      Final energy     : 44.93509 
      Topographic error: 0.0053 

      ANOVA                : 

        Degrees of freedom :  24 

             F     pvalue significativity
pah      1.343 0.12156074                
pss      1.300 0.14868987                
ms      16.401 0.00000000             ***
cc       1.695 0.01827619               *
ec      17.853 0.00000000             ***

找到最佳簇数：

plot(superClass(fit))

fit1 <- superClass(fit, k = 4)

summary(fit1)
#part of output
SOM Super Classes
     Initial number of clusters :  25 
     Number of super clusters   :  4 


  Frequency table
1 2 3 4 
6 9 4 6 

  Clustering
 1  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 
 1  1  2  2  2  1  1  2  2  2  1  1  2  2  2  3  3  4  4  4  3  3  4  4  4 


  ANOVA

        Degrees of freedom :  3 

              F     pvalue significativity
pah       1.393 0.24277933                
pss       3.071 0.02664661               *
ms       19.007 0.00000000             ***
cc        2.906 0.03332672               *
ec       23.103 0.00000000             ***

更多内容在此 vignette

基于数值对生物序列进行聚类

Clustering biological sequences based on numeric values

r

cluster-analysis

bigdata