This is a test version of Biostars. For the public version, visit https://www.biostars.org.
How to determine the average distance between known motif in a list of DNA sequences

Hi everybody !

So there is my problem : I am searching for the average distance between a known motif inside sequence, and extend this to a list of sequences... The first part is done, the second part (extend to a list of sequences) is the problematic one ! So, here the way i am doing the first part :

source("motifOccurrence.R") #https://www.r-bloggers.com/calculate-the-average-distance-between-a-given-dna-motif-within-dna-sequences-in-r/

library("seqinr")

df <- readDNAStringSet("X.fasta")

df2 <- df[[1]]

motif <- c("T", "C", "C", "A")

coord <- coordMotif(df2, motif)

motidist <- computeDistance(coord)

motidist

[1] 152

It's appear that the first sequence of my fasta list have an average distance of 152 nucelotide between two TCCA motif. And, i don't know how automatize this to all my list in df...

Thanks by advance for the help.

Kévin

r sequence motif distance

Thanks, i will check it !

1 answer

The Rho function in the seqinr package calculate the indipendence of occurence of a k-words using the expected word count. We can modify it to be run just for our target:

library("seqinr")
data(ec999)

exp_count<-function (sequence, wordsize = 2, alphabet = s2c("acgt")) 
{
    wordcount <- count(sequence, wordsize, freq = FALSE, alphabet = alphabet)
    uni <- count(sequence, 1, freq = TRUE, alphabet = alphabet)
    expected_wordfreq <- function(wordsize, uni) {
        if (wordsize == 1) 
            return(uni)
        else kronecker(uni, expected_wordfreq(wordsize - 1, uni))
    }
    expected_wordcount <- sum(wordcount) * expected_wordfreq(wordsize, 
        uni)

names(expected_wordcount)<-words(wordsize,alphabet = alphabet)
return(expected_wordcount)
}
#exp_count(ec999[[1]],3)

Log in to answer this question.