diff netboxr_r.R @ 0:785ed8621503 draft default tip

planemo upload commit 13e28726551f8751007f28a3c6cd5d00ba278056

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/netboxr_r.R	Wed Aug 24 08:27:40 2022 +0000
@@ -0,0 +1,117 @@
+# Set up R error handling to go to stderr
+options(show.error.messages = FALSE,
+        error = function() {
+        cat(geterrmessage(), file = stderr())
+        q("no", 1, FALSE)})
+# Avoid crashing Galaxy with an UTF8 error on German LC settings
+loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8")
+# Import required libraries and data
+suppressPackageStartupMessages({
+  library(netboxr)
+  library(igraph)
+  library(RColorBrewer)
+})
+
+data(netbox2010)
+args <- commandArgs(TRUE)
+# Vars
+gene_list <- scan(args[2], what = character(), sep = "\n")
+cutoff <- args[4]
+community <- args[6]
+global_model <- args[8]
+global_iterations <- args[10]
+global_number <- args[12]
+local_model <- args[14]
+local_iterations <- args[16]
+
+network_plot <- args[18]
+plot_width <- args[20]
+output_sif <- args[22]
+neighbor_list <- args[24]
+modmem <- args[26]
+nt <- args[28]
+
+sink("metadata.txt")
+sink(stdout(), type = "message")
+# Network analysis as described in netboxr vignette
+sif_network <- netbox2010$network
+graph_reduced <- networkSimplify(sif_network, directed = FALSE)
+threshold <- cutoff
+results <- print(geneConnector(geneList = gene_list, networkGraph = graph_reduced,
+                         directed = FALSE, pValueAdj = "BH", pValueCutoff = threshold,
+                         communityMethod = community, keepIsolatedNodes = FALSE))
+
+# Check the p-value of the selected linker
+linker_df <- results$neighborData
+linker_df[linker_df$pValueFDR < threshold, ]
+graph_layout <- layout_with_fr(results$netboxGraph)
+
+# Global Network Null Model
+if (global_model) {
+  global_test <- globalNullModel(netboxGraph = results$netboxGraph, networkGraph = graph_reduced,
+                                iterations = global_iterations, numOfGenes = global_number)
+  global_test
+}
+
+# Local Network Null Model
+if (local_model) {
+  local_test <- localNullModel(netboxGraph = results$netboxGraph, iterations = local_iterations)
+  local_test
+}
+
+## Output
+# Plot the edge annotated graph
+if (network_plot) {
+
+  edges <- results$netboxOutput
+  interaction_type <- unique(edges[, 2])
+  interaction_type_color <- brewer.pal(length(interaction_type), name = "Spectral")
+  edge_colors <- data.frame(interaction_type, interaction_type_color, stringsAsFactors = FALSE)
+  colnames(edge_colors) <- c("INTERACTION_TYPE", "COLOR")
+  netbox_graph_annotated <- annotateGraph(netboxResults = results, edgeColors =
+                                          edge_colors, directed = FALSE, linker = TRUE)
+  pdf("network_plot.pdf", width = plot_width)
+  plot(results$netboxCommunity, netbox_graph_annotated, layout = graph_layout,
+       vertex.size = 10, vertex.shape = V(netbox_graph_annotated)$shape, edge.color
+       = E(netbox_graph_annotated)$interactionColor, edge.width = 3)
+
+  # Add interaction type annotations
+  legend(x = -1.8, y = -1, legend = interaction_type, col =
+           interaction_type_color, lty = 1, lwd = 2, bty = "n", cex = 1)
+  dev.off()
+}
+
+# Local Network Null Model
+if (local_model) {
+  pdf("localModel_histogram.pdf")
+  h <- hist(local_test$randomModularityScore, breaks = 35, plot = FALSE)
+  h$density <- h$counts / sum(h$counts)
+  plot(h, freq = FALSE, ylim = c(0, 0.1), xlim = c(0.1, 0.6), col = "lightblue")
+  abline(v = local_test$modularityScoreObs, col = "red")
+  dev.off()
+}
+
+# NetBox algorithm output in SIF format.
+if (output_sif) {
+  write.table(results$netboxOutput, file = "network.sif", sep = "\t", quote = FALSE,
+              col.names = FALSE, row.names = FALSE)
+}
+
+# Save neighbor data
+if (neighbor_list) {
+  write.table(results$neighborData, file = "neighbor_data.txt", sep = "\t",
+              quote = FALSE, col.names = TRUE, row.names = FALSE)
+}
+
+#Save identified pathway module numbers
+if (modmem) {
+ write.table(results$moduleMembership, file = "community.membership.txt", sep = "\t",
+              quote = FALSE, col.names = FALSE, row.names = FALSE)
+}
+
+# Save file that indicates whether the node is a 'linker' or 'candidate'
+if (nt) {
+  write.table(results$nodeType, file = "nodeType.txt", sep = "\t", quote = FALSE, col.names = FALSE,
+              row.names = FALSE)
+}
+sink(NULL)