@@ -1,6 +1,6 @@

 Package: CNEr 

-Version: 1.10.0

-Date: 2016-10-06

+Version: 1.10.1

+Date: 2016-10-25

 Title: CNE Detection and Visualization

 Description: Large-scale identification and advanced visualization 

              of sets of conserved noncoding elements.

@@ -25,7 +25,8 @@ Imports: Biostrings (>= 2.33.4),

          poweRlaw (>= 0.60.3),

          annotate (>= 1.50.0),

          GO.db (>= 3.3.0),

-         R.utils (>= 2.3.0)

+         R.utils (>= 2.3.0),

+         KEGGREST (>= 1.14.0)

 Depends: R (>= 3.2.2)

 Suggests: Gviz (>= 1.7.4),

           BiocStyle,

@@ -66,4 +67,5 @@ Collate:

         WholeGenomeAlignment.R

         Ancora.R

         CNE-methods.R

-        GO.R

\ No newline at end of file

+        GO.R

+        KEGG.R

@@ -62,6 +62,7 @@ importFrom(GO.db, GOCCANCESTOR, GOCCOFFSPRING, GOCCCHILDREN,

                   GOMFANCESTOR, GOMFOFFSPRING, GOMFCHILDREN)

 importFrom("grDevices", "jpeg", "pdf", "png", "postscript")

 importFrom(R.utils, gunzip, gzip)

+importFrom(KEGGREST, keggGet)

 ### -----------------------------------------------------------------

 ### Export S4 classes defined in CNEr

@@ -1,3 +1,9 @@

+CHANGES IN Bioc 3.5

+------------------------

+NEW FEATURES

+    o Add function orgKEGGIds2EntrezIDs to fetch the mapping between KEGG IDs

+      and Entrez IDs

+    

 CHANGES IN Bioc 3.4

 ------------------------

 NEW FEATURES

