Blog

You have a protein sequence and someone needs the mass. So you look up a table of amino acid molecular weights, write a loop, add them all up, and hand over the number. Then the mass spec person tells you your value is 360 Daltons too heavy for a 21 residue peptide. This is the single most common bug when people calculate protein molecular weight in Python, and it comes from the same place every time: the table you found lists free amino acids, not residues.

Why the naive sum is wrong

When two amino acids form a peptide bond, they lose a water molecule between them. The carboxyl group gives up an OH, the amino group gives up an H, and those leave as H2O.

So a chain of n amino acids has n-1 peptide bonds, which means n-1 waters gone. If you sum the free amino acid weights you are counting all that water. For a 21 residue peptide that is 20 extra waters, or 360.3 Da. For a 300 residue protein it is 5386.6 Da, which is a 5.4 kDa error nobody will miss.

There are two ways to fix it. You can sum free amino acid weights and subtract 18.01528 times (n-1). Or you can use a residue mass table, where the water is already removed from every entry, and add one water back at the end for the free N and C termini.

The second way is better. It is one addition instead of a length dependent correction, and it is much harder to get wrong when you later start slicing sequences into peptides.

The residue mass table

Two tables, because there are two kinds of mass and you need both eventually. More on that below.

AVERAGE_RESIDUE_MASS = {
    "G": 57.0519, "A": 71.0788, "S": 87.0782, "P": 97.1167, "V": 99.1326,
    "T": 101.1051, "C": 103.1388, "L": 113.1594, "I": 113.1594, "N": 114.1038,
    "D": 115.0886, "Q": 128.1307, "K": 128.1741, "E": 129.1155, "M": 131.1926,
    "H": 137.1411, "F": 147.1766, "R": 156.1875, "Y": 163.1760, "W": 186.2132,
}

MONOISOTOPIC_RESIDUE_MASS = {
    "G": 57.02146, "A": 71.03711, "S": 87.03203, "P": 97.05276, "V": 99.06841,
    "T": 101.04768, "C": 103.00919, "L": 113.08406, "I": 113.08406, "N": 114.04293,
    "D": 115.02694, "Q": 128.05858, "K": 128.09496, "E": 129.04259, "M": 131.04049,
    "H": 137.05891, "F": 147.06841, "R": 156.10111, "Y": 163.06333, "W": 186.07931,
}

WATER_AVERAGE = 18.01528
WATER_MONOISOTOPIC = 18.010565Code language: Python (python)

Both tables trace back to the standard atomic weights published by IUPAC, which is why different tools agree to about the second decimal and then start drifting.

Leucine and isoleucine are identical. That is not a typo, they are structural isomers with the same formula, and no mass measurement will ever tell them apart.

How to calculate protein molecular weight in Python

def protein_mass(sequence, monoisotopic=False):
    """Mass in Daltons of an unmodified polypeptide with free termini."""
    table = MONOISOTOPIC_RESIDUE_MASS if monoisotopic else AVERAGE_RESIDUE_MASS
    water = WATER_MONOISOTOPIC if monoisotopic else WATER_AVERAGE

    seq = "".join(sequence.split()).upper()
    if not seq:
        raise ValueError("empty sequence")

    total = water
    for position, residue in enumerate(seq, start=1):
        if residue not in table:
            raise ValueError(
                "position {}: {!r} is not a standard amino acid".format(position, residue)
            )
        total += table[residue]
    return totalCode language: Python (python)

The whitespace strip matters more than it looks. Sequences copied out of FASTA files, PDB headers, or a paper’s supplementary table arrive full of newlines, spaces and line numbers, and .split() with no argument handles all of it.

The explicit error is the other half. A silent dict.get(residue, 0) will happily return a mass for a sequence full of X and B, and you will never find out. Fail loudly and tell the caller which position broke.

Try it on the two chains of human insulin:

INSULIN_A = "GIVEQCCTSICSLYQLENYCN"
INSULIN_B = "FVNQHLCGSHLVEALYLVCGERGFFYTPKT"

print("A chain avg  {:.2f}".format(protein_mass(INSULIN_A)))
print("A chain mono {:.4f}".format(protein_mass(INSULIN_A, monoisotopic=True)))
print("B chain avg  {:.2f}".format(protein_mass(INSULIN_B)))
print("B chain mono {:.4f}".format(protein_mass(INSULIN_B, monoisotopic=True)))Code language: Python (python)

