Merge pull request #332 from netZoo/master

Merging build fixes into devel

DESCRIPTION

@@ -29,7 +29,6 @@ Depends: R (>= 4.2.0),
    igraph,
 	 reticulate,
 	 pandaR,
-   matrixcalc,
    Biobase
 Remotes: 
 	 stan-dev/cmdstanr,
@@ -47,7 +46,6 @@ Imports:
     RCy3, 
     viridisLite,
     STRINGdb,
-    Biobase,
     GOstats,
     AnnotationDbi,
     matrixStats,
@@ -88,7 +86,8 @@ Imports:
     preprocessCore,
     readr,
     RColorBrewer,
-    quantro
+    quantro,
+    matrixcalc
 License: GPL-3
 Encoding: UTF-8
 LazyData: false
@@ -97,7 +96,7 @@ Suggests:
     knitr,
     rmarkdown,
     pkgdown,
-    dorothea
+    dorothea,
 VignetteEngine: knitr
 VignetteBuilder: knitr
 RoxygenNote: 7.3.1

R/ALPACA.R

@@ -35,7 +35,7 @@ alpaca <- function(net.table,file.stem,verbose=FALSE)
     if (!(is.null(file.stem))) write.table(ctrl.memb, paste(c(file.stem,"_ALPACA_ctrl_memb.txt"),collapse=""),row.names=TRUE,col.names=FALSE,quote=FALSE,sep="\t")
 
     pos.table <- net.table[intersect(which(net.table[,3]>=0),which(net.table[,4]>=0)),]
