# Fileset

[1-s2.0-S0368204826000216-main.pdf](https://mdr.nims.go.jp/filesets/90a0f26b-69f8-49cc-b7db-b8d9ac72844a/download)

## Creator

Shunichi Yoneda, [Ryo Murakami](https://orcid.org/0000-0001-8585-9268), [Hiroshi Shinotsuka](https://orcid.org/0000-0001-5147-1396), [Hideki Yoshikawa](https://orcid.org/0000-0002-7389-8865), [Shigeo Tanuma](https://orcid.org/0000-0003-2628-9941), Hiromi Tanaka, Hayaru Shouno, [Kenji Nagata](https://orcid.org/0000-0001-9894-4461)

## Rights

[Creative Commons BY-NC-ND Attribution-NonCommercial-NoDerivs 4.0 International](https://creativecommons.org/licenses/by-nc-nd/4.0/)

## Other metadata

[A data-driven surrogate model for X-ray photoelectron spectroscopy based on survey spectrum background features](https://mdr.nims.go.jp/datasets/ed0727c8-9bb4-43c9-89ef-6caa558d53b5)

## Fulltext

A data-driven surrogate model for X-ray photoelectron spectroscopy based on survey spectrum background featuresJournal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 A0n Contents lists available at ScienceDirectJournal of Electron Spectroscopy and Related Phenomenajournal homepage: www.elsevier.com/locate/elspec  A data-driven surrogate model for X-ray photoelectron spectroscopy based on survey spectrum background featuresShunichi Yoneda a ,∗, Ryo Murakami b, Hiroshi Shinotsuka b , Hideki Yoshikawa b , Shigeo Tanuma b , Hiromi Tanaka a, Hayaru Shouno c , Kenji Nagata d ,∗∗a National Institute of Technology (KOSEN), Yonago College, Yonago, Japanb National Institute for Materials Science, Tsukuba, Japanc Graduate School of Informatics and Engineering, The University of Electro-Communications, Chofu, Japand Center for Basic Research on Materials, National Institute for Materials Science, Tsukuba, JapanA R T I C L E  I N F OKeywords:X-ray photoelectron spectroscopySurvey spectrumSurrogate modelTotal electron attenuation lengthThickness estimation A B S T R A C TWe propose a surrogate model that uses the full wide-energy-range (survey) X-ray photoelectron spectroscopy (XPS) spectrum to quantitatively capture spectral shapes and changes in the sample structure (e.g., an In2O3∕Sibilayer) as functions of film thickness. Although it was known that XPS survey spectra could be calculated for each sample structure using the NIST Simulation of Electron Spectra for Surface Analysis (SESSA), version 2.2.0 of Standard Reference Data 100 (SRD 100), exploring solutions for the sample structure required an unacceptably large amount of computational time to obtain accurate and error-assessed results. The necessity to solve this problem led us to propose a surrogate model that generalizes the calculation results of SESSA. A key feature is the explicit decomposition of the survey spectrum into contributions from individual electron orbitals, with their thickness dependence characterized by response functions for which we quantify the attenuation length by defining the total electron attenuation length as the average distance a photoelectron with an initial kinetic energy of 𝐸𝑘 travels before losing target energy (𝛥𝐸) through multiple scatterings and being detected at a reduced kinetic energy of 𝐸′𝑘(= 𝐸𝑘−𝛥𝐸). Unlike analyses that rely solely on peak intensities, our approach uses the full spectrum to enable accurate and efficient thickness estimation.. 1. IntroductionX-ray photoelectron spectroscopy (XPS) is a surface-sensitive tech-nique used to determine the chemical compositions and chemical states of materials. Because XPS probes only the top few nanometers of the surface, it is widely used to characterize surface chemistry. When a sample is irradiated with X-rays, electrons in the material are emitted as photoelectrons via the photoelectric effect. By measuring the kinetic energies and peak intensities in the resulting spectrum and converting them to binding energies (BEs) of core-level electrons, elemental species and their chemical states can be identified in the sample.Additionally, to non-destructively extract depth-dependent struc-tural information from XPS data, changes in peak intensities can be an-alyzed as a function of the emission angle, i.e., angle-resolved XPS [1], or measured peak intensities can be compared with those of reference samples. In the latter case, depth information is obtained by expressing the attenuation of the measured peak intensity relative to the intensity ∗ Correspondence to: Graduate School of Informatics and Engineering, The University of Electro-Communications, Chofu, Japan.∗∗ Corresponding author.E-mail addresses: y2530148@gl.cc.uec.ac.jp (S. Yoneda), nagata.kenji@nims.go.jp (K. Nagata).of the reference peak as a ratio. This attenuation occurs because pho-toelectrons generated within the solid are inelastically and elastically scattered on their way to the surface; consequently, as described in standard XPS theory, the peak intensity decays approximately exponen-tially with increasing path length to the surface [2]. The attenuation behavior depends on the sample’s structure and composition. A com-monly used parameter that quantifies this depth-dependent decay is the effective attenuation length (EAL), which depends on the photoelectron kinetic energy and material [3].This approach estimates the depth profile by comparing two or more spectra acquired at different emission angles or measured spec-tra with those obtained from reference samples. However, analysis based on multiple spectra measured under varying conditions faces two challenges: The accurate evaluation of rough surfaces is difficult, and surface mapping requires prolonged measurements.Alternatively, depth-dependent structural information can be in-ferred from a single spectrum using not only the peaks but also the background accompanying them (hereafter, the ‘‘peak background’’) [4]https://doi.org/10.1016/j.elspec.2026.147612Received 11 January 2026; Received in revised form 19 March 2026; Accepted 22 vailable online 31 March 2026 368-2048/© 2026 The Authors. Published by Elsevier B.V. This is an open access artc-nd/4.0/ ). March 2026icle under the CC BY-NC-ND license ( http://creativecommons.org/licenses/by- https://www.elsevier.com/locate/elspechttps://www.elsevier.com/locate/elspechttps://orcid.org/0009-0005-5996-6226https://orcid.org/0000-0001-5147-1396https://orcid.org/0000-0002-7389-8865https://orcid.org/0000-0003-2628-9941https://orcid.org/0000-0002-2412-0184https://orcid.org/0000-0001-9894-4461mailto:y2530148@gl.cc.uec.ac.jpmailto:nagata.kenji@nims.go.jphttps://doi.org/10.1016/j.elspec.2026.147612https://doi.org/10.1016/j.elspec.2026.147612http://creativecommons.org/licenses/by-nc-nd/4.0/http://creativecommons.org/licenses/by-nc-nd/4.0/S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 The peak background arises as electrons travel through the sample and lose kinetic energy through inelastic scattering with valence and core electrons, which is observed on the lower-kinetic-energy side of the peak. As the path length increases, both the number of inelastic events and cumulative energy loss increase, producing a broad distribution on the low-kinetic-energy side of the peak. Thus, the peak back-ground carries information about the depth at which photoelectrons are generated; however, because this information is encoded in a broad distribution, the extraction of depth-resolved information from the spectral shape is not straightforward. A well-known approach is the Tougaard method [5], which models the build-up of background because of multiple inelastic scatterings to reproduce backgrounds consistent with the sample structure. Nevertheless, because the in-elastic mean free path (IMFP) [2] and differential elastic scattering cross section both depend on photoelectrons’ kinetic energy, the depth sensitivity encoded in the background nonuniformly varies with energy; consequently, the analysis becomes increasingly challenging as the energy range is broadened.An experienced XPS analyst can quickly extract qualitative com-positional and depth-dependent structural information from a survey spectrum acquired over a broad energy range. Although survey spec-tra contain rich information about the sample, the wider the energy window of interest, the more carefully one must account for multiple scattering processes that generate the background, thereby increas-ing the analytical complexity. Therefore, analysis of survey spectra is predominantly qualitative rather than quantitative [6].In this study, we quantitatively analyze survey spectra and de-velop a method for analyzing XPS data without explicitly modeling the complex processes underlying the peak background. Specifically, we employ the XPS simulator simulation of electron spectra for sur-face analysis (SESSA) [7], which generates XPS spectra for samples possessing various structures by accounting for photoelectron gener-ation and subsequent elastic and inelastic scatterings during electron transport through the solid. Thus, without requiring users to explic-itly model numerous physical phenomena that affect photoelectron spectra, including the peak background, the proposed method enables the calculation of spectral changes as the sample’s three-dimensional structure varies. Because SESSA can generate a wide range of spectra corresponding to different sample structures, the structure can also be estimated from photoelectron spectra by searching for the structure that best reproduces the data [1].Although SESSA can compute survey spectra that include the peaks of all the constituent elements and the associated peak backgrounds, when estimating the parameters that define the sample structure from XPS spectra, finer sampling of a broader parameter space increases the number of required simulations, thereby lengthening the total computation time [1]. This issue is exacerbated when the calculation includes the peak background, which spans a wide energy range.Against this backdrop, we use SESSA to generate survey spectra for a wide range of sample structures and train a machine-learning-based surrogate model that captures structural variations indicated in the survey spectra. Using this surrogate, survey spectra from unknown samples can be analyzed without running simulations during the search for candidate structures, thereby greatly accelerating the analysis.In the first step, we construct a surrogate model for a simplified sample structure. We consider a two-layer In2O3∕Si, where the depth of each layer is homogeneous, and use SESSA to generate survey spectra (BE range 0.0 eV–1000.0 eV) for varying thicknesses of the overlayer (In2O3). Then, we train the surrogate model on these spectra to capture spectral changes as a function of the overlayer thickness.The surrogate first decomposes the survey spectrum into spectral components derived from individual electron orbitals and then uses ma-chine learning to model how each component varies with the overlayer thickness. In practical XPS measurements, a measured survey spectrum, including the associated peak background, cannot be decomposed into 2 contributions from individual electron orbitals; however, using SESSA, we can obtain such orbital-resolved components.To apply the surrogate model, we employ Bayesian inference to estimate the overlayer thickness from the XPS survey spectrum of a sample with an unknown thickness [8]. The results show that this data-driven Bayesian approach enables thickness estimation from a single survey spectrum. Moreover, whereas conventional peak-intensity-based methods [2] have a practical limit of approximately 6 nm, our approach achieves high precision, even at 10 nm.Beyond enabling thickness estimation directly from the survey spec-trum, the surrogate model eliminates the need to compute XPS spectra using SESSA while searching for the optimum thickness, thereby en-abling highly efficient thickness estimation. For a given material system and measurement condition, the surrogate model can be constructed offline from a limited number of SESSA simulations and then reused for subsequent inference without running SESSA in the search loop, substantially reducing the computational cost. Thus, the practical ad-vantage of the proposed method does not rely on the universality of the response functions, but on the reuse of a surrogate model constructed for a given material system and measurement condition. Because the method leverages peak-background information over a wide energy range, it retains non-negligible thickness sensitivity at larger depths and under low-count conditions, which is advantageous for mapping applications requiring short acquisition times.2. Construction of the surrogate modelThis section details the construction of a surrogate model that captures changes in the XPS survey spectrum as a function of the film thickness. Section 2.1 introduces the surrogate model by expressing the survey spectrum as a linear combination of photoelectron spectra for individual electron orbitals. Section 2.2 demonstrates how orbital-resolved photoelectron spectra are generated using the XPS simulator SESSA. Section 2.3 defines the response function that governs the thick-ness dependence within the surrogate, describes its approximation, and discusses its validation. Finally, Section 2.4 evaluates the accuracy of the surrogate model. The details of the layered structure of a sample and the angular configuration of the measurements are provided in the Appendix  B.2.1. ModelingWe define a surrogate model — a computationally efficient mathe-matical approximation — that reproduces the measured XPS spectrum as a function of the film thickness (𝑧). The model is constructed assuming that the measured spectrum can be decomposed into spectral components derived from individual electron orbitals of the constituent elements in the sample. The details are as follows.First, consider the measured spectrum 𝐼𝐴(𝐸, 𝜃) of a bulk sample (𝐴). If the spectrum 𝐼𝐴(𝐸, 𝜃), acquired at an emission angle of 𝜃, contains 𝑚𝐴 peaks, it can be written as the sum of its components (𝐼 𝑖𝐴(𝐸, 𝜃)) as follows: 𝐼𝐴(𝐸, 𝜃) =𝑚𝐴∑𝑖=1𝐼 𝑖𝐴(𝐸, 𝜃), (1)where 𝐼 𝑖𝐴(𝐸, 𝜃) denotes the spectral contribution, including the peak background, associated with photoelectron peak 𝑖 of the constituent elements in sample 𝐴, expressed as a function of the BE (𝐸).Consider the measured spectrum 𝐼𝐶 (𝑧, 𝐸, 𝜃) of a bilayer sample (𝐶), where sample 𝐴 is the substrate layer, and sample 𝐵 forms an overlayer having a thickness of 𝑧. If the measured spectrum 𝐼𝐶 (𝑧, 𝐸, 𝜃) consists of 𝑚𝐴 +𝑚𝐵 peaks, 𝐼𝐶 (𝑧, 𝐸, 𝜃) can be expressed in terms of its components (𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃)) as follows: 𝐼𝐶 (𝑧, 𝐸, 𝜃) =𝑚𝐴∑𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) +𝑚𝐵∑𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃), (2)𝑖=1 𝑗=1S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 1. The schematic of the surrogate model. The left-hand side shows the spectrum 𝐼𝐶 (𝑧, 𝐸, 𝜃) obtained from a bilayer sample. The right-hand side illustrates the linear combination of orbital-resolved spectral components. The response functions (𝒇𝑨(𝑧, 𝐸, 𝜃) and 𝒇𝑩(𝑧, 𝐸, 𝜃)) characterize how the spectrum varies with the overlayer thickness (𝑧). 𝑰𝑨(𝐸, 𝜃) and 𝑰𝑩(𝐸, 𝜃) denote the spectra measured for bulk samples 𝐴 and 𝐵, respectively.In Eq. (2), 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) denotes the spectral component derived from the electron orbitals of the constituent elements in the substrate layer, and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) denotes the component derived from the electron or-bitals of the constituent elements in the overlayer; both components depend on the thickness (𝑧) because of elastic and inelastic scatterings. The substrate-derived contribution (𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃)) satisfies 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) =𝐼 𝑖𝐴(𝐸, 𝜃) at 𝑧 = 0, and 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) = 0 as 𝑧 → ∞. Conversely, the overlayer-derived contribution (𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃)) satisfies 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) = 0 at 𝑧 = 0, and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) = 𝐼 𝑗𝐵(𝐸, 𝜃) as 𝑧 → ∞. Accordingly, 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃)and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) can be written as follows:𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) = 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃)𝐼 𝑖𝐴(𝐸, 𝜃), (3)𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) = 𝐼 𝑗𝐵(𝐸, 𝜃) − 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)𝐼 𝑗𝐵(𝐸, 𝜃)= {1 − 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)}𝐼 𝑗𝐵(𝐸, 𝜃), (4)where 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) are response functions that character-ize how the spectrum 𝐼𝐶 (𝑧, 𝐸, 𝜃) varies with the overlayer thickness (𝑧). By definition, the response functions satisfy 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃), 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) = 1at 𝑧 = 0, and 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃), 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) → 0 as 𝑧 → ∞.By substituting Eqs. (3) and (4) into Eq. (2), 𝐼𝐶 (𝑧, 𝐸, 𝜃) can be written as follows:𝐼𝐶 (𝑧, 𝐸, 𝜃) =𝑚𝐴∑𝑖=1𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃)𝐼 𝑖𝐴(𝐸, 𝜃) +𝑚𝐵∑𝑗=1(1 − 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃))𝐼 𝑗𝐵(𝐸, 𝜃)= 𝒇𝑨(𝑧, 𝐸, 𝜃)𝑰𝑨(𝐸, 𝜃) + {𝟏 − 𝒇𝑩(𝑧, 𝐸, 𝜃)}𝑰𝑩(𝐸, 𝜃). (5)A schematic is shown in Fig.  1. 𝒇𝑨(𝑧, 𝐸, 𝜃)𝑰𝑨(𝐸, 𝜃) and(𝟏 − 𝒇𝑩(𝑧, 𝐸, 𝜃))𝑰𝑩(𝐸, 𝜃) are expressed in vector form because the re-sponse functions and the orbital-resolved spectral components cor-respond one to one and can be written as inner products: 𝒇𝑨 =(𝑓 1𝐴, 𝑓2𝐴,… , 𝑓𝑚𝐴𝐴), 𝑰𝑨 =(𝐼1𝐴, 𝐼2𝐴,… , 𝐼𝑚𝐴𝐴)𝑇 ,𝒇𝑩 =(𝑓 1𝐵 , 𝑓2𝐵 ,… , 𝑓𝑚𝐵𝐵), and 𝑰𝑩 =(𝐼1𝐵 , 𝐼2𝐵 ,… , 𝐼𝑚𝐵𝐵)𝑇 .Importantly, conventional XPS measurements only yield 𝐼𝐶 (𝑧, 𝐸, 𝜃), 𝐼𝐴(𝐸, 𝜃), and 𝐼𝐵(𝐸, 𝜃); the orbital-resolved spectral components (𝐼 𝑖𝐴(𝐸, 𝜃)and 𝐼 𝑗𝐵(𝐸, 𝜃)) are not directly measurable. By contrast, as explained in the next subsection, an XPS simulator can generate the orbital-resolved spectra (𝐼 𝑖 (𝐸, 𝜃) and 𝐼 𝑗 (𝐸, 𝜃)).𝐴 𝐵3 2.2. Generation of orbital-resolved spectraIn Section 2.1, as shown in Eq. (5), the measured spectrum (𝐼𝐶 (𝑧, 𝐸,𝜃)) of bilayer sample 𝐶 can be written as the sum of components 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) derived from electron orbitals of the con-stituent elements in the substrate layer (𝐴) and overlayer (𝐵), re-spectively. However, in conventional XPS measurements, because the orbital-resolved bulk spectra (𝐼 𝑖𝐴(𝐸, 𝜃) and 𝐼 𝑗𝐵(𝐸, 𝜃)) cannot be mea-sured directly, we use an XPS simulator to generate them.Fig.  2 shows the simulated spectra for a bilayer sample (𝐶)(In2O3∕Si), where the substrate is 𝐴 (Si) and the overlayer is 𝐵 (In2O3). The left panel shows the survey spectrum 𝐼𝐶 (𝑧, 𝐸, 𝜃) measured for sample 𝐶 at the overlayer thickness (𝑧). The central panel shows the components 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) derived from each electron orbital of the constituent elements in the bilayer. The different colors indicate the different electron orbitals from which each spectrum is derived: Red and blue correspond to Si 2p3∕2 and In 3d5∕2 resolved spectra, respectively. The right panel shows the thickness dependence of the orbital-resolved spectra, indicating that the response functions (𝑓𝐴(𝑧, 𝐸, 𝜃) and 𝑓𝐵(𝑧, 𝐸, 𝜃)) characterize the intensity variation with increasing thickness. 𝐼𝐴(𝐸, 𝜃) and 𝐼𝐵(𝐸, 𝜃) are the survey spectra mea-sured for bulk samples Si and In2O3, respectively. Gray indicates spectra derived from the electron orbitals included in the survey spectrum but not highlighted here; these components can likewise be expanded using the corresponding response functions. Specific parameters (e.g., peak positions and widths) are provided in the Appendix  B.2.3. Response functions2.3.1. Formulation of the response functionsAs described in Section 2.2, the response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) in Eqs. (3) and (4), respectively, can be expressed based on spectra calculated using SESSA. Specifically, the response functions 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) are given by• Substrate: 𝐴𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) =𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃)𝑖 , (6)𝐼𝐴(𝐸, 𝜃)S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 2. Example of spectra simulated using SESSA. Left: the XPS spectrum 𝐼𝐶 (𝑧, 𝐸, 𝜃) measured for sample 𝐶 (In2O3∕Si) at the overlayer thickness (𝑧). Center: orbital-resolved spectra 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃). The different colors indicate the different electron orbitals from which each spectrum is derived: Red, blue indicate Si 2p3∕2, In 3d5∕2, respectively. Right: the thickness dependence of the orbital-resolved spectra. The response functions (𝑓𝐴(𝑧, 𝐸, 𝜃) and 𝑓𝐵(𝑧, 𝐸, 𝜃)) describe spectral changes with increasing thickness (𝑧). 𝐼𝐴(𝐸, 𝜃) and 𝐼𝐵(𝐸, 𝜃) are the spectra measured using XPS for the bulk samples, with 𝐴 as Si and 𝐵 as In2O3. Specific parameters (e.g., peak positions and widths) are provided in the Appendix  B. Gray indicates spectra derived from the electron orbitals included in the survey spectrum but not highlighted here; these components can likewise be expanded using the corresponding response functions.Fig. 3. Dependence of 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃)∕𝐼 𝑖𝐴(𝐸, 𝜃) and 1 − 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃)∕𝐼 𝑗𝐵(𝐸, 𝜃) on the film thickness at binding energy 𝐸. Sample 𝐶 is In2O3∕Si; (a) orbital-resolved spectra for Si 2p3∕2 and (b) for In 3d5∕2 Data points correspond to the response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) output by SESSA. Circles indicate the peak ((a): 𝐸 = 99.5 eV; (b): 𝐸 = 445.5 eV), squares indicate a point offset by 10 eV from the peak ((a): 𝐸 = 109.5 eV; (b): 𝐸 = 455.5 eV), and triangles indicate a point far from the peak ((a): 𝐸 = 950.0 eV; (b): 𝐸 = 950.0 eV).• Overlayer: 𝐵𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) = 1 −𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃)𝐼 𝑗𝐵(𝐸, 𝜃). (7)The right-hand sides of Eqs. (6) and (7) are obtained directly from the simulator. Using these expressions, we obtain the response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) shown in Fig.  3. The sample is In2O3∕Si, as in Fig.  2. Fig.  3 presents, at a fixed BE of 𝐸, the simulator-derived response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) as functions of the film thickness (𝑧). The horizontal axis is 𝑧, and the vertical axis shows the right-hand sides of Eqs. (6) and (7) on a logarithmic scale. Data-point shapes and colors indicate the chosen BE. Panels (a) and (b) correspond to different electron orbitals. From Fig.  3, we confirm that over the 4 tested BEs, including the background region, the response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) vary approximately linearly with thickness 𝑧 on a logarithmic scale.According to these validation results, the response functions (𝑓 𝑖𝐴(𝑧,𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) decay exponentially with increasing film thickness (𝑧). Therefore, we approximate 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) with the following expressions:• Substrate: 𝐴𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) ≈ exp(− 𝑧𝑖), (8)𝐿𝑇 𝐴(𝐸) cos 𝜃S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 • Overlayer: 𝐵𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) ≈ exp(− 𝑧𝐿𝑇𝑗𝐵(𝐸) cos 𝜃). (9)In Eqs. (8) and (9), the quantity 𝐿𝑇 (𝐸) on the right-hand sides denotes the distance over which the electron intensity, including the back-ground component, decays to 1∕𝑒; we refer to this as the total electron attenuation length (TEAL), which is interpreted as the mean path length that a photoelectron with an initial kinetic energy of 𝐸𝑘𝑖, 𝐸𝑘𝑗 travels before losing energy (𝛥𝐸) via multiple elastic and inelastic scatterings and being detected in the spectrum as background on the high–BE (low-kinetic-energy) side. The angle 𝜃 is the photoelectron emission angle, which is treated as fixed throughout this paper.In this study, we determine TEAL directly by fitting the approximate expressions in Eqs. (8) and (9) to simulated spectra.2.3.2. Determination of response function parametersIn this section, we fit the proposed approximate expressions in Eqs. (8) and (9) to simulated spectra and validate TEAL.Using the right-hand sides of Eqs. (6) and (7), we fit the approximate expressions in Eqs. (8) and (9) to the simulator-derived data in Fig. 3; the resulting fits are shown in Fig.  4. Fig.  4 uses the same format and plotting conditions as those in Fig.  3. Although the vertical axis uses a logarithmic scale, we fit the data in the original (linear) domain using nonlinear least squares to preserve accuracy in regions where the response functions vary rapidly with the thickness (𝑧). Nonlinear least squares minimize the sum of squared residuals (SSRs), defined as follows: SSR =𝑁∑𝑘=1(𝑓𝑘SESSA(𝑧, 𝐸, 𝜃) − 𝑓𝑘f it (𝑧, 𝐸, 𝜃))2 , (10)where 𝑓SESSA(𝑧, 𝐸, 𝜃) is the simulator-derived response function at thickness 𝑧, 𝑓f it (𝑧, 𝐸, 𝜃) is the value of the approximate expression, and 𝑁 is the number of data points used for fitting [9]. The solid lines represent the optimized approximate expressions in Eqs. (8) and (9) fitted to the data points. As shown in Fig.  4, the fits agree well with the data, supporting the validity of the approximation.Historically, the exponential attenuation in Eqs. (8) and (9) has been assumed to apply only to the relationship between the photoelectron peak intensity and film thickness (𝑧) [2]. For peak intensity changes in inelastic scattering only, the response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) are given by• Substrate: 𝐴𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) =𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃)𝐼 𝑖𝐴(𝐸, 𝜃)= exp(− 𝑧𝜆(𝐸) cos 𝜃), (11)• Overlayer: 𝐵𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) = 1 −𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃)𝐼 𝑗𝐵(𝐸, 𝜃)= exp(− 𝑧𝜆(𝐸) cos 𝜃), (12)where 𝜆(𝐸) denotes the distance over which the photoelectron peak intensity falls to 1∕𝑒, also known as the IMFP, which depends on the photoelectron kinetic energy and varies with the material and its structure.As is evident from these expressions, 𝐿𝑇 (𝐸) and 𝜆(𝐸) are parameters that characterize the attenuation length. Table  1 compares 𝐿𝑇 (𝐸) and 𝜆(𝐸). We evaluate 𝜆(𝐸) using the TPP-2M formula of Tanuma, Powell, and Penn [10], with an Al K𝛼 X-ray source. The table uses the same BE (𝐸) as that used in Fig.  4 and reports 𝐿𝑇 (𝐸) and 𝜆(𝐸) at those energies; the rightmost column lists their ratio, 𝐿𝑇 (𝐸)∕𝜆(𝐸). From the table, the ratio at the peak is slightly less than 1 because 𝐿𝑇 (𝐸) accounts for both inelastic and elastic scatterings. In the background region, the ratio substantially increases with increasing BE.5 In the preceding analysis, we evaluated the TEAL (𝐿𝑇 (𝐸)) and IMFP (𝜆(𝐸)) at selected BEs (𝐸). Then, we set the energy loss from the peak top, 𝛥𝐸, in the range from 0.0 eV to 1000.0 eV, with a step of 0.5 eV, and evaluated the ratio 𝐿𝑇 (𝐸)∕𝜆(𝐸) for each 𝛥𝐸. The results are shown in Fig.  5. In this figure, the horizontal axis represents the energy loss (𝛥𝐸), and the vertical axis represents the ratio of 𝐿𝑇 (𝐸) to 𝜆(𝐸). The two curves correspond to the ratios 𝐿𝑇 (𝐸)∕𝜆(𝐸) for two different elec-tron orbitals, indicating that the TEAL (𝐿𝑇 (𝐸)) takes different values depending on the orbital.Now, we discuss why such differences arise between 𝐿𝑇 (𝐸) and 𝜆(𝐸). The peak component’s attenuation can be simply described by the IMFP (𝜆(𝐸)), which is a physical quantity for deriving the EAL of generated photoelectrons that travel before one inelastic scattering. In contrast, in the peak background region, two types of EAL must also be considered. One is the EAL that a traveling photoelectron that has lost energy (𝛥𝐸) travels before losing further energy in an additional inelastic scattering. The other is the EAL where the sum of the energy losses because of multiple inelastic scatterings equals 𝛥𝐸. For the change in the photoelectron intensity in the peak background region, we will call the former and latter the consumption and supply terms, respectively. In the peak background region, this supply term makes the effective attenuation with respect to the film thickness (𝑧) more gradual and, consequently, 𝐿𝑇 (𝐸) increases. In this situation, the TEAL (𝐿𝑇 (𝐸)) shown in Fig.  5 can be interpreted as the mean path length traveled by a photoelectron with an initial kinetic energy of 𝐸𝑘𝑖 or 𝐸𝑘𝑗 until, through multiple elastic and inelastic scatterings, it loses energy (𝛥𝐸) and is detected in the spectrum as background on the high–BE side (at 𝐸 + 𝛥𝐸), i.e., on the low-kinetic-energy side.Overall, TEAL is an effective approximate parameter for response functions derived from simulator-generated survey spectra. Particu-larly, TEAL’s ability to describe background spectral variations, in-cluding those in the high-BE region, offers practical utility that differs from that of conventional IMFP-based formulations. However, because TEAL is obtained by empirical fitting to numerical data calculated using SESSA rather than from first principles and because TEAL can depend on measurement conditions, its range of applicability and generalizability are subject to certain limitations.2.4. Generation of spectra using the surrogate modelIn Section 2.3, we showed that the response functions are well approximated by Eqs. (8) and (9), and we obtained the TEAL (𝐿𝑇 (𝐸)) at each BE (𝐸) using nonlinear least squares. Consequently, substi-tuting Eqs. (8) and (9) into Eq. (5) enables the survey spectrum to be generated for a layered sample (𝐶) and a film thickness of 𝑧 as input. We refer to this as the surrogate model. Because TEAL and the response functions are defined and estimated for each orbital-resolved component in Eqs. (6) and (7), the same procedure can be applied to all orbitals originating from either the substrate or the overlayer; here we show Si 2p3∕2 and In 3d5∕2 as representative examples.We present validation results comparing spectra generated using the surrogate model with those produced using the simulator. Fig.  6 shows the validation for the survey spectra (𝐼𝐶 (𝑧, 𝐸, 𝜃)) of bilayer sample 𝐶 , for which we constructed response functions and estimated TEAL for all 15 orbital-resolved components (i.e., all peaks listed in Table  2) in the In2O3∕Si, and Fig.  7 shows the validation for the orbital-resolved com-ponents (𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃)) derived from the same sample. In Fig.  7, panels (a), (c), and (e) show the Si 2p3∕2 orbital-resolved spectra, while panels (b), (d), and (f) show the In 3d5∕2 orbital-resolved spectra. The figures compare the spectra generated using the surrogate model with those generated using SESSA at film thicknesses (𝑧) of 0.10, 1.00, and 10.00 nm. In the upper panels, surrogate- and SESSA-generated spectra are shown in red and black, respectively; the lower panels plot the corresponding residuals. The plotted residual spectral intensities are expanded by two orders of magnitude as follows: The plotted spectral S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 4. Dependence of the response functions 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) on the film thickness at the BE (𝐸). Sample 𝐶 is In2O3∕Si; the panels show (a) Si 2p3∕2and (b) In 3d5∕2 orbital-resolved spectra. Data points correspond to the simulator-derived response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)). Solid lines depict Eqs. (8) and (9); TEAL 𝐿𝑇  was obtained by fitting these lines to data points of the corresponding color. Circles indicate the peak top ((a): 𝐸 = 99.5 eV; (b): 𝐸 = 445.5 eV), squares indicate a point offset by 10 eV from the peak ((a): 𝐸 = 109.5 eV; (b): 𝐸 = 455.5 eV), and triangles indicate a point far from the peak top ((a): 𝐸 = 950.0 eV; (b): 𝐸 = 950.0 eV).Table 1The TEAL (𝐿𝑇 (𝐸)) plotted as a function of the BE obtained by fitting the response functions (𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) and 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)) for the electron orbitals Si 2p3∕2 and In 3d5∕2 using Eqs. (8) and (9), respectively, together with the IMFP (𝜆(𝐸)) derived from the TPP-2M formula, and the ratio 𝐿𝑇 (𝐸)∕𝜆(𝐸). The X-ray source is Al K𝛼 radiation. The table also lists the corresponding kinetic energy for each BE. Electron orbitals BE (eV) Kinetic energy (eV) TEAL (nm) IMFP (nm) 𝐿𝑇 (𝐸)∕𝜆(𝐸)  99.5 1387.1 2.32 2.81 0.83   Si 2p3∕2 109.5 1377.1 4.53 2.79 1.62   950.0 536.6 30.39 1.37 22.18   445.5 1041.1 1.86 2.25 0.83   In 3d5∕2 455.5 1031.1 3.65 2.23 1.63   950.0 536.6 20.91 1.37 15.26  Fig. 5. Ratio of the TEAL (𝐿𝑇 (𝐸)) to the IMFP (𝜆(𝐸)) plotted as a function of the energy loss (𝛥𝐸). TEAL (𝐿𝑇 (𝐸)) values were obtained by fitting Eqs. (8) and (9) to the data generated using SESSA. The IMFP (𝜆(𝐸)) was calculated using the TPP-2M formula. Sample 𝐶 is In2O3∕Si, and the dashed line (–) and dash-dotted line (-  ⋅ -) represent data derived from the Si 2p3∕2and In 3d5∕2orbital-resolved spectra, respectively. The X-ray source is monochro-matic Al K𝛼 with a photon energy of ℎ𝜈 = 1486.6 eV. In the kinetic energy representation, the peak positions of the Si 2p3∕2and In 3d5∕2orbitals are 1387.1 eV and 1041.0 eV, respectively, and 𝛥𝐸 = 0 corresponds to these peak positions.intensities (vertical axis) are on the orders of 10−5 and 10−7 in the upper and lower panels, respectively.As shown in Fig.  6, the surrogate model reproduces the spectrum very well. This also holds for the orbital-resolved spectra (Fig.  7).Accordingly, for any BE (𝐸) within the validated range, Eqs. (8) and (9) remain applicable, indicating that the inverse problem of 6 determining the film thickness can be solved using this surrogate model.3. Film thickness estimation using the surrogate model3.1. Film thickness evaluation using Bayesian inferenceUsing the surrogate model developed in the preceding sections, we quantitatively estimate the film thickness (𝑧) from the survey spectrum. Although this analysis considers a relatively simple sample structure and methods such as maximum likelihood estimation could also be used, we adopt Bayesian inference so that the resulting posterior dis-tribution can be analyzed to assess the reliability of the estimate and support the appropriate interpretation.The sample structure is In2O3∕Si; detailed parameters are provided in the Appendix  B. The simulated intensities are generated to represent practical measurement conditions. Thus far, the BE (𝐸) has been treated as continuous; in this section, we discretize it for numerical compu-tation and denote the grid points by 𝐸𝑛 (𝑛 = 1, 2,… , 𝑁) as follows: 𝑦𝑛(𝐸𝑛, 𝜃) ∼ Poi(ℎ × 𝐼(𝐸𝑛, 𝜃)). (13)In Eq. (13), Poi(ℎ×𝐼(𝐸𝑛, 𝜃)) denotes a Poisson distribution with a mean of ℎ × 𝐼(𝐸𝑛, 𝜃); the equation states that the intensity 𝑦𝑛(𝐸𝑛, 𝜃) at BE 𝐸𝑛follows a Poisson distribution and is observed as count data. The mean intensity (ℎ × 𝐼(𝐸𝑛, 𝜃)) equals the simulator (SESSA) output (𝐼(𝐸𝑛, 𝜃)) multiplied by an intensity-scaling parameter (ℎ). We set the scaling parameter at ℎ = 108 so that the counting intensities are on the order of 𝐼 ∼ 103.S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 6. Comparison of survey spectra for sample 𝐶 (In2O3∕Si) generated using the surrogate model and the simulator. Panels (a)–(c) correspond to film thicknesses 𝑧 = 0.10, 1.00, and 10.00 nm, respectively. In each panel, the upper plot compares spectra generated by the surrogate model and the simulator, and the lower plot shows the residuals.By applying Bayesian inference to the measured spectrum (𝑦) gen-erated under these conditions, we obtain the posterior distribution of the mean thickness (𝑧̂). Details are provided in the Appendix  C.Fig.  8 presents Bayesian estimates of the film thickness (𝑧). The target thicknesses are 1.00, 5.00, and 10.00 nm. The figure layout is as follows: Upper left: measured spectrum 𝑦 (black) and surrogate-generated spectrum (red), overlaid with the Bayesian fit. Lower left: Poisson negative-log-likelihood at each BE. Right: posterior distribution of the mean thickness (𝑧̂) derived from the Bayesian analysis. The right panel’s legend lists the true thickness (𝑧) and posterior mean for 𝑧̂ and the corresponding 95% confidence interval.As shown in Fig.  8, the surrogate model reproduces the measured data (𝑦) well. Comparing the true thickness (𝑧) with the posterior mean and 95% confidence interval of 𝑧̂ obtained using Bayesian inference, we find that (a) 𝑧 = 1.00 nm, 𝑧̂ = 0.84 ± 0.02 nm; (b) 𝑧 = 5.00 nm, 𝑧̂ = 4.75 ± 0.12 nm; (c) 𝑧 = 10.00 nm, 𝑧̂ = 13.22 ± 0.49 nm. These 7 results indicate high accuracy for thin overlayers, whereas for thicker overlayers, the estimates exceed the true values.3.2. Estimation range of the proposed methodFig.  9 shows the accuracy of the film thickness estimation by plotting 𝑧̂ as functions of 𝑧 (the identity line is 𝑧 = 𝑧̂). We set the prescribed film thicknesses at 𝑧 = 0.50, 1.00,… , 10.00 nm (i.e., 𝑧 =0.50 𝑛1 nm with 𝑛1 = 1, 2,… , 20). The figure layout is as follows: Upper panel: prescribed thickness (𝑧) (horizontal) versus estimated thickness (𝑧̂) (vertical); points show the posterior mean of 𝑧̂, and error bars indicate the 95% confidence interval. The lower panel shows the mean estimation cost. To emulate changes in the count level with increasing measurement time, the overall intensity scale is set at (a) 𝐼 ∼ 104, (b) 𝐼 ∼ 103, and (c) 𝐼 ∼ 102 in Fig.  9. The signal-to-noise ratio (SNR) is the highest for (a) 𝐼 ∼ 104 and the lowest for (c) 𝐼 ∼ 102.As shown in Fig.  9, the error bar centers in panels (a)–(c) are nearly identical. At higher 𝑧 values, the estimation bias shows a broadly similar trend across the panels. The error bars are the narrowest in panel (a), with the highest SNR, and the widest in panel (c), with the lowest SNR.3.3. Discussion3.3.1. Effectiveness of the surrogate modelIn Fig.  8(c) (10.00 nm), the estimated thickness is 𝑧̂ = 13.22±0.49 nm, indicating overestimation. The posterior distribution is also broader than those in panels (a) and (b). This broadening arises because data sensitivity decreases with increasing film thickness (𝑧). As shown in Fig. 4, the peak intensity at 10.00 nm is reduced to approximately 1∕1000of its maximum value; a similar reduction occurs at other BE points. Consequently, the estimation sensitivity is reduced, thereby broadening the posterior distribution.The overestimation is because of the underestimation of the expo-nential fitting near 10.00 nm for the response functions in the back-ground regions in Fig.  4. In the background region, the estimated response functions (square and triangle markers) consistently shift in the same direction as the underestimation in Fig.  8.As shown in Fig.  9, the estimates in panels (a)–(c) follow similar trends, although the posterior width varies with the noise level, indi-cating that the use of the broad energy range and numerous data points in Bayesian inference with the surrogate model enables high-accuracy thickness estimation that is relatively robust to noise, provided that the systematic overestimation of the film thickness is considered for thick films. High-precision and high-accuracy thickness estimations are also achievable for thin overlayers.The present surrogate model assumes a laterally uniform overlayer thickness. If the sample exhibits lateral thickness variations within the analysis area, the measured spectrum represents an average over regions with different local thicknesses. In principle, such effects can be incorporated by constructing the surrogate model from SESSA sim-ulations that include the assumed thickness distribution. For mea-surements with sufficient spatial resolution, the present framework may also be applied point by point to map the thickness distribution. Because XPS is inherently an area-averaging technique, the measured spectrum should be interpreted as the averaged response of the sample within the analysis area.3.3.2. Comparison of the proposed and conventional methodsIn this section, we demonstrate the effectiveness of the proposed method using both peaks and backgrounds by comparing it with a conventional method using background-free peaks. The conventional method determines the film thickness by evaluating the ratio of the peak areas of the In2O3∕Si sample of interest to those of the reference bulk In2O3 sample or bulk Si sample. In this work, the peak areas of Si 2p and In 3d  were evaluated after removing the background, which 5∕2S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 7. Comparison of orbital-resolved spectra for sample 𝐶 (In2O3∕Si) generated using the surrogate model and the simulator. Panels (a), (c), and (e) show the Si 2p3∕2 resolved spectra, and panels (b), (d), and (f) show the In 3d5∕2 resolved spectra. Film thicknesses are 𝑧 = 0.10 nm for panels (a) and (b), 1.00 nm for panels (c) and (d), and 10.00 nm for panels (e) and (f). In each panel, the upper plot compares spectra generated by the surrogate model and the simulator, and the lower plot shows the residuals.was assumed to be linear within a narrow energy window around each peak.The resulting peak areas exhibit exponential attenuation with in-creasing film thickness, analogous to Eqs. (8) and (9). Accordingly, the response functions (𝑓 𝑖𝐴(𝑧, 𝜃) and 𝑓𝑗𝐵(𝑧, 𝜃)) for the ratios of peak areas can be written as follows:• Substrate: 𝐴𝑓 𝑖𝐴(𝑧, 𝜃) =𝐼 𝑖𝐶𝐴(𝑧, 𝜃)𝐼 𝑖𝐴(𝜃)≈ exp(− 𝑧𝐿 cos 𝜃), (14)• Overlayer: 𝐵𝑓 𝑗𝐵(𝑧, 𝜃) = 1 −𝐼 𝑗𝐶𝐵(𝑧, 𝜃)𝐼 𝑗𝐵(𝜃)≈ exp(− 𝑧𝐿 cos 𝜃). (15)In Eqs. (14) and (15), the middle expressions denote ratios of peak area intensities: 𝐼 𝑖𝐶𝐴(𝑧, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝜃) are the peak area intensities of the specified electron orbitals of the constituent elements in samples having a thickness of 𝑧, whereas 𝐼 𝑖 (𝜃) and 𝐼 𝑗 (𝜃) are the reference (bulk) 𝐴 𝐵8 peak area intensities for the corresponding orbitals. The right-hand sides indicate that these intensity ratios follow an exponential trend; the photoelectron emission angle in the exponent is fixed at 𝜃 = 45◦. In the denominator of the exponent, the parameter 𝐿 is the EAL, an empirical value obtained from the attenuation of measured intensities [3]. 𝐿 was determined by averaging the values derived from 𝐼 𝑖𝐴(𝜃) with 𝐼 𝑖𝐶𝐴(𝑧 =1.0 nm, 𝜃) and with 𝐼 𝑖𝐶𝐴(𝑧 = 3.0 nm, 𝜃) and similarly from 𝐼 𝑗𝐵(𝜃) with 𝐼 𝑗𝐶𝐵(𝑧 = 1.0 nm, 𝜃) and with 𝐼 𝑗𝐶𝐵(𝑧 = 3.0 nm, 𝜃). The resulting values are Si 2p: 2.25 < 𝐿 < 2.35 nm and In 3d5∕2: 1.80 < 𝐿 < 1.90 nm, approximately 0.8 times the IMFP values set in SESSA (Si 2p: 2.87 nm; In 3d5∕2: 2.30 nm), indicating that the obtained 𝐿 values are reasonable.Next, to estimate the film thickness (𝑧) from the measured peak areas (𝐼 𝑖𝐶𝐴(𝜃) and 𝐼 𝑗𝐶𝐵(𝜃)), for which the corresponding 𝑧 values of 𝐼 𝑖𝐶𝐴(𝑧, 𝜃) and 𝐼 𝑗𝐶𝐵(𝑧, 𝜃) are unknown, we solve Eqs. (14) and (15) for 𝑧 to obtain• Substrate: 𝐴𝑧 ≈ − 𝐿 cos 𝜃 log𝐼 𝑖𝐶𝐴(𝜃)𝑖 , (16)𝐼𝐴(𝜃)S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 8. Bayesian estimation results for film thicknesses 𝑧 = 1.00, 5.00, and 10.00 nm at an overall intensity scale of 𝐼 ∼ 103. Panels (a)–(c) correspond to the prescribed thicknesses 𝑧 = 1.00, 5.00, and 10.00 nm, respectively. In each panel, the upper-left plot shows the fit between the measured and surrogate-generated spectra, the lower-left plot shows the Poisson negative log-likelihood at each binding-energy point, and the right plot shows the posterior distribution obtained by Bayesian inference. The posterior mean 𝑧̂ and the 95% confidence interval are as follows: (a) 𝑧̂ = 0.84±0.02 nm, (b) 𝑧̂ = 4.75±0.12 nm, and (c) 𝑧̂ = 13.22±0.49 nm.• Overlayer: 𝐵𝑧 ≈ − 𝐿 cos 𝜃 log⎛⎜⎜⎝1 −𝐼 𝑗𝐶𝐵(𝜃)𝐼 𝑗𝐵(𝜃)⎞⎟⎟⎠. (17)Fig.  10 compares the thickness estimates (𝑧̂) calculated using these equations with those obtained using the proposed method.We add Poisson noise to the data used for estimation; thus, the observed spectral intensities (count statistics) exhibit stochastic fluc-tuations even under identical sample and instrumental conditions. To assess the impact of these fluctuations, we generate 𝑁 = 100 noisy 9 realizations for each prescribed thickness (𝑧) and estimate the thickness for each realization.Because the peak areas 𝐼 𝑖𝐶𝐴(𝜃) and 𝐼 𝑗𝐶𝐵(𝜃) are obtained after back-ground subtraction, Poisson fluctuations can yield non-physical values. In such cases, Eqs. (16) and (17) do not return a real thickness because the arguments of the logarithms become non-positive. We therefore evaluate the failure probabilities of the conventional inversion as𝑝fail,𝐴(𝑧) =1𝑁∑I[ 𝐼 𝑖𝐶𝐴 ,𝑛(𝜃)𝑖 ≤ 0], (18)𝑁 𝑛=1 𝐼𝐴(𝜃)S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 9. Estimated film thickness (𝑧̂) obtained using the surrogate model with Bayesian inference, plotted as a function of the designated thickness (𝑧), for spectra recorded at different XPS acquisition times. Data points: posterior mean of 𝑧̂; error bars: 95% confidence intervals. Overall intensity scale: (a) 𝐼 ∼ 104; (b) 𝐼 ∼ 103; (c) 𝐼 ∼ 102. Lower panels: mean estimation cost.𝑝fail,𝐵(𝑧) =1𝑁𝑁∑𝑛=1I⎡⎢⎢⎣1 −𝐼 𝑗𝐶𝐵 ,𝑛(𝜃)𝐼 𝑗𝐵(𝜃)≤ 0⎤⎥⎥⎦, (19)where I[⋅] is the indicator function. We then define the quantification limits as the largest thicknesses satisfying 𝑝fail(𝑧) < 0.05: 𝑧QL,𝐴 = max{𝑧 ∣ 𝑝fail,𝐴(𝑧) < 0.05}, 𝑧QL,𝐵 = max{𝑧 ∣ 𝑝fail,𝐵(𝑧) < 0.05}.(20)Realizations counted as failures are excluded from the plotted esti-mates.The validation range matches that in Section 3.1, i.e., 𝑧 = 0.50 𝑛1 nmwith 𝑛1 = 1, 2,… , 20. In conventional analysis, narrow spectra are extracted around each peak. For both the Si 2p and In 3d5∕2 peaks, we use a discrete BE grid with a constant step of 𝛥𝐸 = 0.05 eV, identical to the sampling interval used for the survey spectra.In Fig.  10, panels (a), (d), and (g) show estimates obtained using the proposed surrogate model; panels (b), (e), and (h) show estimates based on attenuation of the Si 2p peak area using the conventional method, and panels (c), (f), and (i) show estimates based on attenuation of the In 3d5∕2 peak area using the conventional method. The intensity scale (proportional to the data acquisition time) is set at 𝐼 ∼ 104, 𝐼 ∼ 103, and 𝐼 ∼ 102 in the top, middle, and bottom rows, respectively. Black dots denote individual thickness estimates (𝑧̂) for each dataset, and shading indicates the density of estimates, with darker shading representing higher densities.For reference, we indicate a vertical lines corresponding to the information depths for 95% of the signal, 𝐼𝐷95%, in Fig.  10. For an exponential attenuation model in the peak-area method, the detected peak intensity scales with exp[−𝑧∕(𝜆 cos 𝜃)]; therefore, the depth from which 95% of the analytical signal originates is𝐼𝐷95% = 𝜆 cos 𝜃 ln[ 11 − 0.95]≈ 3.0 𝜆 cos 𝜃.Using 𝜃 = 45◦, this 𝐼𝐷95% line provides a practical upper bound for quantitative peak-area analysis, beyond which the peak becomes strongly noise-dominated and the inversion becomes unstable.The film thickness evaluated using the proposed method is much less scattered than that evaluated using the conventional method, and this advantage is the most pronounced for short acquisition times (𝐼 ∼102) and at higher film thicknesses (𝑧). We attribute this performance to the use of more data points from a wide-energy-range (survey) XPS spectrum and the explicit use of the peak background component in the estimation. In addition, within our validation range up to 𝑧 = 10.0 nm, 10 the proposed method returned valid thickness estimates for all noisy realizations, i.e., no non-estimable cases were observed, whereas the conventional peak-area inversion can fail when the logarithm argu-ments become non-positive. Consistent with this hypothesis, Table  1 shows that whereas at the top of the peak, TEAL is 𝐿𝑇 (𝐸) = 2.32 and 1.86 nm, in the peak background region at 950.0 eV, TEAL increases to 𝐿𝑇 = 30.39 and 20.91 nm, respectively. Thus, the peak-background region provides non-negligible thickness sensitivity even at larger 𝑧(up to 𝑧 = 10.0 nm in this study). By exploiting this background, the proposed approach achieves high-accuracy thickness estimates.3.3.3. Comparison of spectrum generation timesIn this study, the surrogate model enabled rapid spectrum gener-ation. We compared the spectrum generation times of the surrogate model and the NIST Simulation of Electron Spectra for Surface Analysis (SESSA) on a benchmark machine equipped with an Intel(R) Core(TM) i7-10870H @ 2.20 GHz (8 cores) and 32.0 GB of random access memory, using In2O3∕Si as the benchmark sample.In SESSA, the generation time depends on the Monte Carlo con-vergence criterion (convergence rate, CR). We therefore measured the generation time at CR = 10−4, 10−5, and 10−6 and report the results as relative times normalized to unity for the surrogate model (surrogate model = 1). SESSA required 1.5 × 104, 1.07 × 105, and 7.94 × 105times longer than the surrogate model at CR = 10−4, 10−5, and 10−6, respectively. Thus, even under the loosest convergence setting (CR =10−4), the surrogate model accelerates spectrum generation by four orders of magnitude, and under the strictest setting tested (CR = 10−6) by nearly six orders of magnitude, enabling rapid structural inference.4. Summary and discussionIn this study, we proposed a data-driven surrogate model that enables quantitative analysis of wide-energy-range (survey) XPS spectra by explicitly exploiting peak-background features, with an application to film thickness estimation for a bilayer In2O3∕Si system. Using SESSA-generated orbital-resolved spectra, we expressed the measured survey spectrum as a linear combination of orbital-resolved spectral compo-nents and modeled the thickness dependence of each component using response functions (Eq. (5)). We confirmed that, over the tested binding energies including the background region, the response functions vary approximately linearly with thickness on a logarithmic scale (Fig.  3) and therefore can be approximated by exponential forms (Eqs. (8) and (9)). By fitting these exponential forms to simulator-derived response S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 10. Estimated film thickness (𝑧̂) plotted as a function of the designated thickness (𝑧) for spectra recorded at different XPS acquisition times (overall intensity scale: top row 𝐼 ∼ 104; middle row 𝐼 ∼ 103; bottom row 𝐼 ∼ 102). (a), (d), and (g): thickness estimates obtained using the proposed surrogate-model approach. (b), (e), and (h): thickness estimates obtained using the conventional peak-area attenuation method based on the Si 2p peak. (c), (f), and (i): thickness estimates obtained using the conventional peak-area attenuation method based on the In 3d5∕2 peak. The vertical red and blue lines indicate the information depth for 95% of the signal, 𝐼𝐷95% = 𝜆 cos 45◦ ln[1∕(1 − 0.95)]. For the conventional method, realizations for which the logarithm arguments become non-positive are treated as failures (non-estimable cases) and are excluded from the plotted estimates. functions (Fig.  4), we determined the total electron attenuation length (TEAL, 𝐿𝑇 (𝐸)) and compared it with the IMFP 𝜆(𝐸) derived from the TPP-2M formula (Table  1). In the peak-background region at binding energy 𝐸 = 950.0 eV, TEAL increases to 𝐿𝑇 = 30.39 nm for Si 2p3∕2and 𝐿𝑇 = 20.91 nm for In 3d5∕2, and the ratio 𝐿𝑇 (𝐸)∕𝜆(𝐸) substantially increases relative to that at the peak (Table  1 and Fig.  5).We validated the surrogate model by comparing surrogate-generatedspectra with those generated using SESSA and confirmed that the surro-gate reproduces the spectra well for both the survey spectrum and the orbital-resolved components (Figs.  6 and 7). Using a Bayesian inference framework with a Poisson observation model (Eq. (13)), we inferred the posterior distribution of the film thickness from a single survey spectrum (Fig.  8). For prescribed thicknesses of 𝑧 = 1.00, 5.00, and 10.00 nm at an overall intensity scale of 𝐼 ∼ 103, the posterior means and 95% confidence intervals were 𝑧̂ = 0.84 ± 0.02 nm, 4.75 ± 0.12 nm, 11 and 13.22±0.49 nm, respectively (Fig.  8). Across 𝑧 = 0.50–10.00 nm, the estimation trends were broadly similar for different overall intensity scales (𝐼 ∼ 104, 103, 102), while the posterior widths reflected the noise level (Fig.  9). Compared with the conventional peak-area attenuation method (Eqs. (16) and (17)), the proposed method yielded less scat-tered estimates, particularly at short acquisition times (𝐼 ∼ 102) and larger 𝑧, and returned valid thickness estimates for all noisy realizations within the validation range up to 𝑧 = 10.0 nm, whereas the conventional inversion can fail when the logarithm arguments become non-positive (Fig.  10). In addition, the surrogate model substantially reduced the spectrum generation time: relative generation times (surrogate model = 1) were 15,000, 107,000, and 794,000 for SESSA at CR = 10−4, 10−5, and 10−6, respectively, and the surrogate reduced the generation time by approximately a factor of 106 relative to SESSA at CR = 10−6.S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 For thicker overlayers, the Bayesian estimates exceeded the true thickness (e.g., 𝑧 = 10.00 nm and 𝑧̂ = 13.22 ± 0.49 nm; Fig.  8), and the Discussion attributes this overestimation to the underestimation of the exponential fitting near 𝑧 = 10.00 nm for the response functions in the background region (Fig.  4). Accordingly, further refinement of the response-function approximation in the peak-background region is needed to improve the accuracy for thicker films.However, because TEAL is obtained by empirical fitting to nu-merical data calculated using SESSA rather than from first principles, and because it can depend on both the material system and the mea-surement condition, its range of applicability and generalizability are subject to certain limitations. For a given material system and mea-surement condition, the surrogate model can be constructed offline using SESSA-generated spectra and then reused for subsequent infer-ence without running SESSA in the search loop. Because the response functions (and thus TEAL) can be determined in advance from a limited number of SESSA simulations for the specified system, the subsequent analysis can be carried out without repeated SESSA-based fitting, sub-stantially reducing the computational cost. Therefore, even though TEAL depends on the material system and may also depend on the mea-surement condition, the practical advantage of the proposed method remains unchanged. The proposed method does not assume a universal TEAL; rather, it enables rapid inverse analysis for a given material system and measurement condition by replacing repeated SESSA-based calculations with a surrogate model.Analysis beyond film thickness to more complex sample structures will require extending the proposed approach to additional material systems and further advancing the surrogate model.CRediT authorship contribution statementShunichi Yoneda: Writing – original draft, Visualization, Valida-tion, Software, Methodology, Formal analysis. Ryo Murakami: Writing – review & editing, Supervision. Hiroshi Shinotsuka: Writing – re-view & editing, Supervision. Hideki Yoshikawa: Writing – review & editing, Supervision. Shigeo Tanuma: Writing – review & editing, Supervision. Hiromi Tanaka: Writing – review & editing, Supervision. Hayaru Shouno: Writing – review & editing, Supervision. Kenji Na-gata: Writing – review & editing, Supervision, Project administration, Conceptualization.Funding sourceThis work was supported by the National Institute for Materials Sci-ence (NIMS) through the NIMS Collaboration Hub Promotion Scheme and the NIMS Internship Program.Declaration of competing interestThe authors declare that they have no known competing finan-cial interests or personal relationships that could have appeared to influence the work reported in this paper.AcknowledgmentThe authors gratefully acknowledge the National Institute for Mate-rials Science (NIMS) for its support of this work.Appendix A. NotationThis section summarizes the main variables used in this paper and their meanings.12 Sample- and layer-structure-related notation𝐴 : Constituent element of the substrate layer (Si in this study).𝐵 : Constituent element of the surface layer (In2O3 in this study).𝐶 : Bilayer sample 𝐶 = 𝐵∕𝐴 = In2O3∕Si.𝑧 : Film thickness of the surface layer (𝐵) (nm).𝑧̂ : Estimated film thickness obtained using Bayesian inference (e.g.,posterior mean).𝑚𝐴 : Total number of electron orbitals of all the constituent elements in sample 𝐴.𝑚𝐵 : Total number of electron orbitals of all the constituent elements in sample 𝐵.𝑖 : Index assigned to the electron orbitals of the constituent elements in the substrate layer (𝐴) (𝑖 = 1,… , 𝑚𝐴).𝑗 : Index assigned to the electron orbitals of the constituent elements in the surface layer (𝐵) (𝑗 = 1,… , 𝑚𝐵).𝑘 : Index of data points used for response function fitting (index for the sum of squared residuals, SSRs).Energy- and angle-related notation𝐸 : Binding energy.𝐸𝑛 : Discretized binding energy point (𝑛 = 1,… , 𝑁).𝑛 : Index of the discretized binding energy point (𝐸𝑛).𝑁 : Total number of energy/data points.𝛥𝐸 : Energy lost by a photoelectron because of multiple scatterings.𝐸𝑘 : Kinetic energy.𝐸𝑘𝑖 : Kinetic energy of the peak derived from electron orbital 𝑖 of the constituent elements in the substrate layer (𝐴).𝐸𝑘𝑗 : Kinetic energy of the peak derived from electron orbital 𝑗 of the constituent elements in the surface layer (𝐵).𝐸′𝑘 = 𝐸𝑘 − 𝛥𝐸 : Kinetic energy at which the photoelectron is detected after energy loss.BE : Abbreviation for binding energy.𝜃 : Emission angle of the photoelectron.Spectrum-related notationWide-energy-range spectra (energy distributions)𝐼𝐴(𝐸, 𝜃) : Measured spectrum of bulk sample 𝐴.𝐼𝐵(𝐸, 𝜃) : Measured spectrum of bulk sample 𝐵.𝐼𝐶 (𝑧, 𝐸, 𝜃) : Measured spectrum of the bilayer sample (𝐶) having a film thickness of 𝑧.𝐼 𝑖𝐴(𝐸, 𝜃) : Component spectrum derived from electron orbital 𝑖 of the constituent elements in sample 𝐴.𝐼 𝑗𝐵(𝐸, 𝜃) : Component spectrum derived from electron orbital 𝑗 of the constituent elements in sample 𝐵.S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) : Spectrum of sample 𝐶 derived from electron orbital 𝑖 of the constituent elements in the substrate layer (𝐴).𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) : Spectrum of sample 𝐶 derived from electron orbital 𝑗 of the constituent elements in the surface layer (𝐵).𝐼(𝐸𝑛, 𝜃) : Simulated spectrum generated using SESSA; the film thick-ness (𝑧) is unknown.𝐼(𝐸, 𝜃; 𝑧) : Approximate spectrum generated using the surrogate model constructed from the SESSA output (approximation of 𝐼𝐶 ).Peak areas (quantities used in the conventional method)𝐼 𝑖𝐴(𝜃) : Area (energy-integrated intensity) of the photoelectron peak derived from electron orbital 𝑖 of the constituent elements in bulk sample 𝐴.𝐼 𝑗𝐵(𝜃) : Area of the photoelectron peak derived from electron orbital 𝑗of the constituent elements in bulk sample 𝐵.𝐼 𝑖𝐶𝐴(𝑧, 𝜃) : Peak area corresponding to electron orbital 𝑖 of the con-stituent elements in the substrate layer (𝐴) of sample 𝐶.𝐼 𝑗𝐶𝐵(𝑧, 𝜃) : Peak area corresponding to electron orbital 𝑗 of the con-stituent elements in the surface layer (𝐵) of sample 𝐶.Response-function- and attenuation-length-related notationResponse functions𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) : Response function corresponding to electron orbital 𝑖 of the constituent elements in the substrate layer (𝐴).Definition: 𝐼 𝑖𝐶𝐴(𝑧, 𝐸, 𝜃) = 𝑓 𝑖𝐴(𝑧, 𝐸, 𝜃) 𝐼 𝑖𝐴(𝐸, 𝜃).𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃) : Response function corresponding to electron orbital 𝑗 of the constituent elements in the surface layer (𝐵).Definition: 𝐼 𝑗𝐶𝐵(𝑧, 𝐸, 𝜃) = {1 − 𝑓 𝑗𝐵(𝑧, 𝐸, 𝜃)} 𝐼 𝑗𝐵(𝐸, 𝜃).𝒇𝐴(𝑧, 𝐸, 𝜃) : Response function vector for the substrate layer (𝐴), (𝑓 1𝐴(𝑧, 𝐸, 𝜃),… , 𝑓𝑚𝐴𝐴 (𝑧, 𝐸, 𝜃))⊤.𝒇𝐵(𝑧, 𝐸, 𝜃) : Response function vector for the surface layer 𝐵, (𝑓 1𝐵(𝑧, 𝐸,𝜃),… , 𝑓𝑚𝐵𝐵 (𝑧, 𝐸, 𝜃))⊤.Attenuation lengths and length scales𝐿𝑇𝑖𝐴(𝐸) : Total electron attenuation length (TEAL) corresponding to electron orbital 𝑖 of the constituent elements in the substrate layer (𝐴).𝐿𝑇𝑗𝐵(𝐸) : TEAL corresponding to electron orbital 𝑗 of the constituent elements in the surface layer (𝐵).𝐿(𝐸) : Effective attenuation length (EAL); attenuation length used in the conventional peak area decay model.𝜆(𝐸) : Inelastic mean free path (IMFP); the mean free path where the photoelectron peak intensity decays to 1∕𝑒.TEAL : Abbreviation for the total electron attenuation length (used as 𝐿𝑇 (𝐸), etc., in equations).IMFP : Abbreviation for the inelastic mean free path (used as 𝜆(𝐸) in equations).𝐼𝐷95% : Information depth from which 95% of the analytical signal originates,𝐼𝐷95% = 𝜆(𝐸) cos 𝜃 ln[ 11 − 0.95].13 Notation for Bayesian estimation and the statistical model𝑦𝑛 : Observed counts at energy point 𝐸𝑛 (assumed to follow a Poisson distribution).𝒚 = (𝑦1,… , 𝑦𝑁 ) : The sequence of counts over the entire spectrum. = {(𝐸𝑛, 𝑦𝑛)}𝑁𝑛=1 : The set of observed data.𝛩 = {𝑧, ℎ} : The set of parameters (film thickness and intensity scale).ℎ : The intensity correction parameter.𝑓 (𝐸𝑛;𝛩) : The fitting function.𝑓 (𝐸𝑛;𝛩) = ℎ × 𝐼(𝐸𝑛, 𝜃; 𝑧) ≈ ℎ × 𝐼(𝐸𝑛, 𝜃; 𝑧).𝜆𝑛 : The mean of the Poisson distribution (𝜆𝑛 = 𝑓 (𝐸𝑛;𝛩)).𝑃 ( ∣ 𝛩) : Likelihood (the product of Poisson distributions at each energy point).𝑃 (𝛩) : Prior distribution.𝑃 (𝛩 ∣ ) : Posterior distribution.𝑃 () : Marginal likelihood (normalizing constant).𝑎𝑧, 𝑏𝑧 : Hyperparameters of the gamma prior distribution Gam(𝑧 ∣𝑎𝑧, 𝑏𝑧) for the film thickness (𝑧).𝑎ℎ, 𝑏ℎ : Hyperparameters of the gamma prior distribution Gam(ℎ ∣𝑎ℎ, 𝑏ℎ) for the intensity correction parameter (ℎ).S/N ratio : Signal-to-noise ratio (discussed in the figure captions as differences in the overall intensity scale, e.g., 𝐼 ≈ 104, 103, 𝑎𝑛𝑑102).Notation for the conventional method and quantification limits𝑝fail,𝐴(𝑧) : Failure probability of the conventional inversion for the substrate peak-area method.𝑝fail,𝐵(𝑧) : Failure probability of the conventional inversion for the overlayer peak-area method.𝑧QL,𝐴 : Quantification limit defined by 𝑝fail,𝐴(𝑧) < 0.05.𝑧QL,𝐵 : Quantification limit defined by 𝑝fail,𝐵(𝑧) < 0.05.I[⋅] : Indicator function.Experimental conditions and other notation𝜎 : The standard deviation of each Gaussian-shaped peak (the width of the peak profile).CR : The convergence condition of the Monte Carlo calculation in SESSA (the convergence ratio, e.g., 𝐶𝑅 = 10−6).Appendix B. Experimental conditionsWe detail the spectrum generation conditions used to construct the surrogate model. Spectra were generated using SESSA, the sample structure and geometry shown in Fig.  11, and the following parameters. The target sample was In2O3∕Si. Survey spectra were generated over a BE range from 0.0 eV to 1000.0 eV at a step of 0.5 eV.For each electron orbital, the shape of the photoelectron peaks are modeled as a Gaussian distribution, and the peak positions and standard deviations (𝜎) are set at the values in Table  2. For the measurement conditions, the X-ray source is monochromatic Al K𝛼radiation. The angle between the incident X-rays and photoelectron detector is 45◦, and the X-rays are incident normal to the sample surface. Additionally, the detector’s aperture (acceptance) angle is set at ±20◦.S. Yoneda et al. Journal of Electron Spectroscopy and Related Phenomena 286 (2026) 147612 Fig. 11. (a) Sample structure and (b) measurement geometry for In2O3∕Si.Table 2Peak positions and Gaussian standard deviations (𝜎) for each electron orbital. Core level Peak position (eV) Standard deviation (eV)   Si 2s 150.8 0.7   Si 2p1∕2 100.1 0.5   Si 2p3∕2 99.5 0.5   In 3s 828.4 1.7   In 3p1∕2 704.9 1.5   In 3p3∕2 667.0 1.6   In 3d3∕2 453.1 0.6   In 3d5∕2 445.6 0.7   In 4s 124.4 1.2   In 4p1∕2 79.5 4.5   In 4p3∕2 79.5 4.5   In 4d3∕2 19.0 0.9   In 4d5∕2 18.2 0.9   O 1s 532.2 1.2   O KL23L23 976.6 4.1  Appendix C. Application of the surrogate model for film thickness estimationsWe outline the Bayesian inference used in Section 3.1.By constructing the surrogate model (𝐼(𝐸𝑛, 𝜃; 𝑧)), we approximate the SESSA output for any film thickness (𝑧) in a given layer stack. Accordingly, for a sample with a known layer configuration but an unknown thickness (𝑧), we estimate 𝑧 by fitting 𝐼(𝐸𝑛, 𝜃; 𝑧) to the XPS spectral data ( = (𝐸𝑛, 𝑦𝑛(𝐸𝑛, 𝜃))𝑁𝑛=1). To assess the reliability of the estimates, we adopt a Bayesian framework.We adopt the surrogate of the SESSA output as the fitting function (𝑓 (𝐸𝑛;𝛩)) as follows:𝑓 (𝐸𝑛;𝛩) = ℎ × 𝐼(𝐸𝑛, 𝜃; 𝑧)≈ ℎ × 𝐼(𝐸𝑛, 𝜃; 𝑧), (21)where the parameter set is 𝛩 = 𝑧, ℎ, and ℎ is an intensity-scaling parameter.To estimate the film thickness (𝑧) using the surrogate model (𝐼(𝐸𝑛,𝜃; 𝑧)) and Bayesian inference, we posit a spectrum generation model. The spectral measurements (𝑦𝑛(𝐸𝑛, 𝜃)) are modeled as a Poisson distri-bution (Poi(𝑦𝑛(𝐸𝑛, 𝜃) ∣ 𝜆)) as follows:𝑦𝑛(𝐸𝑛, 𝜃) ∼ 𝑃 (𝑦𝑛(𝐸𝑛, 𝜃) ∣ 𝛩) = Poi(𝑦𝑛(𝐸𝑛, 𝜃) ∣ 𝜆 = 𝑓 (𝐸𝑛;𝛩))=𝑓 (𝐸𝑛, 𝛩)𝑦𝑛(𝐸𝑛 ,𝜃) exp(−𝑓 (𝐸𝑛, 𝛩))𝑦𝑛(𝐸𝑛, 𝜃)!. (22)Assuming independence across 𝑛, we have 𝑃 ( ∣ 𝛩) =𝑁∏𝑛=1𝑃 (𝑦𝑛(𝐸𝑛, 𝜃) ∣ 𝛩) =𝑁∏𝑛=1𝑓 (𝐸𝑛, 𝛩)𝑦𝑛(𝐸𝑛 ,𝜃) exp(−𝑓 (𝐸𝑛, 𝛩))𝑦𝑛(𝐸𝑛, 𝜃)!.(23)Our goal is to characterize the posterior distribution (𝑃 (𝛩 ∣ )) of the parameters (𝛩) for the measured data (). From the product rule, 14 Table 3Hyperparameters for the prior distributions. Hyperparameter Value   𝑎𝑧 1.50   𝑏𝑧 0.0167   𝑎ℎ 3.00   𝑏ℎ 0.00000500 𝑃 (, 𝛩) = 𝑃 ( ∣ 𝛩)𝑃 (𝛩) = 𝑃 (𝛩 ∣ )𝑃 (), we obtain 𝑃 (𝛩 ∣ ) =𝑃 ( ∣ 𝛩)𝑃 (𝛩)𝑃 ()∝ 𝑃 ( ∣ 𝛩)𝑃 (𝛩), (24)where 𝑃 () is the evidence (a constant with respect to 𝛩). The prior 𝑃 (𝛩) encodes prior information and factorizes as 𝑃 (𝛩) = 𝑃 (𝑧)𝑃 (ℎ)under the assumed independence of 𝑧 and ℎ. Because both 𝑧 and ℎ are positive, we place gamma priors Gam(𝑎, 𝑏) as follows: 𝑃 (𝛩) = 𝑃 (𝑧)𝑃 (ℎ) = Gam(𝑧 ∣ 𝑎𝑧, 𝑏𝑧)Gam(ℎ ∣ 𝑎ℎ, 𝑏ℎ), (25)where 𝑎 and 𝑏 denote the gamma shape and scale parameters (see Table 3), respectively. The hyperparameters 𝑎𝑧, 𝑏𝑧, 𝑎ℎ, and 𝑏ℎ are listed in Table  3. Then, we sample from 𝑃 (𝛩 ∣ ) ∝ 𝑃 ( ∣ 𝛩)𝑃 (𝛩), using the Metropolis algorithm [11] to draw from the joint posterior and obtain marginal posteriors 𝑃 (𝑧 ∣ ) and 𝑃 (ℎ ∣ ). In this study, the intensity-scaling parameter (ℎ) was also analytically estimated [12].Data availabilityData will be made available on request.References[1] H. Shinotsuka, K. Nagata, M. Siriwardana, H. Yoshikawa, H. Shouno, M. Okada, Sample structure prediction from measured XPS data using Bayesian estimation and SESSA simulator, J. Electron Spectrosc. Relat. Phenom. 267 (2023) 147370, http://dx.doi.org/10.1016/j.elspec.2023.147370.[2] J. Matthew, Surface analysis by auger and X-ray photoelectron spectroscopy, in: D. Briggs, J.T. Grant (Eds.), Surface Analysis By Auger and X-Ray Photoelectron Spectroscopy, IM Publications, Chichester, UK and SurfaceSpectra, 2003.[3] A. Jablonski, C.J. Powell, Effective attenuation lengths for different quantitative applications of X-ray photoelectron spectroscopy, J. Phys. Chem. Ref. Data 49 (2020) 033102, http://dx.doi.org/10.1063/5.0008576.[4] S. Tougaard, Surface nanostructure determination by X-ray photoemission spec-troscopy peak shape analysis, J. Vac. Sci. Technol. A 14 (3) (1996) 1415–1423, http://dx.doi.org/10.1116/1.579963.[5] S. Tougaard, Low energy inelastic electron scattering properties of noble and transition metals, Solid State Commun. 61 (9) (1987) 547–549, http://dx.doi.org/10.1016/0038-1098(87)90166-9.[6] Siegfried Hofmann, Qualitative analysis (principle and spectral interpretation), in: Auger- and X-Ray Photoelectron Spectroscopy in Materials Science: A User-Oriented Guide, Springer Berlin Heidelberg, Berlin, Heidelberg, 2013, pp. 43–76, http://dx.doi.org/10.1007/978-3-642-27381-0_3.[7] W.S.M. Werner, W. Smekal, C.J. Powell, J. Gorham, Simulation of Electron Spectra for Surface Analysis (SESSA) Version 2.2 User’s Guide, Technical Report NIST NSRDS 100, National Institute of Standards and Technology, Gaithersburg, MD, 2021, http://dx.doi.org/10.6028/NIST.NSRDS.100-2021.[8] A. Machida, K. Nagata, R. Murakami, H. Shinotsuka, H. Shouno, H. Yoshikawa, M. Okada, Bayesian inference method utilizing SESSA in quantitative layer structure estimation from XPS data, J. Electron Spectrosc. Relat. Phenom. 273 (2024) 147449, http://dx.doi.org/10.1016/j.elspec.2024.147449.[9] SciPy Community, SciPy documentation: scipy.optimize.curve_fit, 2024, URL https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.html. (Accessed 23 October 2024).[10] S. Tanuma, C.J. Powell, D.R. Penn, Calculations of electron inelastic mean free paths. V. Data for 14 organic compounds over the 50–2000 eV range, Surf. Interface Anal. 21 (3) (1994) 165–176, http://dx.doi.org/10.1002/sia.740210302.[11] W.K. Hastings, Monte Carlo sampling methods using Markov chains and their ap-plication, Biometrika 57 (1) (1970) 97–109, http://dx.doi.org/10.1093/biomet/57.1.97.[12] Persi Diaconis, Donald Ylvisaker, Conjugate priors for exponential families, Ann. Statist. 7 (2) (1979) 269–281, http://dx.doi.org/10.1214/aos/1176344611.http://dx.doi.org/10.1016/j.elspec.2023.147370http://refhub.elsevier.com/S0368-2048(26)00021-6/sb2http://refhub.elsevier.com/S0368-2048(26)00021-6/sb2http://refhub.elsevier.com/S0368-2048(26)00021-6/sb2http://refhub.elsevier.com/S0368-2048(26)00021-6/sb2http://refhub.elsevier.com/S0368-2048(26)00021-6/sb2http://dx.doi.org/10.1063/5.0008576http://dx.doi.org/10.1116/1.579963http://dx.doi.org/10.1016/0038-1098(87)90166-9http://dx.doi.org/10.1016/0038-1098(87)90166-9http://dx.doi.org/10.1016/0038-1098(87)90166-9http://dx.doi.org/10.1007/978-3-642-27381-0_3http://dx.doi.org/10.6028/NIST.NSRDS.100-2021http://dx.doi.org/10.1016/j.elspec.2024.147449https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.htmlhttps://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.htmlhttps://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.curve_fit.htmlhttp://dx.doi.org/10.1002/sia.740210302http://dx.doi.org/10.1002/sia.740210302http://dx.doi.org/10.1002/sia.740210302http://dx.doi.org/10.1093/biomet/57.1.97http://dx.doi.org/10.1093/biomet/57.1.97http://dx.doi.org/10.1093/biomet/57.1.97http://dx.doi.org/10.1214/aos/1176344611 A data-driven surrogate model for X-ray photoelectron spectroscopy based on survey spectrum background features Introduction Construction of the Surrogate Model Modeling Generation of Orbital-Resolved Spectra Response Functions Formulation of the Response Functions Determination of Response Function Parameters Generation of Spectra Using the Surrogate Model Film Thickness Estimation Using the Surrogate Model Film Thickness Evaluation Using Bayesian Inference Estimation Range of the Proposed Method Discussion Effectiveness of the Surrogate Model Comparison of the Proposed and Conventional Methods Comparison of Spectrum Generation Times Summary and Discussion CRediT authorship contribution statement Funding Source Declaration of competing interest Acknowledgment Appendix A. Notation Sample- and layer-structure-related notation Energy- and angle-related notation Spectrum-related notation Wide-energy-range spectra (energy distributions) Peak areas (quantities used in the conventional method) Response-function- and attenuation-length-related notation Response functions Attenuation lengths and length scales Notation for Bayesian estimation and the statistical model Notation for the conventional method and quantification limits Experimental conditions and other notation Appendix B. Experimental Conditions Appendix C. Application of the Surrogate Model for Film Thickness Estimations Data availability References