3

我有一个包含蛋白质序列(200 个序列)的文本文件,如下所示。

>ptn1
AAGHM
>ptn2
MGLKKRR

我需要为序列的每个字符提供以下值,并且必须找到每个序列的平均值。

A= 0.2, G= 0.5, L=0.14, M= 0.70, R= 0.55, C=0.48, H= 1.00 , K=0.4

期望的输出

ptn1  - 0.52
ptn2  - 0.462

如何使用 awk 或 python 做到这一点?

您的建议将不胜感激

4

4 回答 4

6
def avg(sequence):
    v= {'A': 0.2, 'C': 0.48, 'R': 0.55, 'G': 0.5, 'H': 1.0,
        'K': 0.4, 'M': 0.7, 'L': 0.14}
    return sum(v[x] for x in sequence) / len(sequence)

avg("AAGHM")  # => 0.5199999999999999
avg("MGLKKRR" # => 0.46285714285714274
于 2012-07-23T06:51:45.807 回答
5

FS=""
http://www.gnu.org/software/gawk/manual/html_node/Single-Character-Fields.html#Single-Character-Fields需要 gawk
用法:
awk -f foo.awk foo.txt

BEGIN {
    FS=""
    k["A"]=0.2; k["G"]=0.5; k["L"]=0.14; k["M"]=0.70
    k["R"]=0.55; k["C"]=0.48; k["H"]=1.00; k["K"]=0.4
}

/^>/{
    $1=""
    name=$0
    next
}

{
    s=0
    for (i=1; i<=NF; i++) {
      s+=k[$(i)]
    }
    printf "%s - %.3f\n", name, s/NF
}
于 2012-07-23T07:20:14.860 回答
1

您的序列文件似乎是 FASTA 格式。您应该使用专门适用于处理序列的工具,而不是使用awk.

如果您的文件格式发生变化怎么办?如果您想根据格式规范提取有意义的数据怎么办?已经构建的解析器非常适合这些问题。

我喜欢Biopython

from Bio import SeqIO
records = SeqIO.parse("sequences.fasta", "fasta")

# Function borrowed from Tichodrama's answer
def avg(sequence):
    v= {'A': 0.2, 'C': 0.48, 'R': 0.55, 'G': 0.5, 'H': 1.0,
        'K': 0.4, 'M': 0.7, 'L': 0.14}
    return sum(v[x] for x in sequence) / len(sequence)

for record in records:
    print("Score for {}: {:.2f}".format(record.id, avg(record.seq)))

请注意,record上述每个SeqRecord对象都是一个包含有用信息的对象,例如解析的定义信息和序列字母表(例如,您可以区分蛋白质序列和 DNA 序列)。

于 2012-07-27T11:42:02.313 回答
0
mapping = {
    'A': 0.2,
    'G': 0.5,
    'L': 0.14,
    'M': 0.7,
    'R': 0.55,
    'C': 0.48,
    'H': 1.0,
    'K': 0.4}

def process(seqs):
    i = 0
    while i < len(seqs) - 1:
        print seqs[i].strip()[1:], '-', avg(seqs[i+1])
        i += 2

def avg(seq):
    return sum(mapping[v] for v in seq) / len(seq)

# assume the file is data.txt
process(open('data.txt').readlines())
于 2012-07-23T06:58:41.063 回答