脚本模版
shell脚本必备模板
区块介绍
依然是分区块写脚本。
-
定义脚本的使用信息
-
使用
getopt定义长参数与短参数
这里的-f | --force与-o|--output 在该脚本中是个假选项，无任何作用，仅是作为模板展示用。-f | --force其理想作用是在有流程续借的情况下还能重新强制运行该脚本。
- 设置各参数默认值
while case的方式解析选项与参数。- 参数的处理
可能需要根据while case 中获取的参数进行处理，以获取路径、文件和ID等等信息。
- 定义流程接续标志文件。
【该脚本此处未启用】
定义流程接续标志文件，这里将其放在了脚本同目录下。需要根据脚本的位置修改标志文件位置。如果你的脚本是写入了全局环境变量的，这里标志文件位置肯定要改成工作目录。
对于单功能脚本，该功能其实用途不大，当要写多个shell脚本的组合流程时，其可以应对更多复杂的运行情况。
这个写法是跟Trinotate中数据下载的perl脚本学的。（可见阅读别人优秀脚本的重要性！）
- 定义函数
模块化的写法，这肯定是必须要学的内容。
注意shell脚本中函数的参数传递方式吧。我写的时候把自己套迷糊了。
这里判断压缩文件的写法总感觉不太准确和优雅，但是先这样，能用就好。
- 检查与运行
如果输入是文件，可以先检查下是否存在。
这里把echo +date部分注释掉了，其它shell脚本中可以尝试启用记录脚本各步骤运行的起止时间。
- 流程接续
【该脚本未启用本部分，这里仅作为模板示例进行展示】
通过检测流程接续标志文件的有无，判断本次是否运行。
不足
- 有些地方不是很理解，模仿着写的。所以可能有不当的地方，可能需要修改。
- 脚本只能接受单个文件作为输入。
让脚本接受多个文件的写法，尝试了3种写法：逗号“,”分割文件转为数组，空格分割文件转为数组与将参数依次放入一个数组。但是写的都有瑕疵，感觉不是很好，那就还不如单个文件作为输入。
-n指定了默认参数，但是当重新指定-n与数值时，只能使用-n3000或者--read_nums=3000的方式指定。不允许使用空格。不知道如何突破这种特点，总感觉怪怪的。
脚本功能
该脚本为了展示shell模板添加了很多无用注释掉的内容。但是该脚本还是能用的。
该脚本功能是返回fastq文件的碱基质量体系版本。
该功能核心部分来自：FASTQ quality scores - Bioinformatics Workbook ：
https://bioinformaticsworkbook.org/introduction/fastqquality-score-encoding.html#gsc.tab=0
脚本运行
脚本名称run_check_fastq_phred_v2.sh
单个文件
bash run_check_fastq_phred_v2.sh -i SRR9850430_1.fastq.gz
bash run_check_fastq_phred_v2.sh -i SRR9850430_1.fastq.gz -n3000
bash run_check_fastq_phred_v2.sh -i SRR9850430_1.fastq.gz --read_nums=3000
三个文件返回结果均为SRR9850430_1.fastq.gz Phred+33
多个文件
ls *.fastq.gz |xargs -i bash run_check_fastq_phred_v2.sh -i {}
ls *.fastq.gz |xargs -i bash run_check_fastq_phred_v2.sh -i {} -n3000
结果
SRR9850430_1.fastq.gz Phred+33
SRR9850430_2.fastq.gz Phred+33
SRR9850431_1.fastq.gz Phred+33
SRR9850431_2.fastq.gz Phred+33
SRR9850432_1.fastq.gz Phred+33
SRR9850432_2.fastq.gz Phred+33
SRR9850433_1.fastq.gz Phred+33
SRR9850433_2.fastq.gz Phred+33
SRR9850434_1.fastq.gz Phred+33
SRR9850434_2.fastq.gz Phred+33
脚本
#!/usr/bin/bash
## ---------- 定义脚本的使用方法 ----------------------------------
usage() {
cat << EOF
Usage:bash $0 -i <fastq.gz>
bash $0 -i <fastq.gz> -n3000
bash $0 -i <fastq.gz> --read_nums=3000
ls *.fastq.gz |xargs -i bash $0 -i {}
Options:
-i, --input Input file.
-o, --ouput Output file.
-n, --read_nums A positive integer. Number of reads to test [default: 2500].
-f, --force force tu rerun this step.
-h, --help Print the tips.
EOF
exit 1
}
## ---------- 脚本选项设置 ----------------------------------
OPTIONS=$(getopt -o hfi:o::n:: --long help,force,input:,output::,read_nums:: -n "ERROR:$0" -- "$@")
#if [ $? != 0]; then
# echo "Terminating..." >&2 ;
# exit 1;
#fi
if [ $? -ne 0 ]; then
usage && exit 1
fi
# Note the quotes around `$OPTIONS': they are essential!
eval set -- "$OPTIONS"
## ---------- 设置默认值 ----------------------------------
FILE1=""
FILE2=""
READ_NUMS=2500
FORCE=false
## ---------- 处理解析后的命令行参数 ----------------------------------
while true; do
case "$1" in
-h|--help)
usage
;;
-f|--force)
FORCE=true
shift
;;
-i|--input)
FILE1="$2" # 赋值
shift 2
;;
-o|--output)
FILE2="$2" # 赋值
shift 2
;;
-n|--read_nums)
READ_NUMS=$2
shift 2
;;
--)
shift
break
;;
*)
echo "Cannot recognize: $1"
usage
;;
esac
done
## ---------- 参数的处理 ----------------------------------
if [[ -z $FILE1 ]]; then # 字符串长度是否为0
echo "ERROR: Input file is required."
usage && exit 1
fi
if [[ $READ_NUMS -le 0 ]];then # 小于等于则为真
echo "ERROR: optinon '-n' must be a positive integer!"
usage && exit 1
else
num_lines=$(( READ_NUMS * 4 ))
fi
#echo "$FILE1"
#echo "$READ_NUMS"
## ---------- 流程接续用文件 ----------------------------------
<<EOF
this_scirpt=$(readlink -f "$0")
#echo "$this_scirpt" # 测试脚本所在位置。
end_log="${this_scirpt}.end_log"
EOF
## ---------- 定义函数 ----------------------------------
### 处理压缩文件
function get_fq_lines(){ # $1 fastq.gz file; $2 number of lines of fastq
if [[ "$1" =~ \.gz$ ]]; then
zcat "$1" |head -n "$2"
else
cat "$1" |head -n "$2"
fi
}
#get_fq_lines "$FILE1" "$num_lines" # 简单测试下这个函数
### 检查fastq文件质量体系
function check_fq (){ # $1 fastq.gz file; $2 number of lines of fastq
file=$(basename $1) # print filename only
#zcat $1 |\
#head -n 10000 | \
get_fq_lines $1 $2 | \
awk '{if(NR%4==0) printf("%s",$0);}' | od -A n -t u1 | \
awk -v file=$file 'BEGIN{min=100;max=0;} \
{for(i=1;i<=NF;i++) \
{if($i>max) max=$i; \
if($i<min) min=$i;}}END \
{if(max<=74 && min<59) \
print file"\tPhred+33"; \
else \
if(max>73 && min>=64) \
print file"\tPhred+64"; \
else \
if(min>=59 && min<64 && max>73) \
print file"\tSolexa+64"; else print file"\tUnknown score encoding!";}'
}
## ---------- 检查与运行 ----------------------------------
## 检查输入文件是否存在
if [ ! -f "$FILE1" ]; then
echo "ERROR: $i not exits."
usage
exit 1
fi
#echo "---------- step 1 start: `date`----------------------------------"
#echo $FILE1 $num_lines
check_fq $FILE1 $num_lines # check phred score
#echo "---------- step 1 done: `date`----------------------------------"
## ---------- 流程接续 ----------------------------------
<<EOF
#[[ -f $end_log ]] || ( echo "$FILE1" && touch ${end_log} )
if [[ $FORCE ]]; then rm -rf $end_log ;fi
if [[ ! -f $end_log ]] ;then
for i in $(ls $FILE1); do
check_fq $i $num_lines
done # 主体命令
fi
if [ $? -ne 0 ]; then
usage && exit 1
else
touch ${end_log}
fi
EOF
#echo "---------- DONE: `date`----------------------------------"
收费服务
『收费服务』 SSR分析收费版：MISA+Primer3流程
推荐阅读
关注后可发消息
