new file mode 100644

Binary files /dev/null and b/inst/scripts/import/flyFactorSurvey/flyFactorSurvey.RData differ

new file mode 100644

@@ -0,0 +1,287 @@

+# flyFactorSurvey/import.R

+#------------------------------------------------------------------------------------------------------------------------

+library(org.Dm.eg.db)

+library(biomaRt)

+#------------------------------------------------------------------------------------------------------------------------

+bindingDomainXrefSourceFile <- function() {"TFfile2b.tsv"}

+printf <- function(...) print(noquote(sprintf(...)))

+#------------------------------------------------------------------------------------------------------------------------

+run = function (flyFactorSurveyRootDir)

+{

+  filenames = getMatrixFilenames (flyFactorSurveyRootDir)

+

+     # the rootDir has only matrix files, with one exception:

+     # a file which describes the protein binding domains

+     # remove that from the matrix file list

+  

+  removers <- grep (bindingDomainXrefSourceFile(), filenames)

+  if(length(removers) > 0)

+      filenames <- filenames[-removers]

+  

+  full.filenames = file.path(flyFactorSurveyRootDir, filenames)

+

+  list.BD = createXrefBindingDomain (flyFactorSurveyRootDir)

+  xref = createXref ()  # maps flybase ids to uniprot

+  tbl.ref = createExperimentRefTable ()

+  matrices = readAndParse (full.filenames)

+  tbl.md = createMetadata (matrices, tbl.ref, xref, list.BD)

+     # now normalize, not in readAndParse, because we need the counts to go into tbl.md$sequenceCount

+  matrices = normalizeMatrices (matrices)

+  matrices = renameMatrices (matrices, tbl.md)

+  serializedFile <- "flyFactorSurvey.RData"

+  save (matrices, tbl.md, file=serializedFile)

+  printf("saved %d matrices to %s", length(matrices), serializedFile)

+  printf("copy %s to <packageRoot>/MotifDb/inst/extdata, rebuild package", serializedFile)

+

+} # run

+#------------------------------------------------------------------------------------------------------------------------

+# rec'd extra xref info from michael brodsky on (18 jul 2012)

+# protein binding domains are of immediate interest.

+createXrefBindingDomain = function (flyFactorSurveyRootDir)

+{

+  f1 = file.path(flyFactorSurveyRootDir, bindingDomainXrefSourceFile())

+  tbl.raw = read.table (f1, sep='\t', header=TRUE, as.is=TRUE, comment.char='')

+  list.BD = tbl.raw$DNAbindingDomain_Primary

+  names (list.BD) = tbl.raw$FlybaseID

+  invisible (list.BD)

+

+} # createXrefBindingDomain

+#------------------------------------------------------------------------------------------------------------------------

+# a named list, with uniprot id as content, flybase gene names (FBgn) as names

+createXref = function ()

+{

+  tbl.egUNIPROT = toTable (org.Dm.egUNIPROT)

+  tbl.egFLYBASE = toTable (org.Dm.egFLYBASE)

+  tbl.xref = merge (tbl.egFLYBASE, tbl.egUNIPROT)

+  xref = tbl.xref [, 3]

+  names (xref) = tbl.xref [, 2]

+  invisible (xref)

+

+} # createXref

+#------------------------------------------------------------------------------------------------------------------------

+getMatrixFilenames = function (flyFactorSurveyRootDir)

+{

+  all.files = list.files (flyFactorSurveyRootDir)

+  files.to.exclude = c ('go.R', 'RCS', 'rdata')

+  for (file in files.to.exclude) {

+    removers = grep (file, all.files, ignore.case=T)

+    if (length (removers) > 0)

+      all.files = all.files [-removers]

+    } # for

+

+  invisible (sort(all.files))

+

+} # getMatrixFilenames

+#-----------------------------------------------------------------------------------------------------------------------

+createExperimentRefTable = function ()

