Merge pull request #40 from michaelmarty/v3

V3

.gitignore

@@ -135,5 +135,7 @@ __pycache__/
 unidec_bin/Example Data/ADH_unidecfiles/ADH_peak_areas.dat
 unidec_bin/Example Data/POPC_Nanodiscs_unidecfiles/POPC_Nanodiscs_1D_Mass_Defects.txt
 unidec_bin/Example Data/POPC_Nanodiscs_unidecfiles/POPC_Nanodiscs_2D_Mass_Defects.txt
+unidec_modules/thermo_reader/
+unidec_bin/Presets/Marius/
 
 

GUniDec.py

@@ -150,15 +150,14 @@ def on_open_file(self, filename, directory, skipengine=False, **kwargs):
 
         # Plot 1D
         if self.eng.config.batchflag == 0:
-            self.view.plot1.plotrefreshtop(self.eng.data.data2[:, 0], self.eng.data.data2[:, 1], "Data", "m/z",
-                                           "Intensity", "Data", self.eng.config)
+            self.makeplot1(imfit=False)
         # IM Loading and Plotting
         if self.eng.config.imflag == 1 and self.eng.config.batchflag != 1:
             self.view.controls.ctlmindt.SetValue(str(np.amin(self.eng.data.data3[:, 1])))
             self.view.controls.ctlmaxdt.SetValue(str(np.amax(self.eng.data.data3[:, 1])))
-            if self.eng.config.batchflag == 0:
-                self.view.plot1im.contourplot(self.eng.data.rawdata3, self.eng.config, xlab="m/z (Th)",
-                                              ylab="Arrival Time (ms)", title="IM-MS Data")
+            #if self.eng.config.batchflag == 0:
+            #    self.view.plot1im.contourplot(self.eng.data.rawdata3, self.eng.config, xlab="m/z (Th)",
+            #                                  ylab="Arrival Time (ms)", title="IM-MS Data")
         # tstart = time.perf_counter()
         # Load Config to GUI
         self.import_config()
@@ -332,7 +331,10 @@ def on_dataprep_button(self, e=None):
         self.eng.process_data()
         self.eng.get_auto_peak_width(set=False)
         self.import_config()
+        #print("Data Prep Time1: %.2gs" % (time.perf_counter() - tstart))
         self.view.clear_all_plots()
+        self.makeplot1(imfit=False)
+        '''
         self.view.plot1.plotrefreshtop(self.eng.data.data2[:, 0], self.eng.data.data2[:, 1], "Data Sent to UniDec",
                                        "m/z (Th)", "Normalized Intensity", "Data", self.eng.config)
         if self.eng.config.intthresh != 0 and self.eng.config.imflag == 0:
@@ -344,7 +346,7 @@ def on_dataprep_button(self, e=None):
         if self.eng.config.imflag == 1:
             self.view.plot1im.contourplot(self.eng.data.data3, self.eng.config, xlab="m/z (Th)",
                                           ylab="Arrival Time (ms)", title="IM-MS Data")
-        self.view.plot1.repaint()
+        self.view.plot1.repaint()'''
 
         self.view.SetStatusText("Data Length: " + str(len(self.eng.data.data2)), number=2)
         self.view.SetStatusText("R\u00B2 ", number=3)
@@ -425,6 +427,7 @@ def on_pick_peaks(self, e=None):
             self.makeplot4(1)
 
         self.view.SetStatusText("Peak Pick Done", number=5)
+
         self.on_score()
         pass
 
@@ -483,27 +486,68 @@ def on_auto(self, e=None):
     #
     # ...........................................
 
-    def makeplot1(self, e=None):
+    def makeplot1(self, e=None, intthresh=False, imfit=True):
         """
         Plot data and fit in self.view.plot1 and optionally in plot1fit
         :param e: unused event
         :return: None
         """
         if self.eng.config.batchflag == 0:
             tstart = time.perf_counter()
+            leg = False
+
             if self.eng.config.imflag == 1:
