High-throughput XPS spectrum modeling with autonomous background subtraction for 3d5/2 peak mapping of SnS

ABSTRACT We propose a fitting model that automatically conducts the background subtraction during high-throughput peak fitting. The fitting model consists of the pseudo-Voigt mixture model and the ramp-sum background model, and each model represents the peak and background component, respectively. The optimization of the fitting model is performed by the spectrum adapted ECM algorithm that enables us to perform the peak fitting and background subtraction simultaneously through the high-throughput calculation. Application of the proposed model to the synthetic spectral data showed appropriate decomposition of the peak and background component. We also applied the proposed model to 3721 spectral data collected from the SnS sheet by X-ray photoelectron spectroscopy. The spectral data from the SnS sheet were successfully decomposed to the component of Sn $$3{d_{3/2}}$$3d3/2 peak, Sn $$3{d_{5/2}}$$3d5/2 peak and background, respectively. As the spectrum adapted ECM algorithm can efficiently analyze a large amount of spectral data, we can obtain the color map showing spatial distribution of Sn(II) and Sn(IV) using the parameter of Sn $$3{d_{5/2}}$$3d5/2 peak. The proposed model supports us to obtain the insightful spatial distribution of peak component that has been difficult to obtain by conventional peak fitting. GRAPHICAL ABSTRACT


Introduction
Peak fitting is one of the important approaches in spectroscopy for investigating physical and chemical properties of materials and devices [1][2][3][4].In order to perform high-accurate peak fitting, Nagata et al. [5] introduced the Bayesian spectroscopy using the exchange Monte Carlo method [6].The Bayesian spectroscopy enables us to conduct accurate analysis by incorporating the physical constraints for the shape of peak and background component into the probabilistic model, and to estimate the number of peak components based on the value of marginal likelihood.Although large calculation cost is generally required in the exchange Monte Carlo method, Bayesian spectroscopy has been applied to various spectroscopic techniques, such as core-level X-ray photoelectron spectroscopy (XPS) [7,8], an absorption spectrum [9], a photoluminescence spectrum [10], X-ray magnetic circular dichroism spectrum [11], and a coherent phonon signal [12].
Recently, measurement informatics, which is a research area to integrate techniques in measurements and informatics, has been attracting attention for the efficient measurement and analysis of large amount of spectral data in material science [13][14][15][16][17][18][19][20][21].In spectroscopy, performing peak fitting on a large amount of spectral data has following significances: (1) it is possible to extract low-dimensional features from a high-dimensional spectral data set; (2) lowdimensional features are worth descriptors for understanding property and functionality of materials.In particular, the low-dimensional features represented by the model parameters, such as a Lorentzian mixture model, contribute to analyze the statistical properties of large amount of spectral data.For example, a color map of the peak-position (i.e.location parameter of Lorentzian distribution) reveals the spatial characteristics of the peak shift, which is difficult to determine by pinpoint analysis on material surface.However, although it is possible to collect a large amount of spectral data by recent development of spectroscopic devices, it is almost impossible to complete the analysis of over thousands sets of spectral data by the conventional procedure of the peak fitting.Especially, the background subtraction and the parameter estimation of the fitting model is a seriously time-consuming work in the analysis.As researchers generally adopt suitable background component of spectral data and parameters of fitting model according to previous reports and their experiences, conventional peak fitting strongly resorts to the manual trial and error.Such a manual trial-anderror is extremely inefficient when analyzing large amount of spectral data, and it is difficult to systematically reproduce individual result of the analysis.Thus, it is necessary to develop the method that can perform the background subtraction and peak fitting simultaneously and quickly even for a large amount of spectral data.In the case of XPS, the background subtraction is conducted prior to the parameter estimation by using various methods, such as Shirley [22], Tougaard [23] and linear method.As the parameter estimation is affected by the background subtraction, enormous trialand-error is required in the peak fitting.Herra-Gomes et al. [24] proposed the approach for the background subtraction during the peak fitting to improve the efficiency of peak fitting.This approach was used in the XPS analysis [25,26], whereas the approach for high-throughput analysis has not been developed.
In this study, we propose a fitting model that incorporates the background component and is applicable in the high-throughput peak fitting.By applying this model, peak fitting and background subtraction are efficiently conducted at the same calculation.Calculation is conducted by using the framework of the spectrum adapted expectation-maximization (EM) algorithm which is a high-throughput peak fitting method [19].The spectrum adapted EM algorithm has advantages that are monotonous increase and convergence of likelihood in iterative calculation [27] and sufficiently quick calculation speed in practical peak fitting [19].However, the spectrum adapted EM algorithm has the limitation that it can only use the Gaussian mixture model.To overcome a limitation of available fitting model, the spectrum adapted EM algorithm was extended by introducing the expectation-conditional maximization (ECM) algorithm [28] to be the spectrum adapted ECM algorithm [29], which enables us to use various fitting models, such as the pseudo-Voigt mixture model and Doniach-Šunjić-Gauss mixture model.We further extend the framework for high-throughput peak fitting by introducing the fitting model including the background component.The proposed model is applied to the analysis of synthetic and experimental spectral data.In synthetic spectral data analysis, we show that the proposed model is applicable for synthetic spectral data with various noise levels.In experimental spectral data analysis, we obtain the color map indicating a spatial distribution of Sn(II) and Sn(IV) by the analysis of 3721 spectral data collected by surface mapping of the SnS sheet.