+{

+  experiment.names = c ('SANGER', 'SOLEXA', 'NAR', 'Cell', 'FlyReg', 'NBT')

+  pubmedID = c (NA_character_, NA_character_, '1858360', '18332042', NA_character_, NA_character_)

+  experimentType = c ('bacterial 1-hybrid, SANGER sequencing',

+                      'bacterial 1-hybrid, SOLEXA sequencing',

+                      'bacterial 1-hybrid, SANGER sequencing',

+                      'bacterial 1-hybrid, SANGER sequencing',

+                      'bacterial 1-hybrid, SANGER sequencing',

+                      'bacterial 1-hybrid, SANGER sequencing')

+

+  tbl.ref = data.frame (pmid=pubmedID, experimentType=experimentType, stringsAsFactors=FALSE)

+  rownames (tbl.ref) = experiment.names

+  invisible (tbl.ref)

+

+} # createExperimentRefTable

+#------------------------------------------------------------------------------------------------------------------------

+parsePWMfromText = function (lines.of.text)

+{

+  stopifnot (sum (sapply (c ('A', 'C', 'T', 'G'), function (token) length (grep (token, lines.of.text)))) == 4)

+

+

+  for (line in lines.of.text) {

+    tokens = strsplit (line, '\\s*[:\\|]') [[1]]

+    nucleotide = tokens [1]

+    numbers.raw = tokens [2]

+    number.tokens = strsplit (numbers.raw, '\\s+', perl=T)[[1]]

+    while (nchar (number.tokens [1]) == 0)

+      number.tokens = number.tokens [-1]

+    numbers = as.numeric (number.tokens)

+    if (!exists ('pwm.result'))

+      pwm.result = matrix (nrow=4, ncol=length (numbers), byrow=TRUE, dimnames=list (c ('A', 'C', 'G', 'T'), 1:length(numbers)))

+    pwm.result [nucleotide,] = numbers

+    }

+

+  pwm.result

+

+} # parsePWMfromText

+#------------------------------------------------------------------------------------------------------------------------

+# we expect exactly one matrix per file

+readAndParse = function (filenames)

+{

+  matrices = list ()

+

+  for (filename in filenames) {

+    lines.of.text = scan (filename, what=character(0), sep='\n', comment='#', quiet=TRUE)

+    mtx = parsePWMfromText (lines.of.text)

+    #mtx.normalized = apply (mtx, 2, function (colvector) colvector / sum (colvector))

+    matrices [[basename(filename)]] = mtx

+    }

+ 

+   invisible (matrices)

+

+} # readAndParse

+#-----------------------------------------------------------------------------------------------------------------------

+createMetadata = function (matrices, tbl.ref, xref, list.BD)