-                self.view.plot1fit.contourplot(self.eng.data.fitdat2d, self.eng.config, xlab="m/z (Th)",
-                                               ylab="Arrival Time (ms)", title="IM-MS Fit")
-            self.view.plot1.plotrefreshtop(self.eng.data.data2[:, 0], self.eng.data.data2[:, 1], "Data and UniDec Fit",
-                                           "m/z (Th)", "Normalized Intensity", "Data", self.eng.config, nopaint=True)
-            try:
-                self.view.plot1.plotadd(self.eng.data.data2[:, 0], self.eng.data.fitdat, 'red', "Fit Data")
-                pass
-            except:
-                pass
+                try:
+                    self.view.plot1im.contourplot(self.eng.data.data3, self.eng.config, xlab="m/z (Th)",
+                                                  ylab="Arrival Time (ms)", title="IM-MS Data")
+                except:
+                    pass
+                if imfit:
+                    try:
+                        self.view.plot1fit.contourplot(self.eng.data.fitdat2d, self.eng.config, xlab="m/z (Th)",
+                                                       ylab="Arrival Time (ms)", title="IM-MS Fit")
+                    except:
+                        pass
+
+            if self.eng.config.reductionpercent<0:
+                print("Making Dot Plot")
+                data2 = ud.dataprep(self.eng.data.rawdata, self.eng.config, peaks=False, intthresh=False)
+                self.view.plot1.plotrefreshtop(data2[:, 0], data2[:, 1],
+                                               "Data and UniDec Fit",
+                                               "m/z (Th)", "Normalized Intensity", "Data", self.eng.config,
+                                               nopaint=True)
+                self.view.plot1.plotadddot(self.eng.data.data2[:,0], self.eng.data.data2[:,1], 'blue', "o", "Peaks")
+
+                try:
+                    if len(self.eng.data.fitdat) > 0:
+                        self.view.plot1.plotadddot(self.eng.data.data2[:, 0], self.eng.data.fitdat, 'red', "s", "Fit Data")
+                        leg = True
+                    pass
+                except:
+                    pass
+
+            else:
+                self.view.plot1.plotrefreshtop(self.eng.data.data2[:, 0], self.eng.data.data2[:, 1], "Data and UniDec Fit",
+                                               "m/z (Th)", "Normalized Intensity", "Data", self.eng.config, nopaint=True)
+
+                if self.eng.config.intthresh != 0 and self.eng.config.imflag == 0 and intthresh:
+                    self.view.plot1.plotadd(self.eng.data.data2[:, 0],
+                                            np.zeros_like(self.eng.data.data2[:, 1]) + self.eng.config.intthresh, "red",
+                                            "Noise Threshold")
+                    leg=True
+
+                try:
+                    if len(self.eng.data.fitdat) > 0:
+                        self.view.plot1.plotadd(self.eng.data.data2[:, 0], self.eng.data.fitdat, 'red', "Fit Data")
+                        leg = True
+                    pass
+                except:
+                    pass
             if self.eng.config.aggressiveflag != 0 and len(self.eng.data.baseline) == len(self.eng.data.fitdat):
                 self.view.plot1.plotadd(self.eng.data.data2[:, 0], self.eng.data.baseline, 'blue', "Baseline")
-            self.view.plot1.add_legend()
+            if leg:
+                self.view.plot1.add_legend()
+            self.view.plot1.repaint()
             tend = time.perf_counter()
             print("Plot 1: %.2gs" % (tend - tstart))
 

metaunidec/gui_elements/ud_cont_meta.py

@@ -301,7 +301,7 @@ def __init__(self, parent, config, pres, panel, iconfile):
         i += 1
 
         self.ctlpoolflag = wx.RadioBox(panel2b, label="m/z to Mass Transformation",
-                                       choices=["Integration", "Interpolation"])
+                                       choices=["Integrate", "Interpolate", "Smart"])
         gbox2b.Add(self.ctlpoolflag, (i, 0), span=(1, 2), flag=wx.ALIGN_CENTER_VERTICAL)
         i += 1
 
