TFBS search with TFBSTools and a for loop

Question

I am currently looking for TFBS motifs in DNA. The original code I wrote to look for 10 different transcription factor motifs worked beautifully. But now I want to search for more than 100 motifs. Therefore I want to use loops so that I dont have to type in all the motif names manually (and can easily change my list if my interests change). But when I do this my code no longer works. I can print the result of the loops and it looks like the list contains all the correct information, but when I then use this list to compare with my DNA the list seems to only remember the last item of the list.

Here is my original code that works:

initializeJASPARDB("jaspar.sqlite", version="2014")
MA0018.2 <- getMatrixByID(db, ID="MA0018.2")
MA0024.1 <- getMatrixByID(db, ID="MA0024.1")
MA0062.1 <- getMatrixByID(db, ID="MA0062.1")
MA0099.2 <- getMatrixByID(db, ID="MA0099.2")
MA0112.1 <- getMatrixByID(db, ID="MA0112.1")
MA0153.1 <- getMatrixByID(db, ID="MA0153.1")
MA0156.1 <- getMatrixByID(db, ID="MA0156.1")
MA0157.1 <- getMatrixByID(db, ID="MA0157.1")
MA0158.1 <- getMatrixByID(db, ID="MA0158.1")
MA0259.1 <- getMatrixByID(db, ID="MA0259.1")

pwm1 <- toPWM(MA0018.2)
pwm2 <- toPWM(MA0024.1)
pwm3 <- toPWM(MA0062.1)
pwm4 <- toPWM(MA0099.2)
pwm5 <- toPWM(MA0112.1)
pwm6 <- toPWM(MA0153.1)
pwm7 <- toPWM(MA0156.1)
pwm8 <- toPWM(MA0157.1)
pwm9 <- toPWM(MA0158.1)
pwm10 <- toPWM(MA0259.1)

pwmList <- PWMatrixList(pwm1=pwm1, pwm2=pwm2, pwm3=pwm3, pwm4=pwm4, pwm5=pwm5, pwm6=pwm6, pwm7=pwm7, pwm8=pwm8, pwm9=pwm9, pwm10=pwm10)

seq3 <- readDNAStringSet(file = "sequences/HBV/HBV_genotypeEplusFplusGplusHplusRF.fasta")

sitesetList <- searchSeq(pwmList, seq3, min.score="90%")

as(sitesetList, "data.frame")
as(sitesetList, "DataFrame")
as(sitesetList, "GRanges")
#writeGFF3(sitesetList)
Viruses_TFBS <- writeGFF2(sitesetList)

write.csv(Viruses_TFBS, file = "analysis_HBV/HBV_combmanually_10TFBS.csv")

head(Viruses_TFBS, 5)

Here is my new code with the loop:

TFBS_motifs <- read.csv("TFBS_motifs_test.csv")
seqHSV <- readDNAStringSet(file = "sequences/HSV1.fasta")

# for loop objects
TFBS_matrices <- list()
pwms <- list()
TFBS_ID <- TFBS_motifs$ID

for (i in 1:nrow(TFBS_motifs)) {
  
 pwms[[i]] <- toPWM(TFBS_matrices[[i]] <- getMatrixByID(db, ID = TFBS_ID[[i]]))
  
}
 #print(pwms)

pwmList <- PWMatrixList(pwms[[i]])

sitesetList <- searchSeq(pwmList, seqHSV, min.score="90%")


as(sitesetList, "data.frame")
as(sitesetList, "DataFrame")
as(sitesetList, "GRanges")
#writeGFF3(sitesetList)
HSV_TFBS <- writeGFF2(sitesetList)

write.csv(HSV_TFBS, file = "analysis_HBV/HSV_test.csv")

Thank you all so much for any suggestions...

user_12728748 · Answer

You don't really need to run this in a loop; most of these commands can deal with lists. Try something like this:

library(TFBSTools)
library(JASPAR2018) # or use 2014 version if you prefer
library(Biostrings)
db <- file.path(system.file("extdata", package="JASPAR2018"), "JASPAR2018.sqlite")

# provide your motifs here or read in from file
mymotifs <- c("MA0018.2", "MA0024.1", "MA0062.1", "MA0099.2") 
PFMatrixList <- getMatrixByID(db, ID=mymotifs)

# provide your genome fasta file here
HSV1 <- readDNAStringSet("ftp://ftp.ebi.ac.uk/pub/databases/ena/wgs/public/fkj/FKJZ01.fasta.gz") 
sitesetList <- searchSeq(toPWM(PFMatrixList), HSV1, min.score="90%")

# generate GRangesList (and unlist if you prefer a single GRanges object)
unlist(as(sitesetList, "GRangesList"))

TFBS search with TFBSTools and a for loop

1 Answers