+{

+  options ('stringsAsFactors'=FALSE)

+  dataSource = 'FlyFactorSurvey'

+  organism = 'Dmelanogaster'

+  stopifnot (length (grep ('\\.pfm$', names (matrices))) == length (matrices))

+  native.names.raw = gsub ('\\.pfm$', '', names (matrices))

+  #geneSymbols = sapply (strsplit (native.names.raw, '_'), function (tokens) return (tokens [1]))

+  flybase.ids = sapply (strsplit (native.names.raw, '_'), function (tokens) return (tokens [length (tokens)]))

+  tbl.ids <- fbgnToIDs(flybase.ids)

+     # we may have duplicate FBgns in the otherwise unique matrix names,

+     # creating a tbl.ids with fewer rows than we have matrices.

+     # expand tbl.ids to align with the matrix nmes

+  tbl.ids <- tbl.ids[flybase.ids,]

+  stopifnot(nrow(tbl.ids) == length(matrices))

+     #  next up: call fbgnToIDs(flybase.ids), assign geneSymbols

+     #  and other variables from the returned data.frame

+  stopifnot (length (grep ('^FBgn', flybase.ids)) ==  length (flybase.ids))

+  stopifnot (all (nchar (flybase.ids) == 11))

+  native.names.cooked = gsub ('-', '.', native.names.raw)   # so that the dash reliably separates the three parts of full.names

+  full.names = paste (organism, dataSource, native.names.cooked, sep='-')

+

+  tbl.md = data.frame (providerName=native.names.raw)

+  rownames (tbl.md) = full.names

+  tbl.md = cbind (tbl.md, providerId=flybase.ids)

+  tbl.md = cbind (tbl.md, dataSource = rep (dataSource, nrow (tbl.md)))

+

+  tbl.md = cbind (tbl.md, geneSymbol=tbl.ids$symbol)

+  tbl.md = cbind (tbl.md, geneId=tbl.ids$gene_id)

+  #geneIdType <- rep("ENTREZ", nrow(tbl.md))

+     # not all FBgn ids can be mapped (by fbgnToIds) to entrez genes; for those,

+     # the gene_id field is the FBgn id

+  #geneIdIsFBgn <- grep("^FBgn", tbl.ids$gene_id)

+  #if(length(geneIdIsFBgn) > 0)

+  #    geneIdType [geneIdIsFBgn] <- "FLYBASE"

+  tbl.md = cbind (tbl.md, geneIdType=tbl.ids$geneIdType)

+

+  uniprot.ids = as.character (xref [flybase.ids])

+  tbl.md = cbind (tbl.md, proteinId=uniprot.ids)

+  protein.id.types = rep ('UNIPROT', nrow (tbl.md))

+  NA.protein.indices = which (is.na (tbl.md$proteinId))

+  if (length (NA.protein.indices) > 0)

+    protein.id.types [NA.protein.indices] = NA

+  tbl.md = cbind (tbl.md, proteinIdType=protein.id.types)

+  #tbl.md = cbind (tbl.md, proteinIdType=rep ('UNIPROT', nrow (tbl.md)))

+

+  tbl.md = cbind (tbl.md, organism = rep (organism, nrow (tbl.md)))

+

+  counts = as.integer (sapply (matrices, function (mtx) as.integer (mean (round (colSums (mtx))))))

+  tbl.md = cbind (tbl.md, sequenceCount=counts)

+

+  tbl.md = cbind (tbl.md, bindingSequence=rep (NA_character_, nrow (tbl.md)))

+  bindingDomains = rep (NA, nrow (tbl.md))

+  indices = match (tbl.md$providerId, names (list.BD))

+  bindingDomain = as.character (list.BD)[indices]

+  bindingDomain [bindingDomain == 'NA'] = NA

+  bindingDomain [bindingDomain == ''] = NA

+  tbl.md = cbind (tbl.md, bindingDomain)

+

+  experiment.types = rep (NA_character_, nrow (tbl.md))

+  from.Cell = grep ('_Cell_', native.names.raw)

+  from.NAR  = grep ('_NAR_', native.names.raw)

+  from.SOLEXA = grep ('_SOLEXA_', native.names.raw)

+  from.SANGER = grep ('_SANGER_', native.names.raw)

+  from.FlyReg = grep ('_FlyReg_', native.names.raw) 

+

+  experiment.types [from.SOLEXA] = 'bacterial 1-hybrid, SOLEXA sequencing'

+  experiment.types [from.SANGER] = 'bacterial 1-hybrid, SANGER sequencing'

+  experiment.types [from.Cell] = 'bacterial 1-hybrid, SANGER sequencing'

+  experiment.types [from.NAR] = 'bacterial 1-hybrid, SANGER sequencing'

+  experiment.types [from.FlyReg] = 'DNASE I footprinting'

+

+  tbl.md = cbind (tbl.md, tfFamily=rep(NA_character_, nrow (tbl.md)))

+

+  tbl.md = cbind (tbl.md, experimentType=experiment.types)

+

+  pubmedIDs = rep (NA_character_, nrow (tbl.md))

+  pubmedIDs [from.Cell] = '18332042'

+  pubmedIDs [from.NAR] = '1858360'

+

+  tbl.md = cbind (tbl.md, pubmedID=pubmedIDs)

+

+  

+  invisible (tbl.md)

+

+} # createMetadata

+#------------------------------------------------------------------------------------------------------------------------

+normalizeMatrices = function (matrices)

+{

+  mtx.normalized = sapply (matrices, function (mtx) apply (mtx, 2, function (colvector) colvector / sum (colvector)))

+  invisible (mtx.normalized)

+

+} # normalizeMatrices

+#------------------------------------------------------------------------------------------------------------------------

+renameMatrices = function (matrices, tbl.md)

+{

+  stopifnot (length (matrices) == nrow (tbl.md))

+  names (matrices) = rownames (tbl.md)

+  invisible (matrices)

+

+} # renameMatrices

+#-------------------------------------------------------------------------------

+fbgnToIDs <- function(fbgns, useInputForMissingValues=TRUE)

