Speed up prediction and outlier detection.

Using threading for parallel processing speeds up prediction by a factor of 100: most of the time consumed was on passing data to the processes. Switched from the EllipticEnvelope to an IQR-based sigma estimation, followed by a chi^2 test. Since this happens after the PCA, the data is already decorrelated.

The chi^2 test produces the test statistics and uses the number of degrees of freedom to calculate the chi^2 variance and apply a cut at chi^2 mean + sigma * sqrt(chi^2 variance). This implies that the PCA-decorrelated data should be Gaussian, which is not true, since we know it is Poisson. Nevertheless, in the limit on which we have a lot of data (ie: the XGM intensity is above the 500 uJ cut-off), this is probably a good approximation.

pes_to_spec/model.py

@@ -937,13 +937,13 @@ class Model(TransformerMixin, BaseEstimator):
         E = np.fft.fftfreq(len(e), de)
         e_axis = np.linspace(-0.5*(e[-1] - e[0]), 0.5*(e[-1] - e[0]), len(e))
         # generate a gaussian
-        gaussian = np.exp(-0.5*(e_axis)**2/self.high_res_sigma**2)
-        gaussian /= np.sum(gaussian, axis=0, keepdims=True)
-        gaussian = np.clip(gaussian, a_min=1e-6, a_max=None)
-        gaussian_ft = np.fft.fft(gaussian)
+        #gaussian = np.exp(-0.5*(e_axis)**2/self.high_res_sigma**2)
+        #gaussian /= np.sum(gaussian, axis=0, keepdims=True)
+        #gaussian = np.clip(gaussian, a_min=1e-6, a_max=None)
+        #gaussian_ft = np.fft.fft(gaussian)
 
         H = np.mean(Z/D, axis=0)
-        N = np.absolute(gaussian_ft*V)**2
+        N = np.absolute(V)**2
         S = np.mean(np.absolute(D)**2, axis=0)
         H2 = np.absolute(H)**2
         nonzero = np.absolute(H) > 0.2