modified RNA-seq.R in R (525104d3) · Commits · Chaos / chaos.tools

R/RNA-seq.R

+41 −47

Original line number	Diff line number	Diff line
		@@ -297,6 +297,7 @@ de.gene <- function(

		c.volcano <- function(
		.data,
		compare = NULL,
		top_gene_number = 10,
		fold_change = 2,
		p_value = 0.05,
		@@ -323,7 +324,14 @@ c.volcano <- function(
		font = "serif",
		extend_size = 3
		){
		n_compare <- .data$compare \|> unique() \|> length()
		all_compare <- .data$compare \|> unique()
		selected_compare <- compare %>% .[. %in% all_compare]

		if(length(selected_compare) == 0){
		selected_compare <- all_compare
		}

		n_compare <- selected_compare \|> length()

		sub_dt <- rbind(
		.data[gene_name %in% sub_gene],
		@@ -355,6 +363,20 @@ c.volcano <- function(
		regulation == "down" & gene_name %in% sub_gene & gene_id %in% sub_id, "down",
		default = "not_sig"
		)
		][
		compare %in% selected_compare
		][
		,color := sub_regulation
		][
		order(compare,sub_regulation,-abs(log2FC))
		][
		,rank := 1:.N,.(compare,sub_regulation)
		][
		,label_up := fcase(rank %in% 1:top_gene_number & sub_regulation == "up",gene_name)
		][
		,label_down := fcase(rank %in% 1:top_gene_number & sub_regulation == "down",gene_name)
		][
		,label_special := fcase(gene_name %in% special_gene,gene_name)
		]

		stat <- merge(
		@@ -371,36 +393,6 @@ c.volcano <- function(
		)
		]

		top.genes <- c.comb_dt(compare,c("up","down")) \|>
		apply(
		1,
		\(x){
		dt[
		compare == x[1] & sub_regulation == x[2]
		][
		order(-abs(log2FC))
		][
		0:top_gene_number,
		gene_name
		] \|>
		as.data.table() %>%
		.[,.(gene_name = V1,sub_regulation = x[2],compare = x[1])]
		}
		)\|>
		rbindlist()

		dt[
		,color := sub_regulation
		][
		gene_name %in% top.genes$up,
		label_up := gene_name
		][
		gene_name %in% top.genes$down,
		label_down := gene_name
		][
		gene_name %in% special_gene,
		label_special := gene_name
		]

		if(!is.null(special_gene)){
		dt[
		@@ -413,36 +405,40 @@ c.volcano <- function(
		label_all <- c("Down","Not significant","Up")
		color_all <- c(color_down,color_ns,color_up)
		}else{
		label_all <- c(special_label,"Down","Not significant","Up")
		color_all <- c(color_special,color_down,color_ns,color_up)
		label_all <- c("Down","Not significant",special_label,"Up")
		color_all <- c(color_down,color_ns,color_special,color_up)
		}

		plot <- ggplot() +
		geom_point(
		data = dt[color == "not_sig"],
		aes(log2FC,log10P,color = color_ns),
		aes(log2FC,log10P),
		color = color_ns,
		size = 0.8,
		show.legend = show_legend
		) +
		geom_point(
		data = dt[color == "up"],
		aes(log2FC,log10P,color = color_up),
		aes(log2FC,log10P),
		color = color_up,
		size = 0.8,
		show.legend = show_legend
		) +
		geom_point(
		data = dt[color == "down"],
		aes(log2FC,log10P,color = color_down),
		aes(log2FC,log10P),
		color = color_down,
		size = 0.8,
		show.legend = show_legend
		) +
		geom_point(
		data = dt[color == special_label],
		aes(log2FC,log10P,color = color_special),
		aes(log2FC,log10P),
		color = color_special,
		size = 0.8,
		show.legend = show_legend
		) +
		scale_color_manual(values = color_all, label = label_all)
		)
		facet_wrap(~ compare)

		xlimits <- max(abs(min(dt$log2FC)),abs(max(dt$log2FC))) + extend_size

		@@ -589,19 +585,17 @@ c.volcano <- function(
		if(isTRUE(show_gene_number))
		{
		plot <- plot +
		annotate(
		"text",
		geom_text(
		x = xlimits/2,
		y = max(dt$log10P) + 0.5,
		label = stat[regulation == "up",label],
		family = font
		aes(label = label),
		data = stat[regulation == "up"]
		) +
		annotate(
		"text",
		geom_text(
		x = xlimits/2 * -1,
		y = max(dt$log10P) + 0.5,
		label = stat[regulation == "down",label],
		family = font
		aes(label = label),
		data = stat[regulation == "down"]
		)
		}

Admin message