+{

+   fbgns <- unique(fbgns)

+   if(!exists("tbl.names")) {

+       tbl.sym <- toTable(org.Dm.egSYMBOL)

+       tbl.fbgn <- toTable(org.Dm.egFLYBASE)

+       tbl.names <<- merge(tbl.fbgn, tbl.sym)[, c("flybase_id", "symbol", "gene_id")]

+       flyMart <<- useMart(biomart="ensembl", dataset="dmelanogaster_gene_ensembl")

+       }

+   tbl.bm <- getBM(attributes=c("flybase_gene_id", "entrezgene"), filters=c("flybase_gene_id"),

+                   values=fbgns, mart=flyMart)

+   fbgns.notfound <- setdiff(fbgns, tbl.bm$flybase_gene_id)

+   if(length(fbgns.notfound) > 0){

+       new.rows <- data.frame(flybase_gene_id=fbgns.notfound,

+                              entrezgene=NA,

+                              stringsAsFactors=FALSE)

+       tbl.bm <- rbind(tbl.bm, new.rows)

+       }

+

+   tbl.bm <- unique(tbl.bm)

+

+      # some fbgns could be mapped to multiple entrez genes.  remove all but the

+      # one which occurs first

+   

+   duplicates <- which(duplicated(tbl.bm$flybase_gene_id))

+   if(length(duplicates) > 0)

+       tbl.bm <- tbl.bm[-duplicates,]

+   

+   stopifnot(sort(fbgns) == sort(tbl.bm$flybase_gene_id))

+

+       # order the biomart results so that it matches the order of the incoming fgbns

+   tbl.bm <- tbl.bm[match(fbgns, tbl.bm$flybase_gene_id),]

+   hits <- match(tbl.bm$entrezgene, tbl.names$gene_id)

+   result <- tbl.names[hits,]

+   rownames(result) <- fbgns

+

+       # not all FBgn ids have an associated entrez geneID.

+       # use these orphan FBgn ids as their own "geneID"

+   

+   failed.lookups <- which(is.na(result$flybase_id))

+   failed.count <- length(failed.lookups)

+   if(useInputForMissingValues & failed.count > 0) {

+      failed.names <- rownames(result)[failed.lookups]

+      names.matrix <- matrix(rep(failed.names, failed.count), nrow=failed.count)

+      result[failed.lookups,] <- names.matrix

+      }

+

+   geneIdType <- rep("ENTREZ", nrow(result))

+   geneIdType [failed.lookups] <- "FLYBASE"

+   result <- cbind(result, geneIdType, stringsAsFactors=FALSE)

+   result

+

+

+} # fbgnToIDs

+#------------------------------------------------------------------------------------------------------------------------

new file mode 100644

@@ -0,0 +1,300 @@

+# flyFactorSurvey/test.R

+#-------------------------------------------------------------------------------

+library (RUnit)

+#-------------------------------------------------------------------------------

+source("import.R")

+#-------------------------------------------------------------------------------

+run.tests = function (flyFactorSurveyRootDir="~/s/data/public/TFBS/flyFactorSurvey/unpacked")

+{

+    freshStart ()

+    test.fbgnToIDs()

+    x.list.BD <- test.createXrefBindingDomain (flyFactorSurveyRootDir)

+    x.xref <- test.createXref ()

+    x.filenames <- test.getMatrixFilenames (flyFactorSurveyRootDir)

+    x.tbl.ref <- test.createExperimentRefTable ()

+    x.pwm <- test.parsePWMfromText (flyFactorSurveyRootDir)

+    x.mat3 <- test.readAndParse (flyFactorSurveyRootDir)

+    x.tbl.md <- test.createMetadata (x.mat3, x.tbl.ref, x.xref, x.list.BD)

+    x.mat3f <- test.normalizeMatrices (x.mat3)

+    x.mat3fr <- test.renameMatrices (x.mat3f, x.tbl.md [1:3,])

+

+} # run.tests

+#-------------------------------------------------------------------------------

+freshStart = function ()

+{

+    output.files.easy.to.regenerate = grep ('RData$', dir (), value=TRUE)

+    if (length (output.files.easy.to.regenerate) > 0)

+        unlink (output.files.easy.to.regenerate)

+

+} # freshStart

+#-------------------------------------------------------------------------------

+test.createXrefBindingDomain = function (flyFactorSurveyRootDir)

