计算成对的汉明距离并保留较小的距离。

Question

我有一个矩阵：

A = matrix( c(1, 0, 0, 0, 1, 0, 0, 0, 1), nrow=3, ncol=3, byrow = TRUE)
dimnames(A) = list(c("Taxa1", "Taxa2", "Taxa3"), c("A1", "B1", "C1"))
df <- data.frame("ID" = c("A1", "B1", "C1"), "Triplicate" = c("T1", "T1", "T1"))

我希望计算一个分类单元内的A1与B1、A1与C1之间的汉明距离，并将两者中的最小值作为A1的指定值。然后计算B1与A1（我们已经完成了这一步）以及B1与C1之间的汉明距离，保留该值中的最小值作为B1的值，同样适用于C1。注意，样本A1、B1、C1是由T1标识的同一三重复样本的一部分。我会有类似的其他Triplicate样本T2、T3或T4，我想根据Triplicate列的值对样本进行分组，以计算成对的汉明距离。

最终的矩阵应该是：

df$Hamming <- c(2, 2, 2)

由于A1、B1与A1、C1之间的距离为2，因此保留值为2。

PS：汉明距离的简单描述（1分钟）： https://www.youtube.com/watch?v=P02mJhS9qQ4 添加我正在使用的确切数据： https://www.dropbox.com/s/wrlwmdipeyhbcok/Hamming.RData?dl=0

- Manasi Shah

---

您是想要每个分类单元中A1与B1、C1之间的距离，还是所有分类单元中的距离？能否提供您期望的结果呢？ - Lamia

矩阵中的数值表示二进制还是十进制数字？ - Sathish

https://dev59.com/Dl_Va4cB1Zd3GeqPVahu - Sathish

@Lamia添加了结果。我想要所有分类单元的距离总和。如果您有任何问题，请告诉我。@Sathish如果可以的话，矩阵可以是二进制的。您分享的链接没有按照我的triplicate变量所需的方式分组答案。 - Manasi Shah

Answer 1

通过使用 xor() 函数可以识别出位翻转，而从 xor() 函数的结果中累加可得到位翻转的总数。我没有优化 hamm_dist_min() 中的代码。

xor(0,0)
# [1] FALSE
xor(1,1)
# [1] FALSE
xor(0,1)
# [1] TRUE
xor(1,0)
# [1] TRUE

根据OP的要求，计算汉明距离的两个方向。例如：AB、BA、AC、CA、BC和CB，这些形成了三元组ABC、BCA和CAB。如果您只想要一个方向，例如：AB、AC和BC，您可以在hamm_dist_min()函数中使用combn()函数设置列号。

数据：

A = matrix( c(1, 0, 0, 0, 1, 0, 0, 0, 1), nrow=3, ncol=3, byrow = TRUE)
dimnames(A) = list(c("Taxa1", "Taxa2", "Taxa3"), c("A1", "B1", "C1"))
df <- data.frame("ID" = c("A1", "B1", "C1"), "Triplicate" = c("T1", "T1", "T1"))

海明距离

# minimum of hamming distance
hamm_dist_min <- function(data)
{
  # setup combinations of column numbers
  n_col <- ncol(data)
  x <- expand.grid(seq_len(n_col), seq_len(n_col))
  x <- x[ x[, 1] != x[, 2], ]
  x <- x[order(x[, 1]), ]
  x <- split(x, cut(x[, 1], breaks = c(0, seq_len(n_col)), labels = colnames(data) ))

  # minimum of hamming distance
  h_d <- unlist(lapply(x, function(y){
    min( colSums(apply(y, 1, function(z) xor(data[, z[1]], data[, z[2]]))))
  }))
  return(h_d)
}

hamm_dist_min(data = A)
# A1 B1 C1 
# 2  2  2

df$Hamming <- hamm_dist_min(data = A)
df
#   ID Triplicate Hamming
# 1 A1         T1       2
# 2 B1         T1       2
# 3 C1         T1       2

Youtube示例：

df1 <- matrix( c(0,1,1,1,1,1,0,0), ncol = 2, byrow = FALSE)
colnames(df1) <- LETTERS[1:2]
hamm_dist_min(data = df1)
# A B 
# 3 3 

---

编辑： 基于问题中新增的数据集。

注意：如果一个样本类型只有一列，则将0作为汉明距离，因为我们需要至少2列来计算汉明距离。看一下df的三个副本列中的T71。您可以返回NA，表示不可用，而不是值0。

load("Hamming.RData")

# setup unique colnames pattern
col_list <- unique(unlist(lapply( colnames(A), function(x){
  substr(x = x, start = 1, stop = nchar(x) - 1)
} )
))

# get hamming distance
my_results <- lapply( col_list, function(x){
  cols_x <- grep(x, colnames(A) )
  if(length(cols_x) == 1 ){  # return 0 for one column
    return( setNames( object = rep( 0, length(cols_x)), nm = colnames(A)[cols_x]))
  } else{ # return minimum of hamming distance
    return(hamm_dist_min(data = A[, cols_x]))
  }
})