A chain avg  2383.71
A chain mono 2382.0000
B chain avg  3429.96
B chain mono 3427.6845Code language: plaintext (plaintext)

The published reduced A chain mass is 2383.7. If you had summed free amino acid weights instead you would have printed 2744.01.

Average or monoisotopic

Average mass uses the standard atomic weights, which are isotope abundance weighted averages. Monoisotopic mass uses only the lightest stable isotope of each element: carbon 12, hydrogen 1, nitrogen 14, oxygen 16.

Which one you want depends entirely on what you are comparing against.

Below roughly 3 kDa, a decent instrument resolves the individual isotope peaks and the leftmost one is the monoisotopic mass. Use monoisotopic. Above roughly 10 kDa the isotope envelope smears into one broad peak whose centroid is the average mass, and the monoisotopic peak is too rare to see at all. Use average. In between, it depends on your resolution, and you should check which peak you are actually picking.

Getting this backwards is a subtle error because for a small peptide the two values differ by about 1 Da, which looks like a plausible protonation artifact rather than a bug.

Speaking of protonation, mass spec reports m/z, not mass. The bridge is one line:

PROTON = 1.007276

def mz(mass, charge):
    return (mass + charge * PROTON) / chargeCode language: Python (python)

Neutral monoisotopic insulin is 5803.6376, so the 5+ ion lands at 1161.7348.

Sequences from a FASTA file

Nobody types sequences in by hand. Most of the time you calculate protein molecular weight in Python straight off a FASTA file, so here is a parser small enough to paste anywhere:

def read_fasta(path):
    header, chunks = None, []
    with open(path) as handle:
        for line in handle:
            line = line.rstrip()
            if line.startswith(">"):
                if header is not None:
                    yield header, "".join(chunks)
                header, chunks = line[1:], []
            elif line:
                chunks.append(line)
    if header is not None:
        yield header, "".join(chunks)


for header, seq in read_fasta("chains.fasta"):
    name = header.split()[0]
    print("{}\t{}\t{:.2f}\t{:.4f}".format(
        name, len(seq), protein_mass(seq), protein_mass(seq, monoisotopic=True)
    ))Code language: Python (python)

With the two insulin chains in chains.fasta:

insulin_A	21	2383.71	2382.0000
insulin_B	30	3429.96	3427.6845Code language: plaintext (plaintext)

If you already have Biopython in the environment, Bio.SeqUtils.molecular_weight(seq, seq_type="protein") does the same job and handles DNA and RNA too. Its numbers will differ from the ones above in the second decimal place because the underlying tables are rounded differently. For a check outside Python entirely, paste the sequence into Expasy ProtParam and compare. That difference is noise next to everything in the next section.

What the mass table does not know

Your number is the mass of a bare, unmodified, fully reduced chain. Real proteins are rarely that.

Disulfide bonds are the easy one. Each one costs two hydrogens, so subtract 2.01565 monoisotopic or 2.01588 average per bond. Insulin has three, two between the chains and one inside the A chain:

DISULFIDE_MONO = 2.01565
DISULFIDE_AVERAGE = 2.01588

reduced_avg = protein_mass(INSULIN_A) + protein_mass(INSULIN_B)
reduced_mono = (protein_mass(INSULIN_A, monoisotopic=True)
                + protein_mass(INSULIN_B, monoisotopic=True))

print("insulin avg  {:.2f}".format(reduced_avg - 3 * DISULFIDE_AVERAGE))
print("insulin mono {:.4f}".format(reduced_mono - 3 * DISULFIDE_MONO))Code language: Python (python)

insulin avg  5807.62
insulin mono 5803.6376Code language: plaintext (plaintext)

Both land on the published values for intact human insulin, UniProt P01308.

Everything else is a delta you add. The monoisotopic ones worth memorising: phosphorylation +79.96633, oxidised methionine +15.99491, acetylation +42.01057, methylation +14.01565, deamidation +0.98402, and carbamidomethyl cysteine +57.02146 if the sample saw iodoacetamide, which in a proteomics workflow it almost always did. The rest live in Unimod, which is where the search engines get theirs. A cleaved signal peptide or an initiator methionine that got trimmed will move you by hundreds of Daltons, and no amount of decimal precision saves you from starting with the wrong sequence.

So when your calculated mass is off by 80, do not go hunting for a rounding bug in the table. Go look at what the protein is wearing.