@@ -813,7 +813,9 @@ def setup_tool_tips(self):
         self.ctlpoolflag.SetToolTip(wx.ToolTip(
             "Sets type of conversion from m/z to mass.\nIntegration:\n\tEach m/z bin goes to the nearest mass bin"
             "\n\tBest for undersampled masses\nInterpolation:\n\tEach mass value interpolates its value in m/z space"
-            "\n\tBest for oversampled mass data"))
+            "\n\tBest for oversampled mass data"
+            "\nSmart:\n\tDetermines interpolation vs. integration automatically\n\tfrom the data and blends the two."
+            "\n\tUse Smart in most cases (default). "))
         self.ctlnorm2.SetToolTip(wx.ToolTip(
             "Sets normalization of extraction.\nMaximum will normalize so that the maximum value of a timepoint is 100%."
             "\nSum will normalize so that the sum of all peaks of a timepoint is 100%."

readme.md

@@ -194,6 +194,28 @@ The main GUI class is GUniDec.UniDecApp.
 
 ## Change Log
 
+v.4.2.0
+
+**Added Smart Transform**. From the beginning, UniDec has asked users to decide how it should perform the nonlinear transform from m/z to mass (Integration or Interpolation). Now, the Smart Transform decides for the user based on the density of the data. Where the data is sparse relative to the mass sampling, it will use interpolation to fill in the gaps. Where data is dense relative to the mass sampling, it will use integration to ensure that all the data in-between is accounted for. Smart Transform is located under the Additional Deconvolution Parameters tab and has now been made the default. In some preliminary testing, the Smart Transform showed more robust peak areas that were less sensitive to the mass sampling. The peak heights may change some, but the peak areas should be a lot more reliable (remember Ctrl+I to integrate). Note: you may not notice much difference, but using the Smart Transform should avoid some glitches that are hard to spot. Let me know how it goes!
+
+**Major improvements in 2D plotting speed** due to behind-the-scenes changes and optimzation of the plotting code.
+
+**Sparse data support**. Improvements to the algorithm to support sparse data, including just peak centroids. This fundamentally changes the underlying algorithm in minor ways. You may notice minor differences but likely will not. Let me know if it causes any problems.
+
+Better handling of click and zoom at the edges of the plot. Now, it will default to the edge if you go off the plot. 
+
+Added drag and drop for loading the state from a zip file.
+
+Upgraded pymzml version for improved mzML compatibility and made mzML import more error tolerant.
+
+Fixed issue with discrete plot and data that was non-uniformly sampled.
+
+Fix to bug in calculating the error from the weighted standard deviation of charge state masses. 
+
+Fixed bug in Load and Save State with non-text files.
+
+Fixed IM-MS plotting bugs.
+
 v.4.1.2
 
 **Added button in UniDec to turn off data normalization.** Note: the beta values will be automatically scaled to match this. 

unidec.py

@@ -104,7 +104,8 @@ def open_file(self, file_name, file_directory=None, *args, **kwargs):
             newname = os.path.join(os.getcwd(), self.config.outfname + "_imraw.txt")
         if not os.path.isfile(newname):
             try:
-                shutil.copy(file_directory, newname)
+                # shutil.copy(file_directory, newname)
+                np.savetxt(newname, self.data.rawdata)
             except Exception as e:
                 pass
 
@@ -113,7 +114,7 @@ def open_file(self, file_name, file_directory=None, *args, **kwargs):
         else:
             refresh = False
 
-        if os.path.isfile(self.config.infname) and not refresh:
+        if os.path.isfile(self.config.infname) and not refresh and self.config.imflag == 0:
             self.data.data2 = np.loadtxt(self.config.infname)
             self.config.procflag = 1
         else:
@@ -329,7 +330,11 @@ def unidec_imports(self, efficiency=False, everything=False):
         self.pks = peakstructure.Peaks()
         self.data.massdat = np.loadtxt(self.config.outfname + "_mass.txt")
         self.data.ztab = np.arange(self.config.startz, self.config.endz + 1)
-        self.config.massdatnormtop = np.amax(self.data.massdat[:, 1])
+        try:
+            self.config.massdatnormtop = np.amax(self.data.massdat[:, 1])
+        except:
+            self.data.massdat = np.array([self.data.massdat])
+            self.config.massdatnormtop = np.amax(self.data.massdat[:, 1])
         if not efficiency:
             self.data.massgrid = np.fromfile(self.config.outfname + "_massgrid.bin", dtype=float)
 
@@ -399,8 +404,11 @@ def pick_peaks(self):
         self.export_config()
         # Detect Peaks and Normalize
         peaks = ud.peakdetect(self.data.massdat, self.config)
-        # print(peaks)
-        self.setup_peaks(peaks)
+        if len(peaks) > 0:
+            self.setup_peaks(peaks)
+        else:
+            print("No peaks detected", peaks, self.config.peakwindow, self.config.peakthresh)
+            print("Mass Data:", self.data.massdat)
 
     def setup_peaks(self, peaks):
         if self.config.peaknorm == 1:
@@ -422,9 +430,12 @@ def setup_peaks(self, peaks):
         # Generate Intensities of Each Charge State for Each Peak
         mztab = ud.make_peaks_mztab(self.data.mzgrid, self.pks, self.config.adductmass)
         # Calculate errors for peaks with FWHM
-        ud.peaks_error_FWHM(self.pks, self.data.massdat)
         try:
-
+            ud.peaks_error_FWHM(self.pks, self.data.massdat)
+        except Exception as e:
+            print("Error in FWHM calculation:", e)
+            print(self.data.massdat)
+        try:
             ud.peaks_error_mean(self.pks, self.data.massgrid, self.data.ztab, self.data.massdat, self.config)
         except Exception as e:
             print("Error in error calculations:", e)
@@ -1344,6 +1355,7 @@ def pks_fscore(self):
                     fscore2 = self.score_minimum(height, umin)
                     # print("Fscore5", fscore2, umin, height)
                     p.fscore *= fscore2
+
     '''
     def tscore(self):
         try:
@@ -1368,34 +1380,35 @@ def dscore(self, xfwhm=2, pow=2):
 
         :return:
         """
-
-        self.pks_mscore(xfwhm=xfwhm, pow=pow)
-        self.pks_uscore(pow=pow, xfwhm=xfwhm)
-        self.pks_csscore(xfwhm=xfwhm)
-        self.pks_fscore()
-        #self.tscore()
-        tscores = []
-        ints = []
-        for p in self.pks.peaks:
-            scores = np.array([p.mscore, p.uscore, p.fscore, p.cs_score])  # p.rsquared,
-            p.dscore = np.product(scores)
-            p.lscore = np.product(scores[:-1])
-            tscores.append(p.dscore)
-            ints.append(p.height)
-            print("Mass:", p.mass,
-                  "Peak Shape:", round(p.mscore * 100, 2),
-                  "Uniqueness:", round(p.uscore * 100, 2),
-                  # "Fitting R^2", round(p.rsquared * 100, 2),
-                  "Charge:", round(p.cs_score * 100, 2),
-                  "FWHM:", round(p.fscore * 100, 2),
-                  "Combined:", round(p.dscore * 100, 2))
-            # print(p.intervalFWHM)
-        ints = np.array(ints)
-        self.pks.uniscore = ud.weighted_avg(tscores, ints ** pow) * self.config.error
-        print("R Squared:", self.config.error)
-        #print("TScore:", self.data.tscore)
-        print("Average Peaks Score (UniScore):", self.pks.uniscore)
-
+        if len(self.pks.peaks) > 0:
+            self.pks_mscore(xfwhm=xfwhm, pow=pow)
+            self.pks_uscore(pow=pow, xfwhm=xfwhm)
+            self.pks_csscore(xfwhm=xfwhm)
+            self.pks_fscore()
+            # self.tscore()
+            tscores = []
+            ints = []
+            for p in self.pks.peaks:
+                scores = np.array([p.mscore, p.uscore, p.fscore, p.cs_score])  # p.rsquared,
+                p.dscore = np.product(scores)
+                p.lscore = np.product(scores[:-1])
+                tscores.append(p.dscore)
+                ints.append(p.height)
+                print("Mass:", p.mass,
+                      "Peak Shape:", round(p.mscore * 100, 2),
+                      "Uniqueness:", round(p.uscore * 100, 2),
+                      # "Fitting R^2", round(p.rsquared * 100, 2),
+                      "Charge:", round(p.cs_score * 100, 2),
+                      "FWHM:", round(p.fscore * 100, 2),
+                      "Combined:", round(p.dscore * 100, 2))
+                # print(p.intervalFWHM)
+            ints = np.array(ints)
+            self.pks.uniscore = ud.weighted_avg(tscores, ints ** pow) * self.config.error
+            print("R Squared:", self.config.error)
+            # print("TScore:", self.data.tscore)
+            print("Average Peaks Score (UniScore):", self.pks.uniscore)
+        else:
+            print("No peaks present")
 
     def filter_peaks(self, minscore=0.4):
         # w = deepcopy(self.config.peakwindow)
@@ -1497,9 +1510,24 @@ def open_test_spectrum(self, masslist=None, n=1, **kwargs):
         np.savetxt(os.path.join(testdir, fname), data)
 
         self.open_file(fname, testdir, **kwargs)
+        self.config.maxmz, self.config.minmz = "", ""
         self.data.ztab = ztab
         pass
 
+    def pass_data_in(self, data, n=1, **kwargs):
+        testdir = os.path.join(self.config.UniDecDir, "TestSpectra")
+        if not os.path.isdir(testdir):
+            os.mkdir(testdir)
+
+        fname = "test_" + str(n) + ".txt"
+        path = os.path.join(testdir, fname)
+        print("Creating temp file:", path)
+        np.savetxt(path, data)
+
+        self.open_file(fname, testdir, **kwargs)
+        self.config.maxmz, self.config.minmz = "", ""
+        self.config.default_isotopic_res()
+
     def get_spectrum_peaks(self, threshold=0.05, window=None):
         if window is None:
             window = self.config.mzsig

unidec_bin/Example Data/ADH_unidecfiles/ADH_conf.dat

@@ -34,7 +34,7 @@ rawflag 0
 adductmass 1.007276467
 nativezub 100.0
 nativezlb -100.0
-poolflag 1
+poolflag 2
 accvol 0.0
 peakshapeinflate 1.0
 noiseflag 0.0

unidec_bin/Example Data/BSA_unidecfiles/BSA_conf.dat

@@ -8,17 +8,17 @@ startz 10
 zzsig 1.0
 psig 0.0
 beta 0.0
-mzsig 0.0
+mzsig 1.50166
 psfun 0
 discreteplot 0
 massub 150000.0
 masslb 1000.0
 msig 0.0
-molig 0.0
+molig 1.0
 massbins 1.0
 mtabsig 0.0
-minmz 2715.4349487890913
-maxmz 8616.21283935621
+minmz 3617.539730726923
+maxmz 7066.941420450001
 subbuff 0.0
 subtype 2
 smooth 0.0