diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Alignment/PeakAligner.cs b/src/MSDIAL5/MsdialCore/Algorithm/Alignment/PeakAligner.cs index 64ae19fa1..99d3c846f 100644 --- a/src/MSDIAL5/MsdialCore/Algorithm/Alignment/PeakAligner.cs +++ b/src/MSDIAL5/MsdialCore/Algorithm/Alignment/PeakAligner.cs @@ -1,8 +1,10 @@ using CompMs.Common.Components; using CompMs.Common.DataObj; using CompMs.Common.DataObj.Property; +using CompMs.Common.Enum; using CompMs.Common.Extension; using CompMs.Common.Utility; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using CompMs.MsdialCore.DataObj; using CompMs.MsdialCore.Parameter; using CompMs.MsdialCore.Parser; @@ -27,6 +29,12 @@ public class PeakAligner { public IDataProviderFactory ProviderFactory { get; set; } private IProgress? Progress { get; set; } + /// + /// What the lipid RT check did on the last alignment (hook B): its result for the review list, or why it did not + /// run. Null when the check does not apply to this kind of project. + /// + public LipidRtUnitRun? LipidRtRun { get; private set; } + public PeakAligner(AlignmentProcessFactory factory, IProgress? progress) { Accessor = factory.CreateDataAccessor(); Joiner = factory.CreatePeakJoiner(); @@ -55,6 +63,14 @@ public AlignmentResultContainer Alignment( IsotopeTracking.SetIsotopeTrackingID(container, Param, MspDB, null); } + // Hook B of the lipid RT check, shared by the Console and the GUI. It runs after every step that changes the + // spots' annotations (isotope tracking clears those of spots represented by a labelled file) and before the + // result is saved and exported, so its verdicts and tags describe the annotations that are saved. Isotope + // tracking replaces the spot collection, so it is read from the container here. LC-MS only in this version. + if (Param.MachineCategory == MachineCategory.LCMS) { + LipidRtRun = LipidRtCheckHooks.RunOnSpots(container.AlignmentSpotProperties, Param); + } + if (chromPeakInfoSerializer != null) SerializeSpotInfo( FlattenSpots(refined).ToList(), diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtAntichain.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtAntichain.cs new file mode 100644 index 000000000..657de75ca --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtAntichain.cs @@ -0,0 +1,459 @@ +using System; +using System.Collections.Generic; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// + /// The exact maximum-weight consistent subset of a discordance order (weighted Dilworth/Konig by max-flow). + /// + /// + /// Within one class, element i dominates j when C_i >= C_j, DBe_i <= DBe_j and O_i <= O_j (not all equal): + /// i is predicted to elute later. The pair is discordant when RT_i < RT_j - eps. Discordance is a strict + /// partial order, so a consistent subset (no discordant pair inside) is an antichain of it, and the maximum-weight + /// consistent subset is a maximum-weight antichain: W* = sum(w) - maxflow on the network s -> L_i (w_i), + /// L_i -> R_j (infinite) for each discordant pair, R_j -> t (w_j). + /// + /// Membership comes from one maximum flow and residual reachability: i is in EVERY maximum antichain when L_i is + /// reachable from s and R_i reaches t; in SOME when the residual closure of {s, L_i} contains neither t nor R_i; + /// otherwise in NONE. Membership is defined over all maximum antichains, so the result does not depend on the order + /// of the input or on which maximum flow the algorithm happens to find. + /// + /// Margins in weight units: e_i = W* - W*(without i) for an element in every maximum antichain, and + /// d_i = W* - (w_i + W*(without i and its discordant partners)) for an element in none. Each margin is a maximum + /// flow on the network with some elements removed; it starts from the component's maximum flow with the flow + /// through the removed elements cancelled (every path is s -> L_a -> R_b -> t, so cancelling is local), + /// which gives the same exact value as a flow from scratch at a fraction of the cost. + /// + /// Integer capacities keep every status and margin exact. One Dinic network per connected component. + /// + public static class LipidRtAntichain + { + private const long Infinity = 1L << 60; + + public sealed class Result + { + internal Result(LipidRtFeatureStatus[] status, long[] essential, long[] deficit, long weight, int[] component, HashSet[] partners) { + Status = status; + Essential = essential; + Deficit = deficit; + Weight = weight; + Component = component; + Partners = partners; + } + + public IReadOnlyList Status { get; } + + /// e_i for core elements, the own weight for free elements, 0 otherwise. + public IReadOnlyList Essential { get; } + + /// d_i for excluded elements, 0 otherwise. + public IReadOnlyList Deficit { get; } + + /// The maximum consistent weight W* of the whole set. + public long Weight { get; } + + /// The connected component of the discordance graph, or -1 for a free element. + public IReadOnlyList Component { get; } + + /// The discordant partners of each element, in both directions. + public IReadOnlyList> Partners { get; } + } + + /// + /// Analyses elements. are the discordant pairs (i, j) with i the + /// predicted-later element; the relation must be transitively closed, which discordance is. + /// are positive integers. + /// + public static Result Analyze(int n, IEnumerable<(int Later, int Earlier)> pairs, IReadOnlyList weights, bool margins = true) { + var below = new List[n]; + var partners = new HashSet[n]; + for (int i = 0; i < n; i++) { + below[i] = new List(); + partners[i] = new HashSet(); + } + var seenPairs = new HashSet<(int, int)>(); + foreach (var (i, j) in pairs) { + if (!seenPairs.Add((i, j))) { + continue; + } + below[i].Add(j); + partners[i].Add(j); + partners[j].Add(i); + } + // connected components of the comparability graph + var component = Enumerable.Repeat(-1, n).ToArray(); + var componentCount = 0; + for (int i = 0; i < n; i++) { + if (component[i] >= 0 || partners[i].Count == 0) { + continue; + } + var stack = new Stack(); + stack.Push(i); + component[i] = componentCount; + while (stack.Count > 0) { + var u = stack.Pop(); + foreach (var v in partners[u]) { + if (component[v] < 0) { + component[v] = componentCount; + stack.Push(v); + } + } + } + componentCount++; + } + + var status = new LipidRtFeatureStatus[n]; + var essential = new long[n]; + var deficit = new long[n]; + long total = 0; + for (int i = 0; i < n; i++) { + status[i] = LipidRtFeatureStatus.Free; + essential[i] = weights[i]; + if (partners[i].Count == 0) { + total += weights[i]; + } + } + var members = new List[componentCount]; + for (int c = 0; c < componentCount; c++) { + members[c] = new List(); + } + for (int i = 0; i < n; i++) { + if (component[i] >= 0) { + members[component[i]].Add(i); + } + } + foreach (var index in members) { + var network = new ComponentNetwork(index, weights, below); + var w = network.Sum - network.MaxFlow(); + total += w; + var fromSource = network.ReachFrom(network.Source); + var toSink = network.CoReach(network.Sink); + var m = index.Count; + for (int k = 0; k < m; k++) { + var x = index[k]; + int l = k, r = m + k; + if (fromSource[l] && toSink[r]) { + status[x] = LipidRtFeatureStatus.Core; + continue; + } + bool some; + if (fromSource[l]) { + some = !fromSource[r]; + } + else { + var closure = network.ReachFrom(l); + for (int v = 0; v < closure.Length; v++) { + closure[v] |= fromSource[v]; + } + some = !closure[network.Sink] && !closure[r]; + } + status[x] = some ? LipidRtFeatureStatus.Ambiguous : LipidRtFeatureStatus.Excluded; + } + if (margins) { + network.Snapshot(); + for (int k = 0; k < m; k++) { + var x = index[k]; + if (status[x] == LipidRtFeatureStatus.Core) { + essential[x] = w - network.WeightWithout(new[] { k }); + } + else { + essential[x] = 0; + } + if (status[x] == LipidRtFeatureStatus.Excluded) { + var removed = new List { k }; + for (int j = 0; j < m; j++) { + if (partners[x].Contains(index[j])) { + removed.Add(j); + } + } + deficit[x] = w - (weights[x] + network.WeightWithout(removed)); + } + } + } + } + return new Result(status, essential, deficit, total, component, partners); + } + + /// + /// The network of one component: s -> L_k (w), L_k -> R_j (infinite) per discordant pair, R_k -> t (w), for + /// the component's elements k = 0..m-1. L_k is node k, R_k is node m + k. + /// + private sealed class ComponentNetwork + { + private readonly int _m; + private readonly long[] _weights; + private readonly int[] _sourceEdge; // s -> L_k + private readonly int[] _sinkEdge; // R_k -> t + private readonly Dinic _dinic; + private long[] _snapshot; + + public ComponentNetwork(List index, IReadOnlyList weights, List[] below) { + _m = index.Count; + _weights = index.Select(x => weights[x]).ToArray(); + var position = new Dictionary(_m); + for (int k = 0; k < _m; k++) { + position[index[k]] = k; + } + _dinic = new Dinic(2 * _m + 2); + _sourceEdge = new int[_m]; + _sinkEdge = new int[_m]; + for (int k = 0; k < _m; k++) { + _sourceEdge[k] = _dinic.Add(Source, k, _weights[k]); + _sinkEdge[k] = _dinic.Add(_m + k, Sink, _weights[k]); + Sum += _weights[k]; + } + for (int k = 0; k < _m; k++) { + foreach (var j in below[index[k]]) { + _dinic.Add(k, _m + position[j], Infinity); + } + } + _dinic.Seal(); + } + + public int Source => 2 * _m; + public int Sink => 2 * _m + 1; + public long Sum { get; } + + public long MaxFlow() => _dinic.MaxFlow(Source, Sink); + + public bool[] ReachFrom(int start) => _dinic.ReachFrom(start); + + public bool[] CoReach(int target) => _dinic.CoReach(target); + + /// Keeps the residual of the component's maximum flow as the start of every margin flow. + public void Snapshot() { + _snapshot = _dinic.Capacities(); + } + + /// The maximum antichain weight of the component without the elements . + public long WeightWithout(IReadOnlyCollection removed) { + _dinic.Restore(_snapshot); + var gone = new bool[_m]; + foreach (var k in removed) { + gone[k] = true; + } + // Cancel the flow of every path s -> L_a -> R_b -> t through a removed element, keeping the flow of the + // remaining paths, then cut the removed elements off. + for (int a = 0; a < _m; a++) { + foreach (var e in _dinic.Edges(a)) { + var b = _dinic.Head(e) - _m; + if (b < 0 || b >= _m || !_dinic.IsForward(e)) { + continue; // not an L_a -> R_b edge + } + var flow = _dinic.Flow(e); + if (flow == 0 || (!gone[a] && !gone[b])) { + continue; + } + if (!gone[a]) { + _dinic.Push(_sourceEdge[a], -flow); + } + if (!gone[b]) { + _dinic.Push(_sinkEdge[b], -flow); + } + _dinic.Push(e, -flow); + } + } + long sum = 0; + for (int k = 0; k < _m; k++) { + if (gone[k]) { + _dinic.Isolate(k); + _dinic.Isolate(_m + k); + } + else { + sum += _weights[k]; + } + } + _dinic.MaxFlow(Source, Sink); + long flowValue = 0; + for (int k = 0; k < _m; k++) { + if (!gone[k]) { + flowValue += _dinic.Flow(_sourceEdge[k]); + } + } + return sum - flowValue; + } + } + + /// Dinic's maximum flow on integer capacities. Edge e and its reverse are e and e ^ 1. + private sealed class Dinic + { + private readonly int _n; + private readonly List[] _building; + private readonly List _toList = new List(); + private readonly List _capList = new List(); + private int[][] _graph; + private int[] _to; + private long[] _capacity; + private long[] _original; + + public Dinic(int n) { + _n = n; + _building = new List[n]; + for (int i = 0; i < n; i++) { + _building[i] = new List(); + } + } + + public int Add(int u, int v, long capacity) { + var e = _toList.Count; + _building[u].Add(e); + _toList.Add(v); + _capList.Add(capacity); + _building[v].Add(e + 1); + _toList.Add(u); + _capList.Add(0); + return e; + } + + public void Seal() { + _graph = _building.Select(l => l.ToArray()).ToArray(); + _to = _toList.ToArray(); + _capacity = _capList.ToArray(); + _original = _capList.ToArray(); + } + + public IReadOnlyList Edges(int u) => _graph[u]; + public int Head(int e) => _to[e]; + public bool IsForward(int e) => (e & 1) == 0; + + /// The flow on a forward edge. + public long Flow(int e) => _original[e] - _capacity[e]; + + /// Changes the flow on forward edge by . + public void Push(int e, long amount) { + _capacity[e] -= amount; + _capacity[e ^ 1] += amount; + } + + /// Cuts every edge at , in both directions. + public void Isolate(int u) { + foreach (var e in _graph[u]) { + _capacity[e] = 0; + _capacity[e ^ 1] = 0; + } + } + + public long[] Capacities() => (long[])_capacity.Clone(); + + public void Restore(long[] capacities) => Array.Copy(capacities, _capacity, capacities.Length); + + public long MaxFlow(int s, int t) { + long flow = 0; + var level = new int[_n]; + var iterator = new int[_n]; + var queue = new Queue(); + var path = new List(); + while (true) { + for (int i = 0; i < _n; i++) { + level[i] = -1; + } + level[s] = 0; + queue.Enqueue(s); + while (queue.Count > 0) { + var u = queue.Dequeue(); + foreach (var e in _graph[u]) { + if (_capacity[e] > 0 && level[_to[e]] < 0) { + level[_to[e]] = level[u] + 1; + queue.Enqueue(_to[e]); + } + } + } + if (level[t] < 0) { + return flow; + } + Array.Clear(iterator, 0, _n); + while (true) { + var f = Augment(s, t, level, iterator, path); + if (f == 0) { + break; + } + flow += f; + } + } + } + + // Iterative depth-first search along the level graph; returns the bottleneck pushed, 0 if none. + private long Augment(int s, int t, int[] level, int[] iterator, List path) { + path.Clear(); + var u = s; + while (true) { + if (u == t) { + var f = Infinity; + foreach (var e in path) { + if (_capacity[e] < f) { + f = _capacity[e]; + } + } + foreach (var e in path) { + _capacity[e] -= f; + _capacity[e ^ 1] += f; + } + return f; + } + var advanced = false; + var edges = _graph[u]; + while (iterator[u] < edges.Length) { + var e = edges[iterator[u]]; + if (_capacity[e] > 0 && level[_to[e]] == level[u] + 1) { + path.Add(e); + u = _to[e]; + advanced = true; + break; + } + iterator[u]++; + } + if (!advanced) { + level[u] = -1; // dead end + if (path.Count == 0) { + return 0; + } + var last = path[path.Count - 1]; + path.RemoveAt(path.Count - 1); + u = _to[last ^ 1]; + iterator[u]++; + } + } + } + + /// Nodes reachable from in the residual graph. + public bool[] ReachFrom(int start) { + var seen = new bool[_n]; + seen[start] = true; + var queue = new Queue(); + queue.Enqueue(start); + while (queue.Count > 0) { + var u = queue.Dequeue(); + foreach (var e in _graph[u]) { + var v = _to[e]; + if (_capacity[e] > 0 && !seen[v]) { + seen[v] = true; + queue.Enqueue(v); + } + } + } + return seen; + } + + /// Nodes that can reach in the residual graph. + public bool[] CoReach(int target) { + var seen = new bool[_n]; + seen[target] = true; + var queue = new Queue(); + queue.Enqueue(target); + while (queue.Count > 0) { + var v = queue.Dequeue(); + foreach (var e in _graph[v]) { + var u = _to[e]; + // the edge u -> v is the reverse of e + if (!seen[u] && _capacity[e ^ 1] > 0) { + seen[u] = true; + queue.Enqueue(u); + } + } + } + return seen; + } + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCandidate.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCandidate.cs new file mode 100644 index 000000000..6370299f4 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCandidate.cs @@ -0,0 +1,182 @@ +using CompMs.Common.Enum; +using System; +using System.Collections.Generic; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// + /// One annotated peak or alignment spot of a unit, as the lipid retention-time order check reads it (SPEC 1.2). + /// + /// + /// Every annotated object of the unit is a candidate and receives an outcome; only the eligible ones + /// (, SPEC 2.1) are evaluated, the others are NotEvaluable(Ineligible). An annotated object + /// that is not eligible still takes part in the MS1 artefact layer as a possible brighter partner, with its + /// composition when its name is a lipid. + /// + /// is a key and the last tie-breaker, never evidence. + /// + public sealed class LipidRtCandidate + { + /// MasterPeakID or MasterAlignmentID. + public int Id { get; set; } + + /// Retention time (min) on the unit's own axis. + public double Rt { get; set; } + + /// Precursor m/z. + public double Mz { get; set; } + + /// Peak height (top, or the alignment average). Only ratios and the class median of log10 are used. + public double Height { get; set; } + + /// The adduct, e.g. "[M+NH4]+". + public string Adduct { get; set; } = string.Empty; + + /// The species-level name (left of "|"), e.g. "PC 34:1". + public string Name { get; set; } = string.Empty; + + /// The chain-level name, or empty when the chains are not resolved. + public string Molecular { get; set; } = string.Empty; + + /// The Ontology (LbmClass name). + public string Ontology { get; set; } = string.Empty; + + /// The Ontology, except "BMP" for BMP-named rows exported with Ontology PG. + public string ClassKey { get; set; } = string.Empty; + + /// Species total carbons. + public int Carbon { get; set; } + + /// Species total double bonds, P- counted as O- plus one. + public int DoubleBondEquivalent { get; set; } + + /// Species total extra oxygens. + public int Oxygen { get; set; } + + /// True when the name is a lipid of a lipid Ontology and has a parsed composition with carbons. + public bool IsLipid { get; set; } + + /// An isotope-labelled name: an internal standard. + public bool IsLabelled { get; set; } + + /// The representative is an MS/MS match (annotation tag 3xx/4xx). + public bool IsMsmsAnnotated { get; set; } + + public LipidRtTagLevel TagLevel { get; set; } + + /// The representative's weighted dot product, 0..1. + public double WeightedDot { get; set; } + + /// The chains with roles, for the plausibility profile. + public IReadOnlyList Chains { get; set; } = Array.Empty(); + + public bool IsEligible { get; set; } + + /// Why the candidate is not eligible, e.g. "not_ref_matched,tag_520"; empty when eligible. + public string IneligibleReason { get; set; } = string.Empty; + + /// The composition key (class key, C, DBe, O), or null when the candidate has none. + public LipidRtComposition Composition => IsLipid && Carbon > 0 + ? new LipidRtComposition(string.IsNullOrEmpty(ClassKey) ? Ontology : ClassKey, Carbon, DoubleBondEquivalent, Oxygen) + : null; + } + + /// The species-level composition: class key, total carbons, double-bond equivalents and extra oxygens. + public sealed class LipidRtComposition : IEquatable, IComparable + { + public LipidRtComposition(string classKey, int carbon, int doubleBondEquivalent, int oxygen) { + ClassKey = classKey ?? string.Empty; + Carbon = carbon; + DoubleBondEquivalent = doubleBondEquivalent; + Oxygen = oxygen; + } + + public string ClassKey { get; } + public int Carbon { get; } + public int DoubleBondEquivalent { get; } + public int Oxygen { get; } + + public int CompareTo(LipidRtComposition other) { + var c = string.CompareOrdinal(ClassKey, other.ClassKey); + if (c != 0) return c; + c = Carbon.CompareTo(other.Carbon); + if (c != 0) return c; + c = DoubleBondEquivalent.CompareTo(other.DoubleBondEquivalent); + if (c != 0) return c; + return Oxygen.CompareTo(other.Oxygen); + } + + public bool Equals(LipidRtComposition other) { + return other is not null && ClassKey == other.ClassKey && Carbon == other.Carbon + && DoubleBondEquivalent == other.DoubleBondEquivalent && Oxygen == other.Oxygen; + } + + public override bool Equals(object obj) => Equals(obj as LipidRtComposition); + + public override int GetHashCode() { + unchecked { + var hash = ClassKey.GetHashCode(); + hash = hash * 31 + Carbon; + hash = hash * 31 + DoubleBondEquivalent; + return hash * 31 + Oxygen; + } + } + + /// "PC 34:1" or "OxPC 36:3;O2". + public override string ToString() { + return Oxygen != 0 + ? $"{ClassKey} {Carbon}:{DoubleBondEquivalent};O{Oxygen}" + : $"{ClassKey} {Carbon}:{DoubleBondEquivalent}"; + } + } + + /// + /// One peak or spot of the unit for the MS1 artefact layer: every peak, annotated or not (SPEC 1.2, "context"). + /// + /// + /// A candidate with the same overrides the context peak, so the context may hold every peak of the + /// unit including the annotated ones. + /// + public sealed class LipidRtContextPeak + { + public int Id { get; set; } + public double Rt { get; set; } + public double Mz { get; set; } + public double Height { get; set; } + public string Adduct { get; set; } = string.Empty; + + /// The representative is an MS/MS match (annotation tag 3xx/4xx). + public bool IsMsmsAnnotated { get; set; } + } + + /// One unit (an analysis file, or an alignment) to check. + public sealed class LipidRtRunInput + { + public LipidRtRunInput(IReadOnlyList candidates, IReadOnlyList context, double rtMin, double rtMax, IonMode polarity) { + Candidates = candidates ?? throw new ArgumentNullException(nameof(candidates)); + Context = context ?? Array.Empty(); + RtMin = rtMin; + RtMax = rtMax; + Polarity = polarity; + } + + /// Every annotated object of the unit; each receives an outcome. + public IReadOnlyList Candidates { get; } + + /// Every peak of the unit for the MS1 artefact layer, annotated or not. + public IReadOnlyList Context { get; } + + /// The smallest RT over all peaks of the unit, not only lipids. + public double RtMin { get; } + + /// The largest RT over all peaks of the unit. + public double RtMax { get; } + + public IonMode Polarity { get; } + + /// The RT span the windows scale with: max(RtMax - RtMin, 1 min); a missing bound counts as 0. + public double Span => Math.Max(Finite(RtMax) - Finite(RtMin), 1d); + + private static double Finite(double value) => double.IsNaN(value) || double.IsInfinity(value) ? 0d : value; + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtChain.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtChain.cs new file mode 100644 index 000000000..dc568bd04 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtChain.cs @@ -0,0 +1,108 @@ +using System; +using System.Collections.Generic; +using System.Globalization; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// One chain of a lipid name with its role, as the plausibility profile needs it. + /// + /// is the chain's own count as written: a plasmenyl "P-18:0" chain has 0 here. The + /// species-level double-bond equivalent that counts P- as O- plus one double bond lives on the composition. + /// + public sealed class LipidRtChain : IEquatable + { + public LipidRtChain(LipidRtChainRole role, int carbon, int doubleBond, int oxygen, string ether) { + Role = role; + Carbon = carbon; + DoubleBond = doubleBond; + Oxygen = oxygen; + Ether = ether ?? string.Empty; + } + + public LipidRtChainRole Role { get; } + public int Carbon { get; } + public int DoubleBond { get; } + public int Oxygen { get; } + + /// "O", "P" or empty. + public string Ether { get; } + + public bool IsFatty => Role == LipidRtChainRole.Acyl || Role == LipidRtChainRole.Alkyl || Role == LipidRtChainRole.NAcyl + || Role == LipidRtChainRole.EsterAcyl || Role == LipidRtChainRole.OAcyl; + + public bool IsUnresolved => Role == LipidRtChainRole.SpeciesCore || Role == LipidRtChainRole.SphingoidCore; + + public bool Equals(LipidRtChain other) { + return other is not null && Role == other.Role && Carbon == other.Carbon && DoubleBond == other.DoubleBond + && Oxygen == other.Oxygen && Ether == other.Ether; + } + + public override bool Equals(object obj) => Equals(obj as LipidRtChain); + + public override int GetHashCode() { + unchecked { + var hash = (int)Role; + hash = hash * 31 + Carbon; + hash = hash * 31 + DoubleBond; + hash = hash * 31 + Oxygen; + return hash * 31 + Ether.GetHashCode(); + } + } + + /// "acyl/18/1/0/" as the reference harness writes it. + public override string ToString() { + return string.Format(CultureInfo.InvariantCulture, "{0}/{1}/{2}/{3}/{4}", RoleToken(Role), Carbon, DoubleBond, Oxygen, Ether); + } + + /// Reads the "role/C/DB/O/ether|..." serialization of the reference harness. + public static IReadOnlyList ParseList(string text) { + var chains = new List(); + if (string.IsNullOrEmpty(text)) { + return chains; + } + foreach (var part in text.Split('|')) { + var fields = part.Split('/'); + if (fields.Length != 5) { + throw new FormatException($"Not a chain: '{part}'"); + } + chains.Add(new LipidRtChain( + ParseRole(fields[0]), + int.Parse(fields[1], CultureInfo.InvariantCulture), + int.Parse(fields[2], CultureInfo.InvariantCulture), + int.Parse(fields[3], CultureInfo.InvariantCulture), + fields[4])); + } + return chains; + } + + public static string RoleToken(LipidRtChainRole role) { + switch (role) { + case LipidRtChainRole.Acyl: return "acyl"; + case LipidRtChainRole.Alkyl: return "alkyl"; + case LipidRtChainRole.SphingoidBase: return "sph_base"; + case LipidRtChainRole.NAcyl: return "n_acyl"; + case LipidRtChainRole.EsterAcyl: return "ester_acyl"; + case LipidRtChainRole.OAcyl: return "o_acyl"; + case LipidRtChainRole.SterolCore: return "sterol_core"; + case LipidRtChainRole.SpeciesCore: return "species_core"; + case LipidRtChainRole.SphingoidCore: return "sph_core"; + default: throw new ArgumentOutOfRangeException(nameof(role)); + } + } + + public static LipidRtChainRole ParseRole(string token) { + switch (token) { + case "acyl": return LipidRtChainRole.Acyl; + case "alkyl": return LipidRtChainRole.Alkyl; + case "sph_base": return LipidRtChainRole.SphingoidBase; + case "n_acyl": return LipidRtChainRole.NAcyl; + case "ester_acyl": return LipidRtChainRole.EsterAcyl; + case "o_acyl": return LipidRtChainRole.OAcyl; + case "sterol_core": return LipidRtChainRole.SterolCore; + case "species_core": return LipidRtChainRole.SpeciesCore; + case "sph_core": return LipidRtChainRole.SphingoidCore; + default: throw new FormatException($"Unknown chain role '{token}'"); + } + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCheck.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCheck.cs new file mode 100644 index 000000000..52e375fc6 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCheck.cs @@ -0,0 +1,1030 @@ +using CompMs.Common.Enum; +using CompMs.MsdialCore.Parameter; +using System; +using System.Collections.Generic; +using System.Globalization; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// + /// The lipid retention-time order check (specification v1.0 with the addendum of 2026-10-15): a pure, deterministic + /// evaluation of one unit (an analysis file or an alignment). + /// + /// + /// A Confirmed means: in this run, the species-level composition of the annotation (class, total C, total DB with + /// P- counted as O- plus one, total O) is consistent with the reversed-phase rules among the run's other MS/MS-backed + /// annotations of the same class (more carbons elute later; more double bonds and more oxygens elute earlier), it is + /// bracketed by such neighbours on two axes, and it is not explained as an MS1 artefact of a brighter co-eluting ion. + /// It does not certify the chain assignment, the choice among isomers of one composition, or the absolute RT against + /// a library; the library RT is not used. + /// + /// The steps mirror the frozen Python reference (lipidrt_final.decide_full, poset.analyze, xpol.corroborate): + /// eligibility, plausibility, MS1 artefacts, features, the run-level RP gate, the exact maximum consistent subset per + /// class, adjacent directional support with ties and the class gate, the magnitude prediction, the feature outcome, + /// propagation to members, and the corroboration step. + /// + /// No randomness and no dependence on input order: every ordering is explicit, the max-flow capacities are + /// integers, medians follow numpy (mean of the two middle values) and weights round half to even. + /// + public sealed class LipidRtCheck + { + private const int Scale = 100; + + private static readonly string[] SecondaryTokens = { "Na", "K]", "K+", "Li", "2-", "2+", "3-", "3+" }; + + private static readonly Dictionary> PrimaryByClass = new Dictionary>(StringComparer.Ordinal) { + ["CL"] = new HashSet(StringComparer.Ordinal) { "[M-2H]2-", "[M-H]-" }, + ["DLCL"] = new HashSet(StringComparer.Ordinal) { "[M-H]-", "[M-2H]2-" }, + ["MLCL"] = new HashSet(StringComparer.Ordinal) { "[M-H]-", "[M-2H]2-" }, + }; + + private static readonly HashSet OxidisedClasses = new HashSet(StringComparer.Ordinal) { + "OxTG", "OxPC", "OxPE", "OxPG", "OxPI", "OxPS", "OxFA", "EtherOxPC", "EtherOxPE", "DMEDOxFA", + }; + + private sealed class Offset + { + public Offset(string kind, LipidRtArtefactKinds flag, double value, string partnerIon, string spotIon) { + Kind = kind; Flag = flag; Value = value; PartnerIon = partnerIon; SpotIon = spotIon; + } + public string Kind { get; } + public LipidRtArtefactKinds Flag { get; } + public double Value { get; } + /// The ion an annotated partner must carry for the relation to hold. + public string PartnerIon { get; } + /// A spot annotated as this ion IS that ion: a duplicate, not an artefact. + public string SpotIon { get; } + } + + // m/z(spot) - m/z(brighter co-eluting partner) + private static readonly Offset[] PositiveOffsets = { + new Offset("iso1", LipidRtArtefactKinds.Iso1, 1.00335, null, null), + new Offset("iso2", LipidRtArtefactKinds.Iso2, 2.00671, null, null), + new Offset("naH", LipidRtArtefactKinds.NaH, 21.98194, "[M+H]+", "[M+Na]+"), + new Offset("naNH4", LipidRtArtefactKinds.NaNH4, 4.95539, "[M+NH4]+", "[M+Na]+"), + new Offset("kH", LipidRtArtefactKinds.KH, 37.95588, "[M+H]+", "[M+K]+"), + new Offset("kNa", LipidRtArtefactKinds.KNa, 15.97394, "[M+Na]+", "[M+K]+"), + new Offset("nh4H", LipidRtArtefactKinds.NH4H, 17.02655, "[M+H]+", "[M+NH4]+"), + new Offset("acn", LipidRtArtefactKinds.Acn, 41.02655, "[M+H]+", "[M+ACN+H]+"), + new Offset("water", LipidRtArtefactKinds.Water, -18.01056, null, "[M+H-H2O]+"), + }; + + private static readonly Offset[] NegativeOffsets = { + new Offset("iso1", LipidRtArtefactKinds.Iso1, 1.00335, null, null), + new Offset("iso2", LipidRtArtefactKinds.Iso2, 2.00671, null, null), + new Offset("acetate", LipidRtArtefactKinds.Acetate, 60.02113, "[M-H]-", "[M+CH3COO]-"), + new Offset("formate", LipidRtArtefactKinds.Formate, 46.00548, "[M-H]-", "[M+HCOO]-"), + new Offset("chloride", LipidRtArtefactKinds.Chloride, 35.97668, "[M-H]-", "[M+Cl]-"), + new Offset("demethyl", LipidRtArtefactKinds.Demethyl, -74.03678, "[M+CH3COO]-", null), + new Offset("demethyl", LipidRtArtefactKinds.Demethyl, -60.02113, "[M+HCOO]-", null), + new Offset("water", LipidRtArtefactKinds.Water, -18.01056, null, "[M-H-H2O]-"), + }; + + /// Fragment class -> parent classes: the parent is the fragment plus one realistic chain, at the same RT. + private static readonly Dictionary InSourceChainLoss = new Dictionary(StringComparer.Ordinal) { + ["LPC"] = new[] { "PC" }, ["LPE"] = new[] { "PE" }, ["LPI"] = new[] { "PI" }, ["LPS"] = new[] { "PS" }, ["LPG"] = new[] { "PG" }, ["LPA"] = new[] { "PA" }, + ["EtherLPC"] = new[] { "EtherPC" }, ["EtherLPE"] = new[] { "EtherPE" }, ["MG"] = new[] { "DG" }, ["DG"] = new[] { "TG" }, ["EtherDG"] = new[] { "EtherTG" }, + ["FA"] = new[] { "PC", "PE", "PI", "PS", "PG", "PA", "TG", "DG", "LPC", "LPE", "EtherPC", "EtherPE" }, + }; + + /// Fragment class -> parent classes with the identical composition, at the same RT. + private static readonly Dictionary InSourceHeadLoss = new Dictionary(StringComparer.Ordinal) { + ["Cer_NS"] = new[] { "SM", "HexCer_NS", "Hex2Cer", "CerP" }, ["Cer_NDS"] = new[] { "SM", "HexCer_NDS" }, ["Cer_AS"] = new[] { "HexCer_AS" }, + ["Cer_NP"] = new[] { "HexCer_AP" }, ["DG"] = new[] { "PC", "PE", "PI", "PS", "PG", "PA", "MGDG", "DGDG", "SQDG" }, + ["EtherDG"] = new[] { "EtherPC", "EtherPE" }, ["LPA"] = new[] { "LPC", "LPE", "LPI", "LPS", "LPG" }, ["PA"] = new[] { "PC", "PE", "PI", "PS", "PG" }, + ["MG"] = new[] { "LPC", "LPE", "LPI", "LPS", "LPG" }, ["Sph"] = new[] { "SM", "Cer_NS" }, ["DHSph"] = new[] { "SM", "Cer_NDS" }, + ["DMPE"] = new[] { "PC" }, // NEG: PC [M-CH3]- has the formula of DMPE [M-H]- + }; + + /// Oxidised class -> its base class. + private static readonly Dictionary OxidisedBase = new Dictionary(StringComparer.Ordinal) { + ["OxTG"] = "TG", ["OxPC"] = "PC", ["OxPE"] = "PE", ["OxPI"] = "PI", ["OxPG"] = "PG", ["OxPS"] = "PS", ["OxFA"] = "FA", + ["EtherOxPC"] = "EtherPC", ["EtherOxPE"] = "EtherPE", ["DMEDOxFA"] = "DMEDFA", + }; + + private readonly LipidRtCheckParameter _parameter; + private readonly ILipidRtPlausibilityProfile _profile; + + /// A check with the plausibility profile the parameter names (). + public LipidRtCheck(LipidRtCheckParameter parameter) + : this(parameter, LipidRtChainRuleProfile.FromParameter((parameter ?? throw new ArgumentNullException(nameof(parameter))).PlausibilityProfile)) { + } + + /// A check with an explicit plausibility profile; null means no profile, whatever the parameter says. + public LipidRtCheck(LipidRtCheckParameter parameter, ILipidRtPlausibilityProfile profile) { + _parameter = parameter ?? throw new ArgumentNullException(nameof(parameter)); + _profile = profile; + } + + public LipidRtCheckParameter Parameter => _parameter; + + public ILipidRtPlausibilityProfile Profile => _profile; + + /// + /// Evaluates one unit and applies the corroboration available within it: the second-adduct flag, the class-level + /// policy and the Strict operating point (step 11a and 11c). This is what a single-polarity Console or GUI run + /// produces. adds the other polarity later. + /// + public LipidRtRunResult Evaluate(LipidRtRunInput input) { + var single = EvaluateSingleRun(input); + return LipidRtCorroboration.Apply(new[] { single }, _parameter)[0]; + } + + /// Steps 1-10 only: the unit's own evidence, before any corroboration. + public LipidRtRunResult EvaluateSingleRun(LipidRtRunInput input) { + if (input is null) { + throw new ArgumentNullException(nameof(input)); + } + var p = _parameter; + var rtMin = Finite(input.RtMin); + var rtMax = Finite(input.RtMax); + var span = Math.Max(rtMax - rtMin, 1d); + var eps = p.DiscordanceFraction * span; + var tie = p.TieFraction * span; + var win = p.CollapseFraction * span; + var tailStart = rtMax - p.TailFraction * span; + var negative = IsNegative(input); + + var results = new Dictionary(); + var candidates = input.Candidates; + var ids = new HashSet(); + foreach (var c in candidates) { + if (!ids.Add(c.Id)) { + throw new ArgumentException($"Candidate Ids must be unique within a unit; {c.Id} is repeated.", nameof(input)); + } + } + foreach (var c in candidates) { + if (!IsEvaluable(c)) { + results[c.Id] = Single(c, LipidRtOutcome.NotEvaluable, LipidRtReason.Ineligible, null, new LipidRtEvidence()); + } + } + var eligible = candidates.Where(IsEvaluable).ToList(); + if (eligible.Count == 0) { + return new LipidRtRunResult(input, candidates.Select(c => results[c.Id]).ToList(), null); + } + + // Step 2: plausibility + var plausibility = new Dictionary(); + if (_profile != null) { + foreach (var c in eligible) { + plausibility[c.Id] = _profile.Evaluate(c.Ontology, c.Chains); + } + } + + // Step 3: MS1 artefact explanations + var pool = BuildPool(candidates, input.Context); + var artefacts = ExplainArtefacts(eligible, pool, negative, span); + + // Step 4: features + var features = BuildFeatures(eligible, win, plausibility, artefacts); + + // Step 5: RP gate + var gate = RunGate(features, span, tie); + if (!gate.Passed) { + foreach (var c in eligible) { + results[c.Id] = Single(c, LipidRtOutcome.NotEvaluable, LipidRtReason.Gate, null, new LipidRtEvidence { Gate = gate }); + } + return new LipidRtRunResult(input, candidates.Select(c => results[c.Id]).ToList(), gate); + } + + // Steps 6-7 per class, then the run-level DB/C step ratio, then step 8 + var byClass = features.GroupBy(f => f.ClassKey, StringComparer.Ordinal) + .OrderBy(g => g.Key, StringComparer.Ordinal) + .Select(g => g.ToList()) + .ToList(); + var classInfos = new List(); + foreach (var fs in byClass) { + classInfos.Add(fs.Count >= 2 ? AnalyzeClass(fs, eps, tie) : ClassInfo.Single()); + } + var ratios = new List(); + foreach (var info in classInfos) { + if (info.StepD is double d && info.StepC is double cs && cs > 0 && d > 0) { + ratios.Add(d / cs); + } + } + double? kdRun = ratios.Count >= 2 ? Median(ratios) : (double?)null; + for (int k = 0; k < byClass.Count; k++) { + Magnitude(byClass[k], classInfos[k], kdRun); + } + + // Steps 9-10 + for (int k = 0; k < byClass.Count; k++) { + foreach (var f in byClass[k]) { + DecideFeature(f, byClass[k], classInfos[k], plausibility, artefacts, tailStart, results); + } + } + return new LipidRtRunResult(input, candidates.Select(c => results[c.Id]).ToList(), gate); + } + + private static double Finite(double value) => double.IsNaN(value) || double.IsInfinity(value) ? 0d : value; + + /// An eligible candidate with a composition; anything else is NotEvaluable(Ineligible). + private static bool IsEvaluable(LipidRtCandidate c) => c.IsEligible && c.Composition != null && !double.IsNaN(c.Rt); + + private static bool IsNegative(LipidRtRunInput input) { + switch (input.Polarity) { + case IonMode.Negative: + return true; + case IonMode.Positive: + return false; + default: + return input.Candidates.Any(c => (c.Adduct ?? string.Empty).EndsWith("-", StringComparison.Ordinal)); + } + } + + private static LipidRtObjectResult Single(LipidRtCandidate c, LipidRtOutcome outcome, LipidRtReason reason, double? score, LipidRtEvidence evidence) { + return new LipidRtObjectResult + { + Id = c.Id, + Candidate = c, + Outcome = outcome, + Reason = reason, + Score = score, + SingleRunOutcome = outcome, + SingleRunReason = reason, + SingleRunScore = score, + Evidence = evidence, + }; + } + + internal static bool IsSecondary(string adduct, string classKey) { + if (adduct != null && PrimaryByClass.TryGetValue(classKey ?? string.Empty, out var primary) && primary.Contains(adduct)) { + return false; + } + var a = adduct ?? string.Empty; + return SecondaryTokens.Any(t => a.Contains(t)); + } + + internal static bool IsOxidisedClass(string classKey) => OxidisedClasses.Contains(classKey ?? string.Empty); + + internal static string FormatFeature(string classKey, int c, int d, int o, double rt) { + return string.Format(CultureInfo.InvariantCulture, "{0} {1}:{2}{3}@{4:0.00}", classKey, c, d, o != 0 ? ";O" + o.ToString(CultureInfo.InvariantCulture) : string.Empty, rt); + } + + private static double NonNegativeHeight(double height) => double.IsNaN(height) ? 0d : height; + + private static double Dot(double dot) => double.IsNaN(dot) ? 0d : dot; + + // ------------------------------------------------------------------------------------------ MS1 artefact layer + + private sealed class PoolRow + { + public int Id; + public double Rt; + public double Mz; + public double Height; + public string Adduct = string.Empty; + public bool IsMsmsAnnotated; + public LipidRtComposition Composition; + public string ClassKey = string.Empty; + public int C, D, O; + public bool IsLabelled; + public bool IsLipid; + public string Name = string.Empty; + } + + private static List BuildPool(IReadOnlyList candidates, IReadOnlyList context) { + var rows = new Dictionary(); + foreach (var c in context) { + rows[c.Id] = new PoolRow + { + Id = c.Id, Rt = c.Rt, Mz = c.Mz, Height = NonNegativeHeight(c.Height), Adduct = c.Adduct ?? string.Empty, + IsMsmsAnnotated = c.IsMsmsAnnotated, + }; + } + foreach (var s in candidates) { + var composition = s.Composition; + rows[s.Id] = new PoolRow + { + Id = s.Id, Rt = s.Rt, Mz = s.Mz, Height = NonNegativeHeight(s.Height), Adduct = s.Adduct ?? string.Empty, + IsMsmsAnnotated = s.IsMsmsAnnotated, + Composition = composition, + ClassKey = string.IsNullOrEmpty(s.ClassKey) ? s.Ontology ?? string.Empty : s.ClassKey, + C = s.Carbon, D = s.DoubleBondEquivalent, O = s.Oxygen, + IsLabelled = s.IsLabelled, + IsLipid = composition != null, + Name = s.Name ?? string.Empty, + }; + } + return rows.Values.Where(r => !double.IsNaN(r.Rt)) + .OrderBy(r => r.Rt).ThenBy(r => r.Id) + .ToList(); + } + + private static int LowerBound(List pool, double value) { + int lo = 0, hi = pool.Count; + while (lo < hi) { + var mid = (lo + hi) >> 1; + if (pool[mid].Rt < value) lo = mid + 1; else hi = mid; + } + return lo; + } + + private static int UpperBound(List pool, double value) { + int lo = 0, hi = pool.Count; + while (lo < hi) { + var mid = (lo + hi) >> 1; + if (pool[mid].Rt <= value) lo = mid + 1; else hi = mid; + } + return lo; + } + + private Dictionary> ExplainArtefacts(List eligible, List pool, bool negative, double span) { + var p = _parameter; + var window = p.CoelutionFraction * span; + var kinds = p.ArtefactKinds; + var offsets = (negative ? NegativeOffsets : PositiveOffsets).Where(o => (kinds & o.Flag) != 0).ToArray(); + var result = new Dictionary>(); + foreach (var s in eligible) { + var mz = s.Mz; + var hs = NonNegativeHeight(s.Height); + var rt = s.Rt; + var mine = s.Composition; + var flags = new List(); + var lo = LowerBound(pool, rt - window); + var hi = UpperBound(pool, rt + window); + if (!double.IsNaN(mz) && offsets.Length > 0) { + var tol = p.MzToleranceDa + p.MzTolerancePpm * 1e-6 * Math.Abs(mz); + for (int k = lo; k < hi; k++) { + var c = pool[k]; + if (c.Id == s.Id || double.IsNaN(c.Mz) || Equals(c.Composition, mine)) { + continue; + } + var d = mz - c.Mz; + foreach (var offset in offsets) { + if (Math.Abs(d - offset.Value) > tol) { + continue; + } + var ratio = offset.Kind == "iso1" ? p.Iso1Ratio : offset.Kind == "iso2" ? p.Iso2Ratio : p.ParentRatio; + if (c.Height < ratio * hs) { + continue; + } + if (offset.SpotIon != null && s.Adduct == offset.SpotIon) { + continue; // the spot IS annotated as that ion: a duplicate, not an artefact + } + if (offset.PartnerIon != null && c.IsMsmsAnnotated && c.Adduct != offset.PartnerIon) { + continue; // an annotated partner must carry the ion the relation assumes + } + var partner = c.Name.Length > 0 ? c.Name : string.Format(CultureInfo.InvariantCulture, "m/z {0:0.0000}", c.Mz); + var fold = hs > 0 ? (c.Height / hs).ToString("0.0", CultureInfo.InvariantCulture) : "inf"; + flags.Add(new LipidRtArtefact(offset.Kind, string.Format(CultureInfo.InvariantCulture, + "{0}: {1:0.0000} is {2} {3}{4:0.0000} of a {5}x brighter co-eluting {6}@{7:0.00}", + offset.Kind, mz, partner, offset.Value >= 0 ? "+" : "-", Math.Abs(offset.Value), fold, c.Adduct, c.Rt), c.Id)); + break; + } + } + } + var cls = mine?.ClassKey ?? string.Empty; + var chainParents = (kinds & LipidRtArtefactKinds.IsfChain) != 0 && InSourceChainLoss.TryGetValue(cls, out var pc) ? pc : Array.Empty(); + var headParents = (kinds & LipidRtArtefactKinds.IsfHead) != 0 && InSourceHeadLoss.TryGetValue(cls, out var ph) ? ph : Array.Empty(); + var oxBase = (kinds & LipidRtArtefactKinds.Ox) != 0 && OxidisedBase.TryGetValue(cls, out var ob) ? ob : null; + if (chainParents.Length > 0 || headParents.Length > 0 || oxBase != null) { + for (int k = lo; k < hi; k++) { + var a = pool[k]; + if (!a.IsLipid || a.IsLabelled || a.Id == s.Id || a.Height < p.ParentRatio * hs) { + continue; + } + int dC = a.C - s.Carbon, dD = a.D - s.DoubleBondEquivalent, dO = a.O - s.Oxygen; + if (chainParents.Contains(a.ClassKey) && 12 <= dC && dC <= 26 && 0 <= dD && dD <= 6 && dC - 2 * dD >= 8 && dO >= 0) { + flags.Add(new LipidRtArtefact("isf_chain", string.Format(CultureInfo.InvariantCulture, + "isf: chain-loss fragment of brighter {0} {1}@{2:0.00}", a.ClassKey, a.Name, a.Rt), a.Id)); + break; + } + if (headParents.Contains(a.ClassKey) && dC == 0 && dD == 0 && dO == 0) { + flags.Add(new LipidRtArtefact("isf_head", string.Format(CultureInfo.InvariantCulture, + "isf: head-group loss of brighter {0} {1}@{2:0.00}", a.ClassKey, a.Name, a.Rt), a.Id)); + break; + } + if (oxBase != null && a.ClassKey == oxBase && a.D == s.DoubleBondEquivalent && a.O == 0 && (a.C == s.Carbon || a.C == s.Carbon - 1)) { + flags.Add(new LipidRtArtefact("ox", string.Format(CultureInfo.InvariantCulture, + "ox: co-elutes with brighter base species {0} {1}@{2:0.00} (more O should elute earlier)", a.ClassKey, a.Name, a.Rt), a.Id)); + break; + } + } + } + if (flags.Count > 0) { + result[s.Id] = flags; + } + } + return result; + } + + // -------------------------------------------------------------------------------------------------- features + + private sealed class Feature + { + public string ClassKey; + public int C, D, O; + public double Rt; + public List Members; + public LipidRtCandidate Best; + public double W0; + public double Multiplier; + public double MaxHeight; + public long W = 1; + public LipidRtTagLevel Level; + public bool Hard; + public bool Artefact; + public bool SecondaryOnly; + public int AdductCount; + public string Name; + public string Molecular; + public double HeightMedian; + public LipidRtFeatureStatus Status = LipidRtFeatureStatus.Free; + public double E; + public double Deficit; + public readonly SortedSet Directions = new SortedSet(StringComparer.Ordinal); + public readonly List<(string Label, int Index, double Gap)> Support = new List<(string, int, double)>(); + public List Discordant = new List(); + public int Ties; + public readonly List TiePartners = new List(); + public double? Predicted, Residual, Delta; + public Feature MinorDuplicateOf; + + public string Format() => FormatFeature(ClassKey, C, D, O, Rt); + } + + private static (bool Primary, int Rank, double Dot, double Height, double NegRt) BestKey(LipidRtCandidate s, string classKey) { + return (!IsSecondary(s.Adduct, classKey), LevelRank(s.TagLevel), Dot(s.WeightedDot), NonNegativeHeight(s.Height), -s.Rt); + } + + private static int LevelRank(LipidRtTagLevel level) { + switch (level) { + case LipidRtTagLevel.Chain: return 3; + case LipidRtTagLevel.Msms: return 2; + case LipidRtTagLevel.Class: return 1; + default: return 0; + } + } + + private double LevelWeight(LipidRtTagLevel level) { + switch (level) { + case LipidRtTagLevel.Chain: return _parameter.LevelWeightChain; + case LipidRtTagLevel.Msms: return _parameter.LevelWeightMsms; + case LipidRtTagLevel.Class: return _parameter.LevelWeightClass; + default: return 0.5; + } + } + + private static int CompareBestKeys((bool Primary, int Rank, double Dot, double Height, double NegRt) a, (bool Primary, int Rank, double Dot, double Height, double NegRt) b) { + var c = a.Primary.CompareTo(b.Primary); + if (c != 0) return c; + c = a.Rank.CompareTo(b.Rank); + if (c != 0) return c; + c = a.Dot.CompareTo(b.Dot); + if (c != 0) return c; + c = a.Height.CompareTo(b.Height); + if (c != 0) return c; + return a.NegRt.CompareTo(b.NegRt); + } + + private List BuildFeatures(List eligible, double win, Dictionary plausibility, Dictionary> artefacts) { + var p = _parameter; + var groups = eligible.GroupBy(s => s.Composition).OrderBy(g => g.Key).ToList(); + var features = new List(); + foreach (var group in groups) { + var key = group.Key; + var rs = group.OrderBy(s => s.Rt).ThenBy(s => s.Id).ToList(); + var chains = new List>(); + var current = new List { rs[0] }; + for (int i = 1; i < rs.Count; i++) { + if (rs[i].Rt - current[current.Count - 1].Rt <= win) { + current.Add(rs[i]); + } + else { + chains.Add(current); + current = new List { rs[i] }; + } + } + chains.Add(current); + foreach (var members in chains) { + // the best member: first maximum in (Rt, Id) order + var best = members[0]; + var bestKey = BestKey(best, key.ClassKey); + for (int i = 1; i < members.Count; i++) { + var k = BestKey(members[i], key.ClassKey); + if (CompareBestKeys(k, bestKey) > 0) { + best = members[i]; + bestKey = k; + } + } + var adducts = members.Select(m => m.Adduct ?? string.Empty).Distinct(StringComparer.Ordinal).Count(); + var w = LevelWeight(best.TagLevel) + p.DotWeight * Dot(best.WeightedDot) + p.DuplicateWeight * Math.Min(adducts - 1, p.DuplicateCap); + var hard = plausibility.Count > 0 && plausibility[best.Id].IsHard; + var artefact = members.All(m => artefacts.ContainsKey(m.Id)); + var multiplier = (hard ? p.ImplausibleWeight : 1d) * (artefact ? p.ArtefactWeight : 1d); + features.Add(new Feature + { + ClassKey = key.ClassKey, C = key.Carbon, D = key.DoubleBondEquivalent, O = key.Oxygen, + Rt = best.Rt, Members = members, Best = best, + W0 = Math.Max(w, 0.1) * multiplier, Multiplier = multiplier, + MaxHeight = members.Max(m => NonNegativeHeight(m.Height)), + Level = best.TagLevel, Hard = hard, Artefact = artefact, + SecondaryOnly = members.All(m => IsSecondary(m.Adduct, key.ClassKey)), + AdductCount = adducts, Name = best.Name ?? string.Empty, Molecular = best.Molecular ?? string.Empty, + }); + } + } + foreach (var fs in features.GroupBy(f => f.ClassKey, StringComparer.Ordinal)) { + var list = fs.ToList(); + var logs = list.Where(f => f.MaxHeight > 0).Select(f => Math.Log10(f.MaxHeight)).ToList(); + var median = logs.Count > 0 ? Median(logs) : 0d; + foreach (var f in list) { + f.HeightMedian = median; + var hb = f.MaxHeight > 0 ? Math.Min(Math.Max(Math.Log10(f.MaxHeight) - median, 0d), 3d) : 0d; + f.W = Math.Max(1L, (long)Math.Round(Scale * (f.W0 + p.HeightWeight * hb * f.Multiplier), MidpointRounding.ToEven)); + } + // minor duplicate: same species name, same (or unresolved) molecular name, at another RT, less intense + foreach (var byName in list.GroupBy(f => f.Name, StringComparer.Ordinal)) { + var gs = byName.ToList(); + if (gs.Count < 2) { + continue; + } + foreach (var f in gs) { + foreach (var g in gs) { + if (ReferenceEquals(g, f) || (f.Molecular.Length > 0 && g.Molecular.Length > 0 && f.Molecular != g.Molecular)) { + continue; + } + if (g.MaxHeight > f.MaxHeight || (g.MaxHeight == f.MaxHeight && -g.Rt > -f.Rt)) { + f.MinorDuplicateOf = g; + break; + } + } + } + } + } + return features; + } + + // ------------------------------------------------------------------------------------------------------ gate + + private static string AdjacentKind(int dC, int dD, int dO) { + int a = Math.Abs(dC), b = Math.Abs(dD), c = Math.Abs(dO); + if (b == 0 && c == 0 && 1 <= a && a <= 2) return "C"; + if (a == 0 && c == 0 && b == 1) return "D"; + if (a == 0 && b == 0 && c == 1) return "O"; + return null; + } + + private LipidRtGateStatistics RunGate(List features, double span, double tie) { + var p = _parameter; + if (p.Chromatography == LipidRtChromatography.Other) { + return new LipidRtGateStatistics { Passed = false, ChromatographyIsOther = true }; + } + var slopes = new List(); + int allConcordant = 0, allTotal = 0, trustedConcordant = 0, trustedTotal = 0; + foreach (var group in features.Where(f => !f.Artefact && !f.Hard).GroupBy(f => f.ClassKey, StringComparer.Ordinal)) { + var fs = group.ToList(); + if (fs.Count >= 5 && fs.Select(f => f.C).Distinct().Count() >= 3) { + var sl = new List(); + foreach (var a in fs) { + foreach (var b in fs) { + if (a.D == b.D && a.O == b.O && b.C > a.C) { + sl.Add((b.Rt - a.Rt) / (b.C - a.C)); + } + } + } + if (sl.Count > 0) { + slopes.Add(Median(sl)); + } + } + for (int i = 0; i < fs.Count; i++) { + var a = fs[i]; + for (int j = i + 1; j < fs.Count; j++) { + var b = fs[j]; + if (AdjacentKind(a.C - b.C, a.D - b.D, a.O - b.O) != "C") { + continue; + } + var (hi, lo) = a.C > b.C ? (a, b) : (b, a); + var ok = (hi.Rt - lo.Rt) > tie; + allConcordant += ok ? 1 : 0; + allTotal++; + if (IsTrusted(a) && IsTrusted(b)) { + trustedConcordant += ok ? 1 : 0; + trustedTotal++; + } + } + } + } + var median = slopes.Count > 0 ? Median(slopes) : 0d; + var useTrusted = trustedTotal >= p.GateMinPairs; + var concordant = useTrusted ? trustedConcordant : allConcordant; + var total = useTrusted ? trustedTotal : allTotal; + double? fraction = total > 0 ? (double)concordant / total : (double?)null; + var slopeOk = slopes.Count > 0 && median > p.SlopeFraction * span; + var concordanceOk = total >= p.GateMinPairs && fraction >= p.GateConcordance; + return new LipidRtGateStatistics + { + Passed = slopeOk && concordanceOk, + MedianCarbonSlope = median, + SlopeThreshold = p.SlopeFraction * span, + ClassCount = slopes.Count, + Concordance = fraction, + PairCount = total, + ConcordanceOn = useTrusted ? "trusted" : "all", + ConcordanceAllPairs = allTotal > 0 ? (double)allConcordant / allTotal : (double?)null, + }; + } + + private static bool IsTrusted(Feature f) { + return f.Level == LipidRtTagLevel.Chain && f.MaxHeight > 0 && Math.Log10(f.MaxHeight) >= f.HeightMedian; + } + + // -------------------------------------------------------------------------------------- class order analysis + + private sealed class ClassInfo + { + public int AdjacentTotal; + public double? AdjacentConcordance; + public double? StepC, StepD, StepO; + public bool GateOk = true; + public Dictionary> Neighbours = new Dictionary>(); + + public static ClassInfo Single() => new ClassInfo(); + + public double? Step(string kind) => kind == "C" ? StepC : kind == "D" ? StepD : StepO; + } + + private static bool Dominates(Feature i, Feature j) { + return i.C >= j.C && i.D <= j.D && i.O <= j.O && !(i.C == j.C && i.D == j.D && i.O == j.O); + } + + private ClassInfo AnalyzeClass(List fs, double eps, double tie) { + var p = _parameter; + var n = fs.Count; + var discordant = new List<(int, int)>(); + var adjacent = new List<(int I, int J)>(); + for (int i = 0; i < n; i++) { + for (int j = 0; j < n; j++) { + if (i == j || !Dominates(fs[i], fs[j])) { + continue; + } + var gap = fs[i].Rt - fs[j].Rt; + if (gap < -eps) { + discordant.Add((i, j)); + } + if (AdjacentKind(fs[i].C - fs[j].C, fs[i].D - fs[j].D, fs[i].O - fs[j].O) != null) { + adjacent.Add((i, j)); + } + } + } + var analysis = LipidRtAntichain.Analyze(n, discordant, fs.Select(f => f.W).ToList(), margins: true); + for (int k = 0; k < n; k++) { + var f = fs[k]; + f.Status = analysis.Status[k]; + f.E = (double)analysis.Essential[k] / Scale; + f.Deficit = (double)analysis.Deficit[k] / Scale; + f.Discordant = analysis.Partners[k].OrderByDescending(j => fs[j].W).ThenBy(j => fs[j].Rt).ThenBy(j => j).ToList(); + } + + var steps = new Dictionary> { ["C"] = new List(), ["D"] = new List(), ["O"] = new List() }; + int concordant = 0, total = 0; + foreach (var (i, j) in adjacent) { + var kind = AdjacentKind(fs[i].C - fs[j].C, fs[i].D - fs[j].D, fs[i].O - fs[j].O); + var g = fs[i].Rt - fs[j].Rt; + if (fs[i].Level == LipidRtTagLevel.Chain && fs[j].Level == LipidRtTagLevel.Chain) { + total++; + concordant += g > tie ? 1 : 0; + } + if (fs[i].Status != LipidRtFeatureStatus.Excluded && fs[j].Status != LipidRtFeatureStatus.Excluded) { + steps[kind].Add(g / (kind == "C" ? Math.Abs(fs[i].C - fs[j].C) : 1)); + } + } + var info = new ClassInfo + { + AdjacentTotal = total, + AdjacentConcordance = total > 0 ? (double)concordant / total : (double?)null, + StepC = steps["C"].Count >= 3 ? Median(steps["C"]) : (double?)null, + StepD = steps["D"].Count >= 3 ? Median(steps["D"]) : (double?)null, + StepO = steps["O"].Count >= 3 ? Median(steps["O"]) : (double?)null, + }; + if (total >= p.ClassGateMinPairs && (double)concordant / total < p.ClassGateConcordance) { + info.GateOk = false; + } + if (info.StepC is double stepC && stepC <= 0) { + info.GateOk = false; + } + for (int k = 0; k < n; k++) { + info.Neighbours[k] = new List(); + } + foreach (var (i, j) in adjacent) { + var kind = AdjacentKind(fs[i].C - fs[j].C, fs[i].D - fs[j].D, fs[i].O - fs[j].O); + var g = fs[i].Rt - fs[j].Rt; // expected > 0: i is predicted later + var step = info.Step(kind) is double st && st > 0 ? st : 0d; + var threshold = Math.Max(tie, p.StepFraction * step * (kind == "C" ? Math.Abs(fs[i].C - fs[j].C) : 1)); + var fi = fs[i]; + var fj = fs[j]; + string labelI, labelJ; + switch (kind) { + case "C": labelI = "C<"; labelJ = "C>"; break; + case "D": labelI = "D>"; labelJ = "D<"; break; + default: labelI = "O>"; labelJ = "O<"; break; + } + var usableI = IsUsable(fi); + var usableJ = IsUsable(fj); + if (usableJ) { + info.Neighbours[i].Add(j); + } + if (usableI) { + info.Neighbours[j].Add(i); + } + if (g > threshold) { + if (usableJ) { + fi.Directions.Add(labelI); + fi.Support.Add((labelI, j, g)); + } + if (usableI) { + fj.Directions.Add(labelJ); + fj.Support.Add((labelJ, i, -g)); + } + } + else if (g >= -eps) { + fi.Ties++; + fj.Ties++; + fi.TiePartners.Add(j); + fj.TiePartners.Add(i); + } + } + return info; + } + + private static bool IsUsable(Feature f) { + return (f.Status == LipidRtFeatureStatus.Core || f.Status == LipidRtFeatureStatus.Free) && !f.Hard && !f.Artefact; + } + + private static void Magnitude(List fs, ClassInfo info, double? kdRun) { + var sC = info.StepC; + var sD = info.StepD; + var sO = info.StepO; + if (sD is null && sC is double c && c > 0 && kdRun is double kd) { + sD = kd * c; + } + if (sO is null && sD.HasValue) { + sO = sD; + } + for (int k = 0; k < fs.Count; k++) { + var f = fs[k]; + var predictions = new List(); + if (info.Neighbours.TryGetValue(k, out var neighbours)) { + foreach (var j in neighbours) { + var g = fs[j]; + int dc = f.C - g.C, dd = f.D - g.D, dO = f.O - g.O; + if ((dc != 0 && sC is null) || (dd != 0 && sD is null) || (dO != 0 && sO is null)) { + continue; + } + predictions.Add(g.Rt + dc * (sC ?? 0d) - dd * (sD ?? 0d) - dO * (sO ?? 0d)); + } + } + if (predictions.Count > 0 && sD is double delta && delta > 0) { + f.Predicted = Median(predictions); + f.Residual = f.Rt - f.Predicted; + f.Delta = delta; + } + } + } + + // ------------------------------------------------------------------------------------------------- decision + + private bool SupportOk(int axes, int directions, bool lower, bool upper) { + int minAxes, minDirections; + switch (_parameter.OperatingPoint) { + case LipidRtOperatingPoint.Balanced: + minAxes = 1; minDirections = 1; + break; + default: + minAxes = 2; minDirections = 2; + break; + } + var axesOk = axes >= minAxes && directions >= minDirections; + var bracketOk = lower && upper && directions >= minDirections; + switch (_parameter.OperatingPoint) { + case LipidRtOperatingPoint.Balanced: + return axesOk; + case LipidRtOperatingPoint.Mid: + return axesOk || bracketOk; + default: // Strict and HighPrecision + return axesOk && bracketOk; + } + } + + private void DecideFeature(Feature f, List fs, ClassInfo info, Dictionary plausibility, + Dictionary> artefacts, double tailStart, Dictionary results) { + var p = _parameter; + var axes = f.Directions.Select(d => d[0]).Distinct().Count(); + var nd = f.Directions.Count; + double? rel = null; + if (f.Residual.HasValue && f.Delta is double delta && delta != 0) { + rel = f.Residual.Value / delta; + } + var magnitudeOk = !(p.MagnitudeK > 0 && rel.HasValue && Math.Abs(rel.Value) > p.MagnitudeK); + + LipidRtOutcome outcome; + LipidRtReason reason; + double score; + Feature strongTie = null; + if (f.Status == LipidRtFeatureStatus.Excluded) { + score = -1d - f.Deficit; + if (f.Deficit >= p.InconsistentMinDeficit) { + if (IsOxidisedClass(f.ClassKey) && !p.FlagOxidisedClasses) { + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.OxNoFlag); + } + else if (f.Rt > tailStart && !p.FlagGradientTail) { + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.TailNoFlag); + } + else { + (outcome, reason) = (LipidRtOutcome.Inconsistent, LipidRtReason.InconsistentOrder); + } + } + else { + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.SmallDeficit); + } + } + else if (f.Status == LipidRtFeatureStatus.Ambiguous) { + score = 0d; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Ambiguous); + } + else if (!info.GateOk) { + score = -0.5; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.ClassGate); + } + else { + score = nd + 0.5 * axes - 0.25 * f.Ties + (f.Status == LipidRtFeatureStatus.Core ? 0.1 * Math.Min(f.E, 5d) : 0.5); + if (rel.HasValue) { + score -= 0.5 * Math.Min(Math.Abs(rel.Value), 2d); + } + var lower = f.Directions.Contains("C<") || f.Directions.Contains("D>") || f.Directions.Contains("O>"); + var upper = f.Directions.Contains("C>") || f.Directions.Contains("D<") || f.Directions.Contains("O<"); + foreach (var j in f.TiePartners) { + if (fs[j].Status != LipidRtFeatureStatus.Excluded && fs[j].W >= p.TieRatio * f.W) { + strongTie = fs[j]; + break; + } + } + if (!SupportOk(axes, nd, lower, upper)) { + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Support); + } + else if (p.TieBlock && strongTie != null) { + score -= p.BlockPenalty; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Tie); + } + else if (f.Artefact) { + score -= p.BlockPenalty; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Artefact); + } + else if (!magnitudeOk) { + score -= p.BlockPenalty; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Magnitude); + } + else if (f.Status == LipidRtFeatureStatus.Core && f.E < p.CoreMargin) { + score -= p.BlockPenalty; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Margin); + } + else if (f.SecondaryOnly && p.SecondaryOnlyPolicy == LipidRtSecondaryOnlyPolicy.Block) { + score -= p.PolicyPenalty; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.PolicySecondary); + } + else if (f.MinorDuplicateOf != null && p.MinorDuplicatePolicy == LipidRtMinorDuplicatePolicy.Block) { + score -= p.PolicyPenalty; + (outcome, reason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.MinorDuplicate); + } + else { + (outcome, reason) = (LipidRtOutcome.Confirmed, LipidRtReason.ConfirmedRt); + } + } + + var max = Math.Max(p.MaxEvidenceItems, 0); + var support = SelectSupport(f.Support, max) + .Select(s => new LipidRtNeighbour(s.Label, fs[s.Index].Format(), fs[s.Index].Rt, s.Gap)) + .ToList(); + var discordantPartners = f.Discordant.Take(max) + .Select(j => new LipidRtDiscordantPartner(fs[j].Format(), fs[j].Rt, fs[j].Status)) + .ToList(); + var memberIds = f.Members.Select(m => m.Id).ToList(); + LipidRtArtefact featureArtefact = f.Artefact && artefacts.TryGetValue(f.Best.Id, out var bestArtefacts) ? bestArtefacts[0] : null; + + foreach (var s in f.Members) { + var memberOutcome = outcome; + var memberReason = reason; + var penalty = 0d; + var softOverridden = false; + if (memberOutcome == LipidRtOutcome.Confirmed) { + if (p.ClassLevelPolicy == LipidRtClassLevelPolicy.Never && s.TagLevel == LipidRtTagLevel.Class) { + penalty += p.PolicyPenalty; + (memberOutcome, memberReason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.PolicyClass); + } + else if (artefacts.ContainsKey(s.Id)) { + penalty += p.BlockPenalty; + (memberOutcome, memberReason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.Artefact); + } + } + var verdict = LipidRtPlausibilityVerdict.NotEvaluated; + if (plausibility.Count > 0) { + verdict = plausibility[s.Id]; + if (verdict.IsImplausible) { + if (verdict.Severity == LipidRtPlausibilitySeverity.Hard) { + penalty += 2 * p.PolicyPenalty; + (memberOutcome, memberReason) = (LipidRtOutcome.Inconsistent, LipidRtReason.InconsistentPlausibility); + } + else if (memberOutcome == LipidRtOutcome.Confirmed) { + if (!(axes >= p.SoftOverrideAxes && Dot(s.WeightedDot) >= p.SoftOverrideDot)) { + penalty += p.PolicyPenalty; + (memberOutcome, memberReason) = (LipidRtOutcome.NotEvaluable, LipidRtReason.PlausibilitySoft); + } + else { + softOverridden = true; + } + } + else { + penalty += p.PolicyPenalty; + } + } + } + var evidence = new LipidRtEvidence + { + Feature = f.Format(), + Composition = new LipidRtComposition(f.ClassKey, f.C, f.D, f.O), + FeatureRt = f.Rt, + Status = f.Status, + Margin = f.E, + Deficit = f.Deficit, + Weight = f.W, + Directions = string.Concat(f.Directions), + Axes = axes, + Ties = f.Ties, + Support = support, + Discordant = discordantPartners, + DiscordantCount = f.Discordant.Count, + TiePartner = reason == LipidRtReason.Tie && strongTie != null ? strongTie.Format() : string.Empty, + MinorDuplicateOf = f.MinorDuplicateOf is null ? string.Empty + : string.Format(CultureInfo.InvariantCulture, "{0} {1}@{2:0.00}", f.MinorDuplicateOf.Name, f.MinorDuplicateOf.Molecular, f.MinorDuplicateOf.Rt), + Predicted = f.Predicted, + Residual = f.Residual, + ResidualDoubleBondSteps = rel, + ClassGatePassed = info.GateOk, + ClassAdjacentConcordance = info.AdjacentConcordance, + Members = memberIds, + AdductCount = f.AdductCount, + Artefacts = artefacts.TryGetValue(s.Id, out var own) ? own : (IReadOnlyList)Array.Empty(), + FeatureArtefact = featureArtefact, + Plausibility = verdict, + SoftPlausibilityOverridden = softOverridden, + }; + results[s.Id] = Single(s, memberOutcome, memberReason, score - penalty, evidence); + } + } + + /// + /// The supporting neighbours kept in the evidence record: at most , listed by |gap| + /// descending then direction (SPEC section 8). + /// + /// + /// Which ones are kept is not simply the first of that order. Far neighbours, isomers among + /// them, have the largest gaps and can all lie on one side, so the strongest lower bound and the strongest upper + /// bound are kept first, then the strongest neighbour of each other direction, then the rest by |gap|. A Confirmed + /// rests on a bracket, and its sentence then still shows both sides of it. The reference implementation keeps the + /// first ; only this free text differs, never an outcome or a score. + /// + internal static List<(string Label, int Index, double Gap)> SelectSupport(IReadOnlyList<(string Label, int Index, double Gap)> support, int max) { + var ordered = support + .OrderByDescending(s => Math.Abs(s.Gap)).ThenBy(s => s.Label, StringComparer.Ordinal) + .ToList(); + if (ordered.Count <= max) { + return ordered; + } + var priority = new List(ordered.Count); + var lower = ordered.FindIndex(s => LipidRtNeighbour.IsLowerBoundDirection(s.Label)); + var upper = ordered.FindIndex(s => !LipidRtNeighbour.IsLowerBoundDirection(s.Label)); + foreach (var i in new[] { lower, upper }.Where(i => i >= 0).OrderBy(i => i)) { + priority.Add(i); + } + var directions = new HashSet(priority.Select(i => ordered[i].Label), StringComparer.Ordinal); + for (var i = 0; i < ordered.Count; i++) { + if (directions.Add(ordered[i].Label)) { + priority.Add(i); + } + } + var kept = new HashSet(priority); + for (var i = 0; i < ordered.Count && priority.Count < max; i++) { + if (kept.Add(i)) { + priority.Add(i); + } + } + var selected = new HashSet(priority.Take(max)); + return ordered.Where((s, i) => selected.Contains(i)).ToList(); + } + + // ------------------------------------------------------------------------------------------------- numerics + + /// The median as numpy computes it: the mean of the two middle values for an even count. + internal static double Median(IReadOnlyCollection values) { + if (values.Count == 0) { + throw new ArgumentException("No values.", nameof(values)); + } + var sorted = values.ToArray(); + Array.Sort(sorted); + var mid = sorted.Length / 2; + return sorted.Length % 2 == 1 ? sorted[mid] : (sorted[mid - 1] + sorted[mid]) / 2d; + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCheckHooks.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCheckHooks.cs new file mode 100644 index 000000000..4180bc020 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCheckHooks.cs @@ -0,0 +1,122 @@ +using CompMs.MsdialCore.DataObj; +using CompMs.MsdialCore.Parameter; +using System; +using System.Collections.Generic; +using System.Globalization; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// What the check did on one unit: its result, what happened to the tags, or why it did not run. + public sealed class LipidRtUnitRun + { + internal LipidRtUnitRun(LipidRtRunResult result, LipidRtTagSummary tags, string parameterHash, Exception error) { + Result = result; + Tags = tags; + ParameterHash = parameterHash ?? string.Empty; + Error = error; + } + + /// The unit's result; null when the check did not run (switched off, or it failed). + public LipidRtRunResult Result { get; } + + /// What the protocol did to the tags, including the withdrawal when the check did not run. + public LipidRtTagSummary Tags { get; } + + public string ParameterHash { get; } + + /// The exception the check stopped with; its verdicts and tags were then withdrawn, not half-written. + public Exception Error { get; } + + public bool Ran => Result != null; + + /// One line for a log: the outcome counts over eligible objects, the tags, and the gate when it failed. + public string Describe() { + if (Error != null) { + return $"failed ({Error.GetType().Name}: {Error.Message}); no verdicts or Confirmed tags were written by it, and its earlier ones were withdrawn ({Tags.Withdrawn})."; + } + if (Result is null) { + return Tags.Withdrawn > 0 + ? $"not run; {Tags.Withdrawn} Confirmed tags it had set earlier were withdrawn." + : "not run."; + } + var eligible = Result.Results.Where(r => r.Reason != LipidRtReason.Ineligible).ToList(); + var text = string.Format(CultureInfo.InvariantCulture, + "{0} annotated, {1} evaluable: {2} Confirmed, {3} NotEvaluable, {4} Inconsistent; Confirmed tags set {5} (already tagged {6}, Misannotation {7}, suppressed {8}, earlier ones withdrawn {9}); rules {10}, parameter {11}.", + Result.Results.Count, eligible.Count, + eligible.Count(r => r.Outcome == LipidRtOutcome.Confirmed), + eligible.Count(r => r.Outcome == LipidRtOutcome.NotEvaluable), + eligible.Count(r => r.Outcome == LipidRtOutcome.Inconsistent), + Tags.Set, Tags.AlreadySelected, Tags.Misannotation, Tags.Suppressed, Tags.Withdrawn, + LipidRtCheckParameter.RuleVersion, ParameterHash); + if (Result.Gate != null && !Result.Gate.Passed) { + text += " Reversed-phase gate not passed: " + Result.Gate + "."; + } + return text; + } + } + + /// + /// Where processing runs the lipid RT check (SPEC section 11, addendum B3): hook A on the peaks of one analysis file + /// between annotation and saving, and hook B on the spots of an alignment once they are packed and their annotations + /// are final (after isotope tracking, which clears some), before they are saved. Both are shared by the Console and + /// the GUI, and both apply the verdicts through . + /// + /// + /// + /// The check runs only when says so: a lipidomics project with the + /// check enabled (the default). Whether the chromatography is reversed phase is decided inside the check by its + /// run-level gate. v1 runs in LC-MS processing only; the callers decide that. + /// + /// + /// When the check does not run, the unit still goes through the protocol with no result, so a Confirmed tag the + /// check set on an earlier run, and its verdicts, do not outlive a re-annotation they no longer describe. Tags a + /// person or another tool set are never touched. + /// + /// + /// Each call reads only its own unit and the parameter, writes only that unit's objects, and keeps no state, so + /// hook A can run inside the per-file parallel loop. + /// + /// + public static class LipidRtCheckHooks + { + /// Hook A: the peaks of one analysis file, after annotation and before they are saved. + public static LipidRtUnitRun RunOnPeaks(IReadOnlyList peaks, ParameterBase parameter) { + if (peaks is null) throw new ArgumentNullException(nameof(peaks)); + if (parameter is null) throw new ArgumentNullException(nameof(parameter)); + return Run(parameter, + () => LipidRtInputBuilder.FromPeaks(peaks, parameter.IonMode), + (result, hash) => LipidRtTagProtocol.ApplyToPeaks(peaks, result, hash)); + } + + /// + /// Hook B: the spots of one alignment, after they are packed and every step that changes their annotations has run, + /// and before the result is saved or exported. + /// + public static LipidRtUnitRun RunOnSpots(IReadOnlyList spots, ParameterBase parameter) { + if (spots is null) throw new ArgumentNullException(nameof(spots)); + if (parameter is null) throw new ArgumentNullException(nameof(parameter)); + return Run(parameter, + () => LipidRtInputBuilder.FromSpots(spots, parameter.IonMode), + (result, hash) => LipidRtTagProtocol.ApplyToSpots(spots, result, hash)); + } + + private static LipidRtUnitRun Run(ParameterBase parameter, Func input, Func apply) { + var checkParameter = parameter.LipidRtCheckParam; + var hash = checkParameter.ComputeHash(); + if (!checkParameter.RunsFor(parameter.TargetOmics)) { + return new LipidRtUnitRun(null, apply(null, hash), hash, null); + } + LipidRtRunResult result; + try { + result = new LipidRtCheck(checkParameter).Evaluate(input()); + } + catch (Exception ex) when (!(ex is OperationCanceledException)) { + // An optional check must not take the processing down with it, and must not leave half a result behind: + // the unit is withdrawn as if the check were off, and the caller reports the error. + return new LipidRtUnitRun(null, apply(null, hash), hash, ex); + } + return new LipidRtUnitRun(result, apply(result, hash), hash, null); + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCorroboration.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCorroboration.cs new file mode 100644 index 000000000..f81390b89 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtCorroboration.cs @@ -0,0 +1,211 @@ +using CompMs.Common.Enum; +using CompMs.MsdialCore.Parameter; +using System; +using System.Collections.Generic; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// The two units after cross-polarity corroboration, and the RT line that mapped them. + public sealed class LipidRtCrossPolarityResult + { + internal LipidRtCrossPolarityResult(LipidRtRunResult positive, LipidRtRunResult negative, LipidRtCrossPolarityLine line) { + Positive = positive; + Negative = negative; + Line = line; + } + + public LipidRtRunResult Positive { get; } + public LipidRtRunResult Negative { get; } + public LipidRtCrossPolarityLine Line { get; } + } + + /// + /// The corroboration post-step (SPEC step 11): the second adduct, the other polarity, and the policies that depend + /// on them. A pure function of the units' single-run results; it reads no files and never chains (a promotion is + /// never the basis of another promotion). + /// + /// + /// (a) A Confirmed(Rt) member whose feature holds at least two distinct adducts is adduct-corroborated. + /// (b) With one positive and one negative unit of the same sample set and LC: RtNeg = a + b RtPos is fitted by + /// Theil-Sen over compositions with a single RT cluster in both units; the twins of a candidate are the eligible + /// candidates of the other unit with the same composition within the tolerance of the mapped RT. A Confirmed(Rt) + /// with a Confirmed(Rt) twin is polarity-corroborated (score + 1); a NotEvaluable for a promotable reason (Support, + /// ClassGate, Ambiguous) with a Confirmed(Rt) twin becomes Confirmed(PromotedByPolarity) (score + 0.5). Inconsistent, + /// tie, artefact, magnitude, margin, deficit, policy, plausibility and gate outcomes are never promoted. A + /// Confirmed(Rt) with no Confirmed twin and an Inconsistent twin records a polarity conflict. + /// (c) ClassLevelPolicy Corroborated: a class-level Confirmed without (a) or (b) becomes NotEvaluable(PolicyClass). + /// The Strict operating point: a Confirmed(Rt) without (a) or (b) becomes NotEvaluable(Uncorroborated). + /// + public static class LipidRtCorroboration + { + /// + /// Applies the corroboration step to the given units. Cross polarity runs when the parameter allows it and the + /// units are exactly one positive and one negative unit; otherwise each unit gets (a) and (c) only. + /// + /// New results, in the order of . + public static IReadOnlyList Apply(IReadOnlyList units, LipidRtCheckParameter parameter) { + return Apply(units, parameter, out _); + } + + /// + /// Cross-polarity corroboration of a positive and a negative unit of the same sample set and LC, from their + /// single-run results. The returned units replace the inputs; the inputs are not changed. + /// + public static LipidRtCrossPolarityResult CrossPolarity(LipidRtRunResult positive, LipidRtRunResult negative, LipidRtCheckParameter parameter) { + if (positive is null) throw new ArgumentNullException(nameof(positive)); + if (negative is null) throw new ArgumentNullException(nameof(negative)); + if (positive.Polarity != IonMode.Positive || negative.Polarity != IonMode.Negative) { + throw new ArgumentException("Cross polarity needs one positive and one negative unit."); + } + var results = Apply(new[] { positive, negative }, parameter, out var line); + return new LipidRtCrossPolarityResult(results[0], results[1], line); + } + + private static IReadOnlyList Apply(IReadOnlyList units, LipidRtCheckParameter parameter, out LipidRtCrossPolarityLine line) { + if (units is null) throw new ArgumentNullException(nameof(units)); + if (parameter is null) throw new ArgumentNullException(nameof(parameter)); + line = null; + var copies = units.Select(u => u.Results.Select(r => r.CopyForCorroboration()).ToList()).ToList(); + + // (a) adduct corroboration + foreach (var results in copies) { + foreach (var r in results) { + if (r.SingleRunReason == LipidRtReason.ConfirmedRt) { + r.CorroboratedByAdduct = r.Evidence.AdductCount >= 2; + } + } + } + + // (b) cross polarity + var positive = Enumerable.Range(0, units.Count).Where(i => units[i].Polarity == IonMode.Positive).ToList(); + var negative = Enumerable.Range(0, units.Count).Where(i => units[i].Polarity == IonMode.Negative).ToList(); + if (parameter.CrossPolarity && positive.Count == 1 && negative.Count == 1) { + int a = positive[0], b = negative[0]; + line = FitLine(units[a], units[b], parameter); + var promotable = new HashSet(parameter.GetPromotableReasons()); + Corroborate(units[a], copies[a], units[b], line, forward: true, promotable, parameter); + Corroborate(units[b], copies[b], units[a], line, forward: false, promotable, parameter); + } + + // (c) policies that depend on corroboration + foreach (var results in copies) { + foreach (var r in results) { + var corroborated = r.CorroboratedByAdduct || r.CorroboratedByPolarity.Count > 0; + if (parameter.ClassLevelPolicy == LipidRtClassLevelPolicy.Corroborated + && (r.Reason == LipidRtReason.ConfirmedRt || r.Reason == LipidRtReason.ConfirmedPromotedByPolarity) + && r.Candidate.TagLevel == LipidRtTagLevel.Class && !corroborated) { + r.Outcome = LipidRtOutcome.NotEvaluable; + r.Reason = LipidRtReason.PolicyClass; + } + } + } + if (parameter.OperatingPoint == LipidRtOperatingPoint.Strict) { + foreach (var results in copies) { + foreach (var r in results) { + if (r.Reason == LipidRtReason.ConfirmedRt && !r.CorroboratedByAdduct && r.CorroboratedByPolarity.Count == 0) { + r.Outcome = LipidRtOutcome.NotEvaluable; + r.Reason = LipidRtReason.Uncorroborated; + } + } + } + } + return units.Select((u, i) => new LipidRtRunResult(u.Input, copies[i], u.Gate)).ToList(); + } + + private static void Corroborate(LipidRtRunResult unit, List results, LipidRtRunResult other, LipidRtCrossPolarityLine line, + bool forward, HashSet promotable, LipidRtCheckParameter parameter) { + var index = EligibleByComposition(other); + foreach (var r in results) { + var s = r.Candidate; + if (!IsEligible(s)) { + continue; + } + var mapped = forward ? line.Intercept + line.Slope * s.Rt : (s.Rt - line.Intercept) / line.Slope; + var twins = index.TryGetValue(s.Composition, out var list) + ? list.Where(t => Math.Abs(t.Rt - mapped) <= line.Tolerance).ToList() + : new List(); + r.Twins = twins.Select(t => new LipidRtTwin(t.Id, t.Adduct, t.Rt, other[t.Id].SingleRunReason)).ToList(); + var confirmedTwins = twins.Where(t => other[t.Id].SingleRunReason == LipidRtReason.ConfirmedRt).ToList(); + if (r.SingleRunReason == LipidRtReason.ConfirmedRt) { + if (confirmedTwins.Count > 0) { + r.CorroboratedByPolarity = confirmedTwins.Select(t => t.Id).ToList(); + r.Score = r.Score.HasValue ? r.Score + 1d : null; + } + else if (twins.Count > 0 && twins.Any(t => other[t.Id].SingleRunReason == LipidRtReason.InconsistentOrder)) { + r.PolarityConflict = true; + if (parameter.PolarityConflict == LipidRtPolarityConflictPolicy.Demote) { + r.Outcome = LipidRtOutcome.NotEvaluable; + r.Reason = LipidRtReason.PolarityConflict; + } + } + } + else if (promotable.Contains((int)r.SingleRunReason) && confirmedTwins.Count > 0) { + r.Outcome = LipidRtOutcome.Confirmed; + r.Reason = LipidRtReason.ConfirmedPromotedByPolarity; + r.Score = (r.Score ?? 0d) + 0.5; + r.CorroboratedByPolarity = confirmedTwins.Select(t => t.Id).ToList(); + } + } + } + + private static bool IsEligible(LipidRtCandidate c) => c.IsEligible && c.Composition != null; + + private static Dictionary> EligibleByComposition(LipidRtRunResult unit) { + var index = new Dictionary>(); + foreach (var c in unit.Input.Candidates) { + if (!IsEligible(c)) { + continue; + } + var key = c.Composition; + if (!index.TryGetValue(key, out var list)) { + index[key] = list = new List(); + } + list.Add(c); + } + return index; + } + + /// Compositions with a single RT cluster (spread at most the collapse window) and their median RT. + private static Dictionary UniqueRt(LipidRtRunResult unit, double window) { + var result = new Dictionary(); + foreach (var pair in EligibleByComposition(unit)) { + var rts = pair.Value.Select(c => c.Rt).ToList(); + if (rts.Max() - rts.Min() <= window) { + result[pair.Key] = LipidRtCheck.Median(rts); + } + } + return result; + } + + /// RtNeg ~ a + b RtPos by Theil-Sen over compositions with a single RT cluster in both units. + private static LipidRtCrossPolarityLine FitLine(LipidRtRunResult positive, LipidRtRunResult negative, LipidRtCheckParameter parameter) { + var spanA = positive.Span; + var spanB = negative.Span; + var ua = UniqueRt(positive, parameter.CollapseFraction * spanA); + var ub = UniqueRt(negative, parameter.CollapseFraction * spanB); + var keys = ua.Keys.Where(ub.ContainsKey).OrderBy(k => k).ToList(); + var x = keys.Select(k => ua[k]).ToArray(); + var y = keys.Select(k => ub[k]).ToArray(); + var floor = parameter.CrossPolarityToleranceFraction * Math.Max(spanA, spanB); + if (keys.Count < 5) { + return new LipidRtCrossPolarityLine(0d, 1d, null, keys.Count, floor); + } + var slopes = new List(); + for (int i = 0; i < x.Length; i++) { + for (int j = i + 1; j < x.Length; j++) { + if (Math.Abs(x[j] - x[i]) > 0.025 * spanA) { + slopes.Add((y[j] - y[i]) / (x[j] - x[i])); + } + } + } + var b = slopes.Count >= 10 ? LipidRtCheck.Median(slopes) : 1d; + var a = LipidRtCheck.Median(y.Select((yi, k) => yi - b * x[k]).ToList()); + var residuals = y.Select((yi, k) => yi - (a + b * x[k])).ToList(); + var center = LipidRtCheck.Median(residuals); + var sd = 1.4826 * LipidRtCheck.Median(residuals.Select(r => Math.Abs(r - center)).ToList()); + var tolerance = Math.Max(floor, sd != 0 ? parameter.CrossPolarityToleranceSd * sd : 0d); + return new LipidRtCrossPolarityLine(a, b, sd, keys.Count, tolerance); + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtEnums.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtEnums.cs new file mode 100644 index 000000000..394425734 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtEnums.cs @@ -0,0 +1,123 @@ +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// What the lipid retention-time order check concluded about one object. + public enum LipidRtOutcome : byte + { + NotEvaluated = 0, + Confirmed = 1, + NotEvaluable = 2, + Inconsistent = 3, + } + + /// + /// Why the check reached its outcome. The numbers are part of the persisted record and must not change. + /// + public enum LipidRtReason : byte + { + None = 0, + /// Bracketed by the run's other MS/MS-backed annotations of the class, in RT order. + ConfirmedRt = 1, + /// Not decidable in this run, but a twin in the other polarity was Confirmed by its own evidence. + ConfirmedPromotedByPolarity = 2, + Ineligible = 10, + /// The run does not show reversed-phase order, or the analyst said the chromatography is not RP. + Gate = 11, + /// The class's adjacent compositions are not ordered like reversed phase. + ClassGate = 12, + /// In some but not all maximum order-consistent subsets. + Ambiguous = 13, + /// Not bracketed by adjacent compositions. + Support = 14, + /// Co-elutes with an adjacent composition that weighs at least as much. + Tie = 15, + /// Explained as an MS1 artefact of a brighter co-eluting ion. + Artefact = 16, + /// The RT is closer to a neighbouring composition than the class steps predict. + Magnitude = 17, + /// Wins its order conflicts only narrowly. + Margin = 18, + /// In no maximum consistent subset, but by less than the Inconsistent threshold. + SmallDeficit = 19, + /// An order conflict in an oxidised class, which is not flagged. + OxNoFlag = 20, + /// An order conflict in the gradient-end tail, which is not flagged. + TailNoFlag = 21, + /// A class-level match withheld by ClassLevelPolicy. + PolicyClass = 22, + /// A feature of secondary ions only, withheld by SecondaryOnlyPolicy. + PolicySecondary = 23, + /// A minor duplicate, withheld by MinorDuplicatePolicy. + MinorDuplicate = 24, + /// A soft plausibility rule (B2/B3) failed and the RT/spectrum override did not apply. + PlausibilitySoft = 25, + /// Neither adduct- nor polarity-corroborated at the Strict operating point. + Uncorroborated = 26, + /// Demoted because every other-polarity twin is Inconsistent (PolarityConflict = Demote). + PolarityConflict = 27, + /// Excluded from every maximum order-consistent subset by a clear deficit. A review item, never a Misannotation tag. + InconsistentOrder = 40, + /// A hard plausibility rule (B1) failed, whatever the RT says. + InconsistentPlausibility = 41, + } + + /// The structural level of the representative MS/MS match. + public enum LipidRtTagLevel : byte + { + None = 0, + /// Lipid class resolved (annotation tag 320/420). + Class = 1, + /// An MS/MS spectrum match without lipid levels (330/430). + Msms = 2, + /// Lipid chains or positions resolved (310/410, 300/400). + Chain = 3, + } + + /// Membership of a feature in the maximum order-consistent subsets of its class. + public enum LipidRtFeatureStatus : byte + { + /// No discordant partner. + Free = 0, + /// In every maximum consistent subset. + Core = 1, + /// In some but not all. + Ambiguous = 2, + /// In none. + Excluded = 3, + } + + /// The role of one chain in a lipid name. + public enum LipidRtChainRole : byte + { + Acyl, + Alkyl, + SphingoidBase, + NAcyl, + /// The ester-linked acyl of EOS/ASM/TG_EST/NAGly, written "(FA x:y)". + EsterAcyl, + /// The O-acyl on a hexose, written "(O-x:y)". + OAcyl, + SterolCore, + /// A species-level total of several fatty chains: unresolved. + SpeciesCore, + /// A species-level total of a sphingoid base and its N-acyl: unresolved. + SphingoidCore, + } + + /// The overall verdict of a plausibility profile. + public enum LipidRtPlausibility : byte + { + /// No profile, or not evaluated. + None = 0, + Plausible = 1, + Implausible = 2, + /// A rule applies but the chain it concerns is not resolved. + Untestable = 3, + } + + public enum LipidRtPlausibilitySeverity : byte + { + None = 0, + Soft = 1, + Hard = 2, + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtEvidenceSentence.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtEvidenceSentence.cs new file mode 100644 index 000000000..b38ca3184 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtEvidenceSentence.cs @@ -0,0 +1,161 @@ +using System.Collections.Generic; +using System.Globalization; +using System.Linq; +using System.Text; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// + /// The sentence a person reads for one outcome, e.g. "Confirmed: TG 52:1 elutes after TG 50:1 (+0.30 min) and + /// TG 52:2 (+0.25), and before TG 54:1 (-0.30) and TG 52:0 (-0.25); residual +0.00 min (0.00 DB steps)." + /// + /// + /// Free text for the GUI and the review list; nothing parses it. It says what a Confirmed does not certify only in + /// the documentation, to keep the sentence short. + /// + public static class LipidRtEvidenceSentence + { + public static string Format(LipidRtObjectResult result) { + if (result is null) { + return string.Empty; + } + var c = result.Candidate; + var e = result.Evidence; + var subject = e.Composition?.ToString() ?? c?.Composition?.ToString() ?? c?.Name ?? string.Empty; + switch (result.Reason) { + case LipidRtReason.ConfirmedRt: + return "Confirmed: " + OrderClause(subject, e) + Residual(e) + Corroboration(result) + "."; + case LipidRtReason.ConfirmedPromotedByPolarity: + return "Confirmed in the other polarity (" + Twins(result) + "); in this run: " + Describe(result.SingleRunReason, e) + "."; + case LipidRtReason.InconsistentOrder: + return string.Format(CultureInfo.InvariantCulture, + "Inconsistent: {0} is out of order in its class (deficit {1:0.00}); discordant with {2}{3}.", + subject, e.Deficit, Discordant(e), Residual(e)); + case LipidRtReason.InconsistentPlausibility: + return string.Format(CultureInfo.InvariantCulture, "Inconsistent: plausibility {0} ({1}).", e.Plausibility.Rule, e.Plausibility.Detail); + case LipidRtReason.Ineligible: + return "Not evaluated: " + (string.IsNullOrEmpty(c?.IneligibleReason) ? "not eligible" : c.IneligibleReason.Replace('_', ' ')) + "."; + case LipidRtReason.Gate: + return "Not evaluable: the run does not show reversed-phase order (" + (e.Gate?.ToString() ?? string.Empty) + ")."; + case LipidRtReason.PolicyClass: + return "Not evaluable: a class-level match is not confirmed without corroboration (policy); " + OrderClause(subject, e) + "."; + case LipidRtReason.Uncorroborated: + return "Not evaluable: not corroborated by a second adduct or the other polarity (Strict); " + OrderClause(subject, e) + "."; + case LipidRtReason.PolarityConflict: + return "Not evaluable: the other polarity contradicts (" + Twins(result) + "); " + OrderClause(subject, e) + "."; + default: + return "Not evaluable: " + Describe(result.Reason, e) + "."; + } + } + + private static string Describe(LipidRtReason reason, LipidRtEvidence e) { + switch (reason) { + case LipidRtReason.Support: + return "insufficient adjacent support (" + Directions(e) + ")"; + case LipidRtReason.ClassGate: + return string.Format(CultureInfo.InvariantCulture, "the class order is not reversed-phase like (adjacent concordance {0})", + e.ClassAdjacentConcordance.HasValue ? e.ClassAdjacentConcordance.Value.ToString("0.00", CultureInfo.InvariantCulture) : "-"); + case LipidRtReason.Ambiguous: + return "an order conflict that the evidence does not resolve; discordant with " + Discordant(e); + case LipidRtReason.Tie: + return "co-elutes with the adjacent composition " + e.TiePartner + " (isotope, adduct or in-source artefact?)"; + case LipidRtReason.Artefact: + var artefact = e.Artefacts.FirstOrDefault() ?? e.FeatureArtefact; + return "MS1 artefact: " + (artefact?.Text ?? "explained by a brighter co-eluting ion"); + case LipidRtReason.Magnitude: + return "the RT is closer to a neighbouring composition than the class steps predict" + Residual(e); + case LipidRtReason.Margin: + return string.Format(CultureInfo.InvariantCulture, "wins its order conflicts only narrowly (margin {0:0.00}, {1} discordant)", e.Margin, e.DiscordantCount); + case LipidRtReason.SmallDeficit: + return string.Format(CultureInfo.InvariantCulture, "out of order by a small deficit ({0:0.00}); discordant with {1}", e.Deficit, Discordant(e)); + case LipidRtReason.OxNoFlag: + return string.Format(CultureInfo.InvariantCulture, "an order conflict in an oxidised class, not flagged (deficit {0:0.00})", e.Deficit); + case LipidRtReason.TailNoFlag: + return string.Format(CultureInfo.InvariantCulture, "an order conflict in the gradient-end tail, not flagged (deficit {0:0.00})", e.Deficit); + case LipidRtReason.PolicySecondary: + return "secondary ions only (policy)"; + case LipidRtReason.MinorDuplicate: + return "a minor duplicate of " + e.MinorDuplicateOf + " (policy)"; + case LipidRtReason.PlausibilitySoft: + return string.Format(CultureInfo.InvariantCulture, "plausibility {0} soft ({1})", e.Plausibility.Rule, e.Plausibility.Detail); + case LipidRtReason.ConfirmedRt: + return "Confirmed"; + default: + return reason.ToString(); + } + } + + private static string OrderClause(string subject, LipidRtEvidence e) { + var lower = e.Support.Where(n => n.IsLowerBound).ToList(); + var upper = e.Support.Where(n => !n.IsLowerBound).ToList(); + if (lower.Count == 0 && upper.Count == 0) { + return subject + " has no supporting neighbour"; + } + var builder = new StringBuilder(subject); + if (lower.Count > 0) { + builder.Append(" elutes after ").Append(Join(lower)); + } + if (upper.Count > 0) { + builder.Append(lower.Count > 0 ? ", and before " : " elutes before ").Append(Join(upper)); + } + return builder.ToString(); + } + + private static string Join(IReadOnlyList neighbours) { + var parts = neighbours.Select((n, i) => string.Format(CultureInfo.InvariantCulture, "{0} ({1:+0.00;-0.00;+0.00}{2})", + Strip(n.Feature), n.Gap, i == 0 ? " min" : string.Empty)).ToList(); + if (parts.Count == 1) { + return parts[0]; + } + return string.Join(", ", parts.Take(parts.Count - 1)) + " and " + parts[parts.Count - 1]; + } + + /// "TG 50:1@12.04" -> "TG 50:1". + private static string Strip(string feature) { + var at = feature.LastIndexOf('@'); + return at > 0 ? feature.Substring(0, at) : feature; + } + + private static string Residual(LipidRtEvidence e) { + if (!e.Residual.HasValue || !e.ResidualDoubleBondSteps.HasValue) { + return string.Empty; + } + return string.Format(CultureInfo.InvariantCulture, "; residual {0:+0.00;-0.00;+0.00} min ({1:+0.00;-0.00;+0.00} DB steps)", e.Residual.Value, e.ResidualDoubleBondSteps.Value); + } + + private static string Discordant(LipidRtEvidence e) { + if (e.Discordant.Count == 0) { + return "-"; + } + return string.Join(", ", e.Discordant.Select(d => d.Feature + " [" + d.Status + "]")); + } + + private static string Directions(LipidRtEvidence e) { + return string.Format(CultureInfo.InvariantCulture, "directions {0}, {1} axes, {2} ties", string.IsNullOrEmpty(e.Directions) ? "-" : e.Directions, e.Axes, e.Ties); + } + + private static string Twins(LipidRtObjectResult result) { + if (result.Twins.Count == 0) { + return "no twin"; + } + return string.Join(", ", result.Twins.Select(t => string.Format(CultureInfo.InvariantCulture, "{0} {1}@{2:0.00} {3}", t.Id, t.Adduct, t.Rt, t.SingleRunReason))); + } + + private static string Corroboration(LipidRtObjectResult result) { + var parts = new List(); + if (result.CorroboratedByAdduct) { + parts.Add("second adduct"); + } + if (result.CorroboratedByPolarity.Count > 0) { + parts.Add("other polarity"); + } + if (result.PolarityConflict) { + parts.Add("conflicting other-polarity twin"); + } + if (result.Evidence.SoftPlausibilityOverridden) { + parts.Add("plausibility " + result.Evidence.Plausibility.Rule + " soft overridden by RT and spectrum"); + } + return parts.Count == 0 ? string.Empty : "; " + string.Join("; ", parts); + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtInputBuilder.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtInputBuilder.cs new file mode 100644 index 000000000..c53a63bd8 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtInputBuilder.cs @@ -0,0 +1,151 @@ +using CompMs.Common.DataObj.Result; +using CompMs.Common.Enum; +using CompMs.Common.Utility; +using CompMs.MsdialCore.DataObj; +using System; +using System.Collections.Generic; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// + /// Builds the check's input from MS-DIAL's own objects: the peaks of one analysis file (hook A) or the spots of one + /// alignment (hook B). Reads only; writes nothing back. + /// + /// + /// Drift children of LC-IM-MS objects are not evaluated in v1: pass the top-level list and the parents are read. + /// + public static class LipidRtInputBuilder + { + /// The peaks of one analysis file, with RT on that file's own axis. + public static LipidRtRunInput FromPeaks(IReadOnlyList peaks, IonMode polarity) { + if (peaks is null) throw new ArgumentNullException(nameof(peaks)); + return Build(peaks.Select(p => (p.MasterPeakID, p.PeakFeature.ChromXsTop?.RT?.Value ?? double.NaN, p.PeakFeature.Mass, p.PeakFeature.PeakHeightTop, + p.AdductType?.AdductIonName, p.MatchResults?.Representative, p.Ontology, p.Name)), polarity); + } + + /// The spots of one alignment, with RT on the aligned axis and the average height. + public static LipidRtRunInput FromSpots(IReadOnlyList spots, IonMode polarity) { + if (spots is null) throw new ArgumentNullException(nameof(spots)); + return Build(spots.Select(s => (s.MasterAlignmentID, s.TimesCenter?.RT?.Value ?? double.NaN, s.MassCenter, (double)s.HeightAverage, + s.AdductType?.AdductIonName, s.MatchResults?.Representative, s.Ontology, s.Name)), polarity); + } + + private static LipidRtRunInput Build(IEnumerable<(int Id, double Rt, double Mz, double Height, string Adduct, MsScanMatchResult Representative, string Ontology, string Name)> objects, IonMode polarity) { + var candidates = new List(); + var context = new List(); + var rtMin = double.PositiveInfinity; + var rtMax = double.NegativeInfinity; + foreach (var o in objects) { + if (!double.IsNaN(o.Rt) && !double.IsInfinity(o.Rt)) { + rtMin = Math.Min(rtMin, o.Rt); + rtMax = Math.Max(rtMax, o.Rt); + } + context.Add(new LipidRtContextPeak + { + Id = o.Id, Rt = o.Rt, Mz = o.Mz, Height = o.Height, Adduct = o.Adduct ?? string.Empty, + IsMsmsAnnotated = IsAnnotated(o.Representative) && o.Representative.IsSpectrumMatch, + }); + var candidate = CreateCandidate(o.Id, o.Rt, o.Mz, o.Height, o.Adduct, o.Representative, o.Ontology, o.Name); + if (candidate != null) { + candidates.Add(candidate); + } + } + if (double.IsInfinity(rtMin)) { + rtMin = rtMax = 0d; + } + return new LipidRtRunInput(candidates, context, rtMin, rtMax, polarity); + } + + /// + /// True when the representative carries an annotation at all: a precursor m/z or a spectrum match (annotation + /// tag other than 999). + /// + private static bool IsAnnotated(MsScanMatchResult r) { + return r != null && !r.IsUnknown && (r.IsPrecursorMzMatch || r.IsSpectrumMatch); + } + + /// + /// The candidate for one object, or null when the object carries no annotation (it is then context only). + /// Eligibility follows SPEC 2.1 with addendum A6. + /// + /// MasterPeakID or MasterAlignmentID. + /// The representative match result, read directly: on a project restored without its + /// library the evaluator would have no entry for the annotator and would call every match unmatched. + /// The object's Ontology (LbmClass name). + /// The object's own name, which carries the "Unsettled:" prefix a person may have set. + public static LipidRtCandidate CreateCandidate(int id, double rt, double mz, double height, string adduct, + MsScanMatchResult representative, string ontology, string featureName) { + var r = representative; + if (!IsAnnotated(r)) { + return null; + } + var rawName = r.Name ?? string.Empty; + var name = AnnotationName.WithoutPrefix(rawName) ?? string.Empty; + var trimmed = name.Trim(); + if (trimmed.Length == 0 || string.Equals(trimmed, "Unknown", StringComparison.OrdinalIgnoreCase) || string.Equals(trimmed, "null", StringComparison.OrdinalIgnoreCase)) { + return null; + } + ontology = ontology ?? string.Empty; + var parsed = LipidRtNameParser.Parse(name, ontology); + var isLipid = LipidRtNameParser.IsLipidOntology(ontology, name); + var hasComposition = isLipid && parsed.IsParsed && parsed.Carbon > 0; + var tagLevel = r.IsSpectrumMatch + ? (r.IsLipidPositionMatch || r.IsLipidChainsMatch ? LipidRtTagLevel.Chain : r.IsLipidClassMatch ? LipidRtTagLevel.Class : LipidRtTagLevel.Msms) + : LipidRtTagLevel.None; + + var reasons = new List(); + if (!r.IsReferenceMatched) { + reasons.Add("not_ref_matched"); + } + if (!r.IsSpectrumMatch) { + reasons.Add("no_spectrum_match"); // m/z (+RT) only, e.g. tag 520: RT evidence about an RT-matched suggestion would be circular + } + if ((r.Source & SourceType.TextDB) != 0 || (r.Source & (SourceType.MspDB | SourceType.GeneratedLipid)) == 0) { + reasons.Add("not_spectral_library"); + } + if (!isLipid) { + reasons.Add("nonlipid"); + } + else if (parsed.IsNonChain || parsed.Carbon <= 0) { + reasons.Add("nonchain"); + } + if (parsed.IsLabelled) { + reasons.Add("internal_standard"); + } + if (r.IsManuallyModified) { + reasons.Add("manual"); + } + if (parsed.IsUnsettled || IsUnsettledName(featureName)) { + reasons.Add("unsettled"); + } + return new LipidRtCandidate + { + Id = id, + Rt = rt, + Mz = mz, + Height = height, + Adduct = adduct ?? string.Empty, + Name = parsed.Species.Length > 0 ? parsed.Species : name, + Molecular = parsed.Molecular, + Ontology = ontology, + ClassKey = parsed.ClassKey.Length > 0 ? parsed.ClassKey : ontology, + Carbon = parsed.Carbon, + DoubleBondEquivalent = parsed.DoubleBondEquivalent, + Oxygen = parsed.Oxygen, + IsLipid = hasComposition, + IsLabelled = parsed.IsLabelled, + IsMsmsAnnotated = r.IsSpectrumMatch, + TagLevel = tagLevel, + WeightedDot = r.WeightedDotProduct, + Chains = parsed.Chains, + IsEligible = reasons.Count == 0 && hasComposition, + IneligibleReason = string.Join(",", reasons), + }; + } + + private static bool IsUnsettledName(string name) { + return !string.IsNullOrEmpty(name) && name.TrimStart().StartsWith("Unsettled:", StringComparison.OrdinalIgnoreCase); + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtNameParser.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtNameParser.cs new file mode 100644 index 000000000..bb155d89e --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtNameParser.cs @@ -0,0 +1,592 @@ +using CompMs.Common.Enum; +using CompMs.Common.Lipidomics; +using System; +using System.Collections.Generic; +using System.Globalization; +using System.Linq; +using System.Text.RegularExpressions; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// How far a lipid name resolves its chains. + public enum LipidRtNameLevel : byte + { + None = 0, + /// No fatty chain and no sphingoid base (a sterol core only, for instance). + NonChain, + /// A species-level total only, e.g. "PC 34:1". + Species, + /// Some chains resolved, the rest a total. + Partial, + /// A class with one fatty chain, or a name with one chain. + SingleChain, + /// Every chain resolved. + Chain, + } + + /// A lipid name read for the retention-time order check. + public sealed class LipidRtParsedName + { + internal LipidRtParsedName() { } + + /// True when the name parsed and its total carbon count is positive. + public bool IsParsed { get; internal set; } + + /// Why the name did not parse, or an empty string. + public string Note { get; internal set; } = string.Empty; + + /// "Unsettled", "no MS2", "low score", "w/o MS2", or empty. + public string Prefix { get; internal set; } = string.Empty; + + public bool IsUnsettled { get; internal set; } + + /// An isotope-labelled name such as "(d7)", "(13C)" or "(15N)": an internal standard. + public bool IsLabelled { get; internal set; } + + /// The species-level half (left of "|"), labels and double-bond positions removed. + public string Species { get; internal set; } = string.Empty; + + /// The chain-level half, or the name itself when a single name resolves several chains, or empty. + public string Molecular { get; internal set; } = string.Empty; + + /// The class token of the name, e.g. "PC", "Cer", "SE". + public string ClassToken { get; internal set; } = string.Empty; + + /// The class key: the Ontology, except "BMP" for a BMP-named row whose Ontology is PG. + public string ClassKey { get; internal set; } = string.Empty; + + public int Carbon { get; internal set; } + + /// The double-bond count as written in the species total. + public int DoubleBond { get; internal set; } + + /// The double-bond count with a plasmenyl (P-) linkage counted as O- plus one double bond. + public int DoubleBondEquivalent { get; internal set; } + + /// Oxygens beyond the class's own, including sphingoid hydroxyls (";O2" counts 2). + public int Oxygen { get; internal set; } + + /// "O", "P" or empty. + public string Ether { get; internal set; } = string.Empty; + + /// The chains of the most resolved half, with their roles. + public IReadOnlyList Chains { get; internal set; } = Array.Empty(); + + public LipidRtNameLevel Level { get; internal set; } + + /// True when there is no chain at all, or the name did not parse. + public bool IsNonChain => !IsParsed || Level == LipidRtNameLevel.NonChain || Level == LipidRtNameLevel.None; + } + + /// + /// Reads an MS-DIAL lipid name into the species-level composition the retention-time order check works on: class + /// key, total carbons, total double bonds (P- counted as O- plus one) and total extra oxygens (SPEC section 2.2), + /// and into chains with roles for the plausibility profile (SPEC section 2.4). + /// + /// + /// The species totals come from + /// applied to the species half of a "species|molecular" name, as the specification says. Three forms are read by + /// the chain grammar below instead, because VS2 totals them differently from the species string MS-DIAL writes, + /// or cannot total them: + /// + /// sterol-type names (SE, SG, ASG, ST, BA): VS2 returns the acyl only and drops the ring core. The core is + /// the same for every row of one Ontology, so the order is unaffected, but the species total includes it; + /// a name VS2 cannot total (it throws, or returns no carbons), such as the LNAPE species form "LPE-N (FA)36:4"; + /// the ester scale: for a single chain-level name of FAHFA, NAGly, TG_EST, ASM, Cer_EOS and their kin, the + /// species string MS-DIAL writes carries one more double bond (and oxygen) for the ester bond than the chains add + /// up to, so that offset is added and one Ontology never mixes two scales. + /// + /// The grammar is a port of the reference parser of the specification (lipidparse.py), and the unit tests check the + /// totals against the reference values for every eligible row of the evaluation tables. + /// + public static class LipidRtNameParser + { + private static readonly Regex LabelRegex = new Regex(@"\((?:d|D)\d+\)|\(13C\d*\)|\(15N\d*\)", RegexOptions.CultureInvariant); + private static readonly Regex PositionRegex = new Regex(@"\((?:\d+[EZ]?)(?:,\d+[EZ]?)*\)", RegexOptions.CultureInvariant); + private static readonly Regex ChainTokenRegex = new Regex( + @"^(?[OP]-)?(?\d+):(?\d+)(?(?:;[A-Za-z0-9]*|;?\((?:\d*OH,?)+\))*)$", RegexOptions.CultureInvariant); + private static readonly Regex SpeciesRegex = new Regex(@"^(?[A-Za-z][A-Za-z0-9_\-]*?)\s+(?.+)$", RegexOptions.CultureInvariant | RegexOptions.Singleline); + private static readonly Regex HydroxylGroupRegex = new Regex(@"\((?:\d*OH,?)+\)", RegexOptions.CultureInvariant); + private static readonly Regex HydroxylGroupWithSeparatorRegex = new Regex(@";?\((?:\d*OH,?)+\)", RegexOptions.CultureInvariant); + private static readonly Regex OxygenTokenRegex = new Regex(@"^(\d*)O(\d*)$", RegexOptions.CultureInvariant); + private static readonly Regex DashBetweenDigitsRegex = new Regex(@"(?<=\d)-(?=\d)", RegexOptions.CultureInvariant); + private static readonly Regex CoreSplitRegex = new Regex(@"[_/]", RegexOptions.CultureInvariant); + private static readonly Regex SterolCoreRegex = new Regex(@"^(\d+:\d+(?:;[A-Za-z0-9]+)*)(.*)$", RegexOptions.CultureInvariant | RegexOptions.Singleline); + private static readonly Regex SterolAcylLeadRegex = new Regex(@"^(?:;|\s)*(?:FA|PA|LPA)\s+", RegexOptions.CultureInvariant); + private static readonly Regex PrefixGroupRegex = new Regex(@"^\((O-|FA)?([^()]*)\)", RegexOptions.CultureInvariant); + private static readonly Regex TailGroupRegex = new Regex(@"\(FA([^()]*)\)$", RegexOptions.CultureInvariant); + private static readonly Regex CarbonDoubleBondRegex = new Regex(@"\d+:\d+", RegexOptions.CultureInvariant); + private static readonly Regex Ms4OxygenRegex = new Regex(@"\+(\d*)O\b", RegexOptions.CultureInvariant); + + /// The isotope label pattern: an internal standard shares the endogenous Ontology in the lab LBM2. + public static bool IsLabelledName(string name) { + return !string.IsNullOrEmpty(name) && LabelRegex.IsMatch(name); + } + + /// Sphingolipid Ontologies: the first core chain is the sphingoid base, the second the N-acyl. + internal static readonly HashSet Sphingo = new HashSet(StringComparer.Ordinal) { + "Sph", "DHSph", "PhytoSph", "SM", "ASM", "CerP", "Cer_ADS", "Cer_AS", "Cer_BDS", "Cer_BS", "Cer_OS", + "Cer_EODS", "Cer_EOS", "Cer_NDS", "Cer_NS", "Cer_NP", "Cer_AP", "Cer_EBDS", "Cer_HS", "Cer_HDS", + "Cer_NDOS", "Cer_ABP", "Cer_AH", "Cer_NH", "HexCer_NS", "HexCer_NDS", "HexCer_AP", "HexCer_HS", + "HexCer_HDS", "HexCer_EOS", "Hex2Cer", "Hex3Cer", "SHexCer", "GM3", "GM1", "GD1a", "GD1b", "GD2", + "GD3", "GT1b", "GQ1b", "NGcGM3", "AHexCer", "ASHexCer", "PE_Cer", "PI_Cer", "MIPC", "SL", + "Cer_NS_d7", "SM_d9", "Cer_NH_d9", "Cer_AH_d9", + }; + + private static readonly HashSet SphingoTokens = new HashSet(StringComparer.Ordinal) { + "SPB", "Cer", "SM", "HexCer", "Hex2Cer", "Hex3Cer", "SHexCer", "GM3", "GM1", "GD1a", + "GD1b", "GD2", "GD3", "GT1b", "GQ1b", "NGcGM3", "AHexCer", "ASHexCer", "PE-Cer", "PI-Cer", + "MIPC", "SL", "CerP", + }; + + /// Sphingoid bases without an N-acyl. + internal static readonly HashSet SphingoidBaseOnly = new HashSet(StringComparer.Ordinal) { "Sph", "DHSph", "PhytoSph" }; + + private static readonly HashSet SterolTokens = new HashSet(StringComparer.Ordinal) { "SE", "SG", "ASG", "ST", "BA" }; + + private static readonly HashSet SingleChain = new HashSet(StringComparer.Ordinal) { + "FA", "OxFA", "CAR", "NAE", "NAGln", "NALeu", "NAPhe", "NASer", "NAVal", "NAAla", "NAGABA", "NAAnt", + "NA5HT", "NATau", "LPC", "LPE", "LPA", "LPS", "LPI", "LPG", "EtherLPC", "EtherLPE", "EtherLPG", + "EtherLPS", "EtherLPI", "MG", "CE", "VAE", "LDGTS", "LDGCC", "LDGTA", "MGMG", "DMEDFA", "DMEDOxFA", + "GPNAE", "BisMeLPA", "BRSE", "CASE", "SISE", "STSE", "EGSE", "DEGSE", "DSMSE", "DCAE", "GDCAE", "GLCAE", + "TDCAE", "TLCAE", "LCAE", "KLCAE", "KDCAE", "AHexBRS", "AHexCAS", "AHexCS", "AHexSIS", "AHexSTS", + "CSLPHex", "BRSLPHex", "CASLPHex", "SISLPHex", "STSLPHex", "CE_d7", "LPC_d5", "LPE_d5", "LPG_d5", + "LPS_d5", "LPI_d5", "PC_d5", + }; + + /// The number of fatty chains (alkyl, N-acyl and ester acyl included; a sphingoid base counted apart). + private static readonly Dictionary FattyChainCount = BuildFattyChainCount(); + + private static Dictionary BuildFattyChainCount() { + var counts = new Dictionary(StringComparer.Ordinal) { + ["TG"] = 3, ["OxTG"] = 3, ["EtherTG"] = 3, ["TG_EST"] = 4, ["DG"] = 2, ["EtherDG"] = 2, ["PC"] = 2, ["PE"] = 2, ["PS"] = 2, ["PG"] = 2, + ["PI"] = 2, ["PA"] = 2, ["PT"] = 2, ["PMeOH"] = 2, ["PEtOH"] = 2, ["PBtOH"] = 2, ["DMPE"] = 2, ["MMPE"] = 2, ["BMP"] = 2, ["HBMP"] = 3, + ["EtherPC"] = 2, ["EtherPE"] = 2, ["EtherPS"] = 2, ["EtherPG"] = 2, ["EtherPI"] = 2, ["OxPC"] = 2, ["OxPE"] = 2, ["OxPS"] = 2, + ["OxPG"] = 2, ["OxPI"] = 2, ["OxPA"] = 2, ["EtherOxPC"] = 2, ["EtherOxPE"] = 2, ["CL"] = 4, ["DLCL"] = 2, ["MLCL"] = 3, + ["MGDG"] = 2, ["DGDG"] = 2, ["SQDG"] = 2, ["DGTS"] = 2, ["DGTA"] = 2, ["DGGA"] = 2, ["ADGGA"] = 3, ["DGCC"] = 2, ["SMGDG"] = 2, + ["EtherMGDG"] = 2, ["EtherDGDG"] = 2, ["EtherSMGDG"] = 2, ["LNAPE"] = 2, ["LNAPS"] = 2, ["FAHFA"] = 2, ["DMEDFAHFA"] = 2, + ["NAGly"] = 2, ["NAGlySer"] = 2, ["NAOrn"] = 2, ["NATryA"] = 2, ["WE"] = 2, ["DGMG"] = 1, ["Ac2PIM1"] = 2, ["Ac2PIM2"] = 2, + ["Ac3PIM2"] = 3, ["Ac4PIM2"] = 4, ["SM"] = 1, ["Cer_NS"] = 1, ["Cer_NDS"] = 1, ["Cer_NP"] = 1, ["Cer_AS"] = 1, ["Cer_ADS"] = 1, + ["Cer_AP"] = 1, ["Cer_HS"] = 1, ["Cer_HDS"] = 1, ["Cer_BS"] = 1, ["Cer_BDS"] = 1, ["Cer_EOS"] = 2, ["Cer_EODS"] = 2, + ["Cer_EBDS"] = 2, ["HexCer_EOS"] = 2, ["ASM"] = 2, ["AHexCer"] = 2, ["ASHexCer"] = 2, ["SISPHex"] = 2, ["STSPHex"] = 2, + ["BRSPHex"] = 2, ["CASPHex"] = 2, ["CSPHex"] = 2, + }; + foreach (var name in SingleChain) { + if (!counts.ContainsKey(name)) { + counts[name] = 1; + } + } + foreach (var name in Sphingo) { + if (!counts.ContainsKey(name)) { + counts[name] = SphingoidBaseOnly.Contains(name) ? 0 : 1; + } + } + return counts; + } + + /// + /// The (double bond, oxygen) offset MS-DIAL's species string carries for the ester bond, relative to the sum of + /// the chains (SPEC 2.2, verified on the dual names). + /// + internal static readonly Dictionary EsterScale = new Dictionary(StringComparer.Ordinal) { + ["FAHFA"] = (1, 0), ["DMEDFAHFA"] = (1, 0), + ["NAGly"] = (1, 1), ["NAGlySer"] = (1, 1), ["NAOrn"] = (1, 1), ["NATryA"] = (1, 1), + ["TG_EST"] = (1, 1), ["ASM"] = (1, 1), + ["Cer_EOS"] = (1, 1), ["Cer_EODS"] = (1, 1), ["Cer_EBDS"] = (1, 1), ["HexCer_EOS"] = (1, 1), + }; + + private static readonly HashSet NonLipidOntologies = new HashSet(StringComparer.Ordinal) { + "", "null", "NA", "Undefined", "Others", "Unknown", "SPLASH", + }; + + private static readonly string[] Prefixes = { "Unsettled:", "w/o MS2:", "no MS2:", "low score:" }; + + /// + /// True when the Ontology names a lipid class (an member other than the placeholders) and the + /// name is not a RIKEN VS1 record. The non-lipid MSP matches of a lipidomics run fail this. + /// + public static bool IsLipidOntology(string ontology, string name) { + if (ontology is null || NonLipidOntologies.Contains(ontology)) { + return false; + } + if (ontology.Length == 0 || !char.IsLetter(ontology[0])) { + return false; + } + if (!System.Enum.TryParse(ontology, ignoreCase: false, out LbmClass lbmClass) || !System.Enum.IsDefined(typeof(LbmClass), lbmClass)) { + return false; + } + if (lbmClass.ToString() != ontology) { + return false; + } + return name is null || !name.TrimStart().StartsWith("RIKEN ", StringComparison.Ordinal); + } + + /// + /// Parses an MS-DIAL lipid name. is the feature's Ontology (LbmClass name); it decides + /// the chain roles and the scale of the totals. + /// + public static LipidRtParsedName Parse(string name, string ontology) { + ontology = ontology ?? string.Empty; + var result = new LipidRtParsedName(); + if (name is null) { + result.Note = "empty"; + return result; + } + var n = name.Trim().Trim('"').Trim(); + foreach (var prefix in Prefixes) { + if (n.StartsWith(prefix, StringComparison.OrdinalIgnoreCase)) { + result.Prefix = prefix.TrimEnd(':'); + if (prefix == "Unsettled:") { + result.IsUnsettled = true; + } + n = n.Substring(prefix.Length).Trim(); + } + } + if (n.Length == 0 || string.Equals(n, "unknown", StringComparison.OrdinalIgnoreCase) || string.Equals(n, "null", StringComparison.OrdinalIgnoreCase)) { + result.Note = "unknown"; + return result; + } + result.IsLabelled = LabelRegex.IsMatch(n); + var clean = LabelRegex.Replace(n, string.Empty); + clean = PositionRegex.Replace(clean, string.Empty); + clean = clean.Replace(";(", "("); + clean = Ms4OxygenRegex.Replace(clean, m => ";O" + m.Groups[1].Value); + var halves = clean.Split('|').Select(h => h.Trim()).ToArray(); + + var species = ParseOne(halves[0], ontology, out var note); + if (species is null) { + result.Note = note; + return result; + } + ParsedHalf molecular = null; + if (halves.Length > 1 && halves[1].Length > 0) { + molecular = ParseOne(halves[1], ontology, out var note2); + if (molecular is null) { + result.Note = "molecular half: " + note2; + } + } + result.ClassToken = species.ClassToken; + var speciesChains = species.Chains; + var speciesResolved = !speciesChains.Any(c => c.IsUnresolved) + && speciesChains.Any(c => c.IsFatty || c.Role == LipidRtChainRole.SphingoidBase); + IReadOnlyList chains; + if (molecular != null) { + chains = molecular.Chains; + result.Species = halves[0]; + result.Molecular = halves[1]; + } + else { + chains = speciesChains; + result.Species = halves[0]; + result.Molecular = speciesResolved && speciesChains.Count > 1 ? halves[0] : string.Empty; + } + + // Totals of the species half. + var carbon = speciesChains.Sum(c => c.Carbon); + var doubleBond = speciesChains.Sum(c => c.DoubleBond); + var oxygen = speciesChains.Sum(c => c.Oxygen); + var plasmenyl = speciesChains.Any(c => c.Ether == "P" && (c.Role == LipidRtChainRole.Alkyl || c.Role == LipidRtChainRole.SpeciesCore || c.Role == LipidRtChainRole.Acyl)); + if (!SterolTokens.Contains(species.ClassToken) + && TryTotalsByLipidomicsConverter(DashBetweenDigitsRegex.Replace(halves[0], "_"), ontology, out var vs2Carbon, out var vs2DoubleBondEquivalent, out var vs2Oxygen)) { + // VS2 counts a P- chain with one more double bond: that is already the equivalent. + carbon = vs2Carbon; + doubleBond = vs2DoubleBondEquivalent - (plasmenyl ? 1 : 0); + oxygen = vs2Oxygen; + } + var applyEsterScale = molecular is null && speciesChains.Count > 1 + && EsterScale.TryGetValue(ontology, out _) + && (speciesChains.Any(c => c.Role == LipidRtChainRole.EsterAcyl) || ontology == "FAHFA" || ontology == "DMEDFAHFA"); + if (applyEsterScale) { + var offset = EsterScale[ontology]; + doubleBond += offset.DoubleBond; + oxygen += offset.Oxygen; + } + + var ether = string.Empty; + foreach (var chain in speciesChains) { + if (chain.Ether.Length > 0 && (chain.Role == LipidRtChainRole.Alkyl || chain.Role == LipidRtChainRole.SpeciesCore || chain.Role == LipidRtChainRole.Acyl)) { + ether = chain.Ether; + break; + } + } + + var fatty = chains.Count(c => c.IsFatty); + var unresolved = chains.Count(c => c.IsUnresolved); + var bases = chains.Count(c => c.Role == LipidRtChainRole.SphingoidBase); + LipidRtNameLevel level; + if (fatty == 0 && unresolved == 0 && bases == 0) { + level = LipidRtNameLevel.NonChain; + } + else if (unresolved > 0 && (fatty > 0 || bases > 0)) { + level = LipidRtNameLevel.Partial; + } + else if (unresolved > 0) { + level = LipidRtNameLevel.Species; + } + else if (SingleChain.Contains(ontology) || SphingoidBaseOnly.Contains(ontology) || fatty + bases == 1) { + level = LipidRtNameLevel.SingleChain; + } + else { + level = LipidRtNameLevel.Chain; + } + + result.IsParsed = carbon > 0; + result.Carbon = carbon; + result.DoubleBond = doubleBond; + result.DoubleBondEquivalent = doubleBond + (ether == "P" ? 1 : 0); + result.Oxygen = oxygen; + result.Ether = ether; + result.Chains = chains; + result.Level = level; + result.ClassKey = ontology == "PG" && species.ClassToken == "BMP" ? "BMP" : ontology; + if (carbon == 0 && result.Note.Length == 0) { + result.Note = "C=0"; + } + return result; + } + + /// + /// The species totals by MS-DIAL's own lipid name converter. False when it cannot total the name. + /// + private static bool TryTotalsByLipidomicsConverter(string species, string ontology, out int carbon, out int doubleBondEquivalent, out int oxygen) { + carbon = doubleBondEquivalent = oxygen = 0; + try { + var molecule = LipidomicsConverter.ConvertMsdialLipidnameToLipidMoleculeObjectVS2(species, ontology); + if (molecule is null || molecule.TotalCarbonCount <= 0) { + return false; + } + carbon = molecule.TotalCarbonCount; + doubleBondEquivalent = molecule.TotalDoubleBondCount; + oxygen = molecule.TotalOxidizedCount; + return true; + } + catch (Exception ex) when (ex is IndexOutOfRangeException || ex is ArgumentException || ex is FormatException || ex is NullReferenceException) { + return false; + } + } + + private sealed class ParsedHalf + { + public string ClassToken; + public List Chains; + } + + private sealed class RawChain + { + public int Carbon; + public int DoubleBond; + public int Oxygen; + public string Ether; + public string Role; // provisional: "core", or a final role token + } + + private static ParsedHalf ParseOne(string text, string ontology, out string note) { + note = string.Empty; + var m = SpeciesRegex.Match(text.Trim()); + if (!m.Success) { + note = "no 'CLASS C:DB' form"; + return null; + } + var classToken = m.Groups["cls"].Value; + var body = m.Groups["body"].Value.Trim(); + if (!CarbonDoubleBondRegex.IsMatch(body)) { + note = "no C:DB"; + return null; + } + var chains = ParseBody(classToken, body, out note); + if (chains is null) { + return null; + } + return new ParsedHalf { ClassToken = classToken, Chains = AssignRoles(chains, ontology, classToken), }; + } + + private static List ParseBody(string classToken, string body, out string note) { + note = string.Empty; + var chains = new List(); + if (SterolTokens.Contains(classToken)) { + // "29:1;O;Hex;PA 7:0_15:3" / "29:2;O;Hex;FA 11:0" / "28:1/18:3" / "24:1;O3;T/20:1" / "27:1;O;PE" + var m = SterolCoreRegex.Match(body.Trim()); + if (!m.Success) { + note = "sterol core not found"; + return null; + } + var core = ParseChainToken(m.Groups[1].Value); + if (core is null) { + note = "bad sterol core"; + return null; + } + core.Role = "sterol_core"; + chains.Add(core); + var rest = m.Groups[2].Value.Trim(); + rest = SterolAcylLeadRegex.Replace(rest, string.Empty, 1); + rest = rest.TrimStart('/', ';').Trim(); + if (rest.Length > 0) { + foreach (var token in SplitCore(rest.Replace(" ", string.Empty))) { + var chain = ParseChainToken(token); + if (chain is null) { + note = $"bad sterol acyl '{token}'"; + return null; + } + if (chain.Carbon == 0 && chain.DoubleBond == 0) { + continue; + } + chain.Role = "acyl"; + chains.Add(chain); + } + } + return chains; + } + + body = body.Replace(" ", string.Empty); + // prefix groups "(O-12:0)" (O-acyl on a hexose) and "(FA 15:0)" / "(FA)" (N-acyl of LNAPE) + while (body.StartsWith("(", StringComparison.Ordinal)) { + var m = PrefixGroupRegex.Match(body); + if (!m.Success) { + break; + } + var inner = m.Groups[2].Value; + if (inner.Length > 0) { + var chain = ParseChainToken(inner); + if (chain is null) { + note = $"bad prefix group '{m.Value}'"; + return null; + } + chain.Role = m.Groups[1].Value == "O-" ? "o_acyl" : "acyl"; + chain.Ether = string.Empty; + chains.Add(chain); + } + body = body.Substring(m.Index + m.Length); + } + // tail groups "(FA 16:0)" (ester-linked acyl of EOS/ASM/TG_EST/NAGly) + var tails = new List(); + while (true) { + var m = TailGroupRegex.Match(body); + if (!m.Success) { + break; + } + var chain = ParseChainToken(m.Groups[1].Value); + if (chain is null) { + note = $"bad tail '{m.Value}'"; + return null; + } + chain.Role = "ester_acyl"; + tails.Insert(0, chain); + body = body.Substring(0, m.Index); + } + var tokens = SplitCore(body); + if (tokens.Count == 0) { + note = "no core chain"; + return null; + } + var cores = new List(); + foreach (var token in tokens) { + var chain = ParseChainToken(token); + if (chain is null) { + note = $"bad chain '{token}'"; + return null; + } + if (chain.Carbon == 0 && chain.DoubleBond == 0 && chain.Oxygen == 0) { + continue; // "LPC 16:0/0:0" + } + chain.Role = "core"; + cores.Add(chain); + } + if (cores.Count == 0) { + note = "only empty chains"; + return null; + } + chains.AddRange(cores); + chains.AddRange(tails); + return chains; + } + + private static List AssignRoles(List chains, string ontology, string classToken) { + var core = chains.Where(c => c.Role == "core").ToList(); + var sphingo = Sphingo.Contains(ontology) || (NonLipidOntologies.Contains(ontology) && SphingoTokens.Contains(classToken)); + if (sphingo) { + if (SphingoidBaseOnly.Contains(ontology) || classToken == "SPB") { + foreach (var c in core) { + c.Role = "sph_base"; + } + } + else if (core.Count >= 2) { + core[0].Role = "sph_base"; + for (int i = 1; i < core.Count; i++) { + core[i].Role = "n_acyl"; + } + } + else { + foreach (var c in core) { + c.Role = "sph_core"; + } + } + } + else { + var single = SingleChain.Contains(ontology); + // "NAGly 8:0;O(FA 19:1)", "LPE-N (FA 15:0)18:2": the explicit prefix/tail acyls plus the one core chain + // already account for every fatty chain of the class, so the core is resolved. + var known = chains.Count(c => c.Role == "acyl" || c.Role == "ester_acyl" || c.Role == "o_acyl"); + var expected = FattyChainCount.TryGetValue(ontology, out var count) ? count : 99; + var complete = core.Count == 1 && known >= 1 && known + 1 >= expected; + if (core.Count == 1 && !single && !complete) { + core[0].Role = "species_core"; + } + else { + foreach (var c in core) { + c.Role = c.Ether.Length > 0 ? "alkyl" : "acyl"; + } + } + } + return chains.Select(c => new LipidRtChain(LipidRtChain.ParseRole(c.Role), c.Carbon, c.DoubleBond, c.Oxygen, c.Ether)).ToList(); + } + + private static List SplitCore(string core) { + core = DashBetweenDigitsRegex.Replace(core, "_"); // MS-DIAL 4 "16:0-18:1" + return CoreSplitRegex.Split(core).Where(t => t.Length > 0).ToList(); + } + + private static RawChain ParseChainToken(string token) { + var m = ChainTokenRegex.Match(token.Trim()); + if (!m.Success) { + return null; + } + var oxygen = CountOxygen(m.Groups["mods"].Value); + return new RawChain { + Carbon = int.Parse(m.Groups["c"].Value, CultureInfo.InvariantCulture), + DoubleBond = int.Parse(m.Groups["db"].Value, CultureInfo.InvariantCulture), + Oxygen = oxygen, + Ether = m.Groups["eth"].Value.TrimEnd('-'), + }; + } + + /// ";O" is 1, ";O2" and ";2O" are 2, "(2OH)" is 1, "(1OH,3OH)" is 2; other ";X" tokens are ignored. + private static int CountOxygen(string mods) { + var n = 0; + foreach (Match g in HydroxylGroupRegex.Matches(mods)) { + n += CountOccurrences(g.Value, "OH"); + } + mods = HydroxylGroupWithSeparatorRegex.Replace(mods, string.Empty); + foreach (var token in mods.Split(';')) { + if (token.Length == 0) { + continue; + } + var m = OxygenTokenRegex.Match(token); + if (m.Success) { + var digits = m.Groups[1].Value.Length > 0 ? m.Groups[1].Value : m.Groups[2].Value; + n += digits.Length > 0 ? int.Parse(digits, CultureInfo.InvariantCulture) : 1; + } + } + return n; + } + + private static int CountOccurrences(string text, string value) { + var count = 0; + var index = 0; + while ((index = text.IndexOf(value, index, StringComparison.Ordinal)) >= 0) { + count++; + index += value.Length; + } + return count; + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtPlausibilityProfile.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtPlausibilityProfile.cs new file mode 100644 index 000000000..8012811ef --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtPlausibilityProfile.cs @@ -0,0 +1,176 @@ +using CompMs.MsdialCore.Parameter; +using System; +using System.Collections.Generic; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// The verdict of a plausibility profile on one candidate. + public sealed class LipidRtPlausibilityVerdict + { + public static LipidRtPlausibilityVerdict NotEvaluated { get; } = new LipidRtPlausibilityVerdict(LipidRtPlausibility.None, LipidRtPlausibilitySeverity.None, string.Empty, string.Empty); + + public LipidRtPlausibilityVerdict(LipidRtPlausibility verdict, LipidRtPlausibilitySeverity severity, string rule, string detail) { + Verdict = verdict; + Severity = severity; + Rule = rule ?? string.Empty; + Detail = detail ?? string.Empty; + } + + public LipidRtPlausibility Verdict { get; } + public LipidRtPlausibilitySeverity Severity { get; } + + /// The failed rules, e.g. "B1" or "B2,B3"; empty unless Implausible. + public string Rule { get; } + + public string Detail { get; } + + public bool IsImplausible => Verdict == LipidRtPlausibility.Implausible; + public bool IsHard => IsImplausible && Severity == LipidRtPlausibilitySeverity.Hard; + } + + /// + /// A biological plausibility profile: chain-level rules for an organism, applied as a layer separate from the + /// retention-time evidence (SPEC 2.4). + /// + /// + /// A profile states chain-level constraints only. It must never exclude or penalise a lipid class: the author's + /// decision of 2026-10-15 (animals synthesise MGDG and DGDG, and diet and gut bacteria bring other classes). + /// + public interface ILipidRtPlausibilityProfile + { + string Name { get; } + + LipidRtPlausibilityVerdict Evaluate(string ontology, IReadOnlyList chains); + } + + /// + /// The chain-level plausibility rules B1-B3 as parameters: maximum double bonds per fatty chain (B1, hard), + /// sphingolipid N-acyl double bonds (B2, soft), sphingoid base carbons and double bonds (B3, soft). + /// + /// + /// Rule semantics: a rule fails when a resolved chain breaks it; passes when every chain it concerns is resolved + /// and passes; is untestable when it applies but the chain is not resolved (a species-level name); and is not + /// applicable when it does not concern the class. The verdict is Implausible if any enabled rule fails (Hard when + /// B1 is among them), Plausible if every applicable rule passes, and Untestable otherwise. + /// + public sealed class LipidRtChainRuleProfile : ILipidRtPlausibilityProfile + { + /// The author's mammalian rules (2026-10-15): B1 DB <= 6 per chain; B2 N-acyl DB <= 2; B3 base C16-20, DB <= 2. + public static LipidRtChainRuleProfile Mammalian { get; } = new LipidRtChainRuleProfile("Mammalian", maxChainDoubleBond: 6, nAcylMaxDoubleBond: 2, baseMinCarbon: 16, baseMaxCarbon: 20, baseMaxDoubleBond: 2); + + public LipidRtChainRuleProfile(string name, int maxChainDoubleBond, int nAcylMaxDoubleBond, int baseMinCarbon, int baseMaxCarbon, int baseMaxDoubleBond, + bool b1 = true, bool b2 = true, bool b3 = true) { + Name = name ?? throw new ArgumentNullException(nameof(name)); + MaxChainDoubleBond = maxChainDoubleBond; + NAcylMaxDoubleBond = nAcylMaxDoubleBond; + BaseMinCarbon = baseMinCarbon; + BaseMaxCarbon = baseMaxCarbon; + BaseMaxDoubleBond = baseMaxDoubleBond; + IsB1Enabled = b1; + IsB2Enabled = b2; + IsB3Enabled = b3; + } + + public string Name { get; } + public int MaxChainDoubleBond { get; } + public int NAcylMaxDoubleBond { get; } + public int BaseMinCarbon { get; } + public int BaseMaxCarbon { get; } + public int BaseMaxDoubleBond { get; } + public bool IsB1Enabled { get; } + public bool IsB2Enabled { get; } + public bool IsB3Enabled { get; } + + /// The profile a parameter selects, or null for none. + public static ILipidRtPlausibilityProfile FromParameter(LipidRtPlausibilityProfileKind kind) { + switch (kind) { + case LipidRtPlausibilityProfileKind.Mammalian: + return Mammalian; + default: + return null; + } + } + + private enum RuleResult { NotApplicable, Pass, Fail, Untestable } + + public LipidRtPlausibilityVerdict Evaluate(string ontology, IReadOnlyList chains) { + ontology = ontology ?? string.Empty; + chains = chains ?? Array.Empty(); + var results = new List<(string Rule, RuleResult Result, string Detail)>(); + if (IsB1Enabled) { + results.Add(("B1", B1(chains, out var detail), detail)); + } + if (IsB2Enabled) { + results.Add(("B2", B2(chains, ontology, out var detail), detail)); + } + if (IsB3Enabled) { + results.Add(("B3", B3(chains, ontology, out var detail), detail)); + } + var fails = results.Where(r => r.Result == RuleResult.Fail).ToList(); + if (fails.Count > 0) { + var severity = fails.Any(r => r.Rule == "B1") ? LipidRtPlausibilitySeverity.Hard : LipidRtPlausibilitySeverity.Soft; + return new LipidRtPlausibilityVerdict(LipidRtPlausibility.Implausible, severity, + string.Join(",", fails.Select(r => r.Rule)), + string.Join("; ", fails.Select(r => r.Rule + ": " + r.Detail))); + } + var applicable = results.Where(r => r.Result != RuleResult.NotApplicable).ToList(); + if (applicable.Count > 0 && applicable.All(r => r.Result == RuleResult.Pass)) { + return new LipidRtPlausibilityVerdict(LipidRtPlausibility.Plausible, LipidRtPlausibilitySeverity.None, string.Empty, string.Empty); + } + return new LipidRtPlausibilityVerdict(LipidRtPlausibility.Untestable, LipidRtPlausibilitySeverity.None, string.Empty, string.Empty); + } + + private RuleResult B1(IReadOnlyList chains, out string detail) { + detail = string.Empty; + var bad = chains.FirstOrDefault(c => c.IsFatty && c.DoubleBond > MaxChainDoubleBond); + if (bad != null) { + detail = $"chain {bad.Carbon}:{bad.DoubleBond} has DB>={MaxChainDoubleBond + 1}"; + return RuleResult.Fail; + } + var fatty = chains.Any(c => c.IsFatty); + var unresolved = chains.Any(c => c.IsUnresolved); + if (!fatty && !unresolved) { + return RuleResult.NotApplicable; + } + return unresolved ? RuleResult.Untestable : RuleResult.Pass; + } + + private RuleResult B2(IReadOnlyList chains, string ontology, out string detail) { + detail = string.Empty; + if (!LipidRtNameParser.Sphingo.Contains(ontology) || LipidRtNameParser.SphingoidBaseOnly.Contains(ontology)) { + return RuleResult.NotApplicable; + } + var nAcyl = chains.Where(c => c.Role == LipidRtChainRole.NAcyl).ToList(); + if (nAcyl.Count > 0) { + var bad = nAcyl.FirstOrDefault(c => c.DoubleBond > NAcylMaxDoubleBond); + if (bad != null) { + detail = $"N-acyl {bad.Carbon}:{bad.DoubleBond} has DB>={NAcylMaxDoubleBond + 1}"; + return RuleResult.Fail; + } + return RuleResult.Pass; + } + return RuleResult.Untestable; + } + + private RuleResult B3(IReadOnlyList chains, string ontology, out string detail) { + detail = string.Empty; + if (!LipidRtNameParser.Sphingo.Contains(ontology)) { + return RuleResult.NotApplicable; + } + var b = chains.FirstOrDefault(c => c.Role == LipidRtChainRole.SphingoidBase); + if (b is null) { + return RuleResult.Untestable; + } + if (b.Carbon < BaseMinCarbon || b.Carbon > BaseMaxCarbon) { + detail = $"base {b.Carbon}:{b.DoubleBond} carbons outside {BaseMinCarbon}-{BaseMaxCarbon}"; + return RuleResult.Fail; + } + if (b.DoubleBond > BaseMaxDoubleBond) { + detail = $"base {b.Carbon}:{b.DoubleBond} has DB>{BaseMaxDoubleBond}"; + return RuleResult.Fail; + } + return RuleResult.Pass; + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtResult.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtResult.cs new file mode 100644 index 000000000..ef01d9c17 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtResult.cs @@ -0,0 +1,326 @@ +using CompMs.Common.Enum; +using System; +using System.Collections.Generic; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// An adjacent neighbour that supports a feature's position. + public sealed class LipidRtNeighbour + { + public LipidRtNeighbour(string direction, string feature, double rt, double gap) { + Direction = direction; + Feature = feature; + Rt = rt; + Gap = gap; + } + + /// + /// "C<", "D>" or "O>": a lower bound (a smaller, more unsaturated or more oxidised neighbour elutes + /// earlier); "C>", "D<" or "O<": an upper bound. + /// + public string Direction { get; } + + /// The neighbour feature, e.g. "TG 50:1@12.04". + public string Feature { get; } + + public double Rt { get; } + + /// The feature's RT minus the neighbour's (min). + public double Gap { get; } + + public bool IsLowerBound => IsLowerBoundDirection(Direction); + + /// True for "C<", "D>" and "O>", the directions of a lower bound. + public static bool IsLowerBoundDirection(string direction) { + return direction == "C<" || direction == "D>" || direction == "O>"; + } + } + + /// A feature whose order contradicts this one. + public sealed class LipidRtDiscordantPartner + { + public LipidRtDiscordantPartner(string feature, double rt, LipidRtFeatureStatus status) { + Feature = feature; + Rt = rt; + Status = status; + } + + public string Feature { get; } + public double Rt { get; } + public LipidRtFeatureStatus Status { get; } + } + + /// An MS1 artefact explanation: a brighter co-eluting peak this ion can be a by-product of. + public sealed class LipidRtArtefact + { + public LipidRtArtefact(string kind, string text, int partnerId) { + Kind = kind; + Text = text; + PartnerId = partnerId; + } + + /// "iso1", "iso2", "naH", "naNH4", "kH", "kNa", "nh4H", "acetate", "formate", "chloride", "demethyl", "isf_chain", "isf_head" or "ox". + public string Kind { get; } + + public string Text { get; } + + /// The Id of the brighter co-eluting peak. + public int PartnerId { get; } + } + + /// An eligible candidate of the other polarity with the same composition at the same (mapped) RT. + public sealed class LipidRtTwin + { + public LipidRtTwin(int id, string adduct, double rt, LipidRtReason singleRunReason) { + Id = id; + Adduct = adduct; + Rt = rt; + SingleRunReason = singleRunReason; + } + + public int Id { get; } + public string Adduct { get; } + public double Rt { get; } + + /// The twin's own single-run reason; corroboration never chains. + public LipidRtReason SingleRunReason { get; } + } + + /// The run-level reversed-phase gate (SPEC section 7). + public sealed class LipidRtGateStatistics + { + public bool Passed { get; internal set; } + + /// True when the parameter said the chromatography is not reversed phase; nothing was measured. + public bool ChromatographyIsOther { get; internal set; } + + /// The median over classes of the median within-class carbon slope (min per carbon). + public double MedianCarbonSlope { get; internal set; } + + public double SlopeThreshold { get; internal set; } + + /// Classes that contributed a slope (at least 5 features and 3 distinct carbon counts). + public int ClassCount { get; internal set; } + + /// The adjacent-carbon concordance on the pairs used, or null when there were none. + public double? Concordance { get; internal set; } + + /// The number of adjacent-carbon pairs used. + public int PairCount { get; internal set; } + + /// "trusted" (chain-level features at or above the class median height) or "all". + public string ConcordanceOn { get; internal set; } = string.Empty; + + /// The concordance over all adjacent-carbon pairs, or null when there were none. + public double? ConcordanceAllPairs { get; internal set; } + + public override string ToString() { + if (ChromatographyIsOther) { + return "chromatography parameter says not RP"; + } + return $"median carbon slope {MedianCarbonSlope:0.0000} vs threshold {SlopeThreshold:0.0000} over {ClassCount} classes; " + + $"adjacent-carbon concordance {(Concordance.HasValue ? Concordance.Value.ToString("0.000") : "-")} over {PairCount} {ConcordanceOn} pairs" + + $" (all pairs {(ConcordanceAllPairs.HasValue ? ConcordanceAllPairs.Value.ToString("0.000") : "-")})"; + } + } + + /// The machine-readable evidence behind one object's outcome (SPEC section 8). + public sealed class LipidRtEvidence + { + internal LipidRtEvidence() { } + + /// The feature, e.g. "TG 52:1@12.34"; empty when the object was not evaluated as part of a feature. + public string Feature { get; internal set; } = string.Empty; + + public LipidRtComposition Composition { get; internal set; } + + public double FeatureRt { get; internal set; } + + public LipidRtFeatureStatus Status { get; internal set; } + + /// The essential margin e (evidence units) of a Core feature, the weight of a Free one, 0 otherwise. + public double Margin { get; internal set; } + + /// The deficit d (evidence units) of an Excluded feature, 0 otherwise. + public double Deficit { get; internal set; } + + /// The integer weight of the feature (x100 evidence units). + public long Weight { get; internal set; } + + /// The supporting directions, sorted, e.g. "C<C>D<D>". + public string Directions { get; internal set; } = string.Empty; + + public int Axes { get; internal set; } + + public int Ties { get; internal set; } + + /// + /// Up to the configured number of supporting neighbours, listed by |gap| descending then direction; the strongest + /// lower and upper bound, and the strongest of each direction, are kept before the rest. + /// + public IReadOnlyList Support { get; internal set; } = Array.Empty(); + + /// Up to the configured number of discordant partners, by weight descending then RT. + public IReadOnlyList Discordant { get; internal set; } = Array.Empty(); + + public int DiscordantCount { get; internal set; } + + /// The adjacent composition this feature ties with, when the Tie rule decided. + public string TiePartner { get; internal set; } = string.Empty; + + /// The less intense duplicate's dominant twin, when the MinorDuplicate rule decided. + public string MinorDuplicateOf { get; internal set; } = string.Empty; + + public double? Predicted { get; internal set; } + public double? Residual { get; internal set; } + public double? ResidualDoubleBondSteps { get; internal set; } + + public bool ClassGatePassed { get; internal set; } = true; + public double? ClassAdjacentConcordance { get; internal set; } + + /// The members of the feature (object Ids) and the number of distinct adducts. + public IReadOnlyList Members { get; internal set; } = Array.Empty(); + public int AdductCount { get; internal set; } + + /// This object's own MS1 artefact explanations. + public IReadOnlyList Artefacts { get; internal set; } = Array.Empty(); + + /// The artefact explanation of the feature's representative member, when every member has one. + public LipidRtArtefact FeatureArtefact { get; internal set; } + + public LipidRtPlausibilityVerdict Plausibility { get; internal set; } = LipidRtPlausibilityVerdict.NotEvaluated; + + /// True when a soft plausibility rule failed and "RT and spectrum both excellent" kept the Confirmed. + public bool SoftPlausibilityOverridden { get; internal set; } + + /// Set for NotEvaluable(Gate). + public LipidRtGateStatistics Gate { get; internal set; } + } + + /// The outcome of the check for one object. + public sealed class LipidRtObjectResult + { + internal LipidRtObjectResult() { } + + public int Id { get; internal set; } + + public LipidRtCandidate Candidate { get; internal set; } + + /// The outcome after corroboration and the corroboration-dependent policies. + public LipidRtOutcome Outcome { get; internal set; } + public LipidRtReason Reason { get; internal set; } + + /// The graded score; null when the object was not evaluated (ineligible, or the run gate failed). + public double? Score { get; internal set; } + + /// The outcome of this unit's own evidence, steps 1-10, before any corroboration. + public LipidRtOutcome SingleRunOutcome { get; internal set; } + public LipidRtReason SingleRunReason { get; internal set; } + public double? SingleRunScore { get; internal set; } + + /// A Confirmed(Rt) member whose feature holds at least two distinct adducts. + public bool CorroboratedByAdduct { get; internal set; } + + /// The other-polarity twins Confirmed by their own evidence that corroborate or promote this object. + public IReadOnlyList CorroboratedByPolarity { get; internal set; } = Array.Empty(); + + /// The other-polarity twins, when cross-polarity corroboration ran. + public IReadOnlyList Twins { get; internal set; } = Array.Empty(); + + /// Confirmed here while every other-polarity twin is Inconsistent. + public bool PolarityConflict { get; internal set; } + + public LipidRtEvidence Evidence { get; internal set; } = new LipidRtEvidence(); + + /// The sentence a person reads, e.g. "Confirmed: TG 52:1 elutes after TG 50:1 (+0.30 min) ...". + public string Sentence => LipidRtEvidenceSentence.Format(this); + + internal LipidRtObjectResult CopyForCorroboration() { + return new LipidRtObjectResult + { + Id = Id, + Candidate = Candidate, + Outcome = SingleRunOutcome, + Reason = SingleRunReason, + Score = SingleRunScore, + SingleRunOutcome = SingleRunOutcome, + SingleRunReason = SingleRunReason, + SingleRunScore = SingleRunScore, + Evidence = Evidence, + }; + } + + /// The persisted record of this result. + public LipidRtVerdict ToVerdict(string parameterHash) { + return new LipidRtVerdict + { + Outcome = Outcome, + Reason = Reason, + Score = Score, + RuleVersion = Parameter.LipidRtCheckParameter.RuleVersion, + ParameterHash = parameterHash ?? string.Empty, + CorroboratedByAdduct = CorroboratedByAdduct, + CorroboratedByPolarity = CorroboratedByPolarity.ToArray(), + Evidence = Sentence, + }; + } + } + + /// The outcome of the check for one unit. + public sealed class LipidRtRunResult + { + internal LipidRtRunResult(LipidRtRunInput input, IReadOnlyList results, LipidRtGateStatistics gate) { + Input = input; + Results = results; + Gate = gate; + _byId = new Dictionary(); + foreach (var result in results) { + _byId[result.Id] = result; + } + } + + private readonly Dictionary _byId; + + public LipidRtRunInput Input { get; } + + /// One result per candidate, in the candidates' order. + public IReadOnlyList Results { get; } + + /// The run-level gate; null when nothing was eligible. + public LipidRtGateStatistics Gate { get; } + + public IonMode Polarity => Input.Polarity; + + public double Span => Input.Span; + + public bool TryGet(int id, out LipidRtObjectResult result) => _byId.TryGetValue(id, out result); + + public LipidRtObjectResult this[int id] => _byId[id]; + } + + /// The fitted RT line between the two polarities (step 11b). + public sealed class LipidRtCrossPolarityLine + { + internal LipidRtCrossPolarityLine(double intercept, double slope, double? residualSd, int compositions, double tolerance) { + Intercept = intercept; + Slope = slope; + ResidualSd = residualSd; + Compositions = compositions; + Tolerance = tolerance; + } + + /// RtNeg = Intercept + Slope x RtPos. + public double Intercept { get; } + public double Slope { get; } + + /// 1.4826 x MAD of the residuals, or null when fewer than 5 compositions were shared. + public double? ResidualSd { get; } + + /// Compositions with a single RT cluster in both units. + public int Compositions { get; } + + /// The twin tolerance (min). + public double Tolerance { get; } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtReviewList.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtReviewList.cs new file mode 100644 index 000000000..96f4089db --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtReviewList.cs @@ -0,0 +1,116 @@ +using System; +using System.Collections.Generic; +using System.Globalization; +using System.IO; +using System.Linq; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// The groups of the review list, in the order a person reads them. + public enum LipidRtReviewGroup + { + InconsistentOrderStrong = 0, + InconsistentOrderWeak = 1, + InconsistentPlausibility = 2, + Artefact = 3, + PolarityConflict = 4, + } + + /// One row of the review list. + public sealed class LipidRtReviewItem + { + internal LipidRtReviewItem(LipidRtReviewGroup group, LipidRtObjectResult result) { + Group = group; + Result = result; + } + + public LipidRtReviewGroup Group { get; } + public LipidRtObjectResult Result { get; } + } + + /// + /// The review list (SPEC section 10): the objects a person should look at, which are not tags. Inconsistent by RT + /// order (Strong when the deficit is at least 2 and the residual is at least one double-bond step or unknown, + /// otherwise Weak), Inconsistent by the hard plausibility rule, NotEvaluable with an MS1 artefact explanation, and + /// Confirmed with a conflicting other-polarity twin. The artefact group holds every NotEvaluable object with an + /// artefact explanation of its own, whatever rule decided it. + /// + /// Builds rows and formats them; where they are written is the caller's business. + public static class LipidRtReviewList + { + public static IReadOnlyList Build(LipidRtRunResult result) { + if (result is null) { + throw new ArgumentNullException(nameof(result)); + } + var items = new List(); + foreach (var r in result.Results) { + var group = GroupOf(r); + if (group.HasValue) { + items.Add(new LipidRtReviewItem(group.Value, r)); + } + } + return items + .OrderBy(i => i.Group) + .ThenByDescending(i => i.Group <= LipidRtReviewGroup.InconsistentOrderWeak ? i.Result.Evidence.Deficit : 0d) + .ThenBy(i => i.Result.Candidate.Rt) + .ThenBy(i => i.Result.Id) + .ToList(); + } + + private static LipidRtReviewGroup? GroupOf(LipidRtObjectResult r) { + switch (r.Reason) { + case LipidRtReason.InconsistentOrder: + var e = r.Evidence; + var strong = e.Deficit >= 2d && (!e.ResidualDoubleBondSteps.HasValue || Math.Abs(e.ResidualDoubleBondSteps.Value) >= 1d); + return strong ? LipidRtReviewGroup.InconsistentOrderStrong : LipidRtReviewGroup.InconsistentOrderWeak; + case LipidRtReason.InconsistentPlausibility: + return LipidRtReviewGroup.InconsistentPlausibility; + default: + if (r.Outcome == LipidRtOutcome.NotEvaluable && (r.Reason == LipidRtReason.Artefact || r.Evidence.Artefacts.Count > 0)) { + return LipidRtReviewGroup.Artefact; + } + if (r.PolarityConflict) { + return LipidRtReviewGroup.PolarityConflict; + } + return null; + } + } + + public static readonly string[] Columns = { + "ID", "Name", "Adduct", "RT (min)", "Group", "Outcome", "Reason", "Score", "Evidence", "Neighbours", + "Predicted RT (min)", "Residual (DB steps)", "Artefact partner", "Twins", + }; + + /// Writes the rows as tab-separated text with a header line. + public static void Write(TextWriter writer, IEnumerable items) { + writer.WriteLine(string.Join("\t", Columns)); + foreach (var item in items) { + var r = item.Result; + var c = r.Candidate; + var e = r.Evidence; + var artefact = e.Artefacts.FirstOrDefault() ?? e.FeatureArtefact; + var fields = new[] { + r.Id.ToString(CultureInfo.InvariantCulture), + string.IsNullOrEmpty(c.Molecular) ? c.Name : c.Name + "|" + c.Molecular, + c.Adduct, + c.Rt.ToString("0.000", CultureInfo.InvariantCulture), + item.Group.ToString(), + r.Outcome.ToString(), + r.Reason.ToString(), + r.Score.HasValue ? r.Score.Value.ToString("0.000", CultureInfo.InvariantCulture) : string.Empty, + r.Sentence, + string.Join("; ", e.Support.Select(n => string.Format(CultureInfo.InvariantCulture, "{0} {1} {2:+0.000;-0.000;+0.000}", n.Direction, n.Feature, n.Gap))), + e.Predicted.HasValue ? e.Predicted.Value.ToString("0.000", CultureInfo.InvariantCulture) : string.Empty, + e.ResidualDoubleBondSteps.HasValue ? e.ResidualDoubleBondSteps.Value.ToString("0.00", CultureInfo.InvariantCulture) : string.Empty, + artefact is null ? string.Empty : artefact.PartnerId.ToString(CultureInfo.InvariantCulture) + " " + artefact.Text, + string.Join("; ", r.Twins.Select(t => string.Format(CultureInfo.InvariantCulture, "{0} {1}@{2:0.000} {3}", t.Id, t.Adduct, t.Rt, t.SingleRunReason))), + }; + writer.WriteLine(string.Join("\t", fields.Select(Clean))); + } + } + + private static string Clean(string value) { + return (value ?? string.Empty).Replace('\t', ' ').Replace('\r', ' ').Replace('\n', ' '); + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtTagProtocol.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtTagProtocol.cs new file mode 100644 index 000000000..f5d3e9495 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtTagProtocol.cs @@ -0,0 +1,160 @@ +using CompMs.MsdialCore.DataObj; +using CompMs.MsdialCore.Parameter; +using System; +using System.Collections.Generic; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// What applying the check's verdicts to one object's Confirmed tag did. + public enum LipidRtTagAction + { + /// The outcome is not Confirmed; the tag was not touched. + None = 0, + /// The check selected Confirmed with its own source. + Set = 1, + /// Confirmed was already selected by a person or another tool; it stays theirs. + AlreadySelected = 2, + /// A person's Misannotation is selected; the check does not contradict it. + Misannotation = 3, + /// A person removed the check's Confirmed earlier; the check does not set it again. + Suppressed = 4, + } + + /// What applying one unit's verdicts did to its tags. + public sealed class LipidRtTagSummary + { + /// Confirmed tags this check had set before and removed before re-applying (its own only). + public int Withdrawn { get; internal set; } + + /// Confirmed tags the check set. + public int Set { get; internal set; } + + /// Confirmed outcomes where Confirmed was already selected by a person or another tool. + public int AlreadySelected { get; internal set; } + + /// Confirmed outcomes left untagged because a person's Misannotation is selected. + public int Misannotation { get; internal set; } + + /// Confirmed outcomes left untagged because a person removed the check's tag before. + public int Suppressed { get; internal set; } + + /// Objects that received a verdict. + public int Verdicts { get; internal set; } + + internal void Count(LipidRtTagAction action) { + switch (action) { + case LipidRtTagAction.Set: Set++; break; + case LipidRtTagAction.AlreadySelected: AlreadySelected++; break; + case LipidRtTagAction.Misannotation: Misannotation++; break; + case LipidRtTagAction.Suppressed: Suppressed++; break; + } + } + } + + /// + /// The tag protocol of the lipid RT check (SPEC section 9 as re-expressed in addendum B1): the one routine every + /// caller uses to put a unit's verdicts on its objects - hook A on the peaks of an analysis file, hook B on the spots + /// of an alignment, and a re-run on a loaded result. + /// + /// + /// + /// The check only ever SELECTS Confirmed (decision D1). It keeps its own record of what it set in the tag's source + /// ("MsdialLipidRtCheck/<rule version>"), so before re-applying it removes only the Confirmed tags whose source + /// is its own, of any rule version, and never a tag a person or another tool set. It then selects Confirmed only + /// where the outcome is Confirmed, Confirmed is not already selected (a person's stays a person's), Misannotation is + /// not selected (it never contradicts a person) and Confirmed is not suppressed (a person removed the check's tag + /// before). It never writes Misannotation: Inconsistent stays in the verdict and the review list. + /// + /// + /// The verdict record is replaced on every object of the unit, and cleared where the check had nothing to say, so + /// that no verdict of an earlier run outlives the annotation it was about. + /// + /// + public static class LipidRtTagProtocol + { + /// True when the source is this check's, of any rule version. + public static bool IsOwnSource(string source) { + return source != null && source.StartsWith(LipidRtCheckParameter.TagSourcePrefix, StringComparison.Ordinal); + } + + /// Removes the Confirmed tag only if this check set it. Returns true when it did. + public static bool Withdraw(PeakSpotTagCollection tags) { + if (tags is null) { + throw new ArgumentNullException(nameof(tags)); + } + var confirmed = PeakSpotTag.CONFIRMED; + if (tags.IsSelected(confirmed) && IsOwnSource(tags.GetSource(confirmed))) { + tags.Deselect(confirmed); + return true; + } + return false; + } + + /// + /// Applies one Confirmed outcome to one object's tags with the given source, unless a person's tag or decision is + /// in the way. Call first. + /// + public static LipidRtTagAction Confirm(PeakSpotTagCollection tags, string source) { + if (tags is null) { + throw new ArgumentNullException(nameof(tags)); + } + var confirmed = PeakSpotTag.CONFIRMED; + if (tags.IsSelected(confirmed)) { + return LipidRtTagAction.AlreadySelected; + } + if (tags.IsSelected(PeakSpotTag.MISANNOTATION)) { + return LipidRtTagAction.Misannotation; + } + if (tags.IsSuppressed(confirmed)) { + return LipidRtTagAction.Suppressed; + } + tags.Select(confirmed, source); + return LipidRtTagAction.Set; + } + + /// + /// The whole protocol for one unit: withdraw the check's own Confirmed tags, record each object's verdict (null + /// where the check has none), and select Confirmed where the rules above allow it. + /// + /// The unit's result; null withdraws the check's tags and clears every verdict, which is + /// what a unit gets when the check is switched off or could not run. + /// The hash of the parameter that produced . + public static LipidRtTagSummary Apply(IEnumerable objects, Func id, Func tags, + Action setVerdict, LipidRtRunResult result, string parameterHash) { + if (objects is null) throw new ArgumentNullException(nameof(objects)); + if (id is null) throw new ArgumentNullException(nameof(id)); + if (tags is null) throw new ArgumentNullException(nameof(tags)); + if (setVerdict is null) throw new ArgumentNullException(nameof(setVerdict)); + + var summary = new LipidRtTagSummary(); + var source = LipidRtCheckParameter.TagSource; + foreach (var o in objects) { + var collection = tags(o); + if (Withdraw(collection)) { + summary.Withdrawn++; + } + if (result != null && result.TryGet(id(o), out var r)) { + setVerdict(o, r.ToVerdict(parameterHash)); + summary.Verdicts++; + if (r.Outcome == LipidRtOutcome.Confirmed) { + summary.Count(Confirm(collection, source)); + } + } + else { + setVerdict(o, null); + } + } + return summary; + } + + /// on the peaks of one analysis file (hook A), keyed by MasterPeakID. + public static LipidRtTagSummary ApplyToPeaks(IEnumerable peaks, LipidRtRunResult result, string parameterHash) { + return Apply(peaks, p => p.MasterPeakID, p => p.TagCollection, (p, v) => p.LipidRtVerdict = v, result, parameterHash); + } + + /// on the spots of one alignment (hook B), keyed by MasterAlignmentID. + public static LipidRtTagSummary ApplyToSpots(IEnumerable spots, LipidRtRunResult result, string parameterHash) { + return Apply(spots, s => s.MasterAlignmentID, s => s.TagCollection, (s, v) => s.LipidRtVerdict = v, result, parameterHash); + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtVerdict.cs b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtVerdict.cs new file mode 100644 index 000000000..097480738 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Algorithm/Annotation/LipidRt/LipidRtVerdict.cs @@ -0,0 +1,50 @@ +using MessagePack; +using System; + +namespace CompMs.MsdialCore.Algorithm.Annotation.LipidRt +{ + /// + /// The persisted record of the lipid retention-time order check for one peak or alignment spot, kept for the GUI's + /// evidence sentence and the review list (SPEC section 9 with addendum B2). + /// + /// + /// The record does not own the tag. Whether the Confirmed tag was set by this check or by a person is recorded on + /// the tag itself, in the tag files (addendum B1), so this record carries no tag provenance. + /// + /// Every member's "not recorded" state is its CLR default, because MessagePack assigns default(T) to a key absent + /// from a shorter array: an older record must read as "not evaluated", never as a claim. + /// + [MessagePackObject] + public sealed class LipidRtVerdict + { + [Key(0)] + public LipidRtOutcome Outcome { get; set; } + + [Key(1)] + public LipidRtReason Reason { get; set; } + + /// Null when the object was not evaluated. + [Key(2)] + public double? Score { get; set; } + + [Key(3)] + public string RuleVersion { get; set; } + + [Key(4)] + public string ParameterHash { get; set; } + + [Key(5)] + public bool CorroboratedByAdduct { get; set; } + + /// The Ids of the other-polarity twins that corroborate this object. + [Key(6)] + public int[] CorroboratedByPolarity { get; set; } = Array.Empty(); + + /// The evidence sentence. + [Key(7)] + public string Evidence { get; set; } + + [IgnoreMember] + public bool IsConfirmed => Outcome == LipidRtOutcome.Confirmed; + } +} diff --git a/src/MSDIAL5/MsdialCore/DataObj/AlignmentResultContainer.cs b/src/MSDIAL5/MsdialCore/DataObj/AlignmentResultContainer.cs index 555dad308..7f2efa7b1 100644 --- a/src/MSDIAL5/MsdialCore/DataObj/AlignmentResultContainer.cs +++ b/src/MSDIAL5/MsdialCore/DataObj/AlignmentResultContainer.cs @@ -50,24 +50,7 @@ public void Save(AlignmentFileBean file) { } var tagfile = Path.Combine(Path.GetDirectoryName(file.FilePath), Path.GetFileNameWithoutExtension(file.FilePath) + "_tags.xml"); - var defsElement = new XElement("Definitions"); - foreach (var tag in PeakSpotTag.AllTypes()) { - var tagElement = new XElement("Tag", new XElement("Id", tag.Id), new XElement("Label", tag.Label)); - defsElement.Add(tagElement); - } - var peaksElement = new XElement("Peaks"); - foreach (var peak in AlignmentSpotProperties.SelectMany(Flatten)) { - if (!peak.TagCollection.Any()) { - continue; - } - var peakElement = new XElement("Peak"); - peakElement.SetAttributeValue("Id", peak.MasterAlignmentID); - foreach (var tag in peak.TagCollection.Selected) { - peakElement.Add(new XElement("Tag", tag.Id)); - } - peaksElement.Add(peakElement); - } - var doc = new XElement("PeakSpotTags", defsElement, peaksElement); + var doc = PeakSpotTagXml.Create(AlignmentSpotProperties.SelectMany(Flatten), spot => spot.MasterAlignmentID, spot => spot.TagCollection); doc.Save(tagfile); } @@ -103,21 +86,7 @@ public static AlignmentResultContainer Load(AlignmentFileBean file) var tagfile = Path.Combine(Path.GetDirectoryName(file.FilePath), Path.GetFileNameWithoutExtension(file.FilePath) + "_tags.xml"); if (File.Exists(tagfile)) { var doc = XElement.Load(tagfile); - var dic = new Dictionary(); - foreach (var peakElement in doc.Descendants("Peak")) { - if (int.TryParse(peakElement.Attribute("Id").Value, out var id)) { - var tags = peakElement.Elements("Tag") - .Select(elem => int.TryParse(elem.Value, out var v) ? PeakSpotTag.GetById(v) : null) - .Where(v => v != null) - .ToArray(); - dic[id] = tags; - } - } - foreach (var peak in collection.SelectMany(Flatten)) { - if (dic.TryGetValue(peak.MasterAlignmentID, out var tags)) { - peak.TagCollection.Set(tags); - } - } + PeakSpotTagXml.Apply(doc, collection.SelectMany(Flatten), spot => spot.MasterAlignmentID, spot => spot.TagCollection); } return result; } @@ -170,21 +139,7 @@ public static AlignmentResultContainer LoadLazy(AlignmentFileBean file, Cancella var tagfile = Path.Combine(Path.GetDirectoryName(file.FilePath), Path.GetFileNameWithoutExtension(file.FilePath) + "_tags.xml"); if (File.Exists(tagfile)) { var doc = XElement.Load(tagfile); - var dic = new Dictionary(); - foreach (var peakElement in doc.Descendants("Peak")) { - if (int.TryParse(peakElement.Attribute("Id").Value, out var id)) { - var tags = peakElement.Elements("Tag") - .Select(elem => int.TryParse(elem.Value, out var v) ? PeakSpotTag.GetById(v) : null) - .Where(v => v != null) - .ToArray(); - dic[id] = tags; - } - } - foreach (var peak in collection.SelectMany(Flatten)) { - if (dic.TryGetValue(peak.MasterAlignmentID, out var tags)) { - peak.TagCollection.Set(tags); - } - } + PeakSpotTagXml.Apply(doc, collection.SelectMany(Flatten), spot => spot.MasterAlignmentID, spot => spot.TagCollection); } return result; } diff --git a/src/MSDIAL5/MsdialCore/DataObj/AlignmentSpotProperty.cs b/src/MSDIAL5/MsdialCore/DataObj/AlignmentSpotProperty.cs index 0e2169479..dbcf34805 100644 --- a/src/MSDIAL5/MsdialCore/DataObj/AlignmentSpotProperty.cs +++ b/src/MSDIAL5/MsdialCore/DataObj/AlignmentSpotProperty.cs @@ -6,6 +6,7 @@ using CompMs.Common.Extension; using CompMs.Common.Interfaces; using CompMs.MsdialCore.Algorithm.Annotation; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using CompMs.MsdialCore.Normalize; using MessagePack; using System; @@ -324,6 +325,19 @@ public int GetMSDecResultID() { [IgnoreMember] public PeakSpotTagCollection TagCollection { get; } = new PeakSpotTagCollection(); + /// + /// The lipid retention-time order check's verdict on this spot, for the evidence sentence and the exports; null when + /// the check did not evaluate it. + /// + /// + /// The next key after the post-curation keys (58, 62-65). An alignment written before this key reads it as null, + /// because a key absent from the stored array is assigned default(T). The Confirmed tag itself is not here: it is + /// in , with the check as its source. copies neither, because a + /// clone is a new spot whose annotation a person is about to change. + /// + [Key(66)] + public LipidRtVerdict LipidRtVerdict { get; set; } + public AlignmentSpotProperty Clone(ref int masterId, int alignmentId) { var driftSpots = new List(AlignmentDriftSpotFeatures.Count); foreach (var driftSpot in AlignmentDriftSpotFeatures) { diff --git a/src/MSDIAL5/MsdialCore/DataObj/ChromatogramPeakFeature.cs b/src/MSDIAL5/MsdialCore/DataObj/ChromatogramPeakFeature.cs index cfd92d29e..3334a1fd6 100644 --- a/src/MSDIAL5/MsdialCore/DataObj/ChromatogramPeakFeature.cs +++ b/src/MSDIAL5/MsdialCore/DataObj/ChromatogramPeakFeature.cs @@ -6,6 +6,7 @@ using CompMs.Common.Extension; using CompMs.Common.Interfaces; using CompMs.MsdialCore.Algorithm.Annotation; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using CompMs.MsdialCore.Utility; using MessagePack; using System; @@ -113,6 +114,19 @@ public int GetMSDecResultID() { public int MS1AccumulatedMs1RawSpectrumIdLeft { get; set; } // used for LC-IM-MS/MS [Key(47)] public int MS1AccumulatedMs1RawSpectrumIdRight { get; set; } // used for LC-IM-MS/MS + + // Keys 54 and 55 belong to the MS-DIAL 5.6 line (tags MSDIAL-v5.6.*), where they are RawDataIDType and + // AccumulatedDataIDType: a SpectrumIDType enum (0 = Index, 1 = Raw) saying what kind of ID the raw spectrum IDs + // above are. This line has no such enum, and its IDs are always indices. The two keys are kept here as plain + // integers so that a peak file of either line opens in the other: a 5.6 file keeps its values when it is read and + // saved here, and a file written here holds integers (0 = Index, which is what its IDs are) where a 5.6 reader + // reads the enum, not the nil MessagePack writes for a missing key, which that reader cannot read. When the 5.6 + // line is merged, replace these two with its members; never give these keys to anything else. + [Key(54)] + public int ReservedRawDataIDType { get; set; } + [Key(55)] + public int ReservedAccumulatedDataIDType { get; set; } + [Key(18)] public int MS2RawSpectrumID { get; set; } = -1; // representative ID [Key(19)] @@ -342,6 +356,19 @@ public bool AllDriftFeaturesAreAnnotated(IMatchResultEvaluator + /// The lipid retention-time order check's verdict on this peak, from the check on its own analysis file; null when + /// the check did not evaluate it. + /// + /// + /// Key 56 is the next key after the highest one in use on any branch or release tag: 53 here, and 54-55 in the + /// MS-DIAL 5.6 line (). The unused numbers below it are left alone because an + /// older file may hold something else there. A peak file written before this key, by this line or by 5.6, reads it + /// as null. The Confirmed tag itself is in , with the check as its source. + /// + [Key(56)] + public LipidRtVerdict LipidRtVerdict { get; set; } + public static ChromatogramPeakFeature FromPeakDetectionResult(PeakDetectionResult peakDetectionResult, Chromatogram chromatogram, double mz) { if (peakDetectionResult == null) { return null; diff --git a/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagCollection.cs b/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagCollection.cs index 80c7c3590..c6a7dea21 100644 --- a/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagCollection.cs +++ b/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagCollection.cs @@ -4,18 +4,52 @@ namespace CompMs.MsdialCore.DataObj { + /// + /// The tags of one peak or alignment spot, with who set each tag and which tags a person removed after a tool had + /// set them. + /// + /// + /// + /// PROVENANCE. A selected tag has an optional source. Null means a person set it, or that nobody recorded who did: + /// and are the person's path and keep the + /// meaning they always had. A tool sets a tag with and its own source, e.g. + /// "MsdialLipidRtCheck/1.0" or "Agent/<name>/<version>", so that it can later undo only what it set. + /// + /// + /// SUPPRESSION. When a person removes a tag that a tool set, the tag is also suppressed: no tool may set it again on + /// that object. A person turning the tag on again lifts the suppression. Suppression is kept apart from selection, + /// so a suppressed tag is simply one that is not selected and must not be selected by a tool. + /// + /// + /// Both live in the "_tags.xml" files beside the peak and alignment results (), which + /// older MS-DIAL versions read as the same set of selected tags. + /// + /// public sealed class PeakSpotTagCollection { + /// The source recorded on a suppression a person made through MS-DIAL. + public const string PersonSource = "User"; + private readonly List _selected; + private readonly Dictionary _sources; + private readonly List _suppressed; + private readonly Dictionary _suppressionSources; public PeakSpotTagCollection() { _selected = new List(); + _sources = new Dictionary(); + _suppressed = new List(); + _suppressionSources = new Dictionary(); Selected = _selected.AsReadOnly(); + Suppressed = _suppressed.AsReadOnly(); } public ReadOnlyCollection Selected { get; } + /// The tags a person removed after a tool had set them. No tool may select these again. + public ReadOnlyCollection Suppressed { get; } + public bool IsSelected(PeakSpotTag type) { return _selected.Contains(type); } @@ -24,27 +58,159 @@ public bool IsSelected(PeakSpotTagSearchQuery query) { return query.IsMatched(_selected); } + /// + /// Selects a tag as a person, or as a caller that does not say who it is. A tag that is already selected is left + /// as it is, its source included. + /// public void Select(PeakSpotTag type) { if (!_selected.Contains(type)) { _selected.Add(type); } } - internal void Set(PeakSpotTag[] types) { - _selected.Clear(); - _selected.AddRange(types); + /// + /// Selects a tag on behalf of a tool. A tag that is already selected is left as it is, its source included, so a + /// person's tag stays a person's. A null or empty source is the person's path, . + /// + /// + /// This does not consult : honouring a suppression is the tool's part, and + /// a tool that ignores it is the one repeating what a person undid. + /// + public void Select(PeakSpotTag type, string source) { + if (string.IsNullOrEmpty(source)) { + Select(type); + return; + } + if (type is null) { + return; + } + if (!_selected.Contains(type)) { + _selected.Add(type); + _sources[type] = source; + } + } + + /// + /// The source that set a selected tag: null when a person set it, when nobody recorded who did, or when the tag is + /// not selected. + /// + public string GetSource(PeakSpotTag type) { + return type != null && _sources.TryGetValue(type, out var source) ? source : null; } public void Deselect(PeakSpotTag type) { if (_selected.Contains(type)) { _selected.Remove(type); } + if (type != null) { + _sources.Remove(type); + } + } + + public bool IsSuppressed(PeakSpotTag type) { + return _suppressed.Contains(type); + } + + /// Records that a person removed this tag, so that no tool sets it again. + public void Suppress(PeakSpotTag type) { + Suppress(type, PersonSource); + } + + /// Records that removed this tag, so that no tool sets it again. + public void Suppress(PeakSpotTag type, string source) { + if (type is null) { + return; + } + if (!_suppressed.Contains(type)) { + _suppressed.Add(type); + } + if (string.IsNullOrEmpty(source)) { + _suppressionSources.Remove(type); + } + else { + _suppressionSources[type] = source; + } } + public void Unsuppress(PeakSpotTag type) { + if (type is null) { + return; + } + _suppressed.Remove(type); + _suppressionSources.Remove(type); + } + + /// Who suppressed the tag; null when it is not suppressed or nobody recorded who did. + public string GetSuppressionSource(PeakSpotTag type) { + return type != null && _suppressionSources.TryGetValue(type, out var source) ? source : null; + } + + /// + /// A person turns a tag on: it is selected as the person's and any suppression of it is lifted. + /// + /// + /// A tag that is already on because a tool set it becomes the person's: it keeps its place among the selected tags + /// and loses its source, so the tool no longer counts it as its own and will not withdraw it. "Mark all as + /// confirmed" relies on this; the GUI's single toggles never call this for a tag that is already selected. + /// + public void SelectByPerson(PeakSpotTag type) { + if (type is null) { + return; + } + Select(type); + _sources.Remove(type); + Unsuppress(type); + } + + /// + /// A person turns a tag off. When a tool had set it, the tag is also suppressed so that the tool does not set it + /// again; a tag a person set is simply removed. + /// + /// + /// A suppression records that a person removed a tool's tag (addendum B1). A tag the person turned on again is + /// the person's, and turning it off once more removes the person's own tag, which is not suppressed: a tool may + /// then set it again. + /// + public void DeselectByPerson(PeakSpotTag type) { + var setByTool = IsSelected(type) && GetSource(type) != null; + Deselect(type); + if (setByTool) { + Suppress(type); + } + } + + /// True when any tag is selected. Suppressions alone do not count. public bool Any() { return _selected.Count > 0; } + /// True when there is anything for the tag file: a selected tag or a suppression. + public bool HasEntries() { + return _selected.Count > 0 || _suppressed.Count > 0; + } + + internal void Set(PeakSpotTag[] types) { + Set(types.Select(type => (type, (string)null)), Enumerable.Empty<(PeakSpotTag, string)>()); + } + + /// Replaces everything with what a tag file recorded. + internal void Set(IEnumerable<(PeakSpotTag Tag, string Source)> selected, IEnumerable<(PeakSpotTag Tag, string Source)> suppressed) { + _selected.Clear(); + _sources.Clear(); + _suppressed.Clear(); + _suppressionSources.Clear(); + foreach (var (tag, source) in selected) { + if (tag != null) { + Select(tag, source); + } + } + foreach (var (tag, source) in suppressed) { + if (tag != null) { + Suppress(tag, source); + } + } + } + public override string ToString() { return string.Join("|", _selected.Select(t => t.Label)); } diff --git a/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagXml.cs b/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagXml.cs new file mode 100644 index 000000000..255723f4f --- /dev/null +++ b/src/MSDIAL5/MsdialCore/DataObj/PeakSpotTagXml.cs @@ -0,0 +1,129 @@ +using System; +using System.Collections.Generic; +using System.Linq; +using System.Xml.Linq; + +namespace CompMs.MsdialCore.DataObj +{ + /// + /// The "_tags.xml" file written beside a peak result (.pai2), an alignment result and a spectrum feature result: + /// the tags of every peak, who set each tag, and the tags a person removed after a tool had set them. + /// + /// + /// + /// The file is the human-readable entry point to the tags: the GUI reloads it, and external tools read and write + /// tags through it. It looks like this: + /// + /// <PeakSpotTags> + /// <Definitions><Tag><Id>1</Id><Label>Confirmed</Label></Tag>...</Definitions> + /// <Peaks> + /// <Peak Id="123"> + /// <Tag>3</Tag> (a person's tag: no source) + /// <Tag Source="MsdialLipidRtCheck/1.0">1</Tag> (set by a tool) + /// <Suppressed Source="User">1</Suppressed> (a person removed a tool's tag; no tool sets it again) + /// </Peak> + /// </Peaks> + /// </PeakSpotTags> + /// + /// + /// + /// BACKWARD COMPATIBLE BOTH WAYS. MS-DIAL before tag provenance reads only the value of each Tag element of a + /// Peak, so it ignores the Source attribute and the Suppressed elements and reads the same set of + /// selected tags. A file it wrote has neither, and reads here as tags set by a person or by nobody recorded, with + /// nothing suppressed. A tag without a source is written exactly as it always was. + /// + /// + /// A peak with only suppressions is written too, because the suppression is what stops a tool from setting the tag + /// again; an older reader sees it as a peak without tags. + /// + /// + public static class PeakSpotTagXml + { + public const string SourceAttribute = "Source"; + + /// The tag file for the given objects; only objects with a tag or a suppression get a Peak element. + public static XElement Create(IEnumerable items, Func id, Func tags) { + var defsElement = new XElement("Definitions"); + foreach (var tag in PeakSpotTag.AllTypes()) { + var tagElement = new XElement("Tag", new XElement("Id", tag.Id), new XElement("Label", tag.Label)); + defsElement.Add(tagElement); + } + var peaksElement = new XElement("Peaks"); + foreach (var item in items) { + var collection = tags(item); + if (collection is null || !collection.HasEntries()) { + continue; + } + var peakElement = new XElement("Peak"); + peakElement.SetAttributeValue("Id", id(item)); + foreach (var tag in collection.Selected) { + var tagElement = new XElement("Tag", tag.Id); + var source = collection.GetSource(tag); + if (!string.IsNullOrEmpty(source)) { + tagElement.SetAttributeValue(SourceAttribute, source); + } + peakElement.Add(tagElement); + } + foreach (var tag in collection.Suppressed) { + var suppressedElement = new XElement("Suppressed", tag.Id); + var source = collection.GetSuppressionSource(tag); + if (!string.IsNullOrEmpty(source)) { + suppressedElement.SetAttributeValue(SourceAttribute, source); + } + peakElement.Add(suppressedElement); + } + peaksElement.Add(peakElement); + } + return new XElement("PeakSpotTags", defsElement, peaksElement); + } + + /// + /// Sets the tags the file records on the given objects. Objects the file does not name are left as they are, and + /// a tag id this version does not know is skipped, as before. + /// + public static void Apply(XElement document, IEnumerable items, Func id, Func tags) { + var records = Read(document); + if (records.Count == 0) { + return; + } + foreach (var item in items) { + if (records.TryGetValue(id(item), out var record)) { + tags(item).Set(record.Selected, record.Suppressed); + } + } + } + + internal sealed class Record + { + public List<(PeakSpotTag Tag, string Source)> Selected { get; } = new List<(PeakSpotTag, string)>(); + public List<(PeakSpotTag Tag, string Source)> Suppressed { get; } = new List<(PeakSpotTag, string)>(); + } + + internal static Dictionary Read(XElement document) { + var records = new Dictionary(); + foreach (var peakElement in document.Descendants("Peak")) { + if (!int.TryParse(peakElement.Attribute("Id")?.Value, out var id)) { + continue; + } + var record = new Record(); + foreach (var element in peakElement.Elements("Tag")) { + if (TryGetTag(element, out var tag)) { + record.Selected.Add((tag, element.Attribute(SourceAttribute)?.Value)); + } + } + foreach (var element in peakElement.Elements("Suppressed")) { + if (TryGetTag(element, out var tag)) { + record.Suppressed.Add((tag, element.Attribute(SourceAttribute)?.Value)); + } + } + records[id] = record; + } + return records; + } + + private static bool TryGetTag(XElement element, out PeakSpotTag tag) { + tag = int.TryParse(element.Value, out var value) ? PeakSpotTag.GetById(value) : null; + return tag != null; + } + } +} diff --git a/src/MSDIAL5/MsdialCore/DataObj/SpectrumFeatureCollection.cs b/src/MSDIAL5/MsdialCore/DataObj/SpectrumFeatureCollection.cs index 6c790b34d..b9065c309 100644 --- a/src/MSDIAL5/MsdialCore/DataObj/SpectrumFeatureCollection.cs +++ b/src/MSDIAL5/MsdialCore/DataObj/SpectrumFeatureCollection.cs @@ -23,49 +23,21 @@ public void Save(Stream stream) { if (stream is FileStream fs) { var file = fs.Name; var tagfile = Path.Combine(Path.GetDirectoryName(file), Path.GetFileNameWithoutExtension(file) + "_tags.xml"); - var defsElement = new XElement("Definitions"); - foreach (var tag in PeakSpotTag.AllTypes()) { - var tagElement = new XElement("Tag", new XElement("Id", tag.Id), new XElement("Label", tag.Label)); - defsElement.Add(tagElement); - } - var peaksElement = new XElement("Peaks"); - foreach (var feature in _spectrumFeatures) { - if (!feature.TagCollection.Any()) { - continue; - } - var peakElement = new XElement("Peak"); - peakElement.SetAttributeValue("Id", feature.AnnotatedMSDecResult.MSDecResult.ScanID); - foreach (var tag in feature.TagCollection.Selected) { - peakElement.Add(new XElement("Tag", tag.Id)); - } - peaksElement.Add(peakElement); - } - var doc = new XElement("PeakSpotTags", defsElement, peaksElement); - doc.Save(tagfile); + CreateTagDocument().Save(tagfile); } } public void Save(Stream stream, Stream tagStream) { LargeListMessagePack.Serialize(stream, _spectrumFeatures); - var defsElement = new XElement("Definitions"); - foreach (var tag in PeakSpotTag.AllTypes()) { - var tagElement = new XElement("Tag", new XElement("Id", tag.Id), new XElement("Label", tag.Label)); - defsElement.Add(tagElement); - } - var peaksElement = new XElement("Peaks"); - foreach (var feature in _spectrumFeatures) { - if (!feature.TagCollection.Any()) { - continue; - } - var peakElement = new XElement("Peak"); - peakElement.SetAttributeValue("Id", feature.AnnotatedMSDecResult.MSDecResult.ScanID); - foreach (var tag in feature.TagCollection.Selected) { - peakElement.Add(new XElement("Tag", tag.Id)); - } - peaksElement.Add(peakElement); - } - var doc = new XElement("PeakSpotTags", defsElement, peaksElement); - doc.Save(tagStream); + CreateTagDocument().Save(tagStream); + } + + private XElement CreateTagDocument() { + return PeakSpotTagXml.Create(_spectrumFeatures, feature => feature.AnnotatedMSDecResult.MSDecResult.ScanID, feature => feature.TagCollection); + } + + private static void ApplyTagDocument(XElement doc, SpectrumFeatureCollection features) { + PeakSpotTagXml.Apply(doc, features._spectrumFeatures, feature => feature.AnnotatedMSDecResult.MSDecResult.ScanID, feature => feature.TagCollection); } public static SpectrumFeatureCollection Load(Stream stream) { @@ -74,22 +46,7 @@ public static SpectrumFeatureCollection Load(Stream stream) { var file = fs.Name; var tagfile = Path.Combine(Path.GetDirectoryName(file), Path.GetFileNameWithoutExtension(file) + "_tags.xml"); if (File.Exists(tagfile)) { - var doc = XElement.Load(tagfile); - var dic = new Dictionary(); - foreach (var peakElement in doc.Descendants("Peak")) { - if (int.TryParse(peakElement.Attribute("Id").Value, out var id)) { - var tags = peakElement.Elements("Tag") - .Select(elem => int.TryParse(elem.Value, out var v) ? PeakSpotTag.GetById(v) : null) - .Where(v => v != null) - .ToArray(); - dic[id] = tags; - } - } - foreach (var feature in features._spectrumFeatures) { - if (dic.TryGetValue(feature.AnnotatedMSDecResult.MSDecResult.ScanID, out var tags)) { - feature.TagCollection.Set(tags); - } - } + ApplyTagDocument(XElement.Load(tagfile), features); } } return features; @@ -98,22 +55,7 @@ public static SpectrumFeatureCollection Load(Stream stream) { public static SpectrumFeatureCollection Load(Stream stream, Stream tagStream) { var features = new SpectrumFeatureCollection(LargeListMessagePack.Deserialize(stream)); if (!(tagStream is null)) { - var doc = XElement.Load(tagStream); - var dic = new Dictionary(); - foreach (var peakElement in doc.Descendants("Peak")) { - if (int.TryParse(peakElement.Attribute("Id").Value, out var id)) { - var tags = peakElement.Elements("Tag") - .Select(elem => int.TryParse(elem.Value, out var v) ? PeakSpotTag.GetById(v) : null) - .Where(v => v != null) - .ToArray(); - dic[id] = tags; - } - } - foreach (var feature in features._spectrumFeatures) { - if (dic.TryGetValue(feature.AnnotatedMSDecResult.MSDecResult.ScanID, out var tags)) { - feature.TagCollection.Set(tags); - } - } + ApplyTagDocument(XElement.Load(tagStream), features); } return features; } diff --git a/src/MSDIAL5/MsdialCore/Export/MztabFormatExport.cs b/src/MSDIAL5/MsdialCore/Export/MztabFormatExport.cs index 4b2f5994b..145fe48d5 100644 --- a/src/MSDIAL5/MsdialCore/Export/MztabFormatExport.cs +++ b/src/MSDIAL5/MsdialCore/Export/MztabFormatExport.cs @@ -8,6 +8,7 @@ using CompMs.MsdialCore.Parameter; using System; using System.Collections.Generic; +using System.Globalization; using System.IO; using System.Linq; using System.Text; @@ -110,7 +111,8 @@ string outfile //SML Header var hasComment = spots.Any(s => !string.IsNullOrEmpty(s.Comment)); var hasMs2 = spots.Any(n => n.IsMsmsAssigned); - var SmlDataHeader = WriteSmlHeader(sw, meta, RawFileMetadataDic, AnalysisFileClassDic,hasComment,hasMs2); + var hasLipidRt = HasLipidRtVerdicts(spots); + var SmlDataHeader = WriteSmlHeader(sw, meta, RawFileMetadataDic, AnalysisFileClassDic,hasComment,hasMs2, hasLipidRt); //SML data foreach (var spot in spots) @@ -118,13 +120,13 @@ string outfile var metadata = metaAccessor.GetContent(spot, msdecResults[spot.MasterAlignmentID]); WriteSmlDataLine( sw, spot, meta, metadata, quantAccessor, stats, RawFileMetadataDic, AnalysisFileClassDic, - database, SmlDataHeader, internalStandardDic, hasComment,hasMs2 + database, SmlDataHeader, internalStandardDic, hasComment,hasMs2, hasLipidRt ); foreach (var driftSpot in spot.AlignmentDriftSpotFeatures ?? Enumerable.Empty()) { WriteSmlDataLine( sw, driftSpot, meta, metadata, quantAccessor, stats, RawFileMetadataDic, AnalysisFileClassDic, - database, SmlDataHeader, internalStandardDic, hasComment, hasMs2 + database, SmlDataHeader, internalStandardDic, hasComment, hasMs2, hasLipidRt ); } } @@ -217,7 +219,8 @@ string outfile var hasComment = spots.Any(s => !string.IsNullOrEmpty(s.Comment)); var hasMs2 = spots.Any(n => n.IsMsmsAssigned); - var SmlDataHeader = WriteSmlHeader(sw, meta, RawFileMetadataDic, AnalysisFileClassDic, hasComment, hasMs2); + var hasLipidRt = HasLipidRtVerdicts(spots); + var SmlDataHeader = WriteSmlHeader(sw, meta, RawFileMetadataDic, AnalysisFileClassDic, hasComment, hasMs2, hasLipidRt); var smfTemp = Path.GetTempFileName(); var smeTemp = Path.GetTempFileName(); @@ -235,7 +238,7 @@ string outfile var metadata = metaAccessor.GetContent(spot, msdec); WriteSmlDataLine( sw, spot, meta, metadata, quantAccessor, stats, RawFileMetadataDic, AnalysisFileClassDic, - database, SmlDataHeader, internalStandardDic, hasComment, hasMs2 + database, SmlDataHeader, internalStandardDic, hasComment, hasMs2, hasLipidRt ); WriteSmfDataLine( smfWriter, spot, meta, quantAccessor, stats, RawFileMetadataDic, AnalysisFileClassDic, @@ -246,7 +249,7 @@ string outfile { WriteSmlDataLine( sw, driftSpot, meta, metadata, quantAccessor, stats, RawFileMetadataDic, AnalysisFileClassDic, - database, SmlDataHeader, internalStandardDic, hasComment, hasMs2 + database, SmlDataHeader, internalStandardDic, hasComment, hasMs2, hasLipidRt ); WriteSmfDataLine( smfWriter, driftSpot, meta, quantAccessor, stats, RawFileMetadataDic, AnalysisFileClassDic, @@ -424,7 +427,8 @@ private void WriteSmlDataLine( IReadOnlyList SmlDataHeader, IReadOnlyDictionary internalStandardDic, bool hasComment, - bool hasMs2 + bool hasMs2, + bool hasLipidRt ) { var matchResult = spot.MatchResults.Representative; @@ -537,6 +541,10 @@ bool hasMs2 LineData.Add(string.IsNullOrEmpty(spot.Comment) ? "null" : spot.Comment); } LineData.AddRange(SetOntology(spot)); // add 20260127 + if (hasLipidRt) + { + LineData.AddRange(SetLipidRtVerdict(spot)); + } sw.WriteLine(string.Join(Separator, LineMetaData) + Separator + string.Join(Separator, LineData)); @@ -1072,7 +1080,7 @@ IReadOnlyList database } private List WriteSmlHeader(StreamWriter sw, ParameterBase meta, IReadOnlyDictionary RawFileMetadataDic, - Dictionary AnalysisFileClassDic, bool hasComment, bool hasMs2) + Dictionary AnalysisFileClassDic, bool hasComment, bool hasMs2, bool hasLipidRt) { var SmlHeaderMeta = new List() @@ -1124,6 +1132,14 @@ private List WriteSmlHeader(StreamWriter sw, ParameterBase meta, IReadOn SmlDataHeader.Add("opt_global_user_comment"); } SmlDataHeader.Add("opt_global_Ontology"); + // The lipid RT check's verdict, after the ontology so that no existing column moves. Written when the check + // left a verdict on any spot; a spot it did not evaluate reads "null". + if (hasLipidRt) + { + SmlDataHeader.Add("opt_global_lipid_rt_outcome"); + SmlDataHeader.Add("opt_global_lipid_rt_reason"); + SmlDataHeader.Add("opt_global_lipid_rt_score"); + } sw.WriteLine(string.Join(Separator, SmlHeaderMeta) + Separator + string.Join(Separator, SmlDataHeader)); return SmlDataHeader; @@ -1264,6 +1280,32 @@ AlignmentSpotProperty spot { return new List() { ValueOrNull(spot.Ontology) }; } + + private static bool HasLipidRtVerdicts(IReadOnlyList spots) + { + return spots.Any(s => s.LipidRtVerdict != null + || (s.AlignmentDriftSpotFeatures?.Any(d => d.LipidRtVerdict != null) ?? false)); + } + + /// + /// opt_global_lipid_rt_outcome, _reason and _score: the outcome and reason names of the lipid RT check + /// (e.g. "Confirmed", "ConfirmedRt") and its graded score, "null" where the spot was not evaluated. + /// + private static List SetLipidRtVerdict( + AlignmentSpotProperty spot + ) + { + var verdict = spot.LipidRtVerdict; + if (verdict is null) + { + return new List() { "null", "null", "null" }; + } + return new List() { + verdict.Outcome.ToString(), + verdict.Reason.ToString(), + verdict.Score is double score ? score.ToString("R", CultureInfo.InvariantCulture) : "null", + }; + } private static IReadOnlyDictionary SetStandardDic( IReadOnlyList spots ) diff --git a/src/MSDIAL5/MsdialCore/Parameter/LipidRtCheckParameter.cs b/src/MSDIAL5/MsdialCore/Parameter/LipidRtCheckParameter.cs new file mode 100644 index 000000000..cf854b677 --- /dev/null +++ b/src/MSDIAL5/MsdialCore/Parameter/LipidRtCheckParameter.cs @@ -0,0 +1,439 @@ +using CompMs.Common.Enum; +using MessagePack; +using System; +using System.Collections.Generic; +using System.Globalization; +using System.Linq; +using System.Security.Cryptography; +using System.Text; + +namespace CompMs.MsdialCore.Parameter +{ + /// + /// The precision/recall knob of the lipid retention-time order check (decision D2). + /// + public enum LipidRtOperatingPoint + { + /// HighPrecision, and a Confirmed must also be corroborated by a second adduct or the other polarity. + Strict = 0, + /// Bracket support: lower and upper neighbours on at least two axes. The default. + HighPrecision = 1, + /// Bracket support or two axes. + Mid = 2, + /// One supporting direction. + Balanced = 3, + } + + /// + /// What the analyst says about the chromatography. ParameterBase has no reversed-phase/HILIC field, so the + /// check also tests the data itself (the run-level gate) unless the analyst says the run is not reversed phase. + /// + public enum LipidRtChromatography + { + /// The data gate decides. + Auto = 0, + /// The analyst says reversed phase; the data gate is still applied as a safety net. + ReversedPhase = 1, + /// HILIC, SFC, direct infusion or anything else: nothing is evaluated. + Other = 2, + } + + /// + /// The biological plausibility profile, a layer separate from the retention-time evidence. + /// + /// + /// A profile holds chain-level rules only. It never excludes or penalises a lipid class for an organism + /// (the author's decision of 2026-10-15: animals synthesise MGDG and DGDG, and diet and gut bacteria bring + /// other classes). + /// + public enum LipidRtPlausibilityProfileKind + { + None = 0, + /// B1 (no fatty chain with 7 or more double bonds, hard), B2 and B3 (sphingolipid N-acyl and base, soft). + Mammalian = 1, + } + + /// What happens to a Confirmed member whose own match is class level only. + public enum LipidRtClassLevelPolicy + { + /// A class-level member is Confirmed with its feature (D3 literal). + Consistent = 0, + /// Confirmed only with a second adduct of the same feature or a Confirmed other-polarity twin. The default. + Corroborated = 1, + /// Never Confirmed, but still used as order evidence. + Never = 2, + } + + /// A less intense feature with the same species name (and the same or an unresolved chain name) at another RT. + public enum LipidRtMinorDuplicatePolicy + { + Allow = 0, + Block = 1, + } + + /// A feature made only of secondary ions (Na, K, Li, multiply charged). + public enum LipidRtSecondaryOnlyPolicy + { + Evaluate = 0, + Block = 1, + } + + /// A Confirmed whose other-polarity twins are all Inconsistent. + public enum LipidRtPolarityConflictPolicy + { + /// Stays Confirmed; the conflict is recorded. + Flag = 0, + /// Becomes NotEvaluable(PolarityConflict). + Demote = 1, + } + + /// + /// The MS1 artefact relations the check looks for between a candidate and a brighter co-eluting peak. + /// + [Flags] + public enum LipidRtArtefactKinds + { + None = 0, + Iso1 = 1 << 0, + Iso2 = 1 << 1, + NaH = 1 << 2, + NaNH4 = 1 << 3, + KH = 1 << 4, + KNa = 1 << 5, + NH4H = 1 << 6, + Acn = 1 << 7, + Water = 1 << 8, + Acetate = 1 << 9, + Formate = 1 << 10, + Chloride = 1 << 11, + Demethyl = 1 << 12, + IsfChain = 1 << 13, + IsfHead = 1 << 14, + Ox = 1 << 15, + /// Every kind except the water loss and the acetonitrile adduct, which were not specific on the development data. + Default = Iso1 | Iso2 | NaH | NaNH4 | KH | KNa | NH4H | Acetate | Formate | Chloride | Demethyl | IsfChain | IsfHead | Ox, + } + + /// + /// Parameters of the lipid retention-time order check (specification v1.0, sections 5 and 6, with the addendum of + /// 2026-10-15). + /// + /// + /// The first group is user facing; the rest are expert constants that are serialized but not shown. Every + /// window is a fraction of the run's RT span, so the same values work on a 20-minute and a 6.5-minute gradient. + /// + /// Every member's "not set" state must be what an older project should mean, because MessagePack assigns + /// default(T) to a key that is absent from a shorter array. That is why is nullable: + /// null means "the default for the project's target omics", see . The + /// other members are only read once the whole object exists, and an object written by this version carries + /// every key. + /// + [MessagePackObject] + public sealed class LipidRtCheckParameter + { + /// The version of the rules. Part of every verdict and of the tag source string. + public const string RuleVersion = "1.0"; + + /// The source written on a tag this check sets, without the version. + public const string TagSourcePrefix = "MsdialLipidRtCheck/"; + + /// The source written on a tag this check sets. + [IgnoreMember] + public static string TagSource => TagSourcePrefix + RuleVersion; + + // ---- user facing ---- + + /// Null means "on for lipidomics projects, off otherwise". + [Key(0)] + public bool? IsEnabled { get; set; } + + [Key(1)] + public LipidRtOperatingPoint OperatingPoint { get; set; } = LipidRtOperatingPoint.HighPrecision; + + [Key(2)] + public LipidRtChromatography Chromatography { get; set; } = LipidRtChromatography.Auto; + + [Key(3)] + public LipidRtPlausibilityProfileKind PlausibilityProfile { get; set; } = LipidRtPlausibilityProfileKind.None; + + [Key(4)] + public LipidRtClassLevelPolicy ClassLevelPolicy { get; set; } = LipidRtClassLevelPolicy.Corroborated; + + [Key(5)] + public LipidRtMinorDuplicatePolicy MinorDuplicatePolicy { get; set; } = LipidRtMinorDuplicatePolicy.Allow; + + [Key(6)] + public LipidRtSecondaryOnlyPolicy SecondaryOnlyPolicy { get; set; } = LipidRtSecondaryOnlyPolicy.Evaluate; + + /// Promotion and corroboration by the other polarity, when both polarities are available. + [Key(7)] + public bool CrossPolarity { get; set; } = true; + + // ---- run-scaled windows (fractions of the RT span) ---- + + /// Discordance epsilon: 0.05 min on a 20-min run. + [Key(8)] + public double DiscordanceFraction { get; set; } = 0.0025; + + /// An adjacent gap below this is a tie, not support. + [Key(9)] + public double TieFraction { get; set; } = 0.001; + + /// A supporting gap must also exceed this fraction of the class median adjacent step. + [Key(10)] + public double StepFraction { get; set; } = 0.25; + + /// Adduct duplicates within this window are one feature. + [Key(11)] + public double CollapseFraction { get; set; } = 0.0025; + + /// Co-elution window of the MS1 artefact layer. + [Key(12)] + public double CoelutionFraction { get; set; } = 0.0015; + + /// No Inconsistent in the last fraction of the span (gradient-end compression). + [Key(13)] + public double TailFraction { get; set; } = 0.05; + + // ---- feature weights (evidence units) ---- + + [Key(14)] + public double LevelWeightChain { get; set; } = 3d; + + [Key(15)] + public double LevelWeightMsms { get; set; } = 2d; + + [Key(16)] + public double LevelWeightClass { get; set; } = 1d; + + [Key(17)] + public double DotWeight { get; set; } = 1d; + + [Key(18)] + public double DuplicateWeight { get; set; } = 0.5; + + [Key(19)] + public int DuplicateCap { get; set; } = 2; + + [Key(20)] + public double HeightWeight { get; set; } = 1d; + + /// Weight multiplier of a feature whose best member is hard-implausible. + [Key(21)] + public double ImplausibleWeight { get; set; } = 0.25; + + /// Weight multiplier of a feature all of whose members have an MS1 artefact explanation. + [Key(22)] + public double ArtefactWeight { get; set; } = 0.25; + + // ---- run-level reversed-phase gate ---- + + /// The median within-class carbon slope must exceed this fraction of the span per carbon. + [Key(23)] + public double SlopeFraction { get; set; } = 0.002; + + [Key(24)] + public double GateConcordance { get; set; } = 0.6; + + [Key(25)] + public int GateMinPairs { get; set; } = 5; + + // ---- class gate ---- + + [Key(26)] + public double ClassGateConcordance { get; set; } = 0.7; + + [Key(27)] + public int ClassGateMinPairs { get; set; } = 4; + + // ---- outcome rules ---- + + /// Block when a co-eluting adjacent composition weighs at least this multiple of the feature. + [Key(28)] + public double TieRatio { get; set; } = 1d; + + /// A feature that wins order conflicts must win them by at least this many evidence units. + [Key(29)] + public double CoreMargin { get; set; } = 2d; + + /// Magnitude veto: |residual| above this many double-bond steps is not Confirmed (0 = off). + [Key(30)] + public double MagnitudeK { get; set; } = 0.75; + + // ---- MS1 artefact layer ---- + + [Key(31)] + public LipidRtArtefactKinds ArtefactKinds { get; set; } = LipidRtArtefactKinds.Default; + + [Key(32)] + public double ParentRatio { get; set; } = 1d; + + [Key(33)] + public double Iso1Ratio { get; set; } = 1.5; + + [Key(34)] + public double Iso2Ratio { get; set; } = 3d; + + [Key(35)] + public double MzToleranceDa { get; set; } = 0.002; + + [Key(36)] + public double MzTolerancePpm { get; set; } = 6d; + + // ---- outcomes ---- + + /// The minimum deficit for Inconsistent. + [Key(37)] + public double InconsistentMinDeficit { get; set; } = 1d; + + /// "RT and spectrum both excellent" override of a soft plausibility rule: this many axes ... + [Key(38)] + public int SoftOverrideAxes { get; set; } = 3; + + /// ... and at least this weighted dot product. + [Key(39)] + public double SoftOverrideDot { get; set; } = 0.8; + + // ---- cross polarity ---- + + [Key(40)] + public double CrossPolarityToleranceFraction { get; set; } = 0.0025; + + [Key(41)] + public double CrossPolarityToleranceSd { get; set; } = 3d; + + /// + /// Single-run reasons that a Confirmed other-polarity twin may promote. Null or empty means the default set + /// (Support, ClassGate, Ambiguous). + /// + [Key(42)] + public List PromotableReasons { get; set; } + + [Key(43)] + public LipidRtPolarityConflictPolicy PolarityConflict { get; set; } = LipidRtPolarityConflictPolicy.Flag; + + // ---- review-list hygiene and penalties ---- + + /// Allow Inconsistent for oxidised classes (off: their conflicts are NotEvaluable). + [Key(44)] + public bool FlagOxidisedClasses { get; set; } = false; + + /// Allow Inconsistent in the gradient-end tail (off: NotEvaluable). + [Key(45)] + public bool FlagGradientTail { get; set; } = false; + + [Key(46)] + public bool TieBlock { get; set; } = true; + + [Key(47)] + public double PolicyPenalty { get; set; } = 2d; + + [Key(48)] + public double BlockPenalty { get; set; } = 1.5; + + /// Supporting and discordant neighbours kept in the evidence record. + [Key(49)] + public int MaxEvidenceItems { get; set; } = 4; + + /// + /// Whether the check runs for a project. An explicit value wins; otherwise it runs for lipidomics projects. + /// + public bool IsEnabledFor(TargetOmics targetOmics) { + return IsEnabled ?? targetOmics == TargetOmics.Lipidomics; + } + + /// + /// Whether processing runs the check: only on a lipidomics project, and only while it is enabled. An explicit + /// True on another target omics does not make it run, because the rules are lipid retention rules. + /// + public bool RunsFor(TargetOmics targetOmics) { + return targetOmics == TargetOmics.Lipidomics && IsEnabledFor(targetOmics); + } + + /// + /// The parameter-file lines of the check. The setting lines use the keys of the Console method file, so the + /// values can be copied into one; the last three lines record what a run did with them. + /// + public List ParametersAsText(TargetOmics targetOmics) { + return new List + { + "Lipid RT check: " + (IsEnabled is bool enabled ? enabled.ToString() : "Auto"), + "Lipid RT check operating point: " + OperatingPoint, + "Lipid RT check chromatography: " + Chromatography, + "Lipid RT check plausibility profile: " + PlausibilityProfile, + "Lipid RT check class-level matches: " + ClassLevelPolicy, + "Lipid RT check minor duplicates: " + MinorDuplicatePolicy, + "Lipid RT check secondary-only features: " + SecondaryOnlyPolicy, + "Lipid RT check runs for this project: " + RunsFor(targetOmics), + "Lipid RT check rule version: " + RuleVersion, + "Lipid RT check parameter hash: " + ComputeHash(), + }; + } + + /// The single-run reasons a Confirmed other-polarity twin may promote, as reason codes. + /// + /// Stored as integers so that the parameter has no dependency on the algorithm's enum; the default is + /// Support (14), ClassGate (12) and Ambiguous (13). On the development data Support promotions were 51/51 + /// in scope, while Margin promotions were 5/7 and were left out. + /// + public IReadOnlyCollection GetPromotableReasons() { + return PromotableReasons is null || PromotableReasons.Count == 0 + ? DefaultPromotableReasons + : PromotableReasons; + } + + [IgnoreMember] + public static IReadOnlyCollection DefaultPromotableReasons { get; } = new[] { 12, 13, 14 }; + + /// A copy, so that a caller can change one value for one run without touching the project's. + public LipidRtCheckParameter Clone() { + var clone = (LipidRtCheckParameter)MemberwiseClone(); + clone.PromotableReasons = PromotableReasons is null ? null : new List(PromotableReasons); + return clone; + } + + /// + /// A short, stable fingerprint of every value that affects an outcome, recorded with each verdict so a reader + /// can tell which settings produced it. + /// + /// + /// Independent of , which decides whether the check runs, not what it concludes. + /// + public string ComputeHash() { + var text = string.Join(";", CanonicalValues()); + using (var sha = SHA256.Create()) { + var bytes = sha.ComputeHash(Encoding.UTF8.GetBytes(text)); + var builder = new StringBuilder(16); + for (int i = 0; i < 8; i++) { + builder.Append(bytes[i].ToString("x2", CultureInfo.InvariantCulture)); + } + return builder.ToString(); + } + } + + private IEnumerable CanonicalValues() { + yield return "rule=" + RuleVersion; + yield return "op=" + OperatingPoint; + yield return "chrom=" + Chromatography; + yield return "profile=" + PlausibilityProfile; + yield return "class=" + ClassLevelPolicy; + yield return "minordup=" + MinorDuplicatePolicy; + yield return "secondary=" + SecondaryOnlyPolicy; + yield return "xpol=" + CrossPolarity; + foreach (var value in new[] { + DiscordanceFraction, TieFraction, StepFraction, CollapseFraction, CoelutionFraction, TailFraction, + LevelWeightChain, LevelWeightMsms, LevelWeightClass, DotWeight, DuplicateWeight, HeightWeight, + ImplausibleWeight, ArtefactWeight, SlopeFraction, GateConcordance, ClassGateConcordance, TieRatio, + CoreMargin, MagnitudeK, ParentRatio, Iso1Ratio, Iso2Ratio, MzToleranceDa, MzTolerancePpm, + InconsistentMinDeficit, SoftOverrideDot, CrossPolarityToleranceFraction, CrossPolarityToleranceSd, + PolicyPenalty, BlockPenalty, }) { + yield return value.ToString("R", CultureInfo.InvariantCulture); + } + yield return "ints=" + string.Join(",", new[] { DuplicateCap, GateMinPairs, ClassGateMinPairs, SoftOverrideAxes, MaxEvidenceItems }); + yield return "art=" + (int)ArtefactKinds; + yield return "promote=" + string.Join(",", GetPromotableReasons().OrderBy(x => x)); + yield return "conflict=" + PolarityConflict; + yield return "flags=" + FlagOxidisedClasses + "," + FlagGradientTail + "," + TieBlock; + } + } +} diff --git a/src/MSDIAL5/MsdialCore/Parameter/ParameterBase.cs b/src/MSDIAL5/MsdialCore/Parameter/ParameterBase.cs index e07ec14c1..ec3186eef 100644 --- a/src/MSDIAL5/MsdialCore/Parameter/ParameterBase.cs +++ b/src/MSDIAL5/MsdialCore/Parameter/ParameterBase.cs @@ -511,6 +511,30 @@ public ParameterBase(bool isLabUseOnly) { [Key(17)] public PostCuratorParameter PostCurationParameter { get; set; } = new PostCuratorParameter(); + /// + /// The lipid retention-time order check. LC-MS only in this version: hook A per analysis file and hook B after + /// alignment, both when and the target omics is lipidomics. + /// + /// + /// Key 18 is the next free key of this class (MsdialDimsParameter uses 20, the other modes 150 and up). A project + /// written before it has no element here, and MessagePack then assigns null rather than the initializer, so the + /// getter hands out the defaults instead of null: an older lipidomics project runs the check with the defaults. + /// The defaults are put in place once, atomically, because hook A reads this from the per-file parallel loop. + /// + [Key(18)] + public LipidRtCheckParameter LipidRtCheckParam { + get { + var parameter = _lipidRtCheckParam; + if (parameter is null) { + System.Threading.Interlocked.CompareExchange(ref _lipidRtCheckParam, new LipidRtCheckParameter(), null); + parameter = _lipidRtCheckParam; + } + return parameter; + } + set => _lipidRtCheckParam = value; + } + private LipidRtCheckParameter _lipidRtCheckParam = new LipidRtCheckParameter(); + public virtual List ParametersAsText() { var pStrings = new List(); pStrings.Add("# Project information"); diff --git a/src/MSDIAL5/MsdialCore/Parser/MsdialPeakSerializer.cs b/src/MSDIAL5/MsdialCore/Parser/MsdialPeakSerializer.cs index 92968e409..2b26e6ab8 100644 --- a/src/MSDIAL5/MsdialCore/Parser/MsdialPeakSerializer.cs +++ b/src/MSDIAL5/MsdialCore/Parser/MsdialPeakSerializer.cs @@ -12,24 +12,7 @@ public static class MsdialPeakSerializer public static void SaveChromatogramPeakFeatures(string file, List chromPeakFeatures) { MessagePackHandler.SaveToFile>(chromPeakFeatures, file); var tagfile = Path.Combine(Path.GetDirectoryName(file), Path.GetFileNameWithoutExtension(file) + "_tags.xml"); - var defsElement = new XElement("Definitions"); - foreach (var tag in PeakSpotTag.AllTypes()) { - var tagElement = new XElement("Tag", new XElement("Id", tag.Id), new XElement("Label", tag.Label)); - defsElement.Add(tagElement); - } - var peaksElement = new XElement("Peaks"); - foreach (var peak in chromPeakFeatures.SelectMany(Flatten)) { - if (!peak.TagCollection.Any()) { - continue; - } - var peakElement = new XElement("Peak"); - peakElement.SetAttributeValue("Id", peak.MasterPeakID); - foreach (var tag in peak.TagCollection.Selected) { - peakElement.Add(new XElement("Tag", tag.Id)); - } - peaksElement.Add(peakElement); - } - var doc = new XElement("PeakSpotTags", defsElement, peaksElement); + var doc = PeakSpotTagXml.Create(chromPeakFeatures.SelectMany(Flatten), peak => peak.MasterPeakID, peak => peak.TagCollection); doc.Save(tagfile); } @@ -38,21 +21,7 @@ public static List LoadChromatogramPeakFeatures(string var tagfile = Path.Combine(Path.GetDirectoryName(file), Path.GetFileNameWithoutExtension(file) + "_tags.xml"); if (File.Exists(tagfile)) { var doc = XElement.Load(tagfile); - var dic = new Dictionary(); - foreach (var peakElement in doc.Descendants("Peak")) { - if (int.TryParse(peakElement.Attribute("Id").Value, out var id)) { - var tags = peakElement.Elements("Tag") - .Select(elem => int.TryParse(elem.Value, out var v) ? PeakSpotTag.GetById(v) : null) - .Where(v => v != null) - .ToArray(); - dic[id] = tags; - } - } - foreach (var peak in peaks.SelectMany(Flatten)) { - if (dic.TryGetValue(peak.MasterPeakID, out var tags)) { - peak.TagCollection.Set(tags); - } - } + PeakSpotTagXml.Apply(doc, peaks.SelectMany(Flatten), peak => peak.MasterPeakID, peak => peak.TagCollection); } return peaks; } diff --git a/src/MSDIAL5/MsdialGuiApp/Icons.xaml b/src/MSDIAL5/MsdialGuiApp/Icons.xaml index 167021793..d76580682 100644 --- a/src/MSDIAL5/MsdialGuiApp/Icons.xaml +++ b/src/MSDIAL5/MsdialGuiApp/Icons.xaml @@ -2641,6 +2641,39 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/AlignmentSpotPropertyModel.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/AlignmentSpotPropertyModel.cs index f7ce2653e..994b7c985 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/AlignmentSpotPropertyModel.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/AlignmentSpotPropertyModel.cs @@ -7,6 +7,7 @@ using CompMs.Common.Interfaces; using CompMs.CommonMVVM; using CompMs.MsdialCore.Algorithm.Annotation; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using CompMs.MsdialCore.DataObj; using Reactive.Bindings; using Reactive.Bindings.Extensions; @@ -17,7 +18,7 @@ namespace CompMs.App.Msdial.Model.DataObj { - public sealed class AlignmentSpotPropertyModel : DisposableModelBase, IPeakSpotModel, IFilterable, IAnnotatedObject, IChromatogramPeak, IMoleculeProperty, IIonProperty + public sealed class AlignmentSpotPropertyModel : DisposableModelBase, IPeakSpotModel, IFilterable, IAnnotatedObject, IChromatogramPeak, IMoleculeProperty, IIonProperty, IPeakSpotTagEditable, ILipidRtCheckedObject { public int AlignmentID => innerModel.AlignmentID; public int MasterAlignmentID => innerModel.MasterAlignmentID; @@ -275,16 +276,76 @@ private bool SetPeakSpotTag(PeakSpotTag tag, bool value, string propertyname) { if (value == innerModel.TagCollection.IsSelected(tag)) { return false; } + // A person's toggle: turning a tag on makes it the person's and lifts a suppression; turning off a tag a tool + // set (e.g. the lipid RT check's Confirmed) also suppresses it, so the tool does not set it again. if (value) { - innerModel.TagCollection.Select(tag); + innerModel.TagCollection.SelectByPerson(tag); } else { - innerModel.TagCollection.Deselect(tag); + innerModel.TagCollection.DeselectByPerson(tag); } OnPropertyChanged(propertyname); + OnPropertyChanged(nameof(ConfirmedDetail)); return true; } + public void SelectPeakSpotTagByPerson(PeakSpotTag tag) { + innerModel.TagCollection.SelectByPerson(tag); + OnPropertyChanged(PeakSpotTagPropertyNames.Of(tag)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + + public PeakSpotTagState CapturePeakSpotTag(PeakSpotTag tag) { + return PeakSpotTagState.Capture(innerModel.TagCollection, tag); + } + + public void RestorePeakSpotTag(PeakSpotTagState state) { + if (state is null) { + throw new ArgumentNullException(nameof(state)); + } + state.RestoreTo(innerModel.TagCollection); + OnPropertyChanged(PeakSpotTagPropertyNames.Of(state.Tag)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + + public LipidRtCheckState CaptureLipidRtCheck() { + return new LipidRtCheckState(PeakSpotTagState.Capture(innerModel.TagCollection, PeakSpotTag.CONFIRMED), innerModel.LipidRtVerdict); + } + + public void WithdrawLipidRtCheck() { + LipidRtTagProtocol.Withdraw(innerModel.TagCollection); + innerModel.LipidRtVerdict = null; + RaiseLipidRtCheckChanged(); + } + + public void RestoreLipidRtCheck(LipidRtCheckState state) { + if (state is null) { + throw new ArgumentNullException(nameof(state)); + } + state.Confirmed.RestoreTo(innerModel.TagCollection); + innerModel.LipidRtVerdict = state.Verdict; + RaiseLipidRtCheckChanged(); + } + + /// The lipid RT check's verdict on this spot; null when the check has none (addendum B2). + public LipidRtVerdict? LipidRtVerdict => innerModel.LipidRtVerdict; + + /// + /// Who set or removed the Confirmed tag, and the lipid RT check's verdict, for the tag's tooltip; null when there is + /// nothing to add to the tag itself. + /// + public string? ConfirmedDetail => LipidRtCheckDescription.ConfirmedDetail(innerModel.TagCollection, innerModel.LipidRtVerdict); + + /// + /// Tells the views that the lipid RT check changed this spot's Confirmed tag and verdict behind the model, as the + /// check works on the spot itself through the shared tag protocol. + /// + internal void RaiseLipidRtCheckChanged() { + OnPropertyChanged(nameof(Confirmed)); + OnPropertyChanged(nameof(LipidRtVerdict)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + public void SwitchPeakSpotTag(PeakSpotTag tag) { if (tag == PeakSpotTag.CONFIRMED) { Confirmed = !Confirmed; diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ChromatogramPeakFeatureModel.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ChromatogramPeakFeatureModel.cs index 29493e743..2e51b70b1 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ChromatogramPeakFeatureModel.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ChromatogramPeakFeatureModel.cs @@ -6,6 +6,7 @@ using CompMs.Common.Interfaces; using CompMs.CommonMVVM; using CompMs.MsdialCore.Algorithm.Annotation; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using CompMs.MsdialCore.DataObj; using Reactive.Bindings.Extensions; using System; @@ -13,7 +14,7 @@ namespace CompMs.App.Msdial.Model.DataObj { - public sealed class ChromatogramPeakFeatureModel : DisposableModelBase, IPeakSpotModel, IFilterable, IMoleculeProperty, IIonProperty, IChromatogramPeak, IAnnotatedObject { + public sealed class ChromatogramPeakFeatureModel : DisposableModelBase, IPeakSpotModel, IFilterable, IMoleculeProperty, IIonProperty, IChromatogramPeak, IAnnotatedObject, IPeakSpotTagEditable, ILipidRtCheckedObject { #region Property public int MasterPeakID => innerModel.MasterPeakID; public double? ChromXValue => innerModel.ChromXs.Value; @@ -176,16 +177,72 @@ private bool SetPeakSpotTag(PeakSpotTag tag, bool value, string propertyname) { if (value == innerModel.TagCollection.IsSelected(tag)) { return false; } + // A person's toggle: turning a tag on makes it the person's and lifts a suppression; turning off a tag a tool + // set (e.g. the lipid RT check's Confirmed) also suppresses it, so the tool does not set it again. if (value) { - innerModel.TagCollection.Select(tag); + innerModel.TagCollection.SelectByPerson(tag); } else { - innerModel.TagCollection.Deselect(tag); + innerModel.TagCollection.DeselectByPerson(tag); } OnPropertyChanged(propertyname); + OnPropertyChanged(nameof(ConfirmedDetail)); return true; } + public void SelectPeakSpotTagByPerson(PeakSpotTag tag) { + innerModel.TagCollection.SelectByPerson(tag); + OnPropertyChanged(PeakSpotTagPropertyNames.Of(tag)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + + public PeakSpotTagState CapturePeakSpotTag(PeakSpotTag tag) { + return PeakSpotTagState.Capture(innerModel.TagCollection, tag); + } + + public void RestorePeakSpotTag(PeakSpotTagState state) { + if (state is null) { + throw new ArgumentNullException(nameof(state)); + } + state.RestoreTo(innerModel.TagCollection); + OnPropertyChanged(PeakSpotTagPropertyNames.Of(state.Tag)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + + public LipidRtCheckState CaptureLipidRtCheck() { + return new LipidRtCheckState(PeakSpotTagState.Capture(innerModel.TagCollection, PeakSpotTag.CONFIRMED), innerModel.LipidRtVerdict); + } + + public void WithdrawLipidRtCheck() { + LipidRtTagProtocol.Withdraw(innerModel.TagCollection); + innerModel.LipidRtVerdict = null; + RaiseLipidRtCheckChanged(); + } + + public void RestoreLipidRtCheck(LipidRtCheckState state) { + if (state is null) { + throw new ArgumentNullException(nameof(state)); + } + state.Confirmed.RestoreTo(innerModel.TagCollection); + innerModel.LipidRtVerdict = state.Verdict; + RaiseLipidRtCheckChanged(); + } + + /// The lipid RT check's verdict on this peak (hook A); null when the check has none (addendum B2). + public LipidRtVerdict? LipidRtVerdict => innerModel.LipidRtVerdict; + + /// + /// Who set or removed the Confirmed tag, and the lipid RT check's verdict, for the tag's tooltip; null when there is + /// nothing to add to the tag itself. + /// + public string? ConfirmedDetail => LipidRtCheckDescription.ConfirmedDetail(innerModel.TagCollection, innerModel.LipidRtVerdict); + + private void RaiseLipidRtCheckChanged() { + OnPropertyChanged(nameof(Confirmed)); + OnPropertyChanged(nameof(LipidRtVerdict)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + public void SwitchPeakSpotTag(PeakSpotTag tag) { if (tag == PeakSpotTag.CONFIRMED) { Confirmed = !Confirmed; diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ILipidRtCheckedObject.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ILipidRtCheckedObject.cs new file mode 100644 index 000000000..7a627d38b --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/ILipidRtCheckedObject.cs @@ -0,0 +1,40 @@ +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; +using System; + +namespace CompMs.App.Msdial.Model.DataObj +{ + /// + /// A peak or spot the lipid RT check gives a verdict to (LC-MS peaks and alignment spots). + /// + /// + /// The verdict, and a Confirmed tag the check set, are about the annotation the check saw. When a person changes the + /// annotation (compound search, MS-FINDER, "Set unknown"), they no longer describe it, and a manual annotation is + /// neither tagged nor used as evidence by the check (addendum A6). The change therefore withdraws both: the check's + /// own Confirmed (a person's or another tool's Confirmed stays) and the verdict. An undo of the change puts both back. + /// + public interface ILipidRtCheckedObject + { + /// The Confirmed tag and the verdict as they are now. + LipidRtCheckState CaptureLipidRtCheck(); + + /// Removes the check's own Confirmed tag and the verdict, and tells the views. + void WithdrawLipidRtCheck(); + + /// Puts the Confirmed tag and the verdict back as captured, and tells the views. + void RestoreLipidRtCheck(LipidRtCheckState state); + } + + /// What the lipid RT check left on one peak or spot at a moment: its Confirmed tag, exactly, and its verdict. + public sealed class LipidRtCheckState + { + public LipidRtCheckState(PeakSpotTagState confirmed, LipidRtVerdict? verdict) { + Confirmed = confirmed ?? throw new ArgumentNullException(nameof(confirmed)); + Verdict = verdict; + } + + /// The Confirmed tag with its source and suppression. + public PeakSpotTagState Confirmed { get; } + + public LipidRtVerdict? Verdict { get; } + } +} diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/IPeakSpotTagEditable.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/IPeakSpotTagEditable.cs new file mode 100644 index 000000000..7aacdb0ff --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/IPeakSpotTagEditable.cs @@ -0,0 +1,47 @@ +using CompMs.MsdialCore.DataObj; + +namespace CompMs.App.Msdial.Model.DataObj +{ + /// + /// A peak or spot whose tags a person edits in the GUI. Setting a tag property is a person's toggle (addendum B1); + /// capturing and restoring a tag is how an undo puts it back exactly, provenance included. + /// + public interface IPeakSpotTagEditable + { + bool Confirmed { get; set; } + + /// + /// A person turns the tag on, or takes it over where a tool set it: afterwards it is selected as the person's, so + /// the tool will not withdraw it, and not suppressed. Setting a tag property does nothing for a tag that is on. + /// + void SelectPeakSpotTagByPerson(PeakSpotTag tag); + + PeakSpotTagState CapturePeakSpotTag(PeakSpotTag tag); + + void RestorePeakSpotTag(PeakSpotTagState state); + } + + /// The names of the tag properties every tag-editable model exposes, for its change notifications. + internal static class PeakSpotTagPropertyNames + { + public static string Of(PeakSpotTag tag) { + if (tag == PeakSpotTag.CONFIRMED) { + return nameof(IPeakSpotTagEditable.Confirmed); + } + if (tag == PeakSpotTag.LOW_QUALITY_SPECTRUM) { + return nameof(AlignmentSpotPropertyModel.LowQualitySpectrum); + } + if (tag == PeakSpotTag.MISANNOTATION) { + return nameof(AlignmentSpotPropertyModel.Misannotation); + } + if (tag == PeakSpotTag.COELUTION) { + return nameof(AlignmentSpotPropertyModel.Coelution); + } + if (tag == PeakSpotTag.OVERANNOTATION) { + return nameof(AlignmentSpotPropertyModel.Overannotation); + } + // An unknown tag: every property may have changed. + return string.Empty; + } + } +} diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/LipidRtCheckDescription.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/LipidRtCheckDescription.cs new file mode 100644 index 000000000..1c18f4d5f --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/LipidRtCheckDescription.cs @@ -0,0 +1,90 @@ +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; +using CompMs.MsdialCore.DataObj; +using System; +using System.Collections.Generic; +using System.Globalization; +using System.Linq; + +namespace CompMs.App.Msdial.Model.DataObj +{ + /// + /// What a person reads in the GUI about the lipid retention-time order check on one peak or spot: the verdict, and who + /// set or removed its Confirmed tag (addendum B1 and B2). + /// + internal static class LipidRtCheckDescription + { + /// One line for the peak information panel: the evidence sentence, or the outcome and reason. + /// Null when the check has no verdict on the object. + public static string? Line(LipidRtVerdict? verdict) { + if (verdict is null) { + return null; + } + return string.IsNullOrEmpty(verdict.Evidence) ? $"{verdict.Outcome} ({verdict.Reason})" : verdict.Evidence; + } + + /// + /// The text behind a Confirmed tag: who set it (or that a person removed a tool's), and the check's verdict with the + /// rule version and parameter hash that produced it. + /// + /// Null when there is nothing to add to the tag itself: no verdict, and Confirmed neither set by a tool + /// nor suppressed. + public static string? ConfirmedDetail(PeakSpotTagCollection tags, LipidRtVerdict? verdict) { + if (tags is null) { + throw new ArgumentNullException(nameof(tags)); + } + var lines = new List(); + var confirmed = PeakSpotTag.CONFIRMED; + if (tags.IsSelected(confirmed)) { + var source = tags.GetSource(confirmed); + if (source is null) { + if (verdict != null) { + lines.Add("Confirmed: set by a person (no tool recorded as its source)."); + } + } + else if (LipidRtTagProtocol.IsOwnSource(source)) { + lines.Add($"Confirmed: set by the lipid RT check ({source})."); + } + else { + lines.Add($"Confirmed: set by {source}."); + } + } + else if (tags.IsSuppressed(confirmed)) { + lines.Add("Confirmed: removed by a person after a tool had set it; no tool sets it again."); + } + if (verdict != null) { + var score = verdict.Score.HasValue ? ", score " + verdict.Score.Value.ToString("0.###", CultureInfo.InvariantCulture) : string.Empty; + lines.Add($"Lipid RT check: {verdict.Outcome} ({verdict.Reason}){score}."); + if (!string.IsNullOrEmpty(verdict.Evidence)) { + lines.Add(verdict.Evidence); + } + if (!string.IsNullOrEmpty(verdict.RuleVersion) || !string.IsNullOrEmpty(verdict.ParameterHash)) { + lines.Add($"Rules {verdict.RuleVersion}, parameter {verdict.ParameterHash}."); + } + } + return lines.Count == 0 ? null : string.Join(Environment.NewLine, lines); + } + + /// + /// What a person is told after GUI processing when the lipid RT check failed on analysis files or on the alignment. + /// The failure does not stop processing, and those results are saved without the check's verdicts and tags, which + /// without this message would look the same as a check that found nothing. + /// + /// The units the check ran on (an analysis file's name, or the alignment's), with what it did. + /// Null when the check failed nowhere. + public static string? ProcessingFailures(IEnumerable<(string Unit, LipidRtUnitRun? Run)> runs) { + if (runs is null) { + throw new ArgumentNullException(nameof(runs)); + } + var failed = runs.Where(r => r.Run?.Error != null).ToList(); + if (failed.Count == 0) { + return null; + } + var lines = new List + { + "The lipid RT check failed during processing. These results were saved without its verdicts and Confirmed tags:", + }; + lines.AddRange(failed.Select(f => $"{f.Unit}: {f.Run!.Describe()}")); + return string.Join(Environment.NewLine, lines); + } + } +} diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/Ms1BasedSpectrumFeature.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/Ms1BasedSpectrumFeature.cs index 31e118f8b..0ce90e1c3 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/Ms1BasedSpectrumFeature.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/Ms1BasedSpectrumFeature.cs @@ -9,7 +9,7 @@ namespace CompMs.App.Msdial.Model.DataObj { - public sealed class Ms1BasedSpectrumFeature : BindableBase, IDisposable, IAnnotatedObject + public sealed class Ms1BasedSpectrumFeature : BindableBase, IDisposable, IAnnotatedObject, IPeakSpotTagEditable { private SpectrumFeature _spectrumFeature; private bool _disposedValue; @@ -64,16 +64,41 @@ private bool SetPeakSpotTag(PeakSpotTag tag, bool value, string propertyname) { if (value == _spectrumFeature.TagCollection.IsSelected(tag)) { return false; } + // A person's toggle; see AlignmentSpotPropertyModel.SetPeakSpotTag. if (value) { - _spectrumFeature.TagCollection.Select(tag); + _spectrumFeature.TagCollection.SelectByPerson(tag); } else { - _spectrumFeature.TagCollection.Deselect(tag); + _spectrumFeature.TagCollection.DeselectByPerson(tag); } OnPropertyChanged(propertyname); + OnPropertyChanged(nameof(ConfirmedDetail)); return true; } + public void SelectPeakSpotTagByPerson(PeakSpotTag tag) { + _spectrumFeature.TagCollection.SelectByPerson(tag); + OnPropertyChanged(PeakSpotTagPropertyNames.Of(tag)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + + public PeakSpotTagState CapturePeakSpotTag(PeakSpotTag tag) { + return PeakSpotTagState.Capture(_spectrumFeature.TagCollection, tag); + } + + public void RestorePeakSpotTag(PeakSpotTagState state) { + if (state is null) { + throw new ArgumentNullException(nameof(state)); + } + state.RestoreTo(_spectrumFeature.TagCollection); + OnPropertyChanged(PeakSpotTagPropertyNames.Of(state.Tag)); + OnPropertyChanged(nameof(ConfirmedDetail)); + } + + /// Who set or removed the Confirmed tag, for the tag's tooltip; null when a person set it or it is off. + /// GC-MS has no lipid RT check, so there is never a verdict here. + public string? ConfirmedDetail => LipidRtCheckDescription.ConfirmedDetail(_spectrumFeature.TagCollection, null); + public void SwitchPeakSpotTag(PeakSpotTag tag) { if (tag == PeakSpotTag.CONFIRMED) { Confirmed = !Confirmed; diff --git a/src/MSDIAL5/MsdialGuiApp/Model/DataObj/PeakSpotTagState.cs b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/PeakSpotTagState.cs new file mode 100644 index 000000000..56c9aea4b --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/Model/DataObj/PeakSpotTagState.cs @@ -0,0 +1,72 @@ +using CompMs.MsdialCore.DataObj; +using System; + +namespace CompMs.App.Msdial.Model.DataObj +{ + /// + /// One tag of one peak or spot as it was at a moment: whether it was selected, who set it, and whether a person had + /// suppressed it. An undo puts the tag back from this exactly, its provenance included (addendum B1). + /// + /// + /// Undoing through a person's toggle would not: turning a tool's Confirmed back on that way makes it the person's, and + /// turning a tag on lifts its suppression. "Mark all as unconfirmed" followed by Undo would then hand every Confirmed + /// the lipid RT check had set to the person, out of the check's reach, and "Mark all as confirmed" followed by Undo + /// would forget what a person had removed before. + /// + public sealed class PeakSpotTagState + { + private PeakSpotTagState(PeakSpotTag tag, bool isSelected, string? source, bool isSuppressed, string? suppressionSource) { + Tag = tag; + IsSelected = isSelected; + Source = source; + IsSuppressed = isSuppressed; + SuppressionSource = suppressionSource; + } + + public PeakSpotTag Tag { get; } + public bool IsSelected { get; } + + /// Who set the tag; null when a person did or nobody recorded who did. + public string? Source { get; } + public bool IsSuppressed { get; } + public string? SuppressionSource { get; } + + public static PeakSpotTagState Capture(PeakSpotTagCollection tags, PeakSpotTag tag) { + if (tags is null) { + throw new ArgumentNullException(nameof(tags)); + } + if (tag is null) { + throw new ArgumentNullException(nameof(tag)); + } + return new PeakSpotTagState(tag, tags.IsSelected(tag), tags.GetSource(tag), tags.IsSuppressed(tag), tags.GetSuppressionSource(tag)); + } + + /// + /// Puts the tag back as it was captured. A part that is already as it was is left alone, so the order of the + /// selected tags (and with it the "Tag:" text of the exports) changes only where the tag itself changes. + /// + public void RestoreTo(PeakSpotTagCollection tags) { + if (tags is null) { + throw new ArgumentNullException(nameof(tags)); + } + if (tags.IsSelected(Tag) != IsSelected || tags.GetSource(Tag) != Source) { + tags.Deselect(Tag); + if (IsSelected) { + if (Source is null) { + tags.Select(Tag); + } + else { + tags.Select(Tag, Source); + } + } + } + if (tags.IsSuppressed(Tag) != IsSuppressed || tags.GetSuppressionSource(Tag) != SuppressionSource) { + tags.Unsuppress(Tag); + if (IsSuppressed) { + // A null source records a suppression nobody signed, which is what it was. + tags.Suppress(Tag, SuppressionSource!); + } + } + } + } +} diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Gcms/GcmsPeakTableModel.cs b/src/MSDIAL5/MsdialGuiApp/Model/Gcms/GcmsPeakTableModel.cs index 3a1803d5d..4830b6c4a 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Gcms/GcmsPeakTableModel.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Gcms/GcmsPeakTableModel.cs @@ -34,13 +34,13 @@ public GcmsAnalysisPeakTableModel(IReadOnlyList spectru public double RiMax { get; } public override void MarkAllAsConfirmed() { - IDoCommand command = new MarkAllAsCommand(PeakSpots, true); + IDoCommand command = new MarkAllAsConfirmedCommand(PeakSpots, true); command.Do(); _undoManager.Add(command); } public override void MarkAllAsUnconfirmed() { - IDoCommand command = new MarkAllAsCommand(PeakSpots, false); + IDoCommand command = new MarkAllAsConfirmedCommand(PeakSpots, false); command.Do(); _undoManager.Add(command); } @@ -51,34 +51,5 @@ public override void SwitchTag(PeakSpotTag tag) { } _target.Value.SwitchPeakSpotTag(tag); } - - sealed class MarkAllAsCommand : IDoCommand { - private readonly IReadOnlyList _peaks; - private readonly List _confirmeds; - private readonly bool _status; - - public MarkAllAsCommand(IReadOnlyList peaks, bool status) - { - _peaks = peaks; - _status = status; - _confirmeds = new List(); - } - - void IDoCommand.Do() { - _confirmeds.Clear(); - foreach (var peak in _peaks) { - if (peak.Confirmed != _status) { - peak.Confirmed = _status; - _confirmeds.Add(peak); - } - } - } - - void IDoCommand.Undo() { - foreach (var peak in _confirmeds) { - peak.Confirmed = !_status; - } - } - } } } diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Information/PeakInformationModel.cs b/src/MSDIAL5/MsdialGuiApp/Model/Information/PeakInformationModel.cs index 5fc3b35e3..12b657cd5 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Information/PeakInformationModel.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Information/PeakInformationModel.cs @@ -2,6 +2,7 @@ using CompMs.App.Msdial.Utility; using CompMs.Common.DataObj.Property; using CompMs.CommonMVVM; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using Reactive.Bindings; using Reactive.Bindings.Extensions; using System; @@ -21,6 +22,13 @@ internal interface IPeakInformationModel : INotifyPropertyChanged string? Ontology { get; } string? InChIKey { get; } string? Comment { get; } + + /// The lipid RT check's verdict in one line; null when the check has none on this peak. + string? LipidRtCheck { get; } + + /// The verdict in full, with who set or removed the Confirmed tag; null when there is nothing to say. + string? LipidRtCheckDetail { get; } + ReadOnlyObservableCollection? PeakPoints { get; } ReadOnlyObservableCollection? PeakAmounts { get; } } @@ -37,7 +45,9 @@ public PeakInformationAnalysisModel(IObservable m (m?.ObserveProperty(m_ => m_.Formula) ?? Observable.Return(new Formula())).Subscribe(m_ => Formula = m_?.FormulaString), (m?.ObserveProperty(m_ => m_.Ontology) ?? Observable.Return(string.Empty)).Subscribe(m_ => Ontology = m_), (m?.ObserveProperty(m_ => m_.InChIKey) ?? Observable.Return(string.Empty)).Subscribe(m_ => InChIKey = m_), - (m?.ObserveProperty(m_ => m_.Comment) ?? Observable.Return(string.Empty)).Subscribe(m_ => Comment = m_) + (m?.ObserveProperty(m_ => m_.Comment) ?? Observable.Return(string.Empty)).Subscribe(m_ => Comment = m_), + (m?.ObserveProperty(m_ => m_.LipidRtVerdict) ?? Observable.Return(null)).Subscribe(v => LipidRtCheck = LipidRtCheckDescription.Line(v)), + (m?.ObserveProperty(m_ => m_.ConfirmedDetail) ?? Observable.Return(null)).Subscribe(d => LipidRtCheckDetail = d), }; }).DisposePreviousValue() .Subscribe() @@ -96,6 +106,18 @@ public string? Comment { } private string? _comment; + public string? LipidRtCheck { + get => _lipidRtCheck; + private set => SetProperty(ref _lipidRtCheck, value); + } + private string? _lipidRtCheck; + + public string? LipidRtCheckDetail { + get => _lipidRtCheckDetail; + private set => SetProperty(ref _lipidRtCheckDetail, value); + } + private string? _lipidRtCheckDetail; + ReadOnlyObservableCollection? IPeakInformationModel.PeakPoints => _peakPoints; public ReadOnlyReactiveCollection? PeakPoints { get => _peakPoints; @@ -130,6 +152,10 @@ internal sealed class PeakInformationAlignmentModel : DisposableModelBase, IPeak public PeakInformationAlignmentModel(IObservable model) { model.Select(m => { if (m is null) { + // No spot, no verdict: the other fields keep their last values as before, but a verdict must not be + // read as being about nothing. + LipidRtCheck = null; + LipidRtCheckDetail = null; return null; } return new CompositeDisposable @@ -139,7 +165,9 @@ public PeakInformationAlignmentModel(IObservable mo m.ObserveProperty(m_ => m_.Formula).Subscribe(m_ => Formula = m_?.FormulaString), m.ObserveProperty(m_ => m_.Ontology).Subscribe(m_ => Ontology = m_), m.ObserveProperty(m_ => m_.InChIKey).Subscribe(m_ => InChIKey = m_), - m.ObserveProperty(m_ => m_.Comment).Subscribe(m_ => Comment = m_) + m.ObserveProperty(m_ => m_.Comment).Subscribe(m_ => Comment = m_), + m.ObserveProperty(m_ => m_.LipidRtVerdict).Subscribe(v => LipidRtCheck = LipidRtCheckDescription.Line(v)), + m.ObserveProperty(m_ => m_.ConfirmedDetail).Subscribe(d => LipidRtCheckDetail = d), }; }).DisposePreviousValue() .Subscribe() @@ -198,6 +226,18 @@ public string? Comment { } private string? _comment; + public string? LipidRtCheck { + get => _lipidRtCheck; + private set => SetProperty(ref _lipidRtCheck, value); + } + private string? _lipidRtCheck; + + public string? LipidRtCheckDetail { + get => _lipidRtCheckDetail; + private set => SetProperty(ref _lipidRtCheckDetail, value); + } + private string? _lipidRtCheckDetail; + ReadOnlyObservableCollection IPeakInformationModel.PeakPoints => _peakPoints; public ReadOnlyReactiveCollection PeakPoints { get => _peakPoints; @@ -286,6 +326,10 @@ public string Comment { } private string _comment = string.Empty; + /// GC-MS has no lipid RT check. + string? IPeakInformationModel.LipidRtCheck => null; + string? IPeakInformationModel.LipidRtCheckDetail => null; + ReadOnlyObservableCollection? IPeakInformationModel.PeakPoints => _peakPoints; public ReadOnlyReactiveCollection? PeakPoints { get => _peakPoints; diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Information/SetAnnotationUsecase.cs b/src/MSDIAL5/MsdialGuiApp/Model/Information/SetAnnotationUsecase.cs index be4443817..c350d9780 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Information/SetAnnotationUsecase.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Information/SetAnnotationUsecase.cs @@ -8,6 +8,12 @@ namespace CompMs.App.Msdial.Model.Information { + /// A person sets the annotation of a peak or spot: a compound from a search, Unsettled, or Unknown. + /// + /// Each of these also withdraws what the lipid RT check left on the object, its own Confirmed tag and its verdict, + /// which were about the annotation it replaces; a manual annotation is not evaluated by the check (addendum A6). + /// Setting Unknown can be undone, and its undo puts them back. + /// internal sealed class SetAnnotationUsecase { private readonly IMoleculeProperty _molecule; @@ -25,6 +31,7 @@ public void SetConfidence(ICompoundResult compoundResult) { compoundResult.MatchResult.IsReferenceMatched = true; _results.RemoveManuallyResults(); _results.AddResult(compoundResult.MatchResult); + WithdrawLipidRtCheck(); } public void SetUnsettled(ICompoundResult compoundResult) { @@ -32,6 +39,7 @@ public void SetUnsettled(ICompoundResult compoundResult) { compoundResult.MatchResult.IsReferenceMatched = true; _results.RemoveManuallyResults(); _results.AddResult(compoundResult.MatchResult); + WithdrawLipidRtCheck(); } public void SetConfidence(MoleculeMsReference msReference, MsScanMatchResult matchResult) { @@ -39,6 +47,7 @@ public void SetConfidence(MoleculeMsReference msReference, MsScanMatchResult mat matchResult.IsReferenceMatched = true; _results.RemoveManuallyResults(); _results.AddResult(matchResult); + WithdrawLipidRtCheck(); } public void SetUnknown() { @@ -46,5 +55,9 @@ public void SetUnknown() { command.Do(); _undoManager.Add(command); } + + private void WithdrawLipidRtCheck() { + (_molecule as ILipidRtCheckedObject)?.WithdrawLipidRtCheck(); + } } } diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Lcms/LcmsMethodModel.cs b/src/MSDIAL5/MsdialGuiApp/Model/Lcms/LcmsMethodModel.cs index 767870544..bbcffed2f 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Lcms/LcmsMethodModel.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Lcms/LcmsMethodModel.cs @@ -4,6 +4,7 @@ using CompMs.App.Msdial.Model.Export; using CompMs.App.Msdial.Model.Loader; using CompMs.App.Msdial.Model.Search; +using CompMs.App.Msdial.Model.Service; using CompMs.App.Msdial.Model.Setting; using CompMs.App.Msdial.Model.Statistics; using CompMs.Common.Components; @@ -11,6 +12,7 @@ using CompMs.Graphics.UI.ProgressBar; using CompMs.MsdialCore.Algorithm; using CompMs.MsdialCore.Algorithm.Annotation; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; using CompMs.MsdialCore.DataObj; using CompMs.MsdialCore.Export; using CompMs.MsdialCore.MSDec; @@ -263,6 +265,7 @@ public override async Task RunAsync(ProcessOption processOption, CancellationTok var parameter = _storage.Parameter; var starttimestamp = DateTime.Now.ToString("yyyyMMddHHmm"); var stopwatch = Stopwatch.StartNew(); + var lipidRtRuns = new List<(string Unit, LipidRtUnitRun? Run)>(); // Run Identification if (processOption.HasFlag(ProcessOption.Identification)) { @@ -282,6 +285,7 @@ public override async Task RunAsync(ProcessOption processOption, CancellationTok if (!ProcessFiles(_storage.AnalysisFiles, runner, processOption)) { return; } + lipidRtRuns.AddRange(processor.LipidRtFailures.Select(failure => (failure.File, (LipidRtUnitRun?)failure.Run))); } // Run second process @@ -292,13 +296,19 @@ public override async Task RunAsync(ProcessOption processOption, CancellationTok // Run Alignment if (processOption.HasFlag(ProcessOption.Alignment)) { - if (!ProcessAlignment(_storage)) + if (!ProcessAlignment(_storage, lipidRtRuns)) return; } stopwatch.Stop(); var ts = stopwatch.Elapsed; AutoParametersSave(starttimestamp, ts, parameter); + // The lipid RT check reports a failure only on the Console, which the GUI does not show; tell the person once, + // after processing, rather than leave results without verdicts looking like ones where it found nothing. + if (LipidRtCheckDescription.ProcessingFailures(lipidRtRuns) is string lipidRtNotice) { + _broker.Publish(new ShortMessageRequest(lipidRtNotice)); + } + await LoadAnalysisFileAsync(AnalysisFileModelCollection.AnalysisFiles.FirstOrDefault(), token).ConfigureAwait(false); } @@ -347,6 +357,12 @@ await Task.Run(() => proteomicsAnnotator.ExecuteSecondRoundAnnotationProcess( } public bool ProcessAlignment(IMsdialDataStorage storage) { + return ProcessAlignment(storage, null); + } + + /// Receives what the lipid RT check did on the alignment (hook B), when not null. + private bool ProcessAlignment(IMsdialDataStorage storage, List<(string Unit, LipidRtUnitRun? Run)>? lipidRtRuns) { + (string Unit, LipidRtUnitRun? Run)? lipidRtRun = null; var request = new ProgressBarRequest("Process alignment..", isIndeterminate: false, async vm => { @@ -360,6 +376,7 @@ public bool ProcessAlignment(IMsdialDataStorage storage) { var alignmentFileModel = AlignmentFiles.Files.Last(); var result = await Task.Run(() => alignmentFileModel.RunAlignment(aligner, CHROMATOGRAM_SPOT_SERIALIZER)).ConfigureAwait(false); + lipidRtRun = (alignmentFileModel.FileName, aligner.LipidRtRun); if (storage.DataBases.ProteomicsDataBases.Any()) { new ProteomeDataAnnotator().MappingToProteinDatabase( @@ -379,6 +396,9 @@ public bool ProcessAlignment(IMsdialDataStorage storage) { await Task.WhenAll(tasks).ConfigureAwait(false); }); _broker.Publish(request); + if (lipidRtRun.HasValue) { + lipidRtRuns?.Add(lipidRtRun.Value); + } return request.Result ?? false; } @@ -489,5 +509,17 @@ public FragmentQuerySettingModel ShowShowFragmentSearchSettingView() { } return new MscleanrSettingModel(_storage.Parameter, AlignmentModel.Ms1Spots); } + + /// + /// The "Re-run lipid RT check" settings for the loaded alignment result, or null when none is loaded. The run saves + /// the result the way the GUI saves it when another result is opened or the project is saved. + /// + public LipidRtCheckSettingModel? ShowLipidRtCheckSettingView() { + var alignmentModel = AlignmentModel; + if (alignmentModel is null) { + return null; + } + return new LipidRtCheckSettingModel(_storage.Parameter, alignmentModel.Ms1Spots, alignmentModel.SaveAsync, alignmentModel.UndoManager); + } } } diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Service/SetUnknownDoCommand.cs b/src/MSDIAL5/MsdialGuiApp/Model/Service/SetUnknownDoCommand.cs index 41ac7961d..5b4d18d82 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Service/SetUnknownDoCommand.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Service/SetUnknownDoCommand.cs @@ -14,6 +14,7 @@ internal class SetUnknownDoCommand : IDoCommand private List? _manuallyResutls; private MsScanMatchResult? _unknownResult; private IMoleculeProperty? _previousMolecule; + private LipidRtCheckState? _previousLipidRtCheck; public SetUnknownDoCommand(IMoleculeProperty molecule, MsScanMatchResultContainerModel container) { @@ -36,6 +37,12 @@ void IDoCommand.Do() DataAccess.ClearMoleculePropertyInfomation(_molecule); _container.RemoveManuallyResults(); _container.AddResult(_unknownResult); + + // The lipid RT check's verdict and its own Confirmed were about the annotation that is now gone (addendum A6). + if (_molecule is ILipidRtCheckedObject checkedObject) { + _previousLipidRtCheck = checkedObject.CaptureLipidRtCheck(); + checkedObject.WithdrawLipidRtCheck(); + } } void IDoCommand.Undo() @@ -51,6 +58,10 @@ void IDoCommand.Undo() _molecule.InChIKey = _previousMolecule.InChIKey; _container.RemoveResult(_unknownResult); _container.AddResults(_manuallyResutls); + + if (_molecule is ILipidRtCheckedObject checkedObject && _previousLipidRtCheck is not null) { + checkedObject.RestoreLipidRtCheck(_previousLipidRtCheck); + } } } } diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Service/UndoManager.cs b/src/MSDIAL5/MsdialGuiApp/Model/Service/UndoManager.cs index 3a3466106..cca4b0da2 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Service/UndoManager.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Service/UndoManager.cs @@ -28,6 +28,16 @@ public void Undo() { Shrink(_redos, REDO_LIMIT); } + /// + /// Forgets every command to undo or redo. For an operation that changes the objects behind the commands' backs, such + /// as a re-run of a tool over the whole result, after which putting back what a command captured would no longer + /// agree with the objects. + /// + public void Clear() { + _undos.Clear(); + _redos.Clear(); + } + public void Redo() { if (!_redos.Any()) { return; diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Setting/LipidRtCheckSettingModel.cs b/src/MSDIAL5/MsdialGuiApp/Model/Setting/LipidRtCheckSettingModel.cs new file mode 100644 index 000000000..7a93167b7 --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/Model/Setting/LipidRtCheckSettingModel.cs @@ -0,0 +1,170 @@ +using CompMs.App.Msdial.Model.DataObj; +using CompMs.App.Msdial.Model.Service; +using CompMs.Common.Enum; +using CompMs.CommonMVVM; +using CompMs.MsdialCore.Algorithm.Annotation.LipidRt; +using CompMs.MsdialCore.Parameter; +using System; +using System.Collections.Generic; +using System.Linq; +using System.Threading.Tasks; + +namespace CompMs.App.Msdial.Model.Setting +{ + /// + /// "Re-run lipid RT check" on a loaded LC-MS alignment result (addendum B3): the user-facing settings of the check, and + /// a run of the check on the result's spots through the same routine as hook B, followed by a save of the result. + /// + /// + /// + /// The settings are the project's: running writes them into , as the + /// MS-CleanR filter writes its own, so the verdicts' parameter hash is the project's and the next processing uses the + /// same settings. The project file keeps them when the project is saved. + /// + /// + /// The check goes through and with it : it + /// withdraws only the Confirmed tags it set itself, sets Confirmed only where no person's tag or decision is in the + /// way, and replaces every spot's verdict. A person's tags, another tool's tags, a Misannotation and a suppressed + /// Confirmed are left as they are. Cross-polarity corroboration is not applied here (v1.1). + /// + /// + /// A run is not undoable, and it clears the result's undo history: an earlier tag edit's undo would put back tags as + /// they were before the run, such as a Confirmed credited to the check that the new verdicts contradict. + /// + /// + internal sealed class LipidRtCheckSettingModel : BindableBase + { + private readonly ParameterBase _parameter; + private readonly IReadOnlyList _spots; + private readonly Func _saveAsync; + private readonly UndoManager _undoManager; + + /// The project's parameter; its lipid RT check settings are read and, on a run, written. + /// The spots of the loaded alignment result. + /// Saves the loaded alignment result (spots with their verdicts, and the tag file). + /// The undo history of the loaded alignment result's edits, cleared by a run. + public LipidRtCheckSettingModel(ParameterBase parameter, IReadOnlyList spots, Func saveAsync, UndoManager undoManager) { + _parameter = parameter ?? throw new ArgumentNullException(nameof(parameter)); + _spots = spots ?? throw new ArgumentNullException(nameof(spots)); + _saveAsync = saveAsync ?? throw new ArgumentNullException(nameof(saveAsync)); + _undoManager = undoManager ?? throw new ArgumentNullException(nameof(undoManager)); + + var current = parameter.LipidRtCheckParam; + _operatingPoint = current.OperatingPoint; + _chromatography = current.Chromatography; + _plausibilityProfile = current.PlausibilityProfile; + _classLevelPolicy = current.ClassLevelPolicy; + _minorDuplicatePolicy = current.MinorDuplicatePolicy; + } + + public LipidRtOperatingPoint OperatingPoint { + get => _operatingPoint; + set => SetProperty(ref _operatingPoint, value); + } + private LipidRtOperatingPoint _operatingPoint; + + public LipidRtChromatography Chromatography { + get => _chromatography; + set => SetProperty(ref _chromatography, value); + } + private LipidRtChromatography _chromatography; + + public LipidRtPlausibilityProfileKind PlausibilityProfile { + get => _plausibilityProfile; + set => SetProperty(ref _plausibilityProfile, value); + } + private LipidRtPlausibilityProfileKind _plausibilityProfile; + + public LipidRtClassLevelPolicy ClassLevelPolicy { + get => _classLevelPolicy; + set => SetProperty(ref _classLevelPolicy, value); + } + private LipidRtClassLevelPolicy _classLevelPolicy; + + public LipidRtMinorDuplicatePolicy MinorDuplicatePolicy { + get => _minorDuplicatePolicy; + set => SetProperty(ref _minorDuplicatePolicy, value); + } + private LipidRtMinorDuplicatePolicy _minorDuplicatePolicy; + + /// Why the check cannot run on this project, or null when it can. + /// The rules are lipid retention rules, so the check runs only on a lipidomics project, as in processing. + public string? CannotRunReason => _parameter.TargetOmics == TargetOmics.Lipidomics + ? null + : $"The lipid RT check runs on lipidomics projects only; this project's target omics is {_parameter.TargetOmics}."; + + /// + /// True when the project's parameter switches the check off ("Lipid RT check: False" in a Console method file). + /// Running it from here switches it back to Auto, which is on for lipidomics, so that processing and this result agree. + /// + public bool IsSwitchedOff => _parameter.LipidRtCheckParam.IsEnabled == false; + + /// What the last run did, for the person: the outcome counts and the tags, or why it did not run. + public string? Summary { + get => _summary; + private set => SetProperty(ref _summary, value); + } + private string? _summary; + + /// The last run; null before the first. + public LipidRtUnitRun? LastRun { + get => _lastRun; + private set => SetProperty(ref _lastRun, value); + } + private LipidRtUnitRun? _lastRun; + + /// + /// Writes the settings into the project's parameter, runs the check on the spots, clears the undo history, tells the + /// spots' views, and saves the alignment result. A failure of the check is reported in (the + /// hook has then withdrawn the check's tags and verdicts, which is saved too); a failure of the save is reported + /// there as well. + /// + /// The run; null when the check cannot run on this project, in which case nothing is changed. + public async Task RunAsync() { + if (CannotRunReason is string reason) { + Summary = reason; + return null; + } + Commit(); + var spots = _spots.Select(spot => spot.innerModel).ToList(); + var parameter = _parameter; + // The check reads and writes the spots themselves, not their views, so it may run off the UI thread; the views + // are told afterwards, on the caller's thread. + var run = await Task.Run(() => LipidRtCheckHooks.RunOnSpots(spots, parameter)); + _undoManager.Clear(); + foreach (var spot in _spots) { + spot.RaiseLipidRtCheckChanged(); + } + LastRun = run; + + var description = $"Lipid RT check on {spots.Count} spots: {run.Describe()}"; + Summary = description + " Saving the alignment result..."; + string saved; + try { + // The save writes the whole result, which takes a while on a large alignment, so it runs off the UI thread + // as the check does. The window is modal and does not close while this runs, so no one edits the spots. + var saveAsync = _saveAsync; + await Task.Run(() => saveAsync()); + saved = "The alignment result was saved."; + } + catch (Exception ex) { + saved = $"The alignment result could not be saved ({ex.GetType().Name}: {ex.Message}); the tags and verdicts are changed in this session only."; + } + Summary = description + " " + saved; + return run; + } + + private void Commit() { + var target = _parameter.LipidRtCheckParam; + target.OperatingPoint = OperatingPoint; + target.Chromatography = Chromatography; + target.PlausibilityProfile = PlausibilityProfile; + target.ClassLevelPolicy = ClassLevelPolicy; + target.MinorDuplicatePolicy = MinorDuplicatePolicy; + if (target.IsEnabled == false) { + target.IsEnabled = null; + OnPropertyChanged(nameof(IsSwitchedOff)); + } + } + } +} diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Table/AlignmentSpotTableModelBase.cs b/src/MSDIAL5/MsdialGuiApp/Model/Table/AlignmentSpotTableModelBase.cs index ea211fa24..4f0767f9f 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Table/AlignmentSpotTableModelBase.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Table/AlignmentSpotTableModelBase.cs @@ -54,13 +54,13 @@ public async Task ExportMatchedSpectraAsync(string directory) { } public override void MarkAllAsConfirmed() { - IDoCommand command = new MarkAllAsCommand(PeakSpots, true); + IDoCommand command = new MarkAllAsConfirmedCommand(PeakSpots, true); command.Do(); _undoManager.Add(command); } public override void MarkAllAsUnconfirmed() { - IDoCommand command = new MarkAllAsCommand(PeakSpots, false); + IDoCommand command = new MarkAllAsConfirmedCommand(PeakSpots, false); command.Do(); _undoManager.Add(command); } @@ -71,34 +71,5 @@ public override void SwitchTag(PeakSpotTag tag) { } _target.Value.SwitchPeakSpotTag(tag); } - - sealed class MarkAllAsCommand : IDoCommand { - private readonly IReadOnlyList _peaks; - private readonly List _confirmeds; - private readonly bool _status; - - public MarkAllAsCommand(IReadOnlyList peaks, bool status) - { - _peaks = peaks; - _status = status; - _confirmeds = new List(); - } - - void IDoCommand.Do() { - _confirmeds.Clear(); - foreach (var peak in _peaks) { - if (peak.Confirmed != _status) { - peak.Confirmed = _status; - _confirmeds.Add(peak); - } - } - } - - void IDoCommand.Undo() { - foreach (var peak in _confirmeds) { - peak.Confirmed = !_status; - } - } - } } } diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Table/AnalysisPeakSpotTableModelBase.cs b/src/MSDIAL5/MsdialGuiApp/Model/Table/AnalysisPeakSpotTableModelBase.cs index 81db085b7..829b2631a 100644 --- a/src/MSDIAL5/MsdialGuiApp/Model/Table/AnalysisPeakSpotTableModelBase.cs +++ b/src/MSDIAL5/MsdialGuiApp/Model/Table/AnalysisPeakSpotTableModelBase.cs @@ -26,44 +26,15 @@ public override void SwitchTag(PeakSpotTag tag) { } public override void MarkAllAsConfirmed() { - IDoCommand command = new MarkAllAsCommand(PeakSpots, true); + IDoCommand command = new MarkAllAsConfirmedCommand(PeakSpots, true); command.Do(); _undoManager.Add(command); } public override void MarkAllAsUnconfirmed() { - IDoCommand command = new MarkAllAsCommand(PeakSpots, false); + IDoCommand command = new MarkAllAsConfirmedCommand(PeakSpots, false); command.Do(); _undoManager.Add(command); } - - sealed class MarkAllAsCommand : IDoCommand { - private readonly IReadOnlyList _peaks; - private readonly List _confirmeds; - private readonly bool _status; - - public MarkAllAsCommand(IReadOnlyList peaks, bool status) - { - _peaks = peaks; - _status = status; - _confirmeds = new List(); - } - - void IDoCommand.Do() { - _confirmeds.Clear(); - foreach (var peak in _peaks) { - if (peak.Confirmed != _status) { - peak.Confirmed = _status; - _confirmeds.Add(peak); - } - } - } - - void IDoCommand.Undo() { - foreach (var peak in _confirmeds) { - peak.Confirmed = !_status; - } - } - } } } diff --git a/src/MSDIAL5/MsdialGuiApp/Model/Table/MarkAllAsConfirmedCommand.cs b/src/MSDIAL5/MsdialGuiApp/Model/Table/MarkAllAsConfirmedCommand.cs new file mode 100644 index 000000000..cfb7ecf31 --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/Model/Table/MarkAllAsConfirmedCommand.cs @@ -0,0 +1,63 @@ +using CompMs.App.Msdial.Model.DataObj; +using CompMs.App.Msdial.Model.Service; +using CompMs.MsdialCore.DataObj; +using System; +using System.Collections.Generic; + +namespace CompMs.App.Msdial.Model.Table +{ + /// + /// "Mark all as confirmed" and "Mark all as unconfirmed" of the peak tables: a person's decision on Confirmed for every + /// peak or spot of the table. + /// + /// + /// + /// The decision is the person's (addendum B1). "Mark all as confirmed" leaves every row with the person's Confirmed: + /// it turns Confirmed on where it is off, and takes it over where a tool such as the lipid RT check had set it, so a + /// later run of that tool, which withdraws only its own tags, does not take back what the person confirmed. "Mark all + /// as unconfirmed" turns Confirmed off wherever it is on; where a tool had set it, that also suppresses it. + /// + /// + /// Undo does not toggle back; it puts each changed tag back exactly as it was, with its source and its suppression, + /// so the check's own tags are the check's again. + /// + /// + internal sealed class MarkAllAsConfirmedCommand : IDoCommand where T : IPeakSpotTagEditable + { + private readonly IReadOnlyList _peaks; + private readonly bool _status; + private readonly List<(T Peak, PeakSpotTagState Before)> _changed; + + public MarkAllAsConfirmedCommand(IReadOnlyList peaks, bool status) { + _peaks = peaks ?? throw new ArgumentNullException(nameof(peaks)); + _status = status; + _changed = new List<(T, PeakSpotTagState)>(); + } + + public void Do() { + _changed.Clear(); + foreach (var peak in _peaks) { + var before = peak.CapturePeakSpotTag(PeakSpotTag.CONFIRMED); + if (_status) { + if (before.IsSelected && before.Source is null && !before.IsSuppressed) { + continue; // already the person's + } + peak.SelectPeakSpotTagByPerson(PeakSpotTag.CONFIRMED); + } + else { + if (!before.IsSelected) { + continue; + } + peak.Confirmed = false; + } + _changed.Add((peak, before)); + } + } + + public void Undo() { + foreach (var (peak, before) in _changed) { + peak.RestorePeakSpotTag(before); + } + } + } +} diff --git a/src/MSDIAL5/MsdialGuiApp/View/Core/MainWindow.xaml.cs b/src/MSDIAL5/MsdialGuiApp/View/Core/MainWindow.xaml.cs index 5d8a0ea30..8239c3925 100644 --- a/src/MSDIAL5/MsdialGuiApp/View/Core/MainWindow.xaml.cs +++ b/src/MSDIAL5/MsdialGuiApp/View/Core/MainWindow.xaml.cs @@ -101,6 +101,8 @@ public MainWindow() { .Subscribe(ShowChildDialog); broker.ToObservable() .Subscribe(ShowChildDialog); + broker.ToObservable() + .Subscribe(ShowChildDialog); broker.ToObservable() .Subscribe(ShowChildSettingDialog("Select library file", height: 600, width: 400)); broker.ToObservable() diff --git a/src/MSDIAL5/MsdialGuiApp/View/Information/PeakInformationView.xaml b/src/MSDIAL5/MsdialGuiApp/View/Information/PeakInformationView.xaml index 816006d1c..85efe69d1 100644 --- a/src/MSDIAL5/MsdialGuiApp/View/Information/PeakInformationView.xaml +++ b/src/MSDIAL5/MsdialGuiApp/View/Information/PeakInformationView.xaml @@ -78,5 +78,13 @@ + + + + + diff --git a/src/MSDIAL5/MsdialGuiApp/View/MainWindowStyle.xaml b/src/MSDIAL5/MsdialGuiApp/View/MainWindowStyle.xaml index 47eb712e8..38f0eac09 100644 --- a/src/MSDIAL5/MsdialGuiApp/View/MainWindowStyle.xaml +++ b/src/MSDIAL5/MsdialGuiApp/View/MainWindowStyle.xaml @@ -539,7 +539,9 @@ O: Overannotation (Ctrl+5) - + + diff --git a/src/MSDIAL5/MsdialGuiApp/View/RibbonControl/MsdialDataVisualizationTab.xaml b/src/MSDIAL5/MsdialGuiApp/View/RibbonControl/MsdialDataVisualizationTab.xaml index 7bf95352d..b7492892a 100644 --- a/src/MSDIAL5/MsdialGuiApp/View/RibbonControl/MsdialDataVisualizationTab.xaml +++ b/src/MSDIAL5/MsdialGuiApp/View/RibbonControl/MsdialDataVisualizationTab.xaml @@ -21,6 +21,11 @@ + diff --git a/src/MSDIAL5/MsdialGuiApp/View/Setting/LipidRtCheckSettingView.xaml b/src/MSDIAL5/MsdialGuiApp/View/Setting/LipidRtCheckSettingView.xaml new file mode 100644 index 000000000..b4870c08d --- /dev/null +++ b/src/MSDIAL5/MsdialGuiApp/View/Setting/LipidRtCheckSettingView.xaml @@ -0,0 +1,89 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +