アミノ酸配列から分子量を計算

FASTA ファイルは核酸またはアミノ酸配列を保存するためのフォーマットの一つである。配列の名前は > で始まる行に書かれ、次の行から配列が書かれている。一つのファイルに複数の配列を書くことも可能である。

>1ALK:A|PDBID|CHAIN|SEQUENCE
TPEMPVLENRAAQGNITAPGGARRLTGDQTAALRNSLSDKPAKNIILLIGDGMGDSEITAARNYAEGAGGFFKGIDALPL
TGQYTHYALNKKTGKPDYVTDSAASATAWSTGVKTYNGALGVDIHEKDHPTILEMAKAAGLATGNVSTAELQDATPAALV
AHVTSRKCYGPSATSQKCPGNALEKGGKGSITEQLLNARADVTLGGGAKTFAETATAGEWQGKTLREEAEARGYQLVSDA
ASLNSVTEANQQKPLLGLFADGNMPVRWLGPKATYHGNIDKPAVTCTPNPQRNDSVPTLAQMTDKAIELLSKNEKGFFLQ
VEGASIDKQDHAANPCGQIGETVDLDEAVQRALEFAKKEGNTLVIVTADHAHASQIVAPDTKAPGLTQALNTKDGAVMVM
SYGNSEEDSQEHTGSQLRIAAYGPHAANVVGLTDQTDLFYTMKAALGLK
>1ALK:B|PDBID|CHAIN|SEQUENCE
TPEMPVLENRAAQGNITAPGGARRLTGDQTAALRNSLSDKPAKNIILLIGDGMGDSEITAARNYAEGAGGFFKGIDALPL
TGQYTHYALNKKTGKPDYVTDSAASATAWSTGVKTYNGALGVDIHEKDHPTILEMAKAAGLATGNVSTAELQDATPAALV
AHVTSRKCYGPSATSQKCPGNALEKGGKGSITEQLLNARADVTLGGGAKTFAETATAGEWQGKTLREEAEARGYQLVSDA
ASLNSVTEANQQKPLLGLFADGNMPVRWLGPKATYHGNIDKPAVTCTPNPQRNDSVPTLAQMTDKAIELLSKNEKGFFLQ
VEGASIDKQDHAANPCGQIGETVDLDEAVQRALEFAKKEGNTLVIVTADHAHASQIVAPDTKAPGLTQALNTKDGAVMVM
SYGNSEEDSQEHTGSQLRIAAYGPHAANVVGLTDQTDLFYTMKAALGLK

アミノ酸配列が保存されている FASTA ファイルを読み込んで、そのペプチドの分子量などを計算するには、次のようにする。

  1. アミノ酸の名前と分子量が 1 セットなった連想配列を作成する
  2. FASTA ファイルを開いて処理する
    1. > から始まるヘッダー行ならば、配列の名前を記録する。
    2. ヘッダー行でなければ、アミノ酸配列を取得し、その分子量を一つずつ足合わせる。
    3. 一つの配列を処理し終えた後に、配列の名前と分子量を出力する。分子量には、ペプチド結合による脱水される水の分子量が含まれるので、これを除いて出力する。
#include <fstream>
#include <iostream>
#include <string>
#include <map>

int main() {
    std::map<std::string, float> aa_weight;

    aa_weight["A"] =  89.09;
    aa_weight["R"] = 174.20;
    aa_weight["N"] = 132.12;
    aa_weight["D"] = 133.10;
    aa_weight["C"] = 121.16;
    aa_weight["Q"] = 146.15;
    aa_weight["E"] = 147.13;
    aa_weight["G"] =  75.07;
    aa_weight["H"] = 155.15;
    aa_weight["I"] = 131.17;
    aa_weight["L"] = 131.17;
    aa_weight["K"] = 146.19;
    aa_weight["M"] = 149.21;
    aa_weight["F"] = 165.19;
    aa_weight["P"] = 115.13;
    aa_weight["S"] = 105.09;
    aa_weight["T"] = 119.12;
    aa_weight["W"] = 204.23;
    aa_weight["Y"] = 181.19;
    aa_weight["V"] = 117.15;

    std::ifstream ifs("./data/1alk.fa");
    std::string buff;

    if (ifs.fail()) {
        std::cerr << "Cannot open file!" << std::endl;
        return -1;
    }

    std::string seq_name;
    float seq_weight = 0.0;
    int n_aa = 0;

    while (getline(ifs, buff)) {
        if (buff.substr(0, 1) == ">") {
            if (seq_weight > 0) {
                std::cout << seq_name << " => " << seq_weight - (n_aa - 1) * 18 << std::endl;
            }

            seq_name = buff;
            seq_weight = 0.0;
            n_aa = 0;
        } else {
            for (int i = 0; i < buff.size(); i++) {
                std::string aa = buff.substr(i, 1);
                auto itr = aa_weight.find(aa);
                if (itr != aa_weight.end()) {
                    seq_weight += aa_weight[aa];
                    n_aa += 1;
                } else {
                    std::cout << "undefined amino acid [" << aa << "] was found" << std::endl;
                }
            }
        }
    }

    std::cout << seq_name << " => " << seq_weight - (n_aa - 1) * 18 << std::endl;

    return 0;
}