Вычисление попарного расстояния Хэмминга между всеми строками двух целых матриц/кадров данных

Question

У меня есть два кадра данных: df1 с справочными данными и df2 с новыми данными. Для каждой строки в df2 мне нужно найти наилучшую (и вторую) подходящую строку до df1 с точки зрения расстояния от хамминга.

Я использовал e1071 пакет для расчета расстояния от помех. Расстояние Хэмминга между двумя векторами x и y может быть вычислена как, например:

x <- c(356739, 324074, 904133, 1025460, 433677, 110525, 576942, 526518, 299386, 
     92497, 977385, 27563, 429551, 307757, 267970, 181157, 3796, 679012, 711274, 
     24197, 610187, 402471, 157122, 866381, 582868, 878) 

y <- c(356739, 324042, 904133, 959893, 433677, 110269, 576942, 2230, 267130, 
     92496, 960747, 28587, 429551, 438825, 267970, 181157, 36564, 677220, 
     711274, 24485, 610187, 404519, 157122, 866413, 718036, 876) 

xm <- sapply(x, intToBits) 
ym <- sapply(y, intToBits) 

distance <- sum(sapply(1:ncol(xm), function(i) hamming.distance(xm[,i], ym[,i]))) 

и полученное расстояние 25. Тем не менее, мне нужно сделать это для всех рядов df1 и df2. Тривиальный метод принимает двойное петлевое гнездо и выглядит ужасно медленным.

Любые идеи, как это сделать более эффективно? В конце концов, мне нужно добавить к df2:

• столбец с идентификатором строки из df1, что дает самое низкое расстояние;
• столбец с самым низким расстоянием;
• столбец с номером строки от df1, который дает 2-е минимальное расстояние;
• колонка со вторым самым низким расстоянием.

Спасибо.

Answer 1

Быстрое вычисление расстояния Хэмминга между двумя целыми числами векторов одинаковой длины

Как я сказал в своем комментарии, мы можем сделать:

hmd0 <- function(x,y) sum(as.logical(xor(intToBits(x),intToBits(y)))) 

вычислить расстояние Хэмминга между двух целых векторов единая длинаx и y. Это использует только базу R, но более эффективно, чем e1071::hamming.distance, , потому что это векторизация!

Для примера x и y в вашем посте, это дает 25. (мой другой ответ покажет, что мы должны делать, если мы хотим, чтобы парный расстояние Хэмминга.)

---

Fast Хэмминга расстояние между матрицей и вектором

Если мы хотим вычислить расстояние от помех между одним y и несколькими x s, т. е. хамом расстояние между вектором и матрицей, мы можем использовать следующую функцию.

hmd <- function(x,y) { 
    rawx <- intToBits(x) 
    rawy <- intToBits(y) 
    nx <- length(rawx) 
    ny <- length(rawy) 
    if (nx == ny) { 
    ## quick return 
    return (sum(as.logical(xor(rawx,rawy)))) 
    } else if (nx < ny) { 
    ## pivoting 
    tmp <- rawx; rawx <- rawy; rawy <- tmp 
    tmp <- nx; nx <- ny; ny <- tmp 
    } 
    if (nx %% ny) stop("unconformable length!") else { 
    nc <- nx/ny ## number of cycles 
    return(unname(tapply(as.logical(xor(rawx,rawy)), rep(1:nc, each=ny), sum))) 
    } 
    } 

Обратите внимание, что:

1. hmd выполняет вычисление столбцов. Он сконструирован так, чтобы быть Кэш-память процессора.Таким образом, если мы хотим сделать некоторые вычисления по строкам, мы должны перенести матрицу сначала;
2. здесь нет очевидной петли; вместо этого мы используем tapply().

---

Быстрого Хэмминг расстояние вычисление между двумя матрицами/кадрами данных

Это то, что вы хотите. Следующая функция foo принимает два кадра данных или матрицы df1 и df2, вычисляя расстояние между df1 и каждой строкой df2. аргумент p - целое число, показывающее, сколько результатов вы хотите сохранить. p = 3 сохранит наименьшие 3 расстояния с идентификаторами строк в df1.