+{

+    print ('--- test.createXrefBindingDomain')

+    list.BD = createXrefBindingDomain (flyFactorSurveyRootDir)

+    checkEquals (length (list.BD), 777)

+    checkTrue (all (grepl ('FBgn', names (list.BD))))

+        # most abundant domain is zf-C2H2

+    checkEquals(names (sort (table (as.character (list.BD)), decreasing=TRUE)[1]),

+                "zf-C2H2")

+    invisible(list.BD)

+

+} # test.createXrefBindingDomain

+#-------------------------------------------------------------------------------

+test.createXref = function ()

+{

+    print ('--- test.createXref')

+    xref = createXref ()

+    checkEquals (xref [['FBgn0261819']], 'F0J881')

+    #checkEquals (xref [['FBgn0259750']], 'E1JHF4') # gone with org.Dm.eg.db_2.9.0

+    checkEquals (xref [['FBgn0000014']], 'E1JIM9')

+    invisible (xref)

+

+} # test.createXref

+#-------------------------------------------------------------------------------

+test.getMatrixFilenames = function (flyFactorSurveyRootDir)

+{

+    print ('--- test.getMatrixFilenames')

+

+    filenames = getMatrixFilenames (flyFactorSurveyRootDir)

+       # as of (03 apr 2013): 614 matrix files, one binding domains extra

+       # file, TFfile2b.tsv

+    

+    checkEquals (length (filenames), 615)

+    checkEquals(filenames[1:3],

+                c("Abd-A_FlyReg_FBgn0000014.pfm", "Abd-B_FlyReg_FBgn0000015.pfm",

+                  "AbdA_Cell_FBgn0000014.pfm"))

+    invisible (filenames)

+

+} # test.getMatrixFilenames

+#-------------------------------------------------------------------------------

+test.createExperimentRefTable = function ()

+{

+    print ('--- test.createExperimentRefTable')

+    tbl.ref = createExperimentRefTable ()

+    checkEquals (nrow (tbl.ref), 6)

+

+    invisible (tbl.ref)

+

+} # test.createExperimentRefTable

+#-------------------------------------------------------------------------------

+test.parsePWMfromText = function (flyFactorSurveyRootDir)

+{

+    print ('--- test.parsePWMfromText')

+    path <- file.path(flyFactorSurveyRootDir, 'Tup_SOLEXA_10_FBgn0003896.pfm')

+

+    lines.of.text =  scan (path, sep='\n',

+                           what=character(0), comment='#', quiet=TRUE)

+    pwm.tup =  parsePWMfromText (lines.of.text)

+    checkEquals (dim (pwm.tup), c (4, 9))

+    checkEquals (colnames (pwm.tup), as.character (1:9))

+    checkEquals (rownames (pwm.tup), c ('A', 'C', 'G', 'T'))

+

+  invisible (pwm.tup)

+  

+} # test.parsePWMfromText

+#-------------------------------------------------------------------------------

+test.readAndParse = function (flyFactorSurveyRootDir)

+{

+    print ('--- test.readAndParse')

+

+    all.files = file.path(flyFactorSurveyRootDir,

+                          getMatrixFilenames (flyFactorSurveyRootDir))

+    

+    sample.files = all.files [1:3]   # grep ('FBgn0003896', all.files)

+    checkEquals (length (sample.files), 3)

+    checkTrue(all(file.exists(sample.files)))

+

+      # by sorting these filenames, we know the order of the three returned

+      # matrices, and results can easily be checked

+    

+    mtx.test = readAndParse (sample.files)

+

+    checkEquals (length (mtx.test), 3)

+    expected.names = c("Abd-A_FlyReg_FBgn0000014.pfm",

+                       "Abd-B_FlyReg_FBgn0000015.pfm",

+                       "AbdA_Cell_FBgn0000014.pfm")

+    actual.names = sort(names(mtx.test))

+    checkEquals(expected.names, actual.names)

+

+    checkEquals (dim (mtx.test [[1]]), c (4, 8))

+    checkEquals (dim (mtx.test [[2]]), c (4, 8))

+    checkEquals (dim (mtx.test [[3]]), c (4, 7))

+

+    invisible (mtx.test)

+

+} # test.readAndParse

+#-------------------------------------------------------------------------------

+test.createMetadata = function (matrices, tbl.ref, xref, list.BD)

