专栏名称: 生信媛

生信媛，从1人分享，到8人同行。坚持分享生信入门方法与课程，持续记录生信相关的分析pipeline, python和R在生物信息学中的利用。内容涵盖服务器使用、基因组转录组分析以及群体遗传。

写了100多个shell脚本的老司机的翻车之旅

生信媛 · 公众号 · 生物 · 2017-07-26 11:32

正文

在写shell脚本上，我就算不是一名老司机，也算是一名合格司机了。

基本上，大部分重复的任务都被我写成了脚本，比如说下面这个脚本就是我用于重测序找variant用的脚本。

#!/bin/bash
# version 1.1# function:  provide different function to align and vairant calling
# author: xuzhougeng, xuzhougeng@163.com
# parameter: samplePath, index, reference threads

set -e
set -u
set -o pipefail

PATH=/home/xzg/miniconda3/bin:/usr/local/bin:/usr/bin:/bin
echo "you should run this program in the project root path "if [ $# -lt 4 ]
then
  echo -e "$# less than 4, please enter samplePath, index, reference and threads"
  exit 1fi
samplePath=$1index=$2suffix=$3threads=$4reference=${index}
mkdir -p alignment
alignDir=alignment
# alignmentfor sample in `cat $samplePath`
do
    filename=${sample##*/}-${suffix}
    echo "processing $filename with bwa"
    if [ ! -f  $alignDir/${filename}.sam ]
    then
        bwa mem -t 8 -B 2 $index ${sample}_1.fq.gz ${sample}_2.fq.gz >\
        $alignDir/${filename}.sam 2>  $alignDir/${filename}.log
        echo "$filename done"
    else
        echo "$filename exists"
    fi
done

# convert sort and index
# warning total memoery > threads x memeoryfor sample in `cat $samplePath`
do
  echo "processing $filename with samtools"
  output=${sample##*/}-${suffix}
  samtools view -@ $threads -b -o $alignDir/${output}.bam $alignDir/${filename}.sam
  samtools sort -@ $threads -m 1G -o $alignDir/${output}.sorted.bam $alignDir/${output}.bam
  samtools index -@ $threads $alignDir/${output}.sorted.bam
  echo "$filename done"done
# vairant calling with bcftools
mkdir -p variant_bcftools
varDir=variant_bcftoolsfor sample in `cat $samplePath`
do
  output=${sample##*/}-${suffix}
  echo "processing $sample with bcftools"
  samtools mpileup  -vu -t AD,DP -f $reference $alignDir/${output}.sorted.bam | \
  bcftools call -vm -Ov > $varDir/${output%%.*}_raw_variants.vcf && echo "$sample done " &
done

但是今天我遇到了一个非常灵异的现象。

我和往常一样，用在我的ATOM文本编辑器中新建了一个配置文件，提供另一个脚本所需要的参数

# test.cfg
dir=/mnt/f/Data/test/Data/Seq
hawkDir=/home/xzg/biosoft/hawk
jellyfishDir=/home/xzg/biosoft/jellyfish-Hawk/bin
sortDir=/home/xzg/biosoft/coreutils/bin
CORES=2KMERSIZE=31

而在另一个脚本，我用awk提取cfg里面的参数，然后遍历里面的文件夹的内容

#!/bin/bash

set -e
set -u
set -o pipefail

#modified from NIKS script
#echo $PATH
...
config=$1if [ ! -f ${config} ] && [ ! -r ${config} ]
then
    echo -e "$1 is not a file or is not readable \n "fi

#  setting
dirs=$(awk 'BEGIN{FS="="}; $1~/dir/ { print $2}' $config)
hawkDir=$(awk 'BEGIN{FS="="}; $1~/hawkDir/ { print $2}'