Mercurial > repos > bgruening > netboxr
comparison netboxr_r.R @ 0:785ed8621503 draft default tip
planemo upload commit 13e28726551f8751007f28a3c6cd5d00ba278056
author | bgruening |
---|---|
date | Wed, 24 Aug 2022 08:27:40 +0000 |
parents | |
children |
comparison
equal
deleted
inserted
replaced
-1:000000000000 | 0:785ed8621503 |
---|---|
1 # Set up R error handling to go to stderr | |
2 options(show.error.messages = FALSE, | |
3 error = function() { | |
4 cat(geterrmessage(), file = stderr()) | |
5 q("no", 1, FALSE)}) | |
6 # Avoid crashing Galaxy with an UTF8 error on German LC settings | |
7 loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") | |
8 # Import required libraries and data | |
9 suppressPackageStartupMessages({ | |
10 library(netboxr) | |
11 library(igraph) | |
12 library(RColorBrewer) | |
13 }) | |
14 | |
15 data(netbox2010) | |
16 args <- commandArgs(TRUE) | |
17 # Vars | |
18 gene_list <- scan(args[2], what = character(), sep = "\n") | |
19 cutoff <- args[4] | |
20 community <- args[6] | |
21 global_model <- args[8] | |
22 global_iterations <- args[10] | |
23 global_number <- args[12] | |
24 local_model <- args[14] | |
25 local_iterations <- args[16] | |
26 | |
27 network_plot <- args[18] | |
28 plot_width <- args[20] | |
29 output_sif <- args[22] | |
30 neighbor_list <- args[24] | |
31 modmem <- args[26] | |
32 nt <- args[28] | |
33 | |
34 sink("metadata.txt") | |
35 sink(stdout(), type = "message") | |
36 # Network analysis as described in netboxr vignette | |
37 sif_network <- netbox2010$network | |
38 graph_reduced <- networkSimplify(sif_network, directed = FALSE) | |
39 threshold <- cutoff | |
40 results <- print(geneConnector(geneList = gene_list, networkGraph = graph_reduced, | |
41 directed = FALSE, pValueAdj = "BH", pValueCutoff = threshold, | |
42 communityMethod = community, keepIsolatedNodes = FALSE)) | |
43 | |
44 # Check the p-value of the selected linker | |
45 linker_df <- results$neighborData | |
46 linker_df[linker_df$pValueFDR < threshold, ] | |
47 graph_layout <- layout_with_fr(results$netboxGraph) | |
48 | |
49 # Global Network Null Model | |
50 if (global_model) { | |
51 global_test <- globalNullModel(netboxGraph = results$netboxGraph, networkGraph = graph_reduced, | |
52 iterations = global_iterations, numOfGenes = global_number) | |
53 global_test | |
54 } | |
55 | |
56 # Local Network Null Model | |
57 if (local_model) { | |
58 local_test <- localNullModel(netboxGraph = results$netboxGraph, iterations = local_iterations) | |
59 local_test | |
60 } | |
61 | |
62 ## Output | |
63 # Plot the edge annotated graph | |
64 if (network_plot) { | |
65 | |
66 edges <- results$netboxOutput | |
67 interaction_type <- unique(edges[, 2]) | |
68 interaction_type_color <- brewer.pal(length(interaction_type), name = "Spectral") | |
69 edge_colors <- data.frame(interaction_type, interaction_type_color, stringsAsFactors = FALSE) | |
70 colnames(edge_colors) <- c("INTERACTION_TYPE", "COLOR") | |
71 netbox_graph_annotated <- annotateGraph(netboxResults = results, edgeColors = | |
72 edge_colors, directed = FALSE, linker = TRUE) | |
73 pdf("network_plot.pdf", width = plot_width) | |
74 plot(results$netboxCommunity, netbox_graph_annotated, layout = graph_layout, | |
75 vertex.size = 10, vertex.shape = V(netbox_graph_annotated)$shape, edge.color | |
76 = E(netbox_graph_annotated)$interactionColor, edge.width = 3) | |
77 | |
78 # Add interaction type annotations | |
79 legend(x = -1.8, y = -1, legend = interaction_type, col = | |
80 interaction_type_color, lty = 1, lwd = 2, bty = "n", cex = 1) | |
81 dev.off() | |
82 } | |
83 | |
84 # Local Network Null Model | |
85 if (local_model) { | |
86 pdf("localModel_histogram.pdf") | |
87 h <- hist(local_test$randomModularityScore, breaks = 35, plot = FALSE) | |
88 h$density <- h$counts / sum(h$counts) | |
89 plot(h, freq = FALSE, ylim = c(0, 0.1), xlim = c(0.1, 0.6), col = "lightblue") | |
90 abline(v = local_test$modularityScoreObs, col = "red") | |
91 dev.off() | |
92 } | |
93 | |
94 # NetBox algorithm output in SIF format. | |
95 if (output_sif) { | |
96 write.table(results$netboxOutput, file = "network.sif", sep = "\t", quote = FALSE, | |
97 col.names = FALSE, row.names = FALSE) | |
98 } | |
99 | |
100 # Save neighbor data | |
101 if (neighbor_list) { | |
102 write.table(results$neighborData, file = "neighbor_data.txt", sep = "\t", | |
103 quote = FALSE, col.names = TRUE, row.names = FALSE) | |
104 } | |
105 | |
106 #Save identified pathway module numbers | |
107 if (modmem) { | |
108 write.table(results$moduleMembership, file = "community.membership.txt", sep = "\t", | |
109 quote = FALSE, col.names = FALSE, row.names = FALSE) | |
110 } | |
111 | |
112 # Save file that indicates whether the node is a 'linker' or 'candidate' | |
113 if (nt) { | |
114 write.table(results$nodeType, file = "nodeType.txt", sep = "\t", quote = FALSE, col.names = FALSE, | |
115 row.names = FALSE) | |
116 } | |
117 sink(NULL) |