mobiusklein / brainpy

A Python implementation of Baffling Recursive Algorithm for Isotopic distributioN calculations
http://mobiusklein.github.io/brainpy
Apache License 2.0
19 stars 10 forks source link

Not computing lighter isotopes #6

Open abefrandsen opened 11 months ago

abefrandsen commented 11 months ago

Greetings. Thanks for the great work on this project. I wanted to raise this issue in case it's not already known.

For some compounds, the isotopic_variants function doesn't seem to compute certain lighter isotopes. I found this in particular for compounds involving tin Sn. In particular, note the isotope pattern for tin as follows:

import brainpy
sn = brainpy.periodic_table["Sn"]
sn.isotopes
>> 
OrderedDict([(-8,
              Isotope(mass=111.904818, abundance=0.0097, neutrons=112, neutron_shift=-8)),
             (-6,
              Isotope(mass=113.902779, abundance=0.0066, neutrons=114, neutron_shift=-6)),
             (-5,
              Isotope(mass=114.903342, abundance=0.0034, neutrons=115, neutron_shift=-5)),
             (-4,
              Isotope(mass=115.901741, abundance=0.1454, neutrons=116, neutron_shift=-4)),
             (-3,
              Isotope(mass=116.902952, abundance=0.0768, neutrons=117, neutron_shift=-3)),
             (-2,
              Isotope(mass=117.901603, abundance=0.2422, neutrons=118, neutron_shift=-2)),
             (-1,
              Isotope(mass=118.903308, abundance=0.0859, neutrons=119, neutron_shift=-1)),
             (0,
              Isotope(mass=119.9021947, abundance=0.3258, neutrons=120, neutron_shift=0)),
             (2,
              Isotope(mass=121.903439, abundance=0.0463, neutrons=122, neutron_shift=2)),
             (4,
              Isotope(mass=123.9052739, abundance=0.0579, neutrons=124, neutron_shift=4))])

The most abundant isotope has mass 119.9, but the next most abundant isotope is lighter (117.9). When we use isotopic_variants on just this atom, we get

isotopic_variants({"Sn":1})
>>
[Peak(mz=119.902195, intensity=0.492386, charge=0),
 Peak(mz=122.042942, intensity=0.335025, charge=0),
 Peak(mz=123.115011, intensity=0.172589, charge=0)]

The lighter isotopes don't show up, but they should.

I'm not sure where in the code this problem arises. I will note that the BRAIN algorithm says to start with the lightest isotope, but maybe instead, your code starts with the most abundant isotope, and only then produces isotope peaks larger than that?

mobiusklein commented 11 months ago

This implementation starts from the monoisotopic form, and doesn't correctly use the lighter isotopes, which was a consequence of the C++ implementation I used as a reference. The problem partially resides in compute_isotopic_coefficients in the C implementation (what is used by default). The Python implementations are more scattered and hard to locate.

I don't have the time to investigate this further myself. If you're interested in looking into it, it'd be much appreciated.