Pseudo-Voigt mixture model
We introduced the pseudo-Voigt mixture model (PVMM) for representing the line shape of the spectral data.When it is not necessary to consider the asymmetry of peak shape, as the PVMM is easy to calculate and is a popular model for peak fitting, we adopt PVMM in the proposed method.Each model for measurement energy steps x n . is constructed by the pseudo-Voigt distribution (pseudoVoigtðx n jμ; σ; g; ηÞ [30] as follows: where (2) and Normalðx n jμ; σÞ ¼ 1 ffi ffi ffi ffi ffiffi 2π p σ exp À 1 2 Here, μ is is the location parameter specifying the peak of spectrum, g (g > 0) is the scale parameter specifying the half-width at half-maximum (HWHM), σ is standard deviation (σ > 0), and η is mixing parameter Since the number of parameters is reduced, we use Equation ( 4) to be the component of PVMM.
In the peak fitting, as these distributions have a wide tail and the measurement energy steps (x n ) are finite, each distribution (Equation 4) is expressed in a truncated distribution in the range of x n as follows: P PV ðx n jμ; σ; ηÞ ¼ pseudoVoigtðx n jμ; σ; ηÞ P N n¼1 pseudoVoigtðx n jμ; σ; ηÞ ; (5) respectively.The truncated pseudo-Voigt distribution (Equation 5) is defined as normalized probability density in the range of measurement energy interval (x 1 À x N ).This normalized probability density means Equation ( 5) is normalized to 1 within the range.By using the truncated pseudo-Voigt distribution, it is possible to perform appropriate analysis on spectral data in various ranges.The PVMM are respectively defined as a linear superposition of k-th normalized pseudo-Voigt distribution (P PV ðx n jμ; σ; ηÞ) as follows: PVMMðx n jλ; μ; σ; ηÞ ¼ where K is the number of components (K > 0), and λ k is the mixing coefficient of k-th component of the peaks (0 � λ k � 1, and . Also, PVMM becomes GMM and LMM in η k = 0 and 1, respectively.

Linear background model
We propose a linear background model based on the probability density function for the spectral data.Our background model consists of two types of probability density functions: uniform and triangular distribution, representing a linear slope and a baseline, respectively, as shown in Figure 1.The proposed model is equivalent to the conventional linear background subtraction, and it is unnecessary to assume hyperparameters such as the side-edge intensity of the linear background.Moreover, joint optimization of both spectral and background models is achieved in the proposed method.
The linear background model is defined by the uniform distribution (P uni ðx n Þ) for a baseline and triangular distribution (P tri ðx n Þ) for linear slope as follows: where and Our background model only has the parameters of the mixing ratio of λ 1 and λ 2 (Equation 7), which are positive values, and the sum of the mixing ratio is unity.The mixing ratio is optimized analytically in the framework of the EM algorithm with low computational cost.Hence, the background model does not need any sophisticated optimizer and preconditioning process.Then, it is suitable for high-throughput analysis.
As the background component approximates the shape of the background to be a monotonically increasing straight line, it is emphasized that Equation (7) does not explicitly incorporate the physical constraints.This background model (Equation 7) can be extended to the ramp-sum model, which is a monotonically increasing nonlinear function.Details are shown in the supporting information.

Fitting model
Here, we define the pseudo-Voigt mixture model with Background (PVMMB) as a fitting model.When the background component is represented by single interval, it corresponds to the background subtraction by the linear method, which is generally used in the peak fitting.In this case, the background component is represented by summing P uni ðx n Þ and P tri ðx n Þ from x 1 to x N (Figure 1(a)).Hence, the fitting model using linear background is represented as follows: PVMMBðx n jλ; μ; σ; ηÞ ¼ As Equation (10) represents the peak and background component in simple form, we apply Equation ( 10) as the fitting model of the spectrum data in this study.

Parameter estimation of the fitting model
Peak fitting using Equation ( 14) is performed by the spectrum adapted ECM algorithm [29].This is because that the peak fitting by the spectrum adapted EM algorithm [19] is difficult to analytically derive the updating rule of the parameters in the maximization (M)-step for PVMMB.This difficulty also occurs in the case of GMM at η k = 0 in PVMMB because we use a truncated distribution to be the component of mixture model.
At first of the calculation, the parameters of fitting model (λ, μ, σ, and η) are initialized, and the initial value of weighted log-likelihood is calculated as follows: Equation ( 11) assumes that the error function of observed intensity is the Poisson distribution.Details of the derivation of Equation ( 11) are shown in supporting information.In the expectation (E)-step, the responsibilities (γ z nk ,γ z nKþ1 and γ z nKþ2 ) are calculated by using current parameter values (λ old , μ old , σ old , and η old ) as follows: and Here, the responsibilities correspond to posterior probabilities when the measurement steps (x) are observed [32].For example, γ z nk is a latent variable associated with x n .In PVMMB, x n is assumed to be generated from one of the pseudo-Voigt and rampsum background components.γ z nk represents the component that generated x n ; i.e. γ z nk is equal to be 1 when x n is generated from k-th component otherwise γ z nk is equal to be 0. Theoretical details are described in Bishop [32] and McLachlan and Krishnan [33].Then, in the conditional maximization (CM)-step for λ, the parameters of λ k are updated from λ old k to λ new k by using current responsibilities (γ z nk ) and intensities (w) that correspond to the measurement steps (x) as follows: where N wγ z nk ¼ P N n¼1 w n γ z nk .λ l in K+l-th and l-th are also calculated from λ old l to λ new l as follows: and where w n γ z nKþ2 , respectively.
In CM-step for μ, the parameters of μ k are updated from μ old k to μ new k .μ new k are obtained by the maximization of the expectation of the complete-data log-likelihood with the weight of the intensity as follows: Here, λ k , σ k and η k are fixed to be λ old k , σ old k and η old k , respectively.This maximization is conducted by using the Brent's method [34].In CM-step for σ k are updated from σ old k to σ new k by maximizing Equation (22)  and σ new k , respectively.After the CM-steps, when the distance of the weighted log-likelihood values (Equation 15) after the update and that immediately before the update is less than 1 � 10 À 8 , the calculation is returned to the E-step.Otherwise, the calculation is determined to be converged and the parameters at that time are adopted to be the solution.Moreover, the calculation is terminated when the iterative calculation reaches 10,000 times in this study.
The most crucial difference between the proposed and conventional regression methods is the assumption of error functions in statistical formalism.Conventional regression assumes the Gaussian error function.On the other hand, the proposed method assumes the Poisson error function.Both are similar in the case of large intensity data but differ in the case of low-intensity data.Joint optimization of both spectral and background models is also remarkable in the high-throughput analysis.In addition, the numerical stability of the proposed method without any preconditioning and hyperparameter tuning is appropriate for applying the high-throughput analysis.
Because the EM algorithm, on which the proposed method is based, guarantees a monotonic increase of the likelihood during calculation, the proposed method converges to the local optimum solution without complex preprocessing and careful choice of initial values.Therefore, the proposed method is numerically robust for the peak fitting of large amounts of spectral data, and the robustness is the advantage of the proposed method for high-throughput analysis against the conventional method, such as non-linear leastsquare optimization.

Generation of the synthetic data
Synthetic spectral data are generated from three-components PVMM and linear background.The data is represented as follows: where w ðspvbÞ are the intensity of synthetic spectral data corresponding to x textitn and spvb mean the synthetic pseudo-Voigt and background distribution.In the XPS, as the intensity is observed as the count data of the number of photoelectrons at each measurement step, we assume that w ðspvbÞ n are randomly generated from the Poisson distribution (PoissonðXj�Þ) as follows: where X is a discrete random variable that is an integer starting from 0 (0, 1, 2, . . .). � is the parameter representing the mean of Poisson distribution.We defined � ðspvbÞ for w ðspvbÞ n , and � ðspvbÞ was obtained as follows: respectively.d is a constant (d � 0) for the number of intensity related to a measurement time.Also, the w ðspvbÞ n in the noise-free is generated as follows: In D ðspvbÞ , the steps (x n ) were collected from the range [0:100] in steps of 0.5, so that the total number of steps was 201.The parameters were set μ

Initial condition of the synthetic data analysis
Spectrum adapted ECM algorithm was performed to each synthetic spectral data.At the calculation, the number of peaks (K) of the fitting model were K = 3, and the initial values of the parameters (μ, σ, η and λ) were randomly collected from the range [20:80], [1:10], [0:0.5] and [0:1] satisfying P K k¼1 ¼ 1, respectively.Table 1 summarizes the initial values of parameters.Calculation was conducted by using our own source code developed in R (http://cran.r-project.org/).R is an open-source programming language and software environment for the statistical analysis.The computer carrying out the calculations had an Intel(R) Xeon W CPU with ten cores at 3 GHz with 128 GB memory.

Result of the synthetic data analysis
Figure 2 shows the obtained fitting curves, and each value of parameters is shown in Table 2. Fitting curves were well-fitted to the synthetic data, and the  2).The values of σ were also close to the true value at d = 4, 5 and noise-free.In contrast, values of σ deviated about � 10% from the true values at d = 3.The values of η were different from the true value even when the spectral data includes moderate to large noise (d = 4, 3).The values of λ were close to the true value at d = 4, 5 and noise free.In d = 3, estimated values of λ deviated about � 10-40% from the true values.Thus, although the peak position (μ) can be estimated accurately, the accuracy of other parameters (σ, η and λ) would decrease in the analysis of very noisy spectral data (d = 3).Variation of the loglikelihood showed monotonically increasing during the iterative calculation, and occasionally showed a large improvement (Figure S1).The values of each parameter fluctuated greatly at the beginning of the iterative calculation (about 1-20 iterations), and then these values became to be stable after about 200-500 iterative calculations at d = 4, 5 and noise free (Figure S2-S4).In contrast, over 700 iterations were required for the parameters to be stable in d = 3 (Figure S5).

Experimental data
A set of experimental data was collected by the surface mapping of the tin monosulfide (SnS) monolayer sheet [35].Since chemical characterization of the nanoscale microstructure is essential to identify the mechanism for film growth of this new monolayer material, the core-level photoelectron spectral mapping was performed.For convenience, we named the set of experimental data from the surface mapping of the SnS sheet to be SnS-data.SnS-data consists of 3721 spectral data from 61 � 61 measurement points.The  SnS was grown on a mica substrate in a three-zone heating tube furnace by physical vapor deposition.More details on the SnS sheet are given in Kawamoto et al. [35].The spectral data were acquired using a core-level photoelectron spectromicroscopy, called '3D nano-ESCA' [36,37], installed at the University of Tokyo outstation beamline, BL07LSU at SPring-8.The photon energy of the synchrotron radiation beam used for measurements was 1000 eV.
In order to obtain a large amount of spectral data in a short time, analyzed spectral data were collected using the hemispherical analyzer in fixed mode.The background of these spectral data is affected by unevenness in the detection efficiency of photoelectrons derived from the micro-channel plate and the background derived from secondary electrons.Therefore, in this case, the background subtraction method that assumes only the background derived from secondary electrons, such as the Tougaard method and Shirley method, is inappropriate.As the background subtraction for such spectral data is generally conducted by the linear method, it is reasonable to use a proposed method.

Initial condition of the analysis
The analysis of SnS-data was conducted to the Sn 3d core-level spectrum.As it is not necessary to consider the asymmetry of Sn 3d 3=2 and Sn 3d 5=2 peaks, PVMM is reasonable for the peak fitting of SnS-data.The number of peaks (K) of SnS-data was K = 2.The initial values of μ, σ and η were set to be μ 1 = 230, μ 2 = 233, η 1 = 0.5 and η 2 = 0.5, σ 1 = 1 and σ 2 = 1, respectively.The initial values of λ were randomly collected from the range [0:1] satisfying P K k¼1 λ k = 1, respectively.The calculation environment is the same as for the synthetic data analysis.

Result of the SnS-data analysis
Figure 3 shows the typical fitting curves attributed to Sn 2 S 3 , SnS islands, and electrodes corresponding to ad in Figure 4.The fitting curves on the SnS islands (#1087, #1897, and #2465) successfully detected Sn 3d 3=2 and Sn 3d 5=2 peaks.As the spectral data from the electrode (#2907) show no clear peaks, the Sn 3d peaks were not identified, and estimated peak components showed a very broad peak to fit the outline of the spectral data (Figure 3).The background model successfully identified the tendency of the spectral data.The values of L PVMMB monotonically increased during calculation, and occasionally showed a large improvement at early stage of the iteration (Figure S6).The values of parameter fluctuated greatly during about 1-200 iterative calculations, and then these values were stable after about 200-400 iterative calculations (Figures S7-S10).The proposed method can perform the peak fitting without significantly deviating from the physical constraints, such as the ratio of λ for Sn 3d 5=2 and 3d 3=2 peaks to be 1/0.667and η for 3d 5=2 and 3d 3=2 peaks to be the same values (Figure S11 and S12).As the proposed method does not explicitly incorporate the physical constraints, the estimated parameters inevitably include variations due to noise in spectral data.
Figure 4 shows a color map of the total intensity of each spectral data (Figure 4(a)) and the binding energy of Sn 3d 5=2 peak in the range of 485.8-486.6 eV (Figure 4(b)).The color map of binding energy shows that the peaks on central part of the SnS island were about 0.2 eV lower than the peaks on the part of Sn 2 S 3 .In contrast, the peaks on the upper side of the SnS island were about 0.3 eV higher binding energy than those on the Sn 2 S 3 .There are the part indicating similar binding energy peaks to the SnS island on the part between the electrodes, whereas the part showing higher binding energy corresponding to the upper side of the SnS island are not observed.As the Sn 3d 5=2 peaks at approximately 485.9 and 486.5 eV correspond to Sn(II) and Sn(IV), respectively [35,38], the contrast of binding energy on the color map is consistent with the mixture of Sn(II) and Sn(IV) on the Sn 2 S 3 and the predominance of Sn(II) on the SnS island.Moreover, at the upper side of the SnS island, it is interpreted that Sn(IV) is predominant and SnS would not be formed.
Figure 5 shows a histogram of the calculation time in the analysis of SnS-data.Analysis of each spectral data was completed in the range of 2.7 to 132.5 seconds and, the total calculation time was 34.5 hours (mean = 33.4sec, median = 7.0 sec).When clear peaks were observed, for example #1087, #1897, #2465 in Figure (3), the calculations were completed less than 15 seconds.In contrast, the calculation time was required over 100 seconds for the spectral data with no clear peaks, such as #2907 in Figure (3).The number of iterative calculations required about 1500 times at most for the data in which a clear peak is observed.In contrast, over 8000 iterations were required to reach the convergence criterion in the data without clear peaks.For the spectral data without clear peaks, many iterative calculations were required to converge the parameter (Figure S7-S10), and the calculation cost increases accordingly.

Discussion and implication
The analysis efficiency of the peak fitting has been improved by replacing the time-consuming procedures with an optimization problem of the PVMMB consisting of the PVMM and background model.The advantage of the peak fitting by using PVMMB is that the peak assignment and background subtraction can be performed simultaneously and effectively.In addition, the spectrum adapted ECM algorithm can perform the peak fitting not only by using the pseudo-Voigt model but also by using a model that represents other peak shapes, such as Doniach-Šunjić-Gauss mixture model [29].In conventional peak fitting, the background subtraction generally requires time-consuming work.For  example, the intensities of background at the ends of measurement point are heuristically selected by a manual trial-and-error, and shape of peak is affected by this trial-and-error.In contrast, proposed procedure conducts the peak fitting without such manual trialand-error because these intensities are optimized in the calculation by spectral adapted ECM algorithm.Moreover, peak fitting is performed avoiding a physically inappropriate situation where the spectral data includes negative intensity due to inappropriate background subtraction.Therefore, the proposed method not only improves the analysis efficiency, but also enables us to perform reasonable peak fitting.
As the parameters of fitting model are converged to a local optimal solution by monotonically increasing log-likelihood value during the iterative calculation (Figure S1 and S6), the peak fitting is stably conducted without the manual trial-and-error [19].In contrast, recent works have demonstrated the peak fitting by using the Bayesian inference that enable to perform an accurate analysis depending on the task.For example, the exchange Monte Carlo method can perform a more accurate analysis than the spectrum adapted EM algorithm.Moreover, the appropriate number of peaks can be determined by calculating the criteria of model selection such as the marginal likelihood.However, a high computational cost would be required for the calculation of Markov Chain Monte Carlo method, and it is not easy to set a suitable prior distribution and inverse temperature for non-expert.As such a setting is unnecessary to use our approach, the peak fitting can be performed easily to the large number of spectral data set.
In the experimental data analysis, we obtained the spatial distribution of Sn 3d 5=2 peaks and visualized the localization of Sn(II) to Sn(IV).This localization cannot be identified by pinpoint and line scan.The proposed procedure enables us to investigate the chemical characterization on a wide range of microstructure without the time-consuming work for the background subtraction and the parameter estimation.We emphasize that the color map (Figure 4(b)) is obtained by the set of μ that is low-dimensional feature extracted from the 3721 spectral data.Such feature is essential to represent the spatial structure of the chemical characterization and provide strong evidence for further investigation of the chemical shifts depending on the state of the substrate.Therefore, the proposed fitting model and spectral adapted ECM algorithm are powerful tool that not only increases the number of available spectral data, but also has the potential to innovate conventional analytical approaches.
The linear background model was used in the analysis because the spectral data are collected by the hemispherical analyzer in fixed mode.As the peak background models are defined individually, other background models based on the generally used methods, such as the Shirley and Tougaard methods, are available in principle.However, to use such a background model, it is necessary to modify the optimization procedure in the CM-step.Therefore, in further study, we extend the proposed method using the model based on general background methods, such as the Shirley and Tougaard methods.
Calculations can be performed sufficiently effective in practical use.The spectrum adapted ECM algorithm completed the calculation in 34.5 hours (mean = 33.4sec, median = 7.0 sec) for 3721 sets of data.It is effectively impossible to manually analyze these large number of spectral data, and to systematically obtain the spatial distribution of peak positions.However, the analysis of experimental data showed that some spectral data required a large calculation time (Figure S6).Such calculation cost may be improved by the technique to accelerate the convergence of iteration [39,40].Further studies are required to develop a spectrum adapted ECM algorithm that utilizes acceleration techniques.
The proposed method requires the fixed number of peaks as an input.Recently, a method for optimizing the number of peaks using an index, such as Bayesian Information Criterion, has been proposed [20,41].To respond to a strong demand for high-throughput peak-number estimation, we are also developing a method to optimize the number of peaks based on the sparse estimation and the optimization using the approximated Bayesian free energy.We intend to advance the further study of the proposed method including these ideas.

Conclusions
Peak fitting generally requires the background subtraction of spectral data and the parameter estimation of fitting model.As these procedures are seriously timeconsuming, conventional peak fitting with manual trial-and-error is unsuitable for systematic analysis of a large number of spectral data.In this study, in order to improve the efficiency of peak fitting, we introduced the fitting model for automatic background subtraction.The fitting model consists of the pseudo-Voigt mixture model and the ramp-sum background model, and both models represent the peak shape and background components, respectively.The optimization of the fitting model is effectively performed by the spectrum adapted ECM algorithm.Application of the proposed model showed that the synthetic spectral data with various noise level were decomposed into the peak and background component less than dozen seconds.We also demonstrated the analysis of 3721 spectral data from the SnS sheet, and these spectral data were successfully decomposed into the component of Sn 3d 3=2 peak, Sn 3d 5=2 peak and background at 34.5 hours.Using the results of analysis, the spatial distribution of the binding energy of Sn 3d 5=2 peak was obtained and clearly indicated the contrast between Sn(II) and Sn(IV).Such spatial distribution provides the insightful information for investigating the properties of the materials and devices.The proposed model plays an important role to obtain the low-dimensional features extracted from a large amount of spectral data through the high-throughput analysis.

Disclosure statement
No potential conflict of interest was reported by the authors.

Figure 1 .
Figure 1.Schematic illustration of the linear background model when the background is divided into a triangle and uniform distribution.It corresponds to the background subtraction by the linear method.P M ðx n jθÞ means any peak distribution, such as the pseudo-Voigt mixture model.
.5}, respectively.The synthetic spectral data were generated in the case of d = {3, 4, 5} and noise-free.

Figure 2 .
Figure 2. Fitting curves of PVMMB to the synthetic spectral data at d = 3, 4, 5 and noise free.The open circles show the spectral data.Red curve is the estimated fitting curve.Red, orange, and yellow dotted curves are each component of the PVMM (K = 3).Light blue and blue curves are P uni and P tri of the background components, respectively.

Figure 3 .
Figure 3.Typical fitting curves of SnS-data.The symbols of a-d are measurement points, which are shown in Figure 4(a).The number in parentheses is the number of spectral data.Open circles show the observed spectral data.Red curve is the estimated fitting curve.Red, orange and yellow dotted curves are each component of the PVMM (K = 2).Light blue and blue curves are P uni and P tri of the background components, respectively.

Figure 4 .
Figure 4. (a) Color map of the total intensity of each spectral data.The symbols of a-d mean the typical part of the sample.Fitting curves of the spectral data from parts are shown in Figure 3.(b) Color map of the binding energy of Sn 2 S 3 peaks in the range of 485.8-486.6 eV.The peaks above 486.6 eV were unified in dark blue.The peaks below 485.8 eV were unified in light yellow.No clear peaks were observed in the spectral data on the electrode (see also Figure 3), and almost peaks of these spectral data were estimated outside the range of 485.8-486.6 eV.

Figure 5 .
Figure 5. Histogram of the calculation time for each spectral data.Each width of bin is 5 seconds.
in which λ k , σ k and η k are fixed to be λ new k , σ new

Table 1 .
Summary of the initial values of parameters in the synthetic data analysis.
calculation completed 5-12 seconds.The values of μ were very close to the true value regardless the magnitude of d (Table

Table 2 .
Summary of the estimated values of parameters in the synthetic data analysis.