+{

+    print ('--- test.createMetadata')

+

+    tbl.md = createMetadata (matrices, tbl.ref, xref, list.BD)

+    checkEquals (nrow (tbl.md), length (matrices))

+    checkEquals (ncol (tbl.md), 15)

+    expected <- c("providerName", "providerId", "dataSource", "geneSymbol",

+                  "geneId", "geneIdType", "proteinId", "proteinIdType",

+                  "organism", "sequenceCount", "bindingSequence",

+                  "bindingDomain", "tfFamily", "experimentType", "pubmedID")

+

+    checkEquals (colnames (tbl.md), expected)

+    tester =  "Dmelanogaster-FlyFactorSurvey-Abd.A_FlyReg_FBgn0000014"

+    checkTrue (tester %in% rownames (tbl.md))

+    x = as.list (tbl.md [tester,])

+    checkEquals (x$providerName, "Abd-A_FlyReg_FBgn0000014")

+    checkEquals (x$providerId, "FBgn0000014")

+    checkEquals (x$geneSymbol, "abd-A")

+    checkEquals (x$geneId, "42037")

+    checkEquals (x$geneIdType, "ENTREZ")

+    checkEquals (x$proteinId, "E1JIM9")

+    checkEquals (x$proteinIdType, "UNIPROT")

+    checkEquals (x$organism, "Dmelanogaster")

+    checkEquals (x$sequenceCount, 37)

+    checkTrue   (is.na (x$bindingSequence))

+    checkEquals (x$bindingDomain, 'Homeobox')

+    checkTrue   (is.na (x$tfFamily))

+    checkEquals (x$experimentType, "DNASE I footprinting")

+    checkTrue   (is.na (x$pubmedID))

+  

+    invisible (tbl.md)

+

+} # test.createMetadata

+#-------------------------------------------------------------------------------

+test.renameMatrices = function (matrices, tbl.md)

+{

+    print ('--- test.renameMatrices')

+    stopifnot (length (matrices) == 3)

+    stopifnot (nrow (tbl.md) == 3)

+

+    old.names = names (matrices)

+    mtxr = renameMatrices (matrices, tbl.md)

+

+      # make sure the dash in the incoming, eg,

+      #     "Abd-A_FlyReg_FBgn0000014.pfm"

+      # is converted to a '.', reserving dashes for 

+      # the separation of dataSource-species-nativeGeneName

+

+

+    checkEquals (names (mtxr) [1],

+                 "Dmelanogaster-FlyFactorSurvey-Abd.A_FlyReg_FBgn0000014")

+

+    checkEquals (names (mtxr) [2],

+                 "Dmelanogaster-FlyFactorSurvey-Abd.B_FlyReg_FBgn0000015")

+

+    checkEquals (names (mtxr) [3],

+                 "Dmelanogaster-FlyFactorSurvey-AbdA_Cell_FBgn0000014")

+

+    invisible (mtxr)

+

+} # renameMatrices

+#-------------------------------------------------------------------------------

+test.fbgnToIDs <- function()

+{

+    print("--- test.fbgnToIDs")

+       # first, make sure that unmappable FBgn ids are handled properly

+    bogus <- fbgnToIDs("bogus")

+    checkEquals(dim(bogus), c(1,4))

+    checkEquals(rownames(bogus), "bogus")

+    checkEquals(as.character(as.list(bogus[1,],)), c(rep("bogus", 3), "FLYBASE"))

+

+       # disable the translation of NA to rowname

+    bogus <- fbgnToIDs("bogus", useInputForMissingValues=FALSE)

+    checkEquals(dim(bogus), c(1,4))

+    checkEquals(rownames(bogus), "bogus")

+    checkTrue(all(is.na(bogus[1,1:3])))

+

+    fbgns <- c ("FBgn0003267", "FBgn0000611", "FBgn0004394", "FBgn0027364",

+                "FBgn0000413")

+    tbl.ids <- fbgnToIDs(fbgns)

+    checkEquals(tbl.ids$flybase_id, fbgns)

+    checkEquals(tbl.ids$symbol, c("ro", "exd", "pdm2", "Six4", "da"))

+    checkEquals(tbl.ids$geneIdType, rep("ENTREZ", 5))

+

+         # 0259750 is obsolete, replace by 0264442

+    tbl.ids <- fbgnToIDs("FBgn0259750")

+    checkEquals(as.list(tbl.ids),

+                list(flybase_id="FBgn0264442", symbol="ab", gene_id="34560",

+                     geneIdType="ENTREZ"))

+

+       # http://flybase.org/static_pages/downloads/IDConv.html

+       # FBgn0003986 	 FBgn0261930 	FBgn0261930	 vnd 

+    tbl.ids <- fbgnToIDs("FBgn0003986")

+    checkEquals(as.character(as.list(tbl.ids[1,],)),

+                c(rep("FBgn0003986", 3), "FLYBASE"))

+    

+    load("../../../extdata/fbgns.RData")

+    checkEquals(length(fbgns), 326)

+    tbl.ids <- fbgnToIDs(fbgns)

+    checkEquals(dim(tbl.ids), c(length(fbgns), 4))

+

+       # no NA's should survive

+    checkTrue(!any(is.na(tbl.ids)))

+

+       # what percentage of the 326 ids fail to map?  do NOT convert NAs

+    tbl.ids <- fbgnToIDs(fbgns, useInputForMissingValues=FALSE)

+    failure.count <- length(which(is.na(tbl.ids$flybase_id)))

+    checkEquals(failure.count, 18)

+

+} # test.fbgnToIDsp