foo <- function(df1, df2, p) { 
    ## check p 
    if (p > nrow(df2)) p <- nrow(df2) 
    ## transpose for CPU cache friendly code 
    xt <- t(as.matrix(df1)) 
    yt <- t(as.matrix(df2)) 
    ## after transpose, we compute hamming distance column by column 
    ## a for loop is decent; no performance gain from apply family 
    n <- ncol(yt) 
    id <- integer(n * p) 
    d <- numeric(n * p) 
    k <- 1:p 
    for (i in 1:n) { 
    distance <- hmd(xt, yt[,i]) 
    minp <- order(distance)[1:p] 
    id[k] <- minp 
    d[k] <- distance[minp] 
    k <- k + p 
    } 
    ## recode "id" and "d" into data frame and return 
    id <- as.data.frame(matrix(id, ncol = p, byrow = TRUE)) 
    colnames(id) <- paste0("min.", 1:p) 
    d <- as.data.frame(matrix(d, ncol = p, byrow = TRUE)) 
    colnames(d) <- paste0("mindist.", 1:p) 
    list(id = id, d = d) 
    } 

Обратите внимание, что:

1. перестановка осуществляется в начале, по причинам прежде;
2. a for цикл используется здесь. Но это действительно эффективно, потому что на каждой итерации проводятся значительные вычисления. Это также более элегантно, чем использование семьи *apply, так как мы запрашиваем несколько выходных данных (строка id id и расстояние d).

---

Эксперимент

Эта часть использует небольшой набор данных для тестирования/продемонстрировать наши функции.

Некоторые игрушечные данные:

set.seed(0) 
df1 <- as.data.frame(matrix(sample(1:10), ncol = 2)) ## 5 rows 2 cols 
df2 <- as.data.frame(matrix(sample(1:6), ncol = 2)) ## 3 rows 2 cols 

Тестовые hmd первый (требуется перестановка):

hmd(t(as.matrix(df1)), df2[1, ]) ## df1 & first row of df2 
# [1] 2 4 6 2 4 

Тестовые foo:

foo(df1, df2, p = 2) 

# $id 
# min1 min2 
# 1 1 4 
# 2 2 3 
# 3 5 2 

# $d 
# mindist.1 mindist.2 
# 1   2   2 
# 2   1   3 
# 3   1   3 

Если вы хотите добавить некоторые столбцы df2, вы знаете, что делать, не так ли?

Answer 2

Пожалуйста, не удивляйтесь, почему я беру другой раздел. Эта часть дает что-то важное. Это не то, что ОП требует, но может помочь любым читателям.

---

Общее расстояние Хэмминга вычисление

В предыдущем ответе, я начинаю с функцией hmd0, которая вычисляет расстояние Хэмминга между двумя целочисленными векторами одинаковой длины.Это означает, что если у нас есть два целочисленных векторы:

set.seed(0) 
x <- sample(1:100, 6) 
y <- sample(1:100, 6) 

мы в конечном итоге с скаляром:

hmd0(x,y) 
# 13 

Что делать, если мы хотим вычислить парного Хэмминг расстояния двух векторов?

В самом деле, простая модификация нашей функции hmd будет делать:

hamming.distance <- function(x, y, pairwise = TRUE) { 
    nx <- length(x) 
    ny <- length(y) 
    rawx <- intToBits(x) 
    rawy <- intToBits(y) 
    if (nx == 1 && ny == 1) return(sum(as.logical(xor(intToBits(x),intToBits(y))))) 
    if (nx < ny) { 
    ## pivoting 
    tmp <- rawx; rawx <- rawy; rawy <- tmp 
    tmp <- nx; nx <- ny; ny <- tmp 
    } 
    if (nx %% ny) stop("unconformable length!") else { 
    bits <- length(intToBits(0)) ## 32-bit or 64 bit? 
    result <- unname(tapply(as.logical(xor(rawx,rawy)), rep(1:ny, each = bits), sum)) 
    } 
    if (pairwise) result else sum(result) 
    } 

Теперь

hamming.distance(x, y, pairwise = TRUE) 
# [1] 0 3 3 2 5 0 
hamming.distance(x, y, pairwise = FALSE) 
# [1] 13 

---

расстояние Хэмминга матрица

Если мы хотим вычислить матрица расстояния от помех, для exa mple,

