Compute spectrum similarity in meta draw (#2299)

* reverted filtering method for psms passed to FlashLFQ

* works yo

* set display spectral angle default to true

* calculate spectral angle on the fly moved to librarySpectrum class

* legend fixed

* unit test plus robust to empty list

* space

* fix string split that caused crash with accessions having hypens

---------

Co-authored-by: Alex <AlexSolivais@gmail.com>
Co-authored-by: MICHAEL SHORTREED <mrshortreed@wisc.edu>
Co-authored-by: Nic Bollis <nbollis@comcast.net>

MetaMorpheus/EngineLayer/SpectralLibrarySearch/LibrarySpectrum.cs

@@ -3,6 +3,9 @@
 using System.Linq;
 using System.Text;
 using System;
+using Easy.Common.Extensions;
+using MassSpectrometry.MzSpectra;
+using ThermoFisher.CommonCore.Data;
 
 namespace EngineLayer
 {
@@ -39,6 +42,34 @@ public LibrarySpectrum(string sequence, double precursorMz, int chargeState, Lis
             Array.Sort(XArray, YArray);
         }
 
+        /// <summary>
+        /// This function enables the spectrum angle to be computed between an individual experimental spectrum and the loaded library spectrum within MetaDraw
+        /// </summary>
+        /// <param name="librarySpectrum"></param>
+        /// <returns></returns>
+        public string CalculateSpectralAngleOnTheFly(List<MatchedFragmentIon> spectrumMatchFragments)
+        {
+            if (spectrumMatchFragments.IsNullOrEmpty())
+            {
+                return "N/A";
+            }
+
+            if(spectrumMatchFragments.IsNotNullOrEmpty()){}
+            SpectralSimilarity spectraComparison = new SpectralSimilarity(
+                spectrumMatchFragments.Select(f => f.Mz).ToArray(),
+                spectrumMatchFragments.Select(f => f.Intensity).ToArray(),
+                MatchedFragmentIons.Select(f => f.Mz).ToArray(),
+                MatchedFragmentIons.Select(f => f.Intensity).ToArray(),
+                SpectralSimilarity.SpectrumNormalizationScheme.mostAbundantPeak,
+                toleranceInPpm: 20,
+                allPeaks: true);
+            double? spectralContrastAngle = spectraComparison.SpectralContrastAngle();
+
+            return spectralContrastAngle == null
+                ? "N/A"
+                : ((double)spectralContrastAngle).ToString("F4");
+        }
+
         public override string ToString()
         {
             StringBuilder spectrum = new StringBuilder();

MetaMorpheus/GuiFunctions/MetaDraw/MetaDrawLogic.cs

@@ -132,7 +132,14 @@ public void DisplaySpectrumMatch(PlotView plotView, PsmFromTsv psm, ParentChildS
 
             spectraFile.InitiateDynamicConnection();
             MsDataScan scan = spectraFile.GetOneBasedScanFromDynamicConnection(psm.Ms2ScanNumber);
+
             LibrarySpectrum librarySpectrum = null;
+            if (SpectralLibrary != null)
+            {
+                SpectralLibrary.TryGetSpectrum(psm.FullSequence, psm.PrecursorCharge, out var librarySpectrum1);
+                librarySpectrum = librarySpectrum1;
+            }
+
             //if not crosslinked
             if (psm.BetaPeptideBaseSequence == null)
             {

MetaMorpheus/GuiFunctions/MetaDraw/MetaDrawSettings.cs

@@ -200,7 +200,7 @@ private static void InitializeDictionaries()
 
                 // lines to be written on the spectrum
                 SpectrumDescription = SpectrumDescriptors.ToDictionary(p => p, p => true);
-                SpectrumDescription["Spectral Angle: "] = false;
+                SpectrumDescription["Spectral Angle: "] = true;
             }
 
             // offset for annotation on base sequence

MetaMorpheus/GuiFunctions/MetaDraw/SpectrumMatch/PeptideSpectrumMatchPlot.cs

@@ -22,7 +22,14 @@ public PeptideSpectrumMatchPlot(OxyPlot.Wpf.PlotView plotView, PsmFromTsv psm, M
         {
             if (annotateProperties)
             {
-                AnnotateProperties();
+                if (librarySpectrum != null)
+                {
+                    AnnotateProperties(librarySpectrum);
+                }
+                else
+                {
+                    AnnotateProperties();
+                }
             }
 
             ZoomAxes(matchedFragmentIons);

MetaMorpheus/GuiFunctions/MetaDraw/SpectrumMatch/SpectrumMatchPlot.cs

@@ -10,11 +10,13 @@
 using System.Windows.Media;
 using System.Windows.Media.Imaging;
 using Chemistry;
+using Easy.Common.Extensions;
 using EngineLayer;
 using iText.IO.Image;
 using iText.Kernel.Pdf;
 using iText.Layout;
 using MassSpectrometry;
+using MassSpectrometry.MzSpectra;
 using mzPlot;
 using OxyPlot;
 using OxyPlot.Annotations;
@@ -42,7 +44,7 @@ public class SpectrumMatchPlot : Plot
         /// <param name="psm">psm to plot</param>
         /// <param name="scan">spectrum to plot</param>
         /// <param name="matchedIons">glyco ONLY child matched ions</param>
-        public SpectrumMatchPlot(OxyPlot.Wpf.PlotView plotView, PsmFromTsv psm, 
+        public SpectrumMatchPlot(OxyPlot.Wpf.PlotView plotView, PsmFromTsv psm,
             MsDataScan scan, List<MatchedFragmentIon> matchedIons = null) : base(plotView)
         {
             Model.Title = string.Empty;
@@ -115,7 +117,8 @@ protected void DrawSpectrum()
                 double mz = Scan.MassSpectrum.XArray[i];
                 double intensity = Scan.MassSpectrum.YArray[i];
 
-                DrawPeak(mz, intensity, MetaDrawSettings.StrokeThicknessUnannotated, MetaDrawSettings.UnannotatedPeakColor, null);
+                DrawPeak(mz, intensity, MetaDrawSettings.StrokeThicknessUnannotated,
+                    MetaDrawSettings.UnannotatedPeakColor, null);
             }
         }
 
@@ -127,7 +130,8 @@ protected void DrawSpectrum()
         /// <param name="strokeWidth"></param>
         /// <param name="color">Color to draw peak</param>
         /// <param name="annotation">text to display above the peak</param>
-        protected void DrawPeak(double mz, double intensity, double strokeWidth, OxyColor color, TextAnnotation annotation)
+        protected void DrawPeak(double mz, double intensity, double strokeWidth, OxyColor color,
+            TextAnnotation annotation)
         {
             // peak line
             var line = new LineSeries();
@@ -150,7 +154,8 @@ protected void DrawPeak(double mz, double intensity, double strokeWidth, OxyColo
         /// <param name="isBetaPeptide"></param>
         /// <param name="matchedFragmentIons"></param>
         /// <param name="useLiteralPassedValues"></param>
-        protected void AnnotateMatchedIons(bool isBetaPeptide, List<MatchedFragmentIon> matchedFragmentIons, bool useLiteralPassedValues = false)
+        protected void AnnotateMatchedIons(bool isBetaPeptide, List<MatchedFragmentIon> matchedFragmentIons,
+            bool useLiteralPassedValues = false)
         {
             List<MatchedFragmentIon> ionsToDisplay = !MetaDrawSettings.DisplayInternalIons
                 ? matchedFragmentIons.Where(p => p.NeutralTheoreticalProduct.SecondaryProductType == null).ToList()
@@ -168,7 +173,8 @@ protected void AnnotateMatchedIons(bool isBetaPeptide, List<MatchedFragmentIon>
         /// <param name="matchedIon">matched ion to annotate</param>
         /// <param name="isBetaPeptide">is a beta x-linked peptide</param>
         /// <param name="useLiteralPassedValues"></param>
-        protected void AnnotatePeak(MatchedFragmentIon matchedIon, bool isBetaPeptide, bool useLiteralPassedValues = false, OxyColor? ionColorNullable = null)
+        protected void AnnotatePeak(MatchedFragmentIon matchedIon, bool isBetaPeptide,
+            bool useLiteralPassedValues = false, OxyColor? ionColorNullable = null)
         {
             OxyColor ionColor;
             if (ionColorNullable == null)
@@ -196,7 +202,8 @@ protected void AnnotatePeak(MatchedFragmentIon matchedIon, bool isBetaPeptide, b
                 ionColor = (OxyColor)ionColorNullable;
             }
 
-            int i = Scan.MassSpectrum.GetClosestPeakIndex(matchedIon.NeutralTheoreticalProduct.NeutralMass.ToMz(matchedIon.Charge));
+            int i = Scan.MassSpectrum.GetClosestPeakIndex(
+                matchedIon.NeutralTheoreticalProduct.NeutralMass.ToMz(matchedIon.Charge));
             double mz = Scan.MassSpectrum.XArray[i];
             double intensity = Scan.MassSpectrum.YArray[i];
 
@@ -221,12 +228,14 @@ protected void AnnotatePeak(MatchedFragmentIon matchedIon, bool isBetaPeptide, b
                     prefix = "A-";
                 }
             }
+
             var peakAnnotation = new TextAnnotation();
             if (MetaDrawSettings.DisplayIonAnnotations)
             {
                 string peakAnnotationText = prefix + matchedIon.NeutralTheoreticalProduct.Annotation;
 
-                if (matchedIon.NeutralTheoreticalProduct.NeutralLoss != 0 && !peakAnnotationText.Contains("-" + matchedIon.NeutralTheoreticalProduct.NeutralLoss.ToString("F2")))
+                if (matchedIon.NeutralTheoreticalProduct.NeutralLoss != 0 &&
+                    !peakAnnotationText.Contains("-" + matchedIon.NeutralTheoreticalProduct.NeutralLoss.ToString("F2")))
                 {
                     peakAnnotationText += "-" + matchedIon.NeutralTheoreticalProduct.NeutralLoss.ToString("F2");
                 }
@@ -256,7 +265,9 @@ protected void AnnotatePeak(MatchedFragmentIon matchedIon, bool isBetaPeptide, b
             {
                 peakAnnotation.Text = string.Empty;
             }
-            if (matchedIon.NeutralTheoreticalProduct.SecondaryProductType != null && !MetaDrawSettings.DisplayInternalIonAnnotations) //if internal fragment
+
+            if (matchedIon.NeutralTheoreticalProduct.SecondaryProductType != null &&
+                !MetaDrawSettings.DisplayInternalIonAnnotations) //if internal fragment
             {
                 peakAnnotation.Text = string.Empty;
             }
@@ -323,11 +334,13 @@ public void ExportPlot(string path, Bitmap combinedBitmaps, double width = 700,
             switch (MetaDrawSettings.ExportType)
             {
                 case "Pdf":
-                    string tempCombinedPath = System.IO.Path.Combine(System.IO.Path.GetDirectoryName(path), "tempCombined.png");
+                    string tempCombinedPath =
+                        System.IO.Path.Combine(System.IO.Path.GetDirectoryName(path), "tempCombined.png");
                     combinedBitmaps.Save(tempCombinedPath, System.Drawing.Imaging.ImageFormat.Png);
 
                     PdfDocument pdfDoc = new(new PdfWriter(path));
-                    Document document = new(pdfDoc, new iText.Kernel.Geom.PageSize((float)width - 30, (float)height - 30));
+                    Document document = new(pdfDoc,
+                        new iText.Kernel.Geom.PageSize((float)width - 30, (float)height - 30));
 
                     ImageData sequenceAndLegendImageData = ImageDataFactory.Create(tempCombinedPath);
                     iText.Layout.Element.Image sequenceAndPtmLegendImage = new(sequenceAndLegendImageData);
@@ -372,7 +385,8 @@ protected void AnnotateLibraryIons(bool isBetaPeptide, List<MatchedFragmentIon>
             {
                 var matchedIon = SpectrumMatch.MatchedIons.FirstOrDefault(p =>
                     p.NeutralTheoreticalProduct.ProductType == libraryIon.NeutralTheoreticalProduct.ProductType
-                    && p.NeutralTheoreticalProduct.FragmentNumber == libraryIon.NeutralTheoreticalProduct.FragmentNumber);
+                    && p.NeutralTheoreticalProduct.FragmentNumber ==
+                    libraryIon.NeutralTheoreticalProduct.FragmentNumber);
 
                 if (matchedIon == null)
                 {
@@ -391,11 +405,15 @@ protected void AnnotateLibraryIons(bool isBetaPeptide, List<MatchedFragmentIon>
 
             foreach (MatchedFragmentIon libraryIon in libraryIons)
             {
-                var neutralProduct = new Product(libraryIon.NeutralTheoreticalProduct.ProductType, libraryIon.NeutralTheoreticalProduct.Terminus,
-                    libraryIon.NeutralTheoreticalProduct.NeutralMass, libraryIon.NeutralTheoreticalProduct.FragmentNumber,
-                    libraryIon.NeutralTheoreticalProduct.AminoAcidPosition, libraryIon.NeutralTheoreticalProduct.NeutralLoss);
+                var neutralProduct = new Product(libraryIon.NeutralTheoreticalProduct.ProductType,
+                    libraryIon.NeutralTheoreticalProduct.Terminus,
+                    libraryIon.NeutralTheoreticalProduct.NeutralMass,
+                    libraryIon.NeutralTheoreticalProduct.FragmentNumber,
+                    libraryIon.NeutralTheoreticalProduct.AminoAcidPosition,
+                    libraryIon.NeutralTheoreticalProduct.NeutralLoss);
 
-                mirroredLibraryIons.Add(new MatchedFragmentIon(ref neutralProduct, libraryIon.Mz, multiplier * libraryIon.Intensity, libraryIon.Charge));
+                mirroredLibraryIons.Add(new MatchedFragmentIon(ref neutralProduct, libraryIon.Mz,
+                    multiplier * libraryIon.Intensity, libraryIon.Charge));
             }
 
             AnnotateMatchedIons(isBetaPeptide, mirroredLibraryIons, useLiteralPassedValues: true);
@@ -408,7 +426,7 @@ protected void AnnotateLibraryIons(bool isBetaPeptide, List<MatchedFragmentIon>
             Model.Axes[1].LabelFormatter = DrawnSequence.YAxisLabelFormatter;
         }
 
-        protected void AnnotateProperties()
+        protected void AnnotateProperties(LibrarySpectrum librarySpectrum = null)
         {
             StringBuilder text = new StringBuilder();
             if (MetaDrawSettings.SpectrumDescription["Precursor Charge: "])
@@ -417,18 +435,25 @@ protected void AnnotateProperties()
                 text.Append(SpectrumMatch.PrecursorCharge);
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Precursor Mass: "])
             {
                 text.Append("Precursor Mass: ");
                 text.Append(SpectrumMatch.PrecursorMass.ToString("F3"));
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Theoretical Mass: "])
             {
                 text.Append("Theoretical Mass: ");
-                text.Append(double.TryParse(SpectrumMatch.PeptideMonoMass, NumberStyles.Any, CultureInfo.InvariantCulture, out var monoMass) ? monoMass.ToString("F3") : SpectrumMatch.PeptideMonoMass);
+                text.Append(
+                    double.TryParse(SpectrumMatch.PeptideMonoMass, NumberStyles.Any, CultureInfo.InvariantCulture,
+                        out var monoMass)
+                        ? monoMass.ToString("F3")
+                        : SpectrumMatch.PeptideMonoMass);
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Protein Accession: "])
             {
                 text.Append("Protein Accession: ");
@@ -438,8 +463,10 @@ protected void AnnotateProperties()
                 }
                 else
                     text.Append(SpectrumMatch.ProteinAccession);
+
                 text.Append("\r\n");
             }
+
             if (SpectrumMatch.ProteinName != null && MetaDrawSettings.SpectrumDescription["Protein: "])
             {
                 text.Append("Protein: ");
@@ -461,53 +488,67 @@ protected void AnnotateProperties()
                 }
                 else
                     text.Append(SpectrumMatch.ProteinName);
+
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Decoy/Contaminant/Target: "])
             {
                 text.Append("Decoy/Contaminant/Target: ");
                 text.Append(SpectrumMatch.DecoyContamTarget);
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Sequence Length: "])
             {
                 text.Append("Sequence Length: ");
                 text.Append(SpectrumMatch.BaseSeq.Length.ToString("F3").Split('.')[0]);
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Ambiguity Level: "])
             {
                 text.Append("Ambiguity Level: ");
                 text.Append(SpectrumMatch.AmbiguityLevel);
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Spectral Angle: "])
             {
-                text.Append("Spectral Angle: ");
                 if (SpectrumMatch.SpectralAngle != null)
-                    text.Append(SpectrumMatch.SpectralAngle.ToString());
-                else
-                    text.Append("N/A");
-                text.Append("\r\n");
+                {
+                    text.Append("Original Spectral Angle: ");
+                    text.Append(SpectrumMatch.SpectralAngle.ToString() + "\r\n");
+                }
+
+                if (librarySpectrum != null)
+                {
+                    text.Append("Displayed Spectral Angle: ");
+                    text.Append(librarySpectrum.CalculateSpectralAngleOnTheFly(this.matchedFragmentIons) + "\r\n");
+                }
             }
+
             if (MetaDrawSettings.SpectrumDescription["Score: "])
             {
                 text.Append("Score: ");
                 text.Append(SpectrumMatch.Score.ToString("F3"));
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["Q-Value: "])
             {
                 text.Append("Q-Value: ");
                 text.Append(SpectrumMatch.QValue.ToString("F3"));
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["PEP: "])
             {
                 text.Append("PEP: ");
                 text.Append(SpectrumMatch.PEP.ToString("F3"));
                 text.Append("\r\n");
             }
+
             if (MetaDrawSettings.SpectrumDescription["PEP Q-Value: "])
             {
                 text.Append("PEP Q-Value: ");

MetaMorpheus/TaskLayer/XLSearchTask/WriteFile.cs

@@ -489,6 +489,7 @@ public static void WriteProteinGlycoLocalization(Dictionary<(string proteinAcces
             Dictionary<string, HashSet<string>> localizedglycans = new Dictionary<string, HashSet<string>>();
             foreach (var item in glycoProteinParsimony.Where(p=>p.Value.IsLocalized && p.Value.MinQValue <= 0.01))
             {
+
                 var key = item.Key.proteinAccession + "#" + item.Key.proteinPosition;
                 if ( localizedglycans.ContainsKey(key))
                 {

MetaMorpheus/Test/LocalizationTest.cs

@@ -24,7 +24,7 @@ public static void TestNonSpecific()
             DigestionParams digestionParams = new DigestionParams(protease: p.Name, maxMissedCleavages: 8, minPeptideLength: 1, maxPeptideLength: 9, initiatorMethionineBehavior: InitiatorMethionineBehavior.Retain);
             var peps = prot.Digest(digestionParams, new List<Modification>(), new List<Modification>()).ToList();
 
-             Assert.AreEqual(1 + 2 + 3 + 4 + 5 + 6 + 7 + 8 + 9, peps.Count);
+            Assert.AreEqual(1 + 2 + 3 + 4 + 5 + 6 + 7 + 8 + 9, peps.Count);
         }
 
         [Test]

MetaMorpheus/Test/MatchIonsOfAllCharges.cs

@@ -491,6 +491,27 @@ public static void TestLibraryExistAfterGPTMDsearch()
 
             Directory.Delete(thisTaskOutputFolder, true);
         }
-
+
+        [Test]
+        public static void TestLibrarySpectrumCalculateSpectralAngleOnTheFly()
+        {
+
+            var librarySpectrumPath = Path.Combine(TestContext.CurrentContext.TestDirectory, @"TestData\SpectralLibrarySearch\SLSNVIAHEISHSWTGNLVTNK.msp");
+            var testLibrary = new SpectralLibrary(new List<string> { librarySpectrumPath });
+            testLibrary.TryGetSpectrum("SLSNVIAHEISHSWTGNLVTNK", 3, out var spectrum);
+
+
+            string psmsPath = Path.Combine(TestContext.CurrentContext.TestDirectory, @"TestData\SpectralLibrarySearch\SLSNVIAHEISHSWTGNLVTNK.psmtsv");
+            List<PsmFromTsv> psms = PsmTsvReader.ReadTsv(psmsPath, out List<string> warnings).Where(p => p.AmbiguityLevel == "1").ToList();
+
+            var computedSpectralSimilarity = spectrum.CalculateSpectralAngleOnTheFly(psms[0].MatchedIons);
+
+            Assert.AreEqual(1,Convert.ToDouble(computedSpectralSimilarity),0.01);
+
+
+            Assert.AreEqual("N/A", spectrum.CalculateSpectralAngleOnTheFly(new List<MatchedFragmentIon>()));
+        }
+
+
     }
 }