1 Star 0 Fork 0

陈同/s-plot

加入 Gitee
与超过 1200万 开发者一起发现、参与优秀开源项目,私有仓库也完全免费 :)
免费加入
克隆/下载
sp_enrichmentPlot.sh 12.60 KB
一键复制 编辑 原始数据 按行查看 历史
Admin 提交于 2017-07-20 21:12 . enrichmentPlot FOR
#!/bin/bash
#set -x
usage()
{
cat <<EOF
${txtcyn}
***CREATED BY Chen Tong (chentong_biology@163.com)***
Usage:
$0 options${txtrst}
${bldblu}Function${txtrst}:
This script is used to do scatter plot and represents two value
variables by color and size using ggplot2.
It is designed for representing the gene annotation data with enriched
GO terms and -pvalues, # of genes.
The parameters for logical variable are either TRUE or FALSE.
Input file (terms only exist in one or a few samples are suitable):
Description GeneRatio qvalue Count Type
ERBB signaling pathway 7/320 0.001836081 7 EHBIO_up
regulation of ERBB signaling pathway 5/320 0.003886659 5 EHBIO_up
negative regulation of cell cycle G1/S phase transition 4/320 0.016153254 4 EHBIO_up
Wnt signaling pathway 13/320 0.01680096 13 EHBIO_up
cell-cell signaling by wnt 13/320 0.0171473 13 EHBIO_up
negative regulation of cell cycle process 8/320 0.019453085 8 EHBIO_up
extrinsic apoptotic signaling pathway 9/320 0.024164034 9 EHBIO_up
positive regulation of extrinsic apoptotic signaling pathway 4/320 0.025708228 4 EHBIO_up
cell cycle G1/S phase transition 7/320 0.035797856 7 EHBIO_up
negative regulation of apoptotic signaling pathway 8/320 0.038684745 8 EHBIO_up
regulation of Notch signaling pathway 4/320 0.041592045 4 EHBIO_up
regulation of cell cycle G1/S phase transition 5/320 0.047407619 5 EHBIO_up
negative regulation of BMP signaling pathway 3/320 0.049460847 3 EHBIO_up
regulation of ERK1 and ERK2 cascade 14/342 0.000629602 14 Baodian_up
positive regulation of cell adhesion 17/342 0.000827275 17 Baodian_up
ERK1 and ERK2 cascade 14/342 0.001086508 14 Baodian_up
regulation of cell growth 17/342 0.002228511 17 Baodian_up
positive regulation of cytoskeleton organization 10/342 0.004406867 10 Baodian_up
regulation of cell-cell adhesion 15/342 0.005075219 15 Baodian_up
regulation of cytoskeleton organization 15/342 0.019685646 15 Baodian_up
negative regulation of Notch signaling pathway 3/342 0.020578211 3 Baodian_up
neuron apoptotic process 10/342 0.040284925 10 Baodian_up
${txtbld}OPTIONS${txtrst}:
-f Data file (with header line, the first column is not the
rowname, tab seperated)${bldred}[NECESSARY]${txtrst}
-t Title of picture[${txtred}Default empty title${txtrst}]
[Scatter plot of horizontal and vertical variable]
-x xlab of picture[${txtred}Default empty xlab${txtrst}]
[The description for horizontal variable]
-y ylab of picture[${txtred}Default empty ylab${txtrst}]
[The description for vertical variable]
-P Legend position[${txtred}Default right. Accept
top, bottom, left, none, or c(0.08, 0.8).${txtrst}]
-R Rotation angle for x-axis value(anti clockwise)
[Default 0]
-H Hjust when rotation angle for x-axis value is not zero(anti clockwise)
[Default 0.5; angle 45, hjust 0 vjust 0]
-V Vjust when rotation angle for x-axis value is not zero(anti clockwise)
[Default 1; angle 90, hjust 0 vjust 1]
-o The variable for horizontal axis.
${bldred}[NECESSARY, such Sample]${txtrst}
-T Type of horizontal axis variable.
${bldred}[NECESSARY, default string, accept numeric]${txtrst}
-O The order for horizontal axis.
${bldred}[Only woroks if horizontal axis variables are strings.
If horizontaol axis are numbers like GeneRatio, this will be treated as
sample order.
Default alphabetical order, accept a string like
"'K562','hESC','GM12878','HUVEC','NHEK','IMR90','HMEC'"
***When -m is used, this default will be ignored too.*********
${txtrst}]
]${txtrst}
-v The variable for vertical axis.
${bldred}[NECESSARY, such as Term]${txtrst}
-c The variable for point color.
${bldred}[NECESSARY, such as p_value]${txtrst}
-s The variable for point size.
${bldred}[NECESSARY, such as count]${txtrst}
-S The variable for sample classification.
${bldred}[Optional, if given will use different shapes for each sample]${txtrst}
-C The specified colors.
${bldred}[Default '"green", "red"'. ]${txtrst}
-l Get log-transformed data for given variable.
[${txtred}Default nolog, means no log10 transform. Accept a variable
like p_value to get (-1) * log10(p_value).${txtrst}]
-w The width of output picture.[${txtred}Default 20${txtrst}]
-a The height of output picture.[${txtred}Default 20${txtrst}]
-E The type of output figures.[${txtred}Default pdf, accept
eps/ps, tex (pictex), png, jpeg, tiff, bmp, svg and wmf)${txtrst}]
-r The resolution of output picture.[${txtred}Default 300 ppi${txtrst}]
-b The formula for facets.[${bldred}Default no facets,
+facet_grid(level ~ .) means divide by levels of 'level' vertcally.
+facet_grid(. ~ level) means divide by levels of 'level' horizontally.
+facet_grid(lev1 ~ lev2) means divide by lev1 vertically and lev2
horizontally.
+facet_wrap(~level, ncol=2) means wrap horizontally with 2
columns.
Example: +facet_wrap(~Size,ncol=6,scale='free')
${txtrst}]
-d If facet is given, you may want to specifize the order of
variable in your facet, default alphabetically.
[${txtred}Accept sth like
(one level one sentence, separate by';')
data\$size <- factor(data\$size, levels=c("l1",
"l2",...,"l10"), ordered=T) ${txtrst}]
-z Other parameters in ggplot format.[${bldred}selection${txtrst}]
-e Execute or not[${bldred}Default TRUE${txtrst}]
-i Install the required packages[${bldred}Default FALSE${txtrst}]
Examples:
# -f: 指定输入文件,格式如上面描述
# -o: 指定横轴的变量,单个样品一般选择GeneRatio或样品名字
# -T: 指定横轴变量的类似,是字符串还是数值
# -v: 指定Y轴显示的内容,一般为富集条目的描述
# -c: 指定用哪一列设置颜色展示,一般为qvalue或p.adjust
# -s: 指定用哪一列设置点的大小,一般为Count
# -l: 指定某一列进行对数操作,一般选qvalue列;如果已做过对数操作,则忽略
# -a: 设置图片输出高度
# -x, -y: 设置横轴和纵轴标题,注意多个单词需要引号括起来
* sp_enrichmentPlot.sh -f GOenrichement.ehbio.xls -o GeneRatio -T numeric -v Description -c qvalue -s Count -l qvalue -a 12 -x "GeneRatio" -y "GO description"
# -o: 指定横轴的变量,单个样品一般选择GeneRatio或样品名字
# -T: 如果是样品名字,指定为字符串
* sp_enrichmentPlot.sh -f GOenrichement.ehbio.xls -o Type -T string -v Description -c qvalue -s Count -l qvalue -a 12 -x "Sample" -y "GO description"
# 多出来的参数是-S用来指定样品分组,不同类型的基因的富集分析用不同的形状表示
* sp_enrichmentPlot.sh -f GOenrichement.xls -o GeneRatio -T numeric -v Description -c qvalue -s Count -l qvalue -a 12 -x "GeneRatio" -y "GO description" -S Type
# 跟单个样品不显示GeneRatio的命令一样,不同的样品分列展示。
* sp_enrichmentPlot.sh -f GOenrichement.xls -o GeneRatio -T numeric -v Description -c qvalue -s Count -l qvalue -a 12 -x "GeneRatio" -y "GO description" -S Type
EOF
}
file=''
title=''
xlab=''
ylab=''
xval=''
xval_ho=''
xval_type='string'
yval=''
execute='TRUE'
ist='FALSE'
color=''
color_v='"green", "red"'
log='nolog'
width=20
height=20
res=300
ext='pdf'
facet=''
size=''
other=''
facet_o=''
legend_pos='right'
xtics_angle=0
hjust=0.5
vjust=1
sample='CTctCT'
while getopts "hf:t:x:y:o:O:P:R:H:V:v:T:c:C:l:w:a:r:E:s:S:b:d:z:e:i:" OPTION
do
case $OPTION in
h)
usage
exit 1
;;
f)
file=$OPTARG
;;
t)
title=$OPTARG
;;
x)
xlab=$OPTARG
;;
y)
ylab=$OPTARG
;;
P)
legend_pos=$OPTARG
;;
R)
xtics_angle=$OPTARG
;;
H)
hjust=$OPTARG
;;
V)
vjust=$OPTARG
;;
o)
xval=$OPTARG
;;
S)
sample=$OPTARG
;;
O)
xval_ho=$OPTARG
;;
T)
xval_type=$OPTARG
;;
v)
yval=$OPTARG
;;
c)
color=$OPTARG
;;
C)
color_v=$OPTARG
;;
l)
log=$OPTARG
;;
w)
width=$OPTARG
;;
a)
height=$OPTARG
;;
r)
res=$OPTARG
;;
E)
ext=$OPTARG
;;
b)
facet=$OPTARG
;;
d)
facet_o=$OPTARG
;;
s)
size=$OPTARG
;;
z)
other=$OPTARG
;;
e)
execute=$OPTARG
;;
i)
ist=$OPTARG
;;
?)
usage
exit 1
;;
esac
done
mid=".scatterplot.dv"
if [ -z $file ] || [ -z $xval ] || [ -z $yval ] ; then
echo 1>&2 "Please give filename, xval and yval."
usage
exit 1
fi
if test "${xval_type}" == "string"; then
if test "${sample}" == "CTctCT"; then
sample="${xval}"
fi
fi
shape=''
if test "${sample}" != "CTctCT" && test "${sample}" != "${xval}"; then
shape=", shape=${sample}"
fi
if test "${log}" != "nolog"; then
color="negLog10_${color}"
fi
cat <<END >${file}${mid}.r
if ($ist){
install.packages("ggplot2", repo="http://cran.us.r-project.org")
}
library(plyr)
library(stringr)
library(ggplot2)
library(grid)
# Function get from https://stackoverflow.com/questions/10674992/convert-a-character-vector-of-mixed-numbers-fractions-and-integers-to-numeric?rq=1
# With little modifications
mixedToFloat <- function(x){
x <- sapply(x, as.character)
is.integer <- grepl("^-?\\\\d+$", x)
is.fraction <- grepl("^-?\\\\d+\\\\/\\\\d+$", x)
is.float <- grepl("^-?\\\\d+\\\\.\\\\d+$", x)
is.mixed <- grepl("^-?\\\\d+ \\\\d+\\\\/\\\\d+$", x)
stopifnot(all(is.integer | is.fraction | is.float | is.mixed))
numbers <- strsplit(x, "[ /]")
ifelse(is.integer, as.numeric(sapply(numbers, \`[\`, 1)),
ifelse(is.float, as.numeric(sapply(numbers, \`[\`, 1)),
ifelse(is.fraction, as.numeric(sapply(numbers, \`[\`, 1)) /
as.numeric(sapply(numbers, \`[\`, 2)),
as.numeric(sapply(numbers, \`[\`, 1)) +
as.numeric(sapply(numbers, \`[\`, 2)) /
as.numeric(sapply(numbers, \`[\`, 3)))))
}
#mixedToFloat(c('1 1/2', '2 3/4', '2/3', '11 1/4', '1'))
data <- read.table(file="$file", sep="\t", quote="", comment="", header=T)
xval_ho <- c(${xval_ho})
if ("${sample}" != "CTctCT" & length(xval_ho) > 1) {
data\$${sample} <- factor(data\$${sample}, levels=sample_ho, ordered=T)
} #else {
# data\$${xval} <- factor(data\$${xval})
#}
# First order by Term, then order by Sample
if ("${sample}" != "CTctCT") {
data <- data[order(data\$${yval}, data\$${sample}), ]
}
if ("${xval_type}" != "string"){
data\$${xval} = mixedToFloat(data\$${xval})
}
if ("${log}" != "nolog"){
log_name = paste0("negLog10_", "${log}")
col_name_data <- colnames(data)
col_name_data <- c(col_name_data, log_name)
data\$log_name <- log10(data\$${log}) * (-1)
colnames(data) <- col_name_data
}
# Get the count of each unique Term
data_freq <- as.data.frame(table(data\$${yval}))
colnames(data_freq) <- c("${yval}", "IDctct")
data2 <- merge(data, data_freq, by="${yval}")
if ("${sample}" != "CTctCT"){
# Collapse sample for each Term
data_samp <- ddply(data2, "${yval}", summarize,
sam_ct_ct_ct=paste(${sample}, collapse="_"))
data2 <- merge(data2, data_samp, by="${yval}")
#print(data2)
if ("${xval_type}" != "string"){
data3 <- data2[order(data2\$IDctct, data2\$sam_ct_ct_ct, data2\$${sample}, data2\$${xval}, data2\$${color}), ]
} else {
data3 <- data2[order(data2\$IDctct, data2\$sam_ct_ct_ct, data2\$${sample}, data2\$${color}), ]
}
} else{
if ("${xval_type}" != "string"){
data3 <- data2[order(data2\$IDctct, data2\$${xval}, data2\$${color}), ]
} else {
data3 <- data2[order(data2\$IDctct, data2\$${color}), ]
}
}
#print(data3)
term_order <- unique(data3\$${yval})
data\$${yval} <- factor(data\$${yval}, levels=term_order, ordered=T)
#print(data)
rm(data_freq, data2, data3)
$facet_o
color_v <- c(${color_v})
p <- ggplot(data, aes(x=${xval},y=${yval})) \
+ labs(x="$xlab", y="$ylab") + labs(title="$title")
if (("${size}" != "") && ("${color}" != "")) {
p <- p + geom_point(aes(size=${size}, color=${color} ${shape})) + \
scale_colour_gradient(low=color_v[1], high=color_v[2], name="${color}")
} else if ("${size}" != "") {
p <- p + geom_point(aes(size=${size} ${shape}))
} else if ("${color}" != "") {
p <- p + geom_point(aes(color=${color} ${shape})) + \
scale_colour_gradient(low="color_v[1]", high=color_v[2], name="${color}")
}
p <- p ${facet}
p <- p + scale_y_discrete(labels=function(x) str_wrap(x, width=60))
p <- p $other
p <- p + theme_bw() + theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank())
if (${xtics_angle} != 0){
p <- p +
theme(axis.text.x=element_text(angle=${xtics_angle},hjust=${hjust},
vjust=${vjust}))
}
top='top'
bottom='bottom'
left='left'
right='right'
none='none'
legend_pos_par <- ${legend_pos}
p <- p + theme(legend.position=legend_pos_par)
ggsave(p, filename="${file}${mid}.${ext}", dpi=$res, width=${width},
height=${height}, units=c("cm"))
END
if [ "$execute" == "TRUE" ]; then
Rscript ${file}${mid}.r
if [ "$?" == "0" ]; then /bin/rm -f ${file}${mid}.r; fi
fi
if [ ! -z "$log" ]; then
log=', trans=\"'$log'\"'
fi
#convert -density 200 -flatten ${file}${mid}.eps ${first}${mid}.png
马建仓 AI 助手
尝试更多
代码解读
代码找茬
代码优化
Shell
1
https://gitee.com/ct5869/s-plot.git
git@gitee.com:ct5869/s-plot.git
ct5869
s-plot
s-plot
master

搜索帮助

A270a887 8829481 3d7a4017 8829481