set.seed(1) 
x <- sample(1:100, 5) 
y <- sample(1:100, 7) 

Матрица расстояние между x и y является:

outer(x, y, hamming.distance) ## pairwise argument has no effect here 

#  [,1] [,2] [,3] [,4] [,5] [,6] [,7] 
# [1,] 2 3 4 3 4 4 2 
# [2,] 7 6 3 4 3 3 3 
# [3,] 4 5 4 3 6 4 2 
# [4,] 2 3 2 5 6 4 2 
# [5,] 4 3 4 3 2 0 2 

Мы также можем сделать:

outer(x, x, hamming.distance) 

#  [,1] [,2] [,3] [,4] [,5] 
# [1,] 0 5 2 2 4 
# [2,] 5 0 3 5 3 
# [3,] 2 3 0 2 4 
# [4,] 2 5 2 0 4 
# [5,] 4 3 4 4 0 

В последнем случае, мы в конечном итоге с симметричной матрицей с 0 по диагонали. Использование outer здесь неэффективно, но оно еще более эффективно, чем запись R-циклов. Поскольку наш hamming.distance написан в R-коде, я остался бы с использованием outer. В my answer до this question я демонстрирую идею использования скомпилированного кода. Это, конечно, требует написания версии C hamming.distance, но я не буду показывать ее здесь.

Answer 3

Вот альтернативное решение, которое использует только базу R, и должно быть очень быстрым, особенно когда ваши df1 и df2 имеют много строк. Основная причина этого заключается в том, что он не использует любой цикл R-уровня для вычисления расстояний Хэмминга, например, для циклов, while-loops или * apply. Вместо этого он использует matrix multiplication for computing the Hamming distance. В R это намного быстрее, чем любой подход с использованием цикла R-уровня. Также обратите внимание, что использование функции * apply не обязательно сделает ваш код более быстрым, чем использование цикла for. Две другие связанные с эффективностью функции этого подхода: (1) Он использует partial sorting для нахождения наилучших двух совпадений для каждой строки в df2 и (2) Он хранит все побитовое представление df1 в одной матрице (то же самое для df2), и делает это за один шаг, без использования каких-либо петель уровня R.

Функция, которая делает всю работу:

# INPUT:  
# X corresponds to your entire df1, but is a matrix 
# Y corresponds to your entire df2, but is a matrix 
# OUTPUT: 
# Matrix with four columns corresponding to the values 
# that you specified in your question 
fun <- function(X, Y) { 

    # Convert integers to bits 
    X <- intToBits(t(X)) 
    # Reshape into matrix 
    dim(X) <- c(ncols * 32, nrows) 

    # Convert integers to bits 
    Y <- intToBits(t(Y)) 
    # Reshape into matrix 
    dim(Y) <- c(ncols * 32, nrows) 

    # Calculate pairwise hamming distances using matrix 
    # multiplication. 
    # Columns of H index into Y; rows index into X. 
    # The code for the hamming() function was retrieved 
    # from this page: 
    # https://johanndejong.wordpress.com/2015/10/02/faster-hamming-distance-in-r-2/ 
    H <- hamming(X, Y) 

    # Now, for each row in Y, find the two best matches 
    # in X. In other words: for each column in H, find 
    # the two smallest values and their row indices. 
    t(apply(H, 2, function(h) { 
    mindists <- sort(h, partial = 1:2) 
    c(
     ind1 = which(h == mindists[1])[1], 
     val1 = mindists[1], 
     hmd2 = which(h == mindists[2])[1], 
     val2 = mindists[2] 
    ) 
    })) 
} 

Для вызова функции на некоторых случайных данных:

# Generate some random test data with no. of columns 
# corresponding to your data 
nrows <- 1000 
ncols <- 26 

# X corresponds to your df1 
X <- matrix(
    sample(1e6, nrows * ncols, replace = TRUE), 
    nrow = nrows, 
    ncol = ncols 
) 

# Y corresponds to your df2 
Y <- matrix(
    sample(1e6, nrows * ncols, replace = TRUE), 
    nrow = nrows, 
    ncol = ncols 
) 

res <- fun(X, Y) 

В приведенном выше примере с 1000 строк в обоих X (DF1) и Y (df2) занял около 1,1 - 1,2 секунды для работы на моем ноутбуке.