Showing posts with label Computational Mass Spectrometry. Show all posts
Showing posts with label Computational Mass Spectrometry. Show all posts

Wednesday, 30 November 2016

Matching a Spectrum to a Protein

In this problem we are given a set of n protein strings and a multiset of weights corresponding to the complete spectrum of some unknown protein string. We are then asked to find the maximum multiplicity and the protein string that corresponds to it. My solution to the problem can be found on my GitHub. For each protein string, it calculates the weight of each prefix and suffix and checks if they correspond to any of the given weights of the protein spectrum. Each time they do a counter is increased by one for that specific protein string. Finally, the highest value is printed together with the corresponding string.

Sample Dataset
4
GSDMQS
VWICN
IASWMQS
PVSMGAD
445.17838
115.02694
186.07931
314.13789
317.1198
215.09061

Sample Output
3
IASWMQS

I actually didn't manage to replicate the result above, which is the example given in the problem description. However, when I tried my program with a larger dataset and uploaded my answer to Rosalind, I got the correct result, so I must have done something right!

Monday, 19 September 2016

Comparing Spectra with the Spectral Convolution

In this task we are given two multisets, S1 and S2, each representing simplified spectra taken from two peptides. We are asked to find the largest multiplicity of the  Minkowski difference (S1⊖S2) and the absolute value of the number x maximizing (S1⊖S2)(x). It took me a while to figure out what I was meant to find, but ultimately I figured that it was the number that occurs most frequently in the multiset formed by S1⊖S2, and its frequency.

Sample Dataset
186.07931 287.12699 548.20532 580.18077 681.22845 706.27446 782.27613 968.35544 968.35544
101.04768 158.06914 202.09536 318.09979 419.14747 463.17369

Expected Output
3
85.03163

Once I had realized what I was after, The programming itself became very easy. All I had to do was generate the multiset formed by S1⊖S2 and pick out the most frequently occurring number. As it turns out, there is a really useful built in library for this called collections. Using the function called most_common, I got precisely what I wanted. Have a look at the code below or on my GitHub.

from collections import Counter

data = []
with open('rosalind_conv.txt','r') as f:
    for line in f:
        data.append(line.strip())
S1 = [float(x) for x in data[0].split()]
S2 = [float(x) for x in data[1].split()]

#calculate Minkowski difference of S1 and S2
min_diff = []
for i in S1:
    for j in S2:
        diff = round(i - j, 5)
        min_diff.append(diff)

#find most common value and its frequency
x = Counter(min_diff).most_common(1)

print(x[0][1])
print(x[0][0])

Thursday, 15 September 2016

Inferring Protein from Spectrum

In this problem we are given a prefix spectrum and are asked to infer a protein sequence from it using the monoisotopic mass table.

Sample Dataset
3524.8542
3710.9335
3841.974
3970.0326
4057.0646

Expected Output
WMQS

This was a fairly simple problem to solve, and I guess the tricky part was to get the rounding right. In order to compare the results from the numbers in the dataset with the numbers in the mass table, we need to round all numbers to four decimals, otherwise they will differ slightly. The code below can also be found on my Github.

L = [float(line) for line in open('rosalind_spec.txt','r')]

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

aa_masses = []
for i in range(len(L) - 1):
    aa_mass = round(L[i + 1] - L[i], 4)
    aa_masses.append(aa_mass)

rnd_mass_table = {}
for k, v in mass_table.items():
    rnd_mass_table[round(v, 4)] = k

prot = ''
for aa in aa_masses:
    prot += rnd_mass_table[aa]

print(prot)

Tuesday, 12 July 2016

Calculating Protein Mass

In this problem we are asked to calculate the molecular weight of a peptide. The peptide is assumed to come from the middle of a protein so we are to assume that it consists entirely of amino acid residues (meaning we don't have to account for the extra weight of the "water molecule" present when we include the two edges of the protein).

I wrote two programs for this problem, the first one uses Biopython and the second doesn't. The following code is the one using Biopython:

from Bio.Seq import Seq                  
from Bio.Alphabet import generic_protein 
from Bio.SeqUtils import molecular_weight

with open('sampledata.txt', 'r') as f:   
    for line in f:                       
        prot_seq = line.strip('\n')      

print('%0.3f' % (molecular_weight(       
    Seq(prot_seq, generic_protein),      
    monoisotopic=True) - 18.01056))      

This program uses the Biopython function molecular_weight from SeqUtils. The function sets monoisotopic=False by default, but because the problem specified that we should use the monoisotopic weitghts we need to set it to monoisotopic=True. Also, the function includes the extra weight of one water molecule, so we need to remove that manually (hence the - 18.01056). Come to think of it there is an option called circular that states that the sequence has no ends. Maybe that would work as well?

Anyhow, I thought I'd try to make this program without the use of Biopython as well, and the following is what I came up with:

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

with open('sampledata.txt', 'r') as f:
    for line in f:                    
        prot_seq = line.strip('\n')   

weight = 0                            
for aa in prot_seq:                   
    weight += weights[aa]             

print('%0.3f' % weight)