+#-------------------------------------------------------------------------------

+# robert stojnic reported that MotifDb 1.1.6 (and before) had some incorrect

+# gene symbols, as seen here:

+#   library(MotifDb)

+#   mdb <- MotifDb

+#   x <- query(mdb, "ab_SANGER_10_FBgn0259750")

+#   t(as.matrix(mcols(x)))

+#                 [,1]                                   

+#   providerName    "ab_SANGER_10_FBgn0259750"             

+#   providerId      "FBgn0259750"

+#   dataSource      "FlyFactorSurvey"                      

+#   geneSymbol      "Ab"                                   

+#   geneId          "FBgn0259750"                          

+#   geneIdType      "FLYBASE"                              

+#   proteinId       "E1JHF4"                               

+#   proteinIdType   "UNIPROT"                              

+#   organism        "Dmelanogaster"                        

+#   sequenceCount   "20"                                   

+#   bindingSequence NA                                     

+#   bindingDomain   "zf-C2H2"                              

+#   tfFamily        NA                                     

+#   experimentType  "bacterial 1-hybrid, SANGER sequencing"

+#   pubmedID        NA

+#

+# FBgn0259750 is obsolete.  flybase provides a converter:

+#    http://flybase.org/static_pages/downloads/IDConv.html shows this, and provide the current FBgn:

+#

+#           Submitted ID          Current ID    Converted ID    Related record 

+#            FBgn0259750          FBgn0264442    FBgn0264442                ab 

+#

+# make sure we get this right

+# the matrix files are:

+#

+#     flyFactorSurveyRootDir/ab_SANGER_10_FBgn0259750.pfm

+#     flyFactorSurveyRootDir/ab_SOLEXA_5_FBgn0259750.pfm

+#

+test.geneNaming <- function(flyFactorSurveyRootDir)

+{

+

+   flyFactorSurveyRootDir <- "/Users/pshannon/s/data/public/TFBS/flyFactorSurvey/unpacked"

+   filenames = getMatrixFilenames (flyFactorSurveyRootDir)

+

+   removers <- grep (bindingDomainXrefSourceFile(), filenames)

+   if(length(removers) > 0)

+       filenames <- filenames[-removers]

+  

+   keepers <- grep("FBgn0259750", filenames)

+   checkTrue(length(keepers) == 2)

+   filenames <- filenames[keepers]

+   full.filenames = file.path(flyFactorSurveyRootDir, filenames)

+   list.BD = createXrefBindingDomain (flyFactorSurveyRootDir)

+   xref = createXref ()  # maps flybase ids to uniprot

+   tbl.ref = createExperimentRefTable ()

+   matrices = readAndParse (full.filenames)

+   tbl.md = createMetadata (matrices, tbl.ref, xref, list.BD)

+     # now normalize, not in readAndParse, because we need the counts to go into tbl.md$sequenceCount

+   checkEquals(tbl.md$geneSymbol, c("ab", "ab"))

+   matrices = normalizeMatrices (matrices)

+   matrices = renameMatrices (matrices, tbl.md)

+

+} # test.geneNaming

+#-------------------------------------------------------------------------------