Cordes de correspondance et de comptage (K-Mer de l'ADN) dans R

Si vous recherchez une vitesse, la solution évidente est stringi code> a> paquet. Il y a stri_count_fixed code> fonction pour compter les modèles. Et maintenant, vérifiez le code et la référence!

DNAlst<-list("CAAACTGATTTT","GATGAAAGTAAAATACCG","ATTATGC","TGGA","CGCGCATCAA")
dna <- stri_paste(rep(c("A","C","G","T"),each=4),c("A","C","G","T"))
result <- t(sapply(DNAlst, stri_count_fixed,pattern=dna,overlap=TRUE))
colnames(result) <- dna
result
     AA AC AG AT CA CC CG CT GA GC GG GT TA TC TG TT
[1,]  2  1  0  1  1  0  0  1  1  0  0  0  0  0  1  3
[2,]  5  1  1  2  0  1  1  0  2  0  0  1  2  0  1  0
[3,]  0  0  0  2  0  0  0  0  0  1  0  0  1  0  1  1
[4,]  0  0  0  0  0  0  0  0  1  0  1  0  0  0  1  0
[5,]  1  0  0  1  2  0  2  0  0  2  0  0  0  1  0  0



fstri <- function(x){
    t(sapply(x, stri_count_fixed,dna,T))
}
fbio <- function(x){
    t(sapply(x, function(x){x1 <-  DNAString(x); oligonucleotideFrequency(x1,2)}))
}

all(fstri(DNAlst)==fbio(DNAlst)) #results are the same
[1] TRUE

longDNA <- sample(DNAlst,100,T)
microbenchmark(fstri(longDNA),fbio(longDNA))
Unit: microseconds
           expr        min         lq        mean     median         uq        max neval
 fstri(longDNA)    689.378    738.184    825.3014    766.862    793.134   6027.039   100
  fbio(longDNA) 118371.825 125552.401 129543.6585 127245.489 129165.711 359335.294   100
127245.489/766.862
## [1] 165.9301

Ma réponse n'était pas aussi rapide que @bartektartanus. Cependant, il est également assez rapide et j'ai écrit le code ...: D

Le côté plus de mon code par rapport aux autres est: H2>

1. n'a pas besoin d'installer la version non différée de stri_count_fixed code> li>
2. probablement Stringi code> sera vraiment lent pour les gros k-mers car il doit générer toutes les combinaisons possibles pour le modèle et ensuite, Vérifiez leur existence dans les données et comptez combien de fois il apparaît. Li>
3. Cela fonctionne également pour une longue simple forte> et et plusieurs séquences fortes> avec la même sortie vraiment rapide. li>
4. Vous pouvez mettre une valeur pour k code> au lieu de créer une chaîne de motif. li>
5. Si vous exécutez oligonucléotidefrequency code> avec un k code> plus grand que 12 Dans une grande séquence, la fonction gèle l'excès d'utilisation de la mémoire et R est redémarré, tandis que ma fonction est assez rapide. Li> ol>

   mon code h2> xxx pré>

   i Utilisez le package code> uniquement pour compter les bases ... Vous pouvez utiliser d'autres options comme Stringi Code> Pour compter ... Si vous supprimez tout le code ci-dessous k_mers lapply code> et retour (k_mers) code> Il renvoie juste la liste ... de tous les k-mers avec les vecteurs répétés respectifs p> séquence code> Voici une séquence de 1000bp h3> xxx pré>

   tests effectués avec ma fonction: strong> p>

   #Slow 1 sequence 12-mers
   > system.time({oligonucleotideFrequency(DNAString(sequence),12)})
        user    system   elapsed 
      150.11      1.14    151.37 
   
   #Biostrings package freezes for a single sequence of 13-mers
   > system.time({oligonucleotideFrequency(sequence,13)})  
   freezes, used all my 8gb RAM
   

Nous avons récemment publié notre forfait 'Kebabs' dans le cadre du bioconducteur 3.0 Libération. Bien que ce paquet vise à fournir des noyaux de séquence pour la classification, la régression et d'autres tâches telles que la similitude Clustering, le colis inclut une fonctionnalité pour calculer les fréquences K-MER efficacement, aussi:

#for the kebabs documentation please see:
browseVignettes("kebabs")