Description Usage Arguments Details Value Author(s) See Also Examples
View source: R/NucleotideOverlap.R
A function for concisely tabulating where genomic features are connected by syntenic hits.
1 2 3 4 5 | NucleotideOverlap(SyntenyObject,
GeneCalls,
LimitIndex = FALSE,
OutputFormat = "Normal",
Verbose = FALSE)
|
SyntenyObject |
An object of class “Synteny” built from the |
GeneCalls |
A named list of objects of class “DFrame” built from |
LimitIndex |
Logical indicating whether to limit which indices in a synteny object to query. |
OutputFormat |
Character string to designate how much information to return. "Sparse" returns only a filled upper triangle of exactly matched positions. "Normal" returns a matrix with associated match information in both the upper and lower triangle of the returned matrix, while "Comprehensive" will return |
Verbose |
Logical indicating whether or not to display a progress bar and print the time difference upon completion. |
Builds a matrix of lists that contain information about linked pairs of genomic features.
An object of class “LinkedPairs”.
Nicholas Cooley npc19@pitt.edu
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 | DBPATH <- system.file("extdata",
"VignetteSeqs.sqlite",
package = "SynExtend")
# Alternatively, to build a database using DECIPHER:
# DBPATH <- tempfile()
# FNAs <- c("ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/006/740/685/GCA_006740685.1_ASM674068v1/GCA_006740685.1_ASM674068v1_genomic.fna.gz",
# "ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/000/956/175/GCA_000956175.1_ASM95617v1/GCA_000956175.1_ASM95617v1_genomic.fna.gz",
# "ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/000/875/775/GCA_000875775.1_ASM87577v1/GCA_000875775.1_ASM87577v1_genomic.fna.gz")
# for (m1 in seq_along(FNAs)) {
# X <- readDNAStringSet(filepath = FNAs[m1])
# X <- X[order(width(X),
# decreasing = TRUE)]
#
# Seqs2DB(seqs = X,
# type = "XStringSet",
# dbFile = DBPATH,
# identifier = as.character(m1),
# verbose = TRUE)
# }
Syn <- FindSynteny(dbFile = DBPATH)
GeneCalls <- vector(mode = "list",
length = ncol(Syn))
GeneCalls[[1L]] <- gffToDataFrame(GFF = system.file("extdata",
"GCA_006740685.1_ASM674068v1_genomic.gff.gz",
package = "SynExtend"),
Verbose = TRUE)
GeneCalls[[2L]] <- gffToDataFrame(GFF = system.file("extdata",
"GCA_000956175.1_ASM95617v1_genomic.gff.gz",
package = "SynExtend"),
Verbose = TRUE)
GeneCalls[[3L]] <- gffToDataFrame(GFF = system.file("extdata",
"GCA_000875775.1_ASM87577v1_genomic.gff.gz",
package = "SynExtend"),
Verbose = TRUE)
# Alternatively:
# GeneCalls <- vector(mode = "list",
# length = ncol(Syn))
# GeneCalls[[1L]] <- rtracklayer::import(system.file("extdata",
# "GCA_006740685.1_ASM674068v1_genomic.gff.gz",
# package = "SynExtend"))
# GeneCalls[[2L]] <- rtracklayer::import(system.file("extdata",
# "GCA_000956175.1_ASM95617v1_genomic.gff.gz",
# package = "SynExtend"))
# GeneCalls[[3L]] <- rtracklayer::import(system.file("extdata",
# "GCA_000875775.1_ASM87577v1_genomic.gff.gz,
# package = "SynExtend"))
names(GeneCalls) <- seq(length(GeneCalls))
Links <- NucleotideOverlap(SyntenyObject = Syn,
GeneCalls = GeneCalls,
LimitIndex = FALSE,
Verbose = TRUE)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.