Creating Gene List

Question

i need to create a list , containing one vector for each gene.Vectors should be the result of using func on each gene.

Answer 1

Hi and welcome to stack overflow. Since you didn't provide any data to test with I made some dummy data as follows. So these are examples of the dummy genes:

valid_codons <- c("aaa", "aac", "aag", "aat", "aca", "acc", "acg", "act",
                  "aga", "agc", "agg", "agt", "ata", "atc", "atg", "att", "caa", "cac", 
                  "cag", "cat", "cca", "ccc", "ccg", "cct", "cga", "cgc", "cgg", "cgt",
                  "cta", "ctc", "ctg", "ctt", "gaa", "gac", "gag", "gat", "gca", "gcc",
                  "gcg", "gct", "gga", "ggc", "ggg", "ggt", "gta", "gtc", "gtg", "gtt",
                  "taa", "tac", "tag", "tat", "tca", "tcc", "tcg", "tct", "tga", "tgc",
                  "tgg", "tgt", "tta", "ttc", "ttg", "ttt")



genes <- replicate(3800, {
  paste0(sample(valid_codons, sample(5:20, 1), replace = TRUE), collapse = "")
})
print(head(genes, 3))
#> [1] "gggtacaaagtgcat"                        
#> [2] "cggaaaaccggggcgtgtccg"                  
#> [3] "ggaccactattactctcctcgggtatagatacccgaggt"

I'm assuming from the function that the data structure you're working with are character vectors, which I made like this:

genes_chars <- strsplit(genes, "")
print(head(genes_chars, 2))
#> [[1]]
#>  [1] "g" "g" "g" "t" "a" "c" "a" "a" "a" "g" "t" "g" "c" "a" "t"
#> 
#> [[2]]
#>  [1] "c" "g" "g" "a" "a" "a" "a" "c" "c" "g" "g" "g" "g" "c" "g" "t" "g" "t" "c"
#> [20] "c" "g"

Now getting to your actual question, I'm wrapping your provided codon_count() function in a lapply loop to calculate the result.

codon_count <- function(gene) {
  answer <- rep(0, 64)
  names(answer) <- valid_codons
  for(i in seq(from=1, to=length(gene), by=3)) {
    codon <- tolower(paste0(gene[i], gene[i+1], gene[i+2]))
    answer[codon] <- answer[codon] + 1
  }
  return(answer[valid_codons])
}

result <- lapply(genes_chars, codon_count)
print(head(result, 2))
#> [[1]]
#> aaa aac aag aat aca acc acg act aga agc agg agt ata atc atg att caa cac cag cat 
#>   1   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   1 
#> cca ccc ccg cct cga cgc cgg cgt cta ctc ctg ctt gaa gac gag gat gca gcc gcg gct 
#>   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0 
#> gga ggc ggg ggt gta gtc gtg gtt taa tac tag tat tca tcc tcg tct tga tgc tgg tgt 
#>   0   0   1   0   0   0   1   0   0   1   0   0   0   0   0   0   0   0   0   0 
#> tta ttc ttg ttt 
#>   0   0   0   0 
#> 
#> [[2]]
#> aaa aac aag aat aca acc acg act aga agc agg agt ata atc atg att caa cac cag cat 
#>   1   0   0   0   0   1   0   0   0   0   0   0   0   0   0   0   0   0   0   0 
#> cca ccc ccg cct cga cgc cgg cgt cta ctc ctg ctt gaa gac gag gat gca gcc gcg gct 
#>   0   0   1   0   0   0   1   0   0   0   0   0   0   0   0   0   0   0   1   0 
#> gga ggc ggg ggt gta gtc gtg gtt taa tac tag tat tca tcc tcg tct tga tgc tgg tgt 
#>   0   0   1   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   0   1 
#> tta ttc ttg ttt 
#>   0   0   0   0

We can check that the dimensions are correct with length() and lengths().

unique(lengths(result))
#> [1] 64
length(result)
#> [1] 3800

However, I think that the code below is somewhat more efficient.

# Split character vectors into groups of three
# Based on https://stackoverflow.com/questions/11619616/how-to-split-a-string-into-substrings-of-a-given-length
splitgenes <- strsplit(genes, "(?<=.{3})", perl = TRUE)

result2 <- t(vapply(splitgenes, function(gene) {
  table(factor(gene, valid_codons))
}, numeric(length(valid_codons))))

# What are the result2 dimensions and content?
dim(result2)
#> [1] 3800   64
result2[1:5, 1:5]
#>      aaa aac aag aat aca
#> [1,]   1   0   0   0   0
#> [2,]   1   0   0   0   0
#> [3,]   0   0   0   0   0
#> [4,]   0   0   0   0   0
#> [5,]   0   1   0   0   0

EDIT:

This is the for-loop equivalent of the lapply statement:

result <- list()
for (i in seq_along(genes_chars)) {
  result[[i]] <- codon_count(genes_chars[[i]])
}

Note however that this is less efficient.