# get triplicate id
triplicate <- paste0( "T", rep(seq_along(my_results), 
                               lengths(my_results)))

# final data
my_results <- unlist(my_results)
df <- data.frame( SampleID = names( my_results ),
                  Hamming = my_results,
                  Triplicate = triplicate,
                  stringsAsFactors = FALSE )

head(df)
#                        SampleID Hamming Triplicate
# Affy22_MDA_1       Affy22_MDA_1       2         T1
# Affy22_MDA_2       Affy22_MDA_2       2         T1
# Affy22_MDA_3       Affy22_MDA_3       3         T1
# GutRef001_MDA_1 GutRef001_MDA_1       4         T2
# GutRef001_MDA_2 GutRef001_MDA_2       4         T2
# GutRef001_MDA_3 GutRef001_MDA_3       6         T2

- Sathish

---

你好！我喜欢这个答案，它给了我一个每个样本的汉明距离。但是对于任何样本，该解决方案没有给出相对于什么的汉明距离？它与矩阵中的所有列进行比较。然而，我想限制它计算在df中提到的“三重复制品”样本之间的最小汉明距离，即我想要A1、B1、C1属于同一三重复制品T1的最小汉明距离。如果我有另一组属于T2的样本A1、B1、C1，则应分别计算这些样本的最小值。有什么办法可以根据“三重复制品”变量分组并计算您的答案呢？谢谢！ - Manasi Shah

你可以将选定的数据传递给函数。例如：hamm_dist_min(data = A[, 1:3])。您可以灵活地更改列数。要查找列名，请使用colnames(A)。 - Sathish

谢谢，我有85个三重复体，所以我想知道是否有办法将该函数应用于三重复体向量（它按Triplicate号码进行分组），然后将结果rbind在一起？ - Manasi Shah

创建一个包含列名的三重复列表 col_list <- list(c("A1", "B1", "C1"), c("A2", "B2", "C2"), c("A3", "B3", "C3"), c("A4", "B4", "C4")，并使用 lapply(col_list, function(x) {hamm_dist_min(data = df[, x]) } ) 获得结果，并从中 rbind do.call('rbind', myresults)。我认为使用 sapply 可以避免 rbind 步骤，但我没有测试过。 - Sathish

嗨，我已经在问题中添加了我正在使用的原始数据。如果您可以提供一些关于如何按三重组分结果的建议，那将非常有帮助！谢谢！ - Manasi Shah

如果您只想严格使用三重样本，则将条件从“length(cols_x) == 1”更改为“length(cols_x)！= 3”，并返回“0”或“NA”。 - Sathish

Answer 2

以下是计算每对列的汉明距离的方法。我不确定您想如何处理对列之间距离相同的情况。在这里，我只是从给定参考列的所有汉明距离相同的列中按字母顺序选择第一列。

library(e1071)  # For the hamming.distance function
library(tidyverse)

# Get Hamming distance for all pairs of columns in matrix A. 
hd = combn(colnames(A), 2, simplify=FALSE) %>% 
  map_df(function(col) data.frame(col1=col[1], col2=col[2], 
                                  Hamming=hamming.distance(A[,col[1]], A[, col[2]]))) %>% 
  # For a given column, keep only the shortest Hamming distance
  group_by(col1) %>% 
  arrange(Hamming, col2) %>% 
  slice(1) %>%
  ungroup %>% 
  # Add a column to mark which Hamming distance pair we kept
  mutate(pair = paste0(col1, "_", col2))

hd

   col1  col2 Hamming  pair
1    A1    B1       2 A1_B1
2    B1    C1       2 B1_C1

现在将汉明距离值与相应的ID连接起来。首先，我们从hd中堆叠col1和col2并删除重复的ID值，然后将结果连接到df。

df = df %>% left_join(
  bind_rows(hd %>% select(col1, Hamming, pair),
            hd %>% select(col1=col2, Hamming, pair)) %>% 
    filter(!duplicated(col1)),
  by=c("ID"="col1")
)

  ID Triplicate Hamming  pair
1 A1         T1       2 A1_B1
2 B1         T1       2 B1_C1
3 C1         T1       2 B1_C1

- eipi10

---

你好，你的答案对我有用，除了hd中的col1和col2应该属于同一个“Triplicate”。现在我正在从不同的“Triplicate”中获取汉明距离对，因为那是最小距离。我尝试在你的代码中添加这个 group_by(c(df$Triplicate, col1)) %>% 但它没有起作用。有什么建议吗？ - Manasi Shah

你能添加一些示例数据（矩阵A和数据框df），以表示当有多个三次重复值时数据是什么样子吗？ - eipi10

嗨，我刚刚在问题中添加了指向我正在使用的确切“A”和“df”数据的链接。请告诉我您的建议！谢谢！ - Manasi Shah