Mercurial > repos > artbio > pathifier
changeset 1:0960bd1161fa draft default tip
planemo upload for repository https://github.com/ARTbio/tools-artbio/tree/main/tools/pathifier commit f7363f7c391bd501fc4d1c10aaf372ec2bd82732
| author | artbio |
|---|---|
| date | Mon, 12 Feb 2024 23:57:01 +0000 |
| parents | fec313f5c889 |
| children | |
| files | pathifier.R pathifier.xml |
| diffstat | 2 files changed, 268 insertions(+), 224 deletions(-) [+] |
line wrap: on
line diff
--- a/pathifier.R Mon Apr 12 09:55:24 2021 +0000 +++ b/pathifier.R Mon Feb 12 23:57:01 2024 +0000 @@ -4,10 +4,10 @@ ################################################################################################## -options(show.error.messages = F, error = function() { - cat(geterrmessage(), file = stderr()); q("no", 1, F) - } -) +options(show.error.messages = FALSE, error = function() { + cat(geterrmessage(), file = stderr()) + q("no", 1, FALSE) +}) # we need that to not crash galaxy with an UTF8 error on German LC settings. loc <- Sys.setlocale("LC_MESSAGES", "en_US.UTF-8") @@ -18,106 +18,112 @@ library(circlize) option_list <- list( - make_option( - "--exp", - type = "character", - help = "Expression matrix"), - make_option( - "--sep", - type = "character", - default = "\t", - help = "File separator [default : '%default']" - ), - make_option( - "--genes", - type = "character", - help = "Gene sets Pathways : gmt format (one pathway per line : Name, description, genes (one by column), tab separated)"), - make_option( - "--is_normal", - default = F, - type = "logical", - help = "Define normals cells in your data"), - make_option( - "--normals", - type = "character", - help = "A vector of sample status : 1 = Healthy, 0 = Tumor. Must be in the same order as in expression data"), - make_option( - "--logfile", - type = "character", - default = "log.txt", - help = "Log file name [default : '%default']" - ), - make_option( - "--max_stability", - type = "logical", - default = T, - help = "If true, throw away components leading to low stability of sampling noise [default : '%default']" - ), - make_option( - "--attempts", - type = "integer", - default = 10, - help = "Number of runs to determine stability. [default : '%default']" - ), - make_option( - "--min_std", - type = "character", - default = "0.4", - help = "Minimum of standard deviation to filter out low variable genes. + make_option( + "--exp", + type = "character", + help = "Expression matrix" + ), + make_option( + "--sep", + type = "character", + default = "\t", + help = "File separator [default : '%default']" + ), + make_option( + "--genes", + type = "character", + help = "Gene sets Pathways : gmt format (one pathway per line : Name, description, genes (one by column), tab separated)" + ), + make_option( + "--is_normal", + default = FALSE, + type = "logical", + help = "Define normals cells in your data" + ), + make_option( + "--normals", + type = "character", + help = "A vector of sample status : 1 = Healthy, 0 = Tumor. Must be in the same order as in expression data" + ), + make_option( + "--logfile", + type = "character", + default = "log.txt", + help = "Log file name [default : '%default']" + ), + make_option( + "--max_stability", + type = "logical", + default = TRUE, + help = "If true, throw away components leading to low stability of sampling noise [default : '%default']" + ), + make_option( + "--attempts", + type = "integer", + default = 10, + help = "Number of runs to determine stability. [default : '%default']" + ), + make_option( + "--min_std", + type = "character", + default = "0.4", + help = "Minimum of standard deviation to filter out low variable genes. Use --min.std data, to use the minimum std of your data [default : '%default']" - ), - make_option( - "--min_exp", - type = "character", - default = "4", - help = "Minimum of gene expression to filter out low expressed genes. + ), + make_option( + "--min_exp", + type = "character", + default = "4", + help = "Minimum of gene expression to filter out low expressed genes. Use --min.exp data, to use the minimum expression of your data [default : '%default']" - ), - make_option( - "--pds", - type = "character", - default = "PDS.tsv", - help = "Output PDS (Pathway deregulation score) of Pathifier in tabular file [default : '%default']" - ), - make_option( - "--heatmap_cluster_cells", - type = "logical", - default = TRUE, - help = "Cluster columns (cells) in the heatmap [default : '%default']" - ), - make_option( - "--heatmap_cluster_pathways", - type = "logical", - default = TRUE, - help = "Cluster rows (pathways) in the heatmap [default : '%default']" - ), - make_option( - "--heatmap_show_cell_labels", - type = "logical", - default = FALSE, - help = "Print column names (cells) on the heatmap [default : '%default']" - ), - make_option( - "--heatmap_show_pathway_labels", - type = "logical", - default = FALSE, - help = "Print row names (pathways) on the heatmap [default : '%default']" - ), - make_option( - "--plot", - type = "character", - default = "./plot.pdf", - help = "Pathifier visualization [default : '%default']" - ), - make_option( - "--rdata", - type = "character", - default = "./results.rdata", - help = "Pathifier object (S4) [default : '%default']")) + ), + make_option( + "--pds", + type = "character", + default = "PDS.tsv", + help = "Output PDS (Pathway deregulation score) of Pathifier in tabular file [default : '%default']" + ), + make_option( + "--heatmap_cluster_cells", + type = "logical", + default = TRUE, + help = "Cluster columns (cells) in the heatmap [default : '%default']" + ), + make_option( + "--heatmap_cluster_pathways", + type = "logical", + default = TRUE, + help = "Cluster rows (pathways) in the heatmap [default : '%default']" + ), + make_option( + "--heatmap_show_cell_labels", + type = "logical", + default = FALSE, + help = "Print column names (cells) on the heatmap [default : '%default']" + ), + make_option( + "--heatmap_show_pathway_labels", + type = "logical", + default = FALSE, + help = "Print row names (pathways) on the heatmap [default : '%default']" + ), + make_option( + "--plot", + type = "character", + default = "./plot.pdf", + help = "Pathifier visualization [default : '%default']" + ), + make_option( + "--rdata", + type = "character", + default = "./results.rdata", + help = "Pathifier object (S4) [default : '%default']" + ) +) parser <- OptionParser(usage = "%prog [options] file", option_list = option_list) args <- parse_args(parser) if (args$sep == "tab") { - args$sep <- "\t" + args$sep <- "\t" } @@ -125,11 +131,13 @@ set.seed(123) # Load expression data for PATHIFIER -exp_matrix <- as.matrix(read.delim(file = args$exp, - as.is = T, - row.names = 1, - sep = args$sep, - check.names = F)) +exp_matrix <- as.matrix(read.delim( + file = args$exp, + as.is = TRUE, + row.names = 1, + sep = args$sep, + check.names = FALSE +)) # Load Genesets annotation gene_sets_file <- file(args$genes, open = "r") @@ -146,19 +154,23 @@ # Prepare data and parameters ################################################## # Extract information from binary phenotypes. 1 = Normal, 0 = Tumor -if (args$is_normal == T) { - normals <- read.delim(file = args$normals, h = F) - normals <- as.logical(normals[, 2]) - n_exp_matrix <- exp_matrix[, normals] -} else n_exp_matrix <- exp_matrix +if (args$is_normal == TRUE) { + normals <- read.delim(file = args$normals, h = FALSE) + normals <- as.logical(normals[, 2]) + n_exp_matrix <- exp_matrix[, normals] +} else { + n_exp_matrix <- exp_matrix +} # Calculate MIN_STD rsd <- apply(n_exp_matrix, 1, sd) min_std <- quantile(rsd, 0.25) # Calculate MIN_EXP -min_exp <- quantile(as.vector(as.matrix(exp_matrix)), - 0.1) # Percentile 10 of data +min_exp <- quantile( + as.vector(as.matrix(exp_matrix)), + 0.1 +) # Percentile 10 of data # Filter low value genes. At least 10% of samples with values over min_exp # Set expression levels < MIN_EXP to MIN_EXP @@ -169,19 +181,23 @@ exp_matrix_filtered[exp_matrix_filtered < min_exp] <- min_exp # Set maximum 5000 genes with more variance -variable_genes <- names(sort(apply(exp_matrix_filtered, 1, var), decreasing = T))[1:5000] +variable_genes <- names(sort(apply(exp_matrix_filtered, 1, var), decreasing = TRUE))[1:5000] variable_genes <- variable_genes[!is.na(variable_genes)] exp_matrix_filtered <- exp_matrix_filtered[variable_genes, ] allgenes <- as.vector(rownames(exp_matrix_filtered)) if (args$min_std == "data") { - args$min_std <- min_std -} else args$min_std <- as.numeric(args$min_std) + args$min_std <- min_std +} else { + args$min_std <- as.numeric(args$min_std) +} if (args$min_exp == "data") { - args$min_exp <- min_exp -} else args$min_exp <- as.numeric(args$min_exp) + args$min_exp <- min_exp +} else { + args$min_exp <- as.numeric(args$min_exp) +} # Open pdf @@ -191,89 +207,108 @@ col_score_fun <- colorRamp2(c(0, 0.5, 1), c("#4575B4", "#FFFFBF", "#D73027")) # Run Pathifier -if (args$is_normal == T) { - pds <- quantify_pathways_deregulation(exp_matrix_filtered, - allgenes, - gs, - pathwaynames, - normals, - maximize_stability = args$max_stability, - attempts = args$attempts, - logfile = args$logfile, - min_std = args$min_std, - min_exp = args$min_exp) - for (i in pathwaynames) { - df <- data.frame(pds$curves[[i]][, 1:3], - normal = normals, - PDS = as.numeric(pds$scores[[i]]), - curve_order = as.vector(pds$curves_order[[i]])) - ordered <- df[df$curve_order, ] +if (args$is_normal == TRUE) { + pds <- quantify_pathways_deregulation(exp_matrix_filtered, + allgenes, + gs, + pathwaynames, + normals, + maximize_stability = args$max_stability, + attempts = args$attempts, + logfile = args$logfile, + min_std = args$min_std, + min_exp = args$min_exp + ) + for (i in pathwaynames) { + df <- data.frame(pds$curves[[i]][, 1:3], + normal = normals, + PDS = as.numeric(pds$scores[[i]]), + curve_order = as.vector(pds$curves_order[[i]]) + ) + ordered <- df[df$curve_order, ] - layout(cbind(1:2, 1:2), heights = c(7, 1)) - sc3 <- scatterplot3d(ordered[, 1:3], - main = paste("Principal curve of", i), - box = F, pch = 19, type = "l") - sc3$points3d(ordered[, 1:3], box = F, pch = 19, - col = col_score_fun(ordered$PDS)) + layout(cbind(1:2, 1:2), heights = c(7, 1)) + sc3 <- scatterplot3d(ordered[, 1:3], + main = paste("Principal curve of", i), + box = FALSE, pch = 19, type = "l" + ) + sc3$points3d(ordered[, 1:3], + box = FALSE, pch = 19, + col = col_score_fun(ordered$PDS) + ) - # Plot color scale legend - par(mar = c(5, 3, 0, 3)) - plot(seq(min(ordered$PDS), max(ordered$PDS), length = 100), rep(0, 100), pch = 15, - axes = TRUE, yaxt = "n", xlab = "Color scale of PDS", ylab = "", bty = "n", - col = col_score_fun(seq(min(ordered$PDS), max(ordered$PDS), length = 100))) + # Plot color scale legend + par(mar = c(5, 3, 0, 3)) + plot(seq(min(ordered$PDS), max(ordered$PDS), length = 100), rep(0, 100), + pch = 15, + axes = TRUE, yaxt = "n", xlab = "Color scale of PDS", ylab = "", bty = "n", + col = col_score_fun(seq(min(ordered$PDS), max(ordered$PDS), length = 100)) + ) - cols_status <- ifelse(ordered$normal, "blue", "red") - sc3 <- scatterplot3d(ordered[, 1:3], - main = paste("Principal curve of", i), - box = F, pch = "", type = "l") - sc3$points3d(ordered[, 1:3], box = F, - pch = ifelse(ordered$normal, 19, 8), - col = ifelse(ordered$normal, "blue", "red")) - legend("topright", pch = c(19, 8), yjust = 0, - legend = c("normal", "cancer"), - col = c("blue", "red"), cex = 1.1) + cols_status <- ifelse(ordered$normal, "blue", "red") + sc3 <- scatterplot3d(ordered[, 1:3], + main = paste("Principal curve of", i), + box = FALSE, pch = "", type = "l" + ) + sc3$points3d(ordered[, 1:3], + box = FALSE, + pch = ifelse(ordered$normal, 19, 8), + col = ifelse(ordered$normal, "blue", "red") + ) + legend("topright", + pch = c(19, 8), yjust = 0, + legend = c("normal", "cancer"), + col = c("blue", "red"), cex = 1.1 + ) - ## annotation for heatmap - sample_status <- data.frame(Status = factor(ifelse(df$normal, "normal", "tumor"))) - rownames(sample_status) <- colnames(exp_matrix_filtered) - color_status_heatmap <- list(Status = c(normal = "blue", tumor = "red")) - } -} else{ - pds <- quantify_pathways_deregulation(exp_matrix_filtered, - allgenes, - gs, - pathwaynames, - maximize_stability = args$max_stability, - attempts = args$attempts, - logfile = args$logfile, - min_std = args$min_std, - min_exp = args$min_exp) - for (i in pathwaynames) { - df <- data.frame(pds$curves[[i]][, 1:3], - PDS = as.numeric(pds$scores[[i]]), - curve_order = as.vector(pds$curves_order[[i]])) - ordered <- df[df$curve_order, ] + ## annotation for heatmap + sample_status <- data.frame(Status = factor(ifelse(df$normal, "normal", "tumor"))) + rownames(sample_status) <- colnames(exp_matrix_filtered) + color_status_heatmap <- list(Status = c(normal = "blue", tumor = "red")) + } +} else { + pds <- quantify_pathways_deregulation(exp_matrix_filtered, + allgenes, + gs, + pathwaynames, + maximize_stability = args$max_stability, + attempts = args$attempts, + logfile = args$logfile, + min_std = args$min_std, + min_exp = args$min_exp + ) + for (i in pathwaynames) { + df <- data.frame(pds$curves[[i]][, 1:3], + PDS = as.numeric(pds$scores[[i]]), + curve_order = as.vector(pds$curves_order[[i]]) + ) + ordered <- df[df$curve_order, ] - layout(cbind(1:2, 1:2), heights = c(7, 1)) - sc3 <- scatterplot3d(ordered[, 1:3], - main = paste("Principal curve of", i), - box = F, pch = 19, type = "l") - sc3$points3d(ordered[, 1:3], box = F, pch = 19, - col = col_score_fun(ordered$PDS)) + layout(cbind(1:2, 1:2), heights = c(7, 1)) + sc3 <- scatterplot3d(ordered[, 1:3], + main = paste("Principal curve of", i), + box = FALSE, pch = 19, type = "l" + ) + sc3$points3d(ordered[, 1:3], + box = FALSE, pch = 19, + col = col_score_fun(ordered$PDS) + ) - # Plot color scale legend - par(mar = c(5, 3, 0, 3)) - plot(seq(min(ordered$PDS), max(ordered$PDS), length = 100), rep(0, 100), pch = 15, - axes = TRUE, yaxt = "n", xlab = "Color scale of PDS", ylab = "", bty = "n", - col = col_score_fun(seq(min(ordered$PDS), max(ordered$PDS), length = 100))) + # Plot color scale legend + par(mar = c(5, 3, 0, 3)) + plot(seq(min(ordered$PDS), max(ordered$PDS), length = 100), rep(0, 100), + pch = 15, + axes = TRUE, yaxt = "n", xlab = "Color scale of PDS", ylab = "", bty = "n", + col = col_score_fun(seq(min(ordered$PDS), max(ordered$PDS), length = 100)) + ) - ## annotation for heatmap (for the moment none for this situation) - sample_status <- NA - color_status_heatmap <- NA - } + ## annotation for heatmap (for the moment none for this situation) + sample_status <- NA + color_status_heatmap <- NA + } } ## Create dataframe from Pathifier list and round score to 4 digits @@ -283,29 +318,31 @@ ## plot heatmap if (ncol(pds_scores) > 1) { pheatmap(t(pds_scores), - main = "Heatmap of Pathway Deregulation Score", # heat map title - cluster_rows = args$heatmap_cluster_pathways, # apply clustering method - cluster_cols = args$heatmap_cluster_cells, # apply clustering method + main = "Heatmap of Pathway Deregulation Score", # heat map title + cluster_rows = args$heatmap_cluster_pathways, # apply clustering method + cluster_cols = args$heatmap_cluster_cells, # apply clustering method - #Additional Options - ## Color labeled columns - annotation_col = sample_status, - annotation_colors = color_status_heatmap, - show_rownames = args$heatmap_show_pathway_labels, - show_colnames = args$heatmap_show_cell_labels, - border_color = NA, - legend = TRUE) + # Additional Options + ## Color labeled columns + annotation_col = sample_status, + annotation_colors = color_status_heatmap, + show_rownames = args$heatmap_show_pathway_labels, + show_colnames = args$heatmap_show_cell_labels, + border_color = NA, + legend = TRUE + ) } dev.off() ## write table write.table(pds_scores, - args$pds, - row.names = T, - col.names = T, - quote = F, - sep = "\t") + args$pds, + row.names = TRUE, + col.names = TRUE, + quote = FALSE, + sep = "\t" +) ## write S4 pathifier object save(pds, file = args$rdata)
--- a/pathifier.xml Mon Apr 12 09:55:24 2021 +0000 +++ b/pathifier.xml Mon Feb 12 23:57:01 2024 +0000 @@ -1,23 +1,28 @@ -<tool id="pathifier" name="Pathifier" version="1.0.1"> +<tool id="pathifier" name="Pathifier" version="@TOOL_VERSION@+galaxy@VERSION_SUFFIX@" profile="@PROFILE@"> <description>: Quantify deregulation of pathways in cancer</description> + <macros> + <token name="@TOOL_VERSION@">1.40.0</token> + <token name="@VERSION_SUFFIX@">0</token> + <token name="@PROFILE@">23.0</token> + </macros> <requirements> - <requirement type="package" version="1.6.2=r35h6115d3f_0">r-optparse</requirement> - <requirement type="package" version="1.22.0=r351_0">bioconductor-pathifier</requirement> - <requirement type="package" version="1.0.12=r351h6115d3f_0">r-pheatmap</requirement> - <requirement type="package" version="0.3_41=r351h6115d3f_0">r-scatterplot3d</requirement> - <requirement type="package" version="0.4.7=r35h6115d3f_0">r-circlize</requirement> + <requirement type="package" version="@TOOL_VERSION@">bioconductor-pathifier</requirement> + <requirement type="package" version="1.7.4">r-optparse</requirement> + <requirement type="package" version="1.0.12">r-pheatmap</requirement> + <requirement type="package" version="0.3_44">r-scatterplot3d</requirement> + <requirement type="package" version="0.4.15">r-circlize</requirement> </requirements> <stdio> <exit_code range="1:" level="fatal" description="Tool exception" /> </stdio> <command detect_errors="exit_code"><![CDATA[ - Rscript $__tool_directory__/pathifier.R + Rscript '$__tool_directory__/pathifier.R' --exp '$input' --sep '$input_sep' --genes '$genes' - #if $reference.is_normal == "Yes": - --is_normal 'TRUE' + #if str($reference.reference_selector) == 'TRUE': + --is_normal '$reference.reference_selector' --normals '$reference.normals' #end if @@ -35,7 +40,6 @@ --logfile '$logfile' --plot '$plot' --rdata '$rdatafile' - ]]></command> <inputs> <param name="input" type="data" format="txt,tabular" label="expression data"/> @@ -46,9 +50,13 @@ <param name="genes" type="data" format="txt" label="Gene sets Pathways" help="Must be in gmt format (one pathway per line : Name, description, genes (one by column), tab separated)" /> <conditional name="reference"> - <param name="is_normal" label="Do you have non cancer transcriptomes in your data set ?" type="boolean" truevalue="Yes" falsevalue="" checked="false" - help="If set the starting curve depends on the matrix of points with in a certain row order (first 'normal' then 'cancer' samples), otherwise the first principal component is used. See help for more informations"/> - <when value="Yes"> + <param name="reference_selector" label="Do you have non cancer transcriptomes in your data set ?" type="select" + help="Yes, if the starting curve depends on the matrix of points with in a certain row order + (first 'normal' then 'cancer' samples), otherwise (No) the first principal component is used. See help section for more informations"> + <option selected="True" value="TRUE">Yes</option> + <option value="">No</option> + </param> + <when value="TRUE"> <param name="normals" type="data" format="tabular" label="Sample status" help="A two-column data frame, first column contains data labels, second column the levels of sample status : 1 = Healthy, 0 = Tumor (no header)" /> </when> @@ -74,7 +82,6 @@ <option value="no" selected="True">No</option> <option value="yes">Yes</option> </param> - </inputs> <outputs> <data name="pds" format="tabular" label="Pathifier Deregulation Score (PDS) of ${on_string}" /> @@ -83,14 +90,14 @@ </data> <data name="plot" format="pdf" label="Pathifier vizualization of ${on_string}" /> <data name="rdatafile" format="rdata" label="Pathifier S4 object of ${on_string}" > - <filter>rds == 'yes'</filter> + <filter>rdata == 'yes'</filter> </data> </outputs> <tests> - <test> + <test expect_num_outputs="3"> <param name="input" value="sheffer.tsv" ftype="tabular"/> <param name="genes" value="kegg_pathways.gmt" ftype="txt" /> - <param name="is_normal" value="True" /> + <param name="reference_selector" value="TRUE" /> <param name="normals" value="normals.tsv" ftype="tabular" /> <param name="log" value="yes" /> <param name="attempts" value="100" /> @@ -98,10 +105,10 @@ <output name="pds" file="sheffer.kegg.tsv" ftype="tabular"/> <output name="plot" file="plot.pdf" ftype="pdf" compare="sim_size" /> </test> - <test> + <test expect_num_outputs="3"> <param name="input" value="sheffer_noref.tsv" ftype="tabular"/> <param name="genes" value="kegg_pathways.gmt" ftype="txt" /> - <param name="is_normal" value="" /> + <param name="reference_selector" value="" /> <param name="log" value="no" /> <param name="rdata" value="yes" /> <param name="attempts" value="50" />