Найти расстояние Хэмминга между последовательностями строк - PullRequest
1 голос
/ 28 января 2020

У меня есть набор данных из 3156 последовательностей ДНК, каждая из которых имеет 98290 символов (SNP), включая (обычно) 5 символов: A, C, G, T, N (пробел) ,

Каков оптимальный способ найти попарно расстояние Хэмминга между этими последовательностями?

Обратите внимание, что для каждой последовательности я действительно хочу найти обратную величину количество последовательностей (включая себя), где расстояние Хемминга для каждого сайта меньше некоторого порогового значения (в данном примере 0,1).

До сих пор я пытался сделать следующее:

library(doParallel)
registerDoParallel(cores=8)
result <- foreach(i = 1:3156) %dopar% {
 temp <- 1/sum(sapply(snpdat, function(x) sum(x != snpdat[[i]])/98290 < 0.1))
}

snpdat - переменная list, где snpdat[[i]] содержит i последовательность ДНК.

Это займет около 36 минут для запуска на ядре i7 - 4790 с 16 ГБ оперативной памяти. Я также попытался использовать пакет stringdist, для получения которого требуется больше времени.

Любая помощь очень важна!

1 Ответ

0 голосов
/ 05 февраля 2020

Я не уверен, что это наиболее оптимальное решение, но мне удалось сократить время выполнения до 15 минут, используя R cpp. Я напишу здесь код на случай, если кто-нибудь когда-нибудь найдет его полезным ...

Это код C ++ (я использовал здесь Sugar операторы) ...

#include <Rcpp.h>
using namespace Rcpp;

double test5(const List& C, const int& x){
  double HD;
  for(int i = 0; i < 3156; i++) if(sum(CharacterVector(C[x])!=CharacterVector(C[i])) < 9829) HD++;
  return HD;
}

После компиляции:

library(Rcpp)
sourceCpp("hd_code.cpp")

Я просто вызываю эту функцию из R:

library(foreach)
library(doParallel)
registerDoParallel(cores = 8)
t =Sys.time()
bla = foreach(i = 1:3156, .combine = "c") %dopar% test5(snpdat,i-1)
Sys.time() - t

Может кто-нибудь придумать еще более быстрый способ сделать это?

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...