linux - 无法让我的程序运行 - %GC 计算器

Question

几周来我一直在尝试不让这个程序运行。我是编程新手，这绝对是一个挑战。我认为我的问题出在我的 if 语句上。我可以让它将名称附加到新文件，但它只是将整个序列附加到文件而不是计算它。我正在处理一个包含多个序列的 fasta 文件，名称以“>”开头，序列在其下方一行。这是我的代码。请帮助，并在此先感谢您!!

#! /bin/bash

#exit program with error if user does not specify input on command line

if [ $# != 1 ]; then
        echo "Please specify fasta input on command line and rerun"
        exit
        else echo "Beginning count"
fi

#collect input from user each time they run the program
input=`cat $1`

#seperate the sequence from the sequence name
name=`grep '>' $1`
sequence=`grep -v '>' $1`


#if name, if sequence
IFS=$'\n'
set -f
for i in $(cat "$1"); 
do
    if [ $i=">" ]; then
        echo "$i" >> GCcontent.txt
        else 
        #count number of occurence of motif ATGC in fasta sequence
        countG=`echo $i | grep -o "G" | wc -l`
        countC=`echo $i | grep -o "C" | wc -l`
        total=`echo $i | wc -m`
        count=`echo "scale=2" ; ($countG+$countC) | bc`

        #calculate percent over total divided by 3bp
        percent=`echo "scale=2 ; ($count/$total*100)" | bc`

        #print output name and percent to file
        echo "$percent" >> GCcontent.txt
    fi
done

echo "Exiting"

exit

编辑: 输入文件:(> 周围没有引号，没有它们我不会把它放在那里) '>'gi|226451773|gb|FJ846591.1 CATTATAGACTGCGTGGTCCGTATTCCCAAGGAGCAGGGAGTTCTGTCCTTCTGGCGCGGTAACCTGGCCAATGTCATCAGATACTTCCCCACCCAGGCTCTTAACTTCGCCTTCAAAGATAAATACAAGCAGATCTTCCTAGGTGGTGTGGACAAGAGGACCCAGTTTTGGCGCTACTTTGCAGGGAATCTGGCATCAGGTGGTGCCGCAGGGGCCACATCCCTGTGTTTTGTGTACCCTCTTGATTTTGCCCGTACCCGTCTAGCAGCTGATGTGGGTAAAGCTGGAGCTGAAAGGGAATTCCGAGGCCTCGGTGACTGCCTGGTTAAGATCTACAAATCTGATGGGATTAAGGGCCTGTACCAAGGCTTTAACGTGTCTGTGCAGGGTATTATCATCTACCGAGCCGCCTACTTCGGTATCTATGACACTGCAAAGGGTAAGTTTGCTGTGGGCTTTAAAGTTGTGTTCTTAGGAGACAATTTAAAAGAGCGTTGTACCAACCTAACATTCCAAGAGCTAGAGAGTTTTTTTAATTGCTGAAGGAAGCCAAGATCATCCAGTGCGACCCTCATGCACAGATGACATGTTTAGGGGATGTGGGGAAAGGAAGTCAGTAAAACTCTACTTTTTGGTAAAAGCATCTCTTTCCTATTCCCAGGAATGCTTCCGGATCCCAAAAACACTCACATCGTCATCAGCTGGATGATCGCACAGACTGTCACTGCTGTTGCTGGGTTGACTTCCTATCCATTTGA (但文件中有多个)

我希望输出文件具有: '>'gi|226451773|gb|FJ846591.1 'Percent 会去这里'

Answer 1

awk 助您一臂之力!

在这里使用您的输入文件(行)是一个概念证明

$ awk '{s=length($2); 
        g=gsub("G","",$2); c=gsub("C","",$2); t=gsub("T","",$2); a=gsub("A","",$2);
        total=a+c+g+t; 
        print a,c,g,t,total,100*(c+g)/total,total==s}' file

185 175 190 209 759 48.0896 1

最后一个值 1 检查单个匹配项的总匹配项和。

您想将逻辑应用于以 >> 符号为前缀的行

$ awk '/^>/{...; print $1, 100*(c+g)/total}' inputfile > outputfile

从上面复制相关部分以填补空白。

如果你有没有任何匹配的行，总计将为零并且除法将失败。你可以加一个守卫。

完整的脚本可以这样

 $ awk '/^>/{g=gsub("G","",$2); c=gsub("C","",$2); total=length($2);             
             print $1, (total?100*(c+g)/total:"no match found")}' inputfile > outputfile