-    pos.graph <- graph.edgelist(as.matrix(pos.table[,seq_len(2)]),directed=TRUE)
+    pos.graph <- graph_from_edgelist(as.matrix(pos.table[,seq_len(2)]),directed=TRUE)
 
     if (length(setdiff(V(pos.graph)$name,names(ctrl.memb)))>0)
     {

R/COBRA.R

@@ -62,7 +62,7 @@ cobra <- function(X, expressionData, method = "pearson"){
     print(dim(expressionData))
   }
   if(method == "pcorsh"){
-   C <- matrix(as.numeric(pcor.shrink(t(expressionData))), dim(expressionData)[1], dim(expressionData)[1])
+   C <- matrix(as.numeric(corpcor::pcor.shrink(t(expressionData))), dim(expressionData)[1], dim(expressionData)[1])
   }
 
   eigenG <- rARPACK::eigs_sym(C,N)

R/CONDOR.R

@@ -10,7 +10,7 @@
 #' \code{condor.object$edges}
 #' @param cs.method is a string to specify which unipartite community 
 #' structure algorithm should be used for the seed clustering. 
-#' Options are \code{LCS} (\code{\link[igraph]{multilevel.community}}), 
+#' Options are \code{LCS} (\code{\link[igraph]{cluster_louvain}}), 
 #' \code{LEC} (\code{\link[igraph]{leading.eigenvector.community}}), 
 #' \code{FG} (\code{\link[igraph]{fastgreedy.community}}).
 #' @param project Provides options for initial seeding of the bipartite 
@@ -60,13 +60,13 @@ condorCluster <- function(condor.object,cs.method="LCS",project=TRUE,low.memory=
   }
 
   #make sure there's only one connected component
-  g.component.test <- graph.data.frame(elist,directed=FALSE)
-  if(!is.connected(g.component.test)){
+  g.component.test <- graph_from_data_frame(elist,directed=FALSE)
+  if(!is_connected(g.component.test)){
     stop("More than one connected component detected,
          method requires only one connected component")
   }  
 
-  G <- graph.data.frame(elist,directed=FALSE)
+  G <- graph_from_data_frame(elist,directed=FALSE)
 
   project.weights <- weights
   #Use unipartite community structure method for first pass
@@ -88,7 +88,7 @@ condorCluster <- function(condor.object,cs.method="LCS",project=TRUE,low.memory=
     gc()
     colnames(gM) <- blue.names
     rownames(gM) <- blue.names
-    G1 = graph.adjacency(gM,mode="undirected",weighted=TRUE,diag=FALSE);
+    G1 = graph_from_adjacency_matrix(gM,mode="undirected",weighted=TRUE,diag=FALSE);
     #if(clusters(G1)$no > 1){print("Warning more than one component! May cause indexing error")}
     #V(G1)$name <- sort(unique(as.vector(esub[,2])))
     #remove loops and multiple edges
@@ -104,7 +104,7 @@ condorCluster <- function(condor.object,cs.method="LCS",project=TRUE,low.memory=
     #blue.indx <- V(G)$name %in% blue.names
   }
 
-  if(cs.method=="LCS"){cs0 = multilevel.community(gcc.initialize, weights=project.weights)}
+  if(cs.method=="LCS"){cs0 = cluster_louvain(gcc.initialize, weights=project.weights)}
   if(cs.method=="LEC"){cs0 = leading.eigenvector.community(gcc.initialize, weights=project.weights)}
   if(cs.method=="FG"){cs0 = fastgreedy.community(gcc.initialize, weights=project.weights)}
   print(paste("modularity of projected graph",max(cs0$modularity)))
@@ -361,8 +361,8 @@ condorMatrixModularity = function(condor.object,T0=cbind(seq_len(q),rep(1,q)),we
   #Convert the edgelist to a sparseMatrix object
   esub <- condor.object$edges
   #make sure there's only one connected component
-  g.component.test <- graph.data.frame(esub,directed=FALSE)
-  if(!is.connected(g.component.test)){
+  g.component.test <- graph_from_data_frame(esub,directed=FALSE)
+  if(!is_connected(g.component.test)){
     stop("More than one connected component,
          method requires only one connected component")
   }
@@ -567,8 +567,8 @@ condorModularityMax = function(condor.object,T0=cbind(seq_len(q),rep(1,q)),weigh
   #Convert the edgelist to a sparseMatrix object
   esub <- condor.object$edges
   #make sure there's only one connected component
-  g.component.test <- graph.data.frame(esub,directed=FALSE)
-  if(!is.connected(g.component.test)){
+  g.component.test <- graph_from_data_frame(esub,directed=FALSE)
+  if(!is_connected(g.component.test)){
     stop("More than one connected component,
          method requires only one connected component")
   }
@@ -835,20 +835,19 @@ condorPlotCommunities = function(condor.object,color_list,point.size=0.01,
 #'  
 condorPlotHeatmap = function(condor.object, main="", xlab="blues", ylab="reds"){
   bo <- condor.object
-  attach(bo)
   # convert edge lists to adjacency matrices (n reds x m blues)
-  adj = get.adjacency(G, attr="weight", sparse=FALSE)
+  adj = as_adjacency_matrix(bo$G, attr="weight", sparse=FALSE)
   # reorder reds according to community membership
-  reds = as.character(red.memb[order(red.memb[,2]),1])
+  reds = as.character(bo$red.memb[order(bo$red.memb[,2]),1])
   adj = adj[reds,]
   # reorder blues according to community membership
-  blues = as.character(blue.memb[order(blue.memb[,2]),1])
+  blues = as.character(bo$blue.memb[order(bo$blue.memb[,2]),1])
   adj = adj[,blues]
-  rowsep = cumsum(as.vector(table(red.memb[,2])))
-  colsep = cumsum(as.vector(table(blue.memb[,2])))
-  labCol <- as.character(sort(blue.memb[,2]))
+  rowsep = cumsum(as.vector(table(bo$red.memb[,2])))
+  colsep = cumsum(as.vector(table(bo$blue.memb[,2])))
+  labCol <- as.character(sort(bo$blue.memb[,2]))
   labCol[duplicated(labCol)] <- ""
-  labRow <- as.character(sort(red.memb[,2]))
+  labRow <- as.character(sort(bo$red.memb[,2]))
   labRow[duplicated(labRow)] <- ""
   heatmap.2(adj, Rowv=FALSE, Colv=FALSE, dendrogram="none", keysize=1.25,
             col=colorpanel(10, "white", "black"), scale="none",
@@ -857,7 +856,6 @@ condorPlotHeatmap = function(condor.object, main="", xlab="blues", ylab="reds"){
             sepwidth = c(0.025, 0.025), ylab=ylab, xlab=xlab, margins=c(3,3),
             labCol=labCol, labRow=labRow, offsetRow=0, offsetCol=0,
             breaks=sort(c(0.1,seq(0, max(adj),length.out=10))))
-  detach(bo)
 }
 
 
@@ -996,7 +994,7 @@ createCondorObject <- function(edgelist,return.gcc=TRUE){
          same column of 'edgelist'.")
   }
 
-  g <- graph.data.frame(edgelist,directed=FALSE)
+  g <- graph_from_data_frame(edgelist,directed=FALSE)
   blue.indx <- V(g)$name %in% unique(edgelist[, 2])
   V(g)$color <- "red"
   V(g)$color[blue.indx] <- "blue"
@@ -1021,9 +1019,9 @@ createCondorObject <- function(edgelist,return.gcc=TRUE){
 
 max.component = function(g){
   # return largest connected component of the iGraph graph object g
-  g.clust = clusters(g);
+  g.clust = components(g);
   maxclust.id = which(g.clust$csize == max(g.clust$csize))[1];
-  h = induced.subgraph(g, which(g.clust$membership == maxclust.id)); # 1-indexed here
+  h = induced_subgraph(g, which(g.clust$membership == maxclust.id)); # 1-indexed here
   return(h);
 }
 

R/DRAGON.R

@@ -175,13 +175,13 @@ estimatePenaltyParameters = function(X1,X2)
   #          + 4.*np.sum(eSqs12))
 
   # T1_1 = -2.*(np.sum(varS1) - np.trace(varS1) + np.sum(eSqs12))
-  T1_1 = -2*(sum(varS1) - matrix.trace(as.matrix(varS1)) + sum(esqS12))
+  T1_1 = -2*(sum(varS1) - matrixcalc::matrix.trace(as.matrix(varS1)) + sum(esqS12))
   # T1_2 = -2.*(np.sum(varS2) - np.trace(varS2) + np.sum(eSqs12))
-  T1_2 = -2*(sum(varS2) - matrix.trace(as.matrix(varS2)) + sum(esqS12))
+  T1_2 = -2*(sum(varS2) - matrixcalc::matrix.trace(as.matrix(varS2)) + sum(esqS12))
   # T2_1 = np.sum(eSqs1) - np.trace(eSqs1)
-  T2_1 = sum(esqS1) - matrix.trace(as.matrix(esqS1))
+  T2_1 = sum(esqS1) - matrixcalc::matrix.trace(as.matrix(esqS1))
   # T2_2 = np.sum(eSqs2) - np.trace(eSqs2)
-  T2_2 = sum(esqS2) - matrix.trace(as.matrix(esqS2))
+  T2_2 = sum(esqS2) - matrixcalc::matrix.trace(as.matrix(esqS2))
   # T3 = 2.*np.sum(eSqs12)
   T3 = 2*sum(esqS12)
   # T4 = 4.*(np.sum(varS12)-np.sum(eSqs12))

R/EGRET.R

@@ -56,6 +56,9 @@
 #' \donttest{
 #' runEgret(qtl,vcf,qbic,motif,expr,ppi,nameGeneMap,tag)
 #' }
+#' file.remove("my_toy_egret_run_egret.RData")
+#' file.remove("my_toy_egret_run_panda.RData")
+#' file.remove("priors_my_toy_egret_run.txt")
 #' @export
 
 runEgret <- function(b,v,q,m,e,p,g,t){
@@ -92,7 +95,7 @@ runEgret <- function(b,v,q,m,e,p,g,t){
   colnames(vcf) <- c("CHROM",  "POS"  ,   "ID"  ,    "REF"  ,   "ALT"   ,  "QUAL"   , "FILTER" , "INFO"   , "FORMAT", "NA12878")
   snp_ids <- paste0(vcf$CHROM,"_",vcf$POS)
   rownames(vcf) <- snp_ids
-  vcf <- tidyr::separate(vcf, NA12878, c("allele1", "allele2"), "\\|", remove = TRUE)
+  vcf <- tidyr::separate(vcf, col="NA12878", into=c("allele1", "allele2"), sep="\\|", remove = TRUE)
   vcf$alt_allele_count <- as.numeric(vcf$allele1) + as.numeric(vcf$allele2)
   vcf$snp_id <- snp_ids
   qbic_ag$alt_allele_count <- vcf$alt_allele_count[match(qbic_ag$snpID, vcf$snp_id)]

R/LIONESS.R

@@ -41,6 +41,7 @@
 #'     motif_file = motif_file_path, ppi_file = ppi_file_path, 
 #'     modeProcess="union",start_sample=1, end_sample=1, precision="single")
 #' }
+#' unlink("lioness_output", recursive=TRUE)
 #' @import reticulate
 #' @export
 lionessPy <- function(expr_file, motif_file=NULL, ppi_file=NULL, computing="cpu", precision="double", save_tmp=TRUE, modeProcess="union", remove_missing=FALSE, start_sample=1, end_sample="None", save_single_network=FALSE, save_dir="lioness_output", save_fmt='npy'){
@@ -179,7 +180,8 @@ lionessPy <- function(expr_file, motif_file=NULL, ppi_file=NULL, computing="cpu"
 #' "coopNet" is the cooperative network
 #' @examples
 #' data(pandaToyData)
-#' lionessRes <- lioness(expr = pandaToyData$expression[,1:3], motif = pandaToyData$motif, ppi = pandaToyData$ppi,hamming=1,progress=FALSE)
+#' lionessRes <- lioness(expr = pandaToyData$expression[,1:3], motif = pandaToyData$motif, 
+#'   ppi = pandaToyData$ppi,hamming=1,progress=FALSE)
 #' @references
 #' Kuijjer, M.L., Tung, M., Yuan, G., Quackenbush, J. and Glass, K., 2015. 
 #' Estimating sample-specific regulatory networks. arXiv preprint arXiv:1505.06440.

R/MONSTER.R

@@ -104,13 +104,12 @@ monsterPrintMonsterAnalysis <- function(x, ...){
 #' data(yeast)
 #' design <- c(rep(0,20),rep(NA,10),rep(1,20))
 #' yeast$exp.cc[is.na(yeast$exp.cc)] <- mean(as.matrix(yeast$exp.cc),na.rm=TRUE)
-#' #monsterRes <- monster(yeast$exp.cc[1:500,], design, yeast$motif, nullPerms=10, numMaxCores=1)
 #' # Example with provided networks
 #' \donttest{
 #' pandaResult <- panda(pandaToyData$motif, pandaToyData$expression, pandaToyData$ppi)
-#' case=getRegNet(pandaResult)
-#' nelemReg=dim(getRegNet(pandaResult))[1]*dim(getRegNet(pandaResult))[2]
-#' nGenes=length(colnames(getRegNet(pandaResult)))
+#' case=pandaResult@regNet
+#' nelemReg=dim(pandaResult@regNet)[1]*dim(pandaResult@regNet)[2]
+#' nGenes=length(colnames(pandaResult@regNet))
 #' control=matrix(rexp(nelemReg, rate=.1), ncol=nGenes)
 #' colnames(control) = colnames(case)
 #' rownames(control) = rownames(case) 
@@ -623,7 +622,7 @@ monsterTransitionNetworkPlot <- function(monsterObj, numEdges=100, numTopTFs=10,
 
   adj.combined <- adj.combined[
     abs(adj.combined[,4])>=sort(abs(adj.combined[,4]),decreasing=TRUE)[numEdges],]
-  tfNet <- graph.data.frame(adj.combined, directed=TRUE)
+  tfNet <- graph_from_data_frame(adj.combined, directed=TRUE)
   vSize <- -log(dTFI_pVals_All)
   vSize[vSize<0] <- 0
   vSize[vSize>3] <- 3
@@ -1165,9 +1164,9 @@ NULL
 #' pandaResult_control <- panda(pandaToyData$motif, pandaToyData$expression[,26:50], pandaToyData$ppi)
 #'
 #' # function takes both panda objects and matrices, or a mixture
-#' monster_res1 <- domonster(pandaResult_exp, pandaResult_control)
-#' monster_res2 <- domonster(pandaResult_exp@regNet, pandaResult_control@regNet)
-#' monster_res3 <- domonster(pandaResult_exp@regNet, pandaResult_control)
+#' monster_res1 <- domonster(pandaResult_exp, pandaResult_control, numMaxCores = 1)
+#' monster_res2 <- domonster(pandaResult_exp@regNet, pandaResult_control@regNet, numMaxCores = 1)
+#' monster_res3 <- domonster(pandaResult_exp@regNet, pandaResult_control, numMaxCores = 1)
 #' }
 domonster <- function(exp_graph, control_graph, nullPerms = 1000, numMaxCores = 3, ...){
   if('panda' %in% class(exp_graph)){

R/PANDA.R

@@ -169,7 +169,6 @@ pandaPy <- function(expr_file, motif_file=NULL, ppi_file=NULL, computing="cpu",
       colnames(panda_net) <- c("TF","Gene","Score")
     }
 
-
     # in-degree of panda network
     py_run_string(paste("indegree=panda_obj.return_panda_indegree()"))
     indegree_net <- py$indegree
@@ -187,14 +186,15 @@ pandaPy <- function(expr_file, motif_file=NULL, ppi_file=NULL, computing="cpu",
       panda_net$Gene <- paste('tar_', panda_net$Gene, sep='')
       message("Rename the content of first two columns with prefix 'reg_' and 'tar_' as there are some duplicate node names between the first two columns" )
     }
-
     output <- list("panda" = panda_net, "indegree" = indegree_net, "outdegree" = outdegree_net)
 
 
-  } else{ py_run_string("panda_network=panda_obj.panda_network",local = FALSE, convert = TRUE) 
-    panda_net <- py$panda_network
+  } else{ 
+    py_run_string("panda_network=panda_obj.panda_network",local = FALSE, convert = TRUE) 
+    panda_net_py <- py$panda_network
+
     # weighted adjacency matrix of PANDA network 
-    output <- list("WAMpanda" = panda_net)
+    output <- list("WAMpanda" = panda_net_py)
   }
 
   message ("...Finish PANDA...")

R/PUMA.R

@@ -60,7 +60,7 @@ puma <- function(motif,expr=NULL,ppi=NULL,alpha=0.1,mir_file,hamming=0.001,
   if(progress)
     print('Initializing and validating')
 
-  if(class(expr)=="ExpressionSet")
+  if(is(expr, "ExpressionSet"))
     expr <- assayData(expr)[["exprs"]]
 
   if (is.null(expr)){
@@ -431,3 +431,5 @@ prepResult <- function(zScale, output, regulatoryNetwork, geneCoreg, tfCoopNetwo
   pandaObj(regNet=regulatoryNetwork, coregNet=geneCoreg, coopNet=tfCoopNetwork, numGenes=numGenes, numTFs=numTFs, numEdges=numEdges)
 }
 
+# Define a new class called pandaObj.
+pandaObj <- setClass("panda", slots=c("regNet","coregNet","coopNet","numGenes","numTFs","numEdges"))

R/SPIDER.R

@@ -68,7 +68,7 @@ spider <- function(motif,expr=NULL,epifilter=NULL,ppi=NULL,alpha=0.1,hamming=0.0
     stop('Chromatin accessibility data does not match motif data size and order.')
   }
 
-  if(class(expr)=="ExpressionSet")
+  if(is(expr, "ExpressionSet"))
     expr <- assayData(expr)[["exprs"]]
 
   if (is.null(expr)){

R/YARN.R

@@ -20,8 +20,11 @@
 #' @importClassesFrom Biobase ExpressionSet
 #'
 #' @examples
-#' download.file('https://netzoo.s3.us-east-2.amazonaws.com/netZooR/unittest_datasets/yarn/bladder.rdata',destfile='netZooR/data/bladder.rdata')
-#' download.file('https://netzoo.s3.us-east-2.amazonaws.com/netZooR/unittest_datasets/yarn/skin.rdata',destfile='netZooR/data/skin.rdata')
+#' u <- 'https://netzoo.s3.us-east-2.amazonaws.com/netZooR/unittest_datasets/'
+#' bladder <- paste0(u, 'yarn/bladder.rdata')
+#' skin <- paste0(u, 'yarn/skin.rdata')
+#' download.file(bladder, destfile='netZooR/data/bladder.rdata')
+#' download.file(skin, destfile='netZooR/data/skin.rdata')
 #' data(skin)
 #' # subsetting and changing column name just for a silly example
 #' skin <- skin[1:10,]
@@ -235,7 +238,7 @@ downloadGTEx <- function(type = "genes", file = NULL, ...) {
 #' head(netZooR:::extractMatrix(skin,normalized=FALSE,log=FALSE))
 #'
 extractMatrix <- function(obj, normalized = FALSE, log = TRUE) {
-  if (class(obj) == "ExpressionSet") {
+  if (is(obj, "ExpressionSet")) {
     if (!normalized) {
       obj <- exprs(obj)
     } else {
@@ -561,8 +564,8 @@ plotDensity <- function(obj, groups = NULL, normalized = FALSE,
 #' # library(RColorBrewer)
 #' data(skin)
 #' tissues <- pData(skin)$SMTSD
-#' heatmapColColors <- brewer.pal(12,"Set3")[as.integer(factor(tissues))]
-#' heatmapCols <- colorRampPalette(brewer.pal(9, "RdBu"))(50)
+#' heatmapColColors <- RColorBrewer::brewer.pal(12,"Set3")[as.integer(factor(tissues))]
+#' heatmapCols <- colorRampPalette(RColorBrewer::brewer.pal(9, "RdBu"))(50)
 #' plotHeatmap(skin,normalized=FALSE,log=TRUE,trace="none",n=10,
 #'  col = heatmapCols,ColSideColors = heatmapColColors,cexRow = 0.6,cexCol = 0.6)
 #'}
@@ -734,10 +737,10 @@ qstats <- function(exprs, groups, window) {
 #' (\href{http://www.ncbi.nlm.nih.gov/pmc/articles/PMC4547484/}{PubMed})
 #'
 #' @source GTEx Portal
-#' @name Skin_data
+#' @name skin
 #' @examples
 #' \donttest{data(skin);
-#' checkMissAnnotation(skin,"GENDER");}
+#' checkMisAnnotation(skin,"GENDER");}
 system('wget https://netzoo.s3.us-east-2.amazonaws.com/netZooR/unittest_datasets/yarn/skin.rdata')
 system('mv skin.rdata data/')
 "skin"
@@ -762,10 +765,10 @@ system('mv skin.rdata data/')
 #' @source GTEx Portal
 #'
 #' @return ExpressionSet object
-#' @name Bladder_data
+#' @name bladder
 #' @examples
 #' \donttest{data(bladder);
-#' checkMissAnnotation(bladder);}
+#' checkMisAnnotation(bladder, "GENDER");}
 system('wget https://netzoo.s3.us-east-2.amazonaws.com/netZooR/unittest_datasets/yarn/bladder.rdata')
 system('mv bladder.rdata data/')
 "bladder"