@@ -50,8 +50,8 @@ readAncora <- function(fn, assembly=NULL,

                                   strand="*",

                                   name=paste0(cne[[1]], ":", 

                                               (cne[[2]]+1), "-", cne[[3]]),

-                                  seqnames2=cne[[1]],

                                   itemRgb=chr2colour(cne[[1]]),

+                                  seqnames2=cne[[1]],

                                   seqinfo=seqinfoQuery)

                      )

   if(is.null(assembly)){

@@ -144,13 +144,13 @@ makeCNEDensity <- function(x, outputDir=".",

   message("Making bedGraph files...")

   bedFirst <- reduce(bedFirst, ignore.strand=TRUE)

   covFirst <- coverage(bedFirst)

-  densityFirst <- runmean(covFirst, k=windowSizeFirst*1000, 

-                          endrule = "constant") * 100

+  densityFirst <- suppressWarnings(runmean(covFirst, k=windowSizeFirst*1000, 

+                          endrule = "constant") * 100)

   bedSecond <- reduce(bedSecond, ignore.strand=TRUE)

   covSecond <- coverage(bedSecond)

-  densitySecond <- runmean(covSecond, k=windowSizeSecond*1000,

-                           endrule = "constant") * 100

+  densitySecond <- suppressWarnings(runmean(covSecond, k=windowSizeSecond*1000,

+                           endrule = "constant") * 100)

   firstTrackLine <- new("GraphTrackLine",

                         name=paste(genomeFirst, "CNEs density", threshold),

@@ -2,7 +2,8 @@

 ### make GRBs from CNEs

 ### Exported!

 makeGRBs <- function(x, winSize=NULL, genes=NULL, ratio=1, 

-                     background=c("chromosome", "genome")){#,

+                     background=c("chromosome", "genome"),

+                     minCNEs=1L){#,

                      #byChrom=FALSE

   if(!is(x, "GRangesList")){

     stop("The input `x` must be a `GRangesList` object!")

@@ -21,7 +22,7 @@ makeGRBs <- function(x, winSize=NULL, genes=NULL, ratio=1,

     names(xList) <- as.character(1:length(xList))

   }

   x <- unlist(x)

-  x <- reduce(x)

+  x <- reduce(x, ignore.strand=TRUE)

   cov <- coverage(x)

   # Guess winSize from genome size if NULL

@@ -36,7 +37,7 @@ makeGRBs <- function(x, winSize=NULL, genes=NULL, ratio=1,

     stop("The `genes` must be a `GRanges` object!")

   }

-  density <- runmean(cov, k=winSize,  endrule="constant")

+  density <- suppressWarnings(runmean(cov, k=winSize,  endrule="constant"))

   if(background == "genome"){

     # calculate the background percentage of coverage

     totalGenomeSize <- 

@@ -93,11 +94,28 @@ makeGRBs <- function(x, winSize=NULL, genes=NULL, ratio=1,

   for(i in 1:length(xList)){

     hitsCNEs <- findOverlaps(xList[[i]], clusterRanges,

                              ignore.strand=TRUE, type="within")

-    cnes <- sapply(split(queryHits(hitsCNEs), subjectHits(hitsCNEs)), length)

+    

+    ## Add the number of CNEs

+    cnes <- lengths(split(queryHits(hitsCNEs), subjectHits(hitsCNEs)))

     mcols(clusterRanges)[[names(xList)[i]]] <- 0L

     mcols(clusterRanges)[[names(xList)[i]]][as.integer(names(cnes))] <- cnes

+    

+    ## Add the CNE ranges

+    cnes <- split(xList[[i]][queryHits(hitsCNEs)], subjectHits(hitsCNEs))

+    missingCNEs <- setdiff(seq_len(length(clusterRanges)),

+                           subjectHits(hitsCNEs))

+    for(j in missingCNEs){

+      ## This is really a bad implementation. 

+      cnes[[as.character(j)]] <- GRanges()

+    }

+    mcols(clusterRanges)[[paste(names(xList)[i], "CNE")]] <- 

+      cnes[order(as.integer(names(cnes)))]

   }

+  # Filter out the GRBs with few CNEs

+  indexToKeep <- apply(as.data.frame(mcols(clusterRanges)[names(xList)]) >= 

+                         minCNEs, 1, any)

+  clusterRanges <- clusterRanges[indexToKeep]

+  

   return(clusterRanges)

-}

-

+}

\ No newline at end of file

new file mode 100644

@@ -0,0 +1,53 @@

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

+### orgKEGGIds2EntrezIDs: This script is supposed to parse the html page of 

+###              certain species's pathway from KEGG and 

+###              download the associated information for each KEGG pathway ID.

+### Exported!

+orgKEGGIds2EntrezIDs <- function(organism="Homo sapiens"){

+  ## Species mapping

+  organismMapping <- read.table("http://rest.kegg.jp/list/organism",

+                                header=FALSE, sep="\t", quote="",

+                                comment.char="")

+  organismID <- organismMapping[grepl(organism, organismMapping$V3,

+                                      ignore.case=TRUE),

+                                2, drop=TRUE]

+  if(length(organismID) == 0L){

+    stop("The provided organism is not available.",

+         "Please refer to http://rest.kegg.jp/list/organism for available",

+         "organisms")

+  }

+  html <- readLines(paste0("http://www.genome.jp/kegg-bin/show_organism?menu_type=pathway_maps&org=", organismID))

+  html <- grep("^\\d{5}&", html, value=TRUE)

+  ### Hopefully the ID is always 5-digit

+  pathwayIDs <- paste0(organismID, substr(html, 1L, 5L))

+  groups <- sample(rep_len(1L:ceiling(length(pathwayIDs) / 10),

+                           length.out=length(pathwayIDs)))

+  pathwayIDs <- split(pathwayIDs, groups)

+  

+  ## query with KEGG Rest server with 10 entries (maximal) a time, 

+  ## and 200s to 400s gap between each query.

+  query <- lapply(pathwayIDs,

+                  function(x){Sys.sleep(sample(200L:400L, size=1L));keggGet(x)})

+  

+  ## re-organise the query object

+  pathways <- list()

+  for(i in 1:length(query)){

+    for(j in 1:length(query[[i]])){

+      pathways[[query[[i]][[j]]$ENTRY]] <-

+        query[[i]][[j]]

+    }

+  }

+  ## Get the Pathway IDs to Entrez Gene IDs mapping

+  pathwayIDs2GeneIDs <- list()

+  for(i in 1:length(pathways)){

+    genesInfo <- pathways[[i]]$GENE

+    if(is.null(genesInfo)){

+      pathwayIDs2GeneIDs[[pathways[[i]]$ENTRY]] <- NA

+      next

+    }

+    pathwayIDs2GeneIDs[[pathways[[i]]$ENTRY]] <-

+      genesInfo[seq(1, length(genesInfo), by=2)]

+  }

+  pathwayIDs2GeneIDs <- pathwayIDs2GeneIDs[!is.na(pathwayIDs2GeneIDs)]

+  return(pathwayIDs2GeneIDs)

+}

\ No newline at end of file

@@ -49,7 +49,7 @@ CNEDensity <- function(dbName, tableName, chr, start, end,

   # Implement get_cne_ranges_in_region_partitioned_by_other_chr later!!!

   ranges <- reduce(ranges)

   covAll <- coverage(ranges, width=context_end)

-  runMeanAll <- runmean(covAll, k=windowSize, "constant")

+  runMeanAll <- suppressWarnings(runmean(covAll, k=windowSize, "constant"))

   ans <- as(runMeanAll, "GRanges")

   ans$score <- ans$score * 100

   return(ans)

@@ -52,7 +52,7 @@ CNE(assembly1Fn=character(1), assembly2Fn=character(1),

     \item{CNE12}{Object of class \code{"GRangePairs"}:

       The preliminary CNEs from axt file with assembly1 as reference.}

     \item{CNE21}{Object of class \code{"GRangePairs"}:

-      The preliminart CNEs from axt file with assembly2 as reference.}

+      The preliminary CNEs from axt file with assembly2 as reference.}

     \item{CNEMerged}{Object of class \code{"GRangePairs"}:

       The CNEs after merging CNE1 and CNE2.}

     \item{CNEFinal}{Object of class \code{"GRangePairs"}:

@@ -9,7 +9,7 @@

 }

 \usage{

   makeGRBs(x, winSize=NULL, genes=NULL, ratio=1,

-           background=c("chromosome", "genome"))

+           background=c("chromosome", "genome"), minCNEs=1L)

 }

 \arguments{

@@ -33,21 +33,26 @@

   }

   \item{background}{

     \code{character}(1):

-    can be "chromosome" or "genome". When slice the CNE density,

+    can be "chromosome" or "genome". When using \code{slice} for the CNE density,

     the background is calculated on a per-chromosome or whole-genome basis.

   }

+  \item{minCNEs}{

+    \code{integer}(1): the minimal number of CNEs that a GRB needs to have.

+  }

 }

 \details{

-  First we calculated the CNE densities from the CNEs.

+  First we calculate the CNE densities from the CNEs.

   Then we segment the regions according to the values of CNE densities.

   The regions with CNE densities above the expected CNE densities * ratio are

   considered as putative GRBs.

-  As a last step, the putative GRBs that do not encompass any gene

-  are filtered out.

+  Putative GRBs that do not encompass any gene are filtered out.

+  Finally, the GRBs that have fewer than \code{minCNEs} number of CNEs will

+  be filtered out.

 }

 \value{

-  A \code{GRanges} object of GRB coordinates is returned, with the number of CNEs

-  within each GRB.

+  A \code{GRanges} object of GRB coordinates is returned.

+  The numbers of CNEs and the coordinates of CNEs within each GRB are

+  returned as a metadata column.

 }

 \author{

@@ -69,5 +74,5 @@

   makeGRBs(cneList, winSize=200, genes=refGenesDanRer10, ratio=1.2,

            background="genome")

   makeGRBs(cneList, winSize=200, genes=refGenesDanRer10, ratio=1.2,

-           background="chromosome")

+           background="chromosome", minCNEs=3L)

 }

new file mode 100644

@@ -0,0 +1,31 @@

+\name{orgKEGGIds2EntrezIDs}

+\alias{orgKEGGIds2EntrezIDs}

+\title{

+  Fetch mapping from KEGG IDs to Entrez IDs

+}

+\description{

+  Given the desired organism name, fetch the mapping between KEGG IDs and Entrez

+  gene IDs.

+}

+\usage{

+  orgKEGGIds2EntrezIDs(organism="Homo sapiens")

+}

+

+\arguments{

+  \item{organism}{

+    \code{character}(1): the name of prganism to query. It has to be available

+    at http://rest.kegg.jp/list/organism.

+  }

+}

+\value{

+  A \code{list} of Entrez gene IDs with KEGG IDs as names.

+}

+\author{

+  Ge Tan

+}

+

+\examples{

+  \donttest{

+    orgKEGGIds2EntrezIDs(organism="Homo sapiens")

+  }

+}

\ No newline at end of file

@@ -6,8 +6,8 @@

 }

 \description{

   CNE widths can follow heavy tailed distribution that are associated with power-laws.

-  This function plot the reverse cumulative density distribution of CNE widths,

-  and fit a discrete power-law distribution.

+  This function plots the reverse cumulative density distribution of CNE widths,

+  and fits a discrete power-law distribution.

   Goodness of fit can also be evaluated.

 }

 \usage{

@@ -1,10 +1,10 @@

 \name{read.rmskFasta}

 \alias{read.rmskFasta}

 \title{

-  Read a soft reoeat masked fasta

+  Read a soft repeat masked fasta

 }

 \description{

-  Read a soft repeat masked fasta file into a \code{GRanges}.

+  Read a soft repeat masked fasta file into a \code{GRanges} object.

 }

 \usage{

   read.rmskFasta(fn)

@@ -5,7 +5,7 @@

   Read the cne file from Ancora format.

 }

 \description{

-  Read the Ancora CNE file into a \code{GRanges} or \code{GRangePairs}.

+  Read the Ancora CNE file into a \code{GRanges} or \code{GRangePairs} object.

 }

 \usage{

   readAncora(fn, assembly=NULL, tAssemblyFn=NULL, qAssemblyFn=NULL)

@@ -6,7 +6,7 @@

 }

 \description{

 Query the SQLite database based on chromosome, 

-coordinates and some other criterias.

+coordinates and some other criteria.

 Primarily not intended to be used directly.

 For the CNE density plot, \code{fetchCNEDensity} function should be used.

 }

@@ -15,7 +15,7 @@

   \item{distance}{

     It can be "far", "medium" or "close". It defines the scoring matrix used in

     \emph{lastz} aligner.

-    Generally, if two species are close to each other at human and chimp level,

+    Generally, if two species are close to each other, for example human and chimp,

     "close" should be used.

     If two species have a divergence time of 100 MYA, "far" should be used.

     In other cases, "medium" should be used.

@@ -55,7 +55,7 @@ Then we can run `lastz` function to generate the _lav_ files.

 assemblyDir <- "/Users/gtan/OneDrive/Project/CSC/CNEr/2bit"

 axtDir <- "/Users/gtan/OneDrive/Project/CSC/CNEr/axt"

 assemblyTarget <- file.path(system.file("extdata",

-                              package="BSgenome.Drerio.UCSC.danRer10"),

+                            package="BSgenome.Drerio.UCSC.danRer10"),

                             "single_sequences.2bit")

 assemblyQuery <- file.path(system.file("extdata",

                                        package="BSgenome.Hsapiens.UCSC.hg38"),

@@ -156,7 +156,7 @@ netSyntenicFile <- chainNetSyntenic(allPreChain, assemblyTarget, assemblyQuery,

                                                       sub("\\.2bit$", "",

                                                       basename(assemblyQuery),

                                                       ignore.case = TRUE),

-                                              ".noClass.net")),

+                                               ".noClass.net")),

                      binaryChainNet="chainNet", binaryNetSyntenic="netSyntenic")

 ```

@@ -184,5 +184,4 @@ document was compiled:

 ```{r sessionInfo, echo=FALSE}

 sessionInfo()

-```

-

+```

\ No newline at end of file