Alternative difference analysis scheme combining R-space EXAFS fit with global optimization XANES fit for X-ray transient absorption spectroscopy

Time-resolved X-ray absorption spectroscopy (TR-XAS), based on the laser-pump/X-ray-probe method, is powerful in capturing the change of the geometrical and electronic structure of the absorbing atom upon excitation. TR-XAS data analysis is generally performed on the laser-on minus laser-off difference spectrum. Here, a new analysis scheme is presented for the TR-XAS difference fitting in both the extended X-ray absorption fine-structure (EXAFS) and the X-ray absorption near-edge structure (XANES) regions. R-space EXAFS difference fitting could quickly provide the main quantitative structure change of the first shell. The XANES fitting part introduces a global non-derivative optimization algorithm and optimizes the local structure change in a flexible way where both the core XAS calculation package and the search method in the fitting shell are changeable. The scheme was applied to the TR-XAS difference analysis of Fe(phen)₃ spin crossover complex and yielded reliable distance change and excitation population.

1. Introduction

Time-resolved X-ray absorption spectroscopy (TR-XAS), based on the laser-pump/X-ray-probe method, is a powerful probe for addressing the geometric and electronic structure of light-generated transient species (Chen et al., 2014; Chergui, 2016). Its time resolution capability is dependent on the X-ray pulse duration. Generally, it can reach 100 ps temporal resolution using synchrotron radiation sources and less than 100 fs using X-ray free-electron lasers. TR-XAS has been widely applied in photochemistry and photophysics fields. Bressler et al. studied the light-induced spin crossover of an iron complex (Bressler et al., 2009). Zhang et al. obtained a 0.01 Å high-resolution structure change of the excited state of an osmium complex (Zhang et al., 2014). Moonshiram et al. revealed the intermediate molecular structure of Co(II) and Co(I) photocatalyst in real H₂ production (Moonshiram et al., 2016). Wen et al. studied BiFeO₃ nanofilm and found an anisotropic unit cell change of the in-plane contraction and out-of-plane extension (Wen et al., 2015).

To benefit from the well developed X-ray absorption fine-structure (XAFS) spectrum analysis methods developed over the last 40 years, different TR-XAS treatments or packages have been developed in both the EXAFS and XANES regions along with extensive TR-XAS applications. The differential EXAFS fitting method in q-space was applied by Borfecchia et al. to study the photoactive metal complexes cis-[Ru(bpy)₂(py)₂] (Borfecchia et al., 2014; Garino et al., 2014). Zhang et al. adopted a fit of the excited-state EXAFS reconstructed according to different excitation fractions (Zhang et al., 2015). FEFF, a popular XANES computation package based on the self-consistent multiple-scattering theory, is used in the study of the PtPOP system by van der Veen et al., who took Kas's Bayes XANES fit method to search for structure change information when the system is excited (van der Veen et al., 2010). The interpolation approach is effective in structure refinement: Smolentsev et al. used it to invoke FEFF and FDMNES in TR-XAS fitting research on platinum dimer (Lockard et al., 2010). They also used the DFT-MO method to simulate the pre-edge feature of different cobalt species (Moonshiram et al., 2016). MXAN (Benfatto et al., 2003), one of the widely used XANES fit algorithms, uses CONTINUUM to perform the XAFS calculation and adopts the MINUIT optimization package. MXAN has been applied to the characterization of the molecular excitation state of [Ru^II(bpy)₃]²⁺ by Benfatto et al. (2006).

Here we present an alternative data analysis scheme for TR-XAS. It combines the R-space EXAFS fit with the global optimization XANES fit for TR-XAS analysis. The fit of the difference spectrum is performed directly in both the EXAFS and XANES regions. We adopt the R-space fit to avoid one more reverse Fourier transform. The XANES difference spectrum fit is similar to that used in MXAN but it is more flexible, as we can choose different calculation cores, like FEFF (Rehr et al., 2010) or FDMNES (Joly et al., 2009), and different adaptive optimization algorithms as well to control the fit process.

2. Method

2.1. Workflow

We present here a TR-XAS fit scheme in both the EXAFS and XANES regions where the R-space difference fit in the EXAFS gives the average bond length changes of the first shell around the absorber, while the detailed structural change is derived in the XANES difference fit. Our XANES fit scheme is flexible, where both the search algorithm (shell) for variables and the spectrum calculation package (core) can be substituted by any suitable package.

The workflow of our TR-XAS fit scheme is composed of two blocks, as shown in Fig. 1, the EXAFS difference fit (block 1) and the XANES difference fit (block 2): each one starts from finding the ground state structure and is followed by searching for the structural change of the intermediate state. The spectrum of the ground state or intermediate state is first calculated using the general EXAFS formula (Newville, 2013) in block 1 or the selected XAS calculation package in block 2. The theoretical difference spectrum is obtained according to the equation

where α is the fraction of the intermediate state, and and are the theoretical spectra of the ground state and the intermediate state, respectively. The theoretical difference spectrum is compared with the experimental difference spectrum to check their consistency in block 2. They are, however, firstly Fourier transformed into R-space in block 1 before comparison. The structure of the intermediate state and α may be modified by a parameter search shell to perform the difference spectrum calculation again, following a comparison with the experimental spectrum until the consistency is acceptable.

Figure 1 (a) Flowchart of the EXAFS difference fit. GS stands for ground state, IntS for intermediate state, diff for difference, th for theoretical and exp for experimental. Data operations are presented in square frames while input/output operations are presented in elliptical frames. (b) Flowchart of the XANES difference fit. The denotions and interpretation of the frames are the same as those in (a).

We give here the definition of the goodness of fit, the minimum of which in the parameter space is what the difference fit searches for,

where x_j is the parameter to be fit, n is the number of independent parameters, m is the number of data points, ∊_i is the individual error in the experimental data set, and w_i is the statistical weight. When w_i = 1, the function reduces to the statistical χ² function.

The conclusion of the difference fit in the EXAFS region is helpful for the XANES fit, as one could limit the varying range of some of the parameters in order to save the fit time.

TR-XAS data analysis is generally performed on the laser-on minus laser-off difference spectrum. The difference highlights the structural change, since the difference spectrum eliminates the experimental system error. Our framework performs the fit directly to the difference spectrum instead of to the ground state and the excitation state spectrum separately. Moreover, pre-processing of the EXAFS spectrum, such as background subtraction, is not needed for the difference fit as the backgrounds are the same in the ground state and intermediate state. Obviously we can see from equation (1) that the fraction of the intermediate state α in the system after pump is critical in producing the XANES fit. There are several ways to estimate α in TR-XAS analysis. Moonshiram et al. estimated the Co(II) intermediate fraction using the Co(II) reference sample's pre-edge fingerprint (Moonshiram et al., 2016). This would fail when the analogue sample is difficult to synthesize or when it is unstable in the experiment, which occurs in most cases. We can estimate the fraction from an optical transient absorption (OTA) experiment carried out under the same conditions as for the TR-XAS (Benfatto et al., 2006; Borfecchia et al., 2013), but some systems are OTA silent, the intermediate in the system cannot be detected by OTA spectroscopy. If we are unable to apply these two methods, an alternative choice for obtaining the intermediate fraction is by treating it as a variable like other structural parameters in the TR-XAS fit (Lockard et al., 2010; Zhang et al., 2015). In our scheme, the intermediate fraction could be treated as a variable or constant.

We use R-space instead of k-space for the EXAFS fit. It is difficult to achieve a high-quality TR-XAS signal across the full EXAFS range. Moreover, it is clear and straightforward to find the trivial difference between the simulation χ(R) and experiment in R-space. Generally the major structure change in the system comes from the nearest atoms to the absorber, e.g. those in the first shell: we can limit the fitting domain to the first shell in R-space during the fit, and the calculation of χ will be much faster as we do not need to consider scattering from other shells.

In the XANES fit, as described above, the core simulates the spectrum and the shell controls the parameter search. FEFF (Rehr et al., 2010) is used to carry out the XANES calculation by default in the scheme, and can be replaced by FDMNES (Joly et al., 2009) at present. The search method adopted by the shell is also substitutable, such as by the mesh adaptive direct search (MADS), the dividing rectangles method (Gablonsky & Kelley, 2001) or the improved stochastic ranking evolution strategy (Runarsson & Yao, 2000). We take here the MADS method as the default for the difference fit. It evolved from the generalized pattern search which is one of the modern grid searches. MADS has proven to be efficient (Audet et al., 2008, 2011) and is now widely used (Miiller et al., 2012; Eisenhower et al., 2012; Berrocal-Plaza et al., 2014). Moreover, NOMAD, which is one of the packages that uses the MADS method, offers a parallel version (Audet et al., 2008).

2.2. Choice of the optimization algorithm

The optimization algorithm is important in the XANES fit. Traditional structure optimization in computation chemistry adopts and develops a series of local optimization algorithms with the initial Hessian guess method and the Hessian update method. However, in crystal structure prediction, the application of heuristic global optimization algorithms, for example the particle swarm optimization in the CALYPSO code (Wang et al., 2012) and the evolutionary optimization in the USPEX code (Glass et al., 2006), has achieved great success. The choice of optimization algorithm should also be deliberated in the XANES fit. First, we prefer to obtain a global minimum rather than a local minimum. Second, since the objective function of the XANES fit is given as a time-consuming flow without analytic derivative information, instead of a simple formula versus the parameters, we do not know whether it is a convex problem before the fit, so generally we adopt a non-linear optimization algorithm. Moreover, we decide to use non-derivative algorithms, as the XANES fit also cannot provide an analytical derivative, and the cost of a numerical derivative computation of the XANES fit is huge; also, the reliability of the numerical derivative is difficult to guarantee. Fortunately many global optimization algorithms do not need information about the derivative, so a global non-derivative algorithm is our choice. The deterministic algorithm and the stochastic and heuristics algorithm are the two main branches of global optimization. Normally the deterministic algorithm can perform a theoretical convergence analysis under given conditions, while some stochastic and heuristics algorithms can only provide a `probabilistic convergence guarantee' or `remain heuristics' (Pardalos et al., 2000). For the deterministic algorithm, taking MADS as an example, Audet et al. used MADS's hierarchy of convergence analysis (Audet & Dennis Jr, 2006) and proved that the lower triangular instance of MADS can produce a dense set of poll directions. Another example is the DIRECT algorithm. This can prove its ability of convergence to the Karush–Kuhn–Tucker point under mild conditions (Gablonsky & Kelley, 2001). Though some optimization algorithms can perform convergence analysis under given conditions, the theoretical condition is usually difficult to test. Therefore the benchmark of the algorithm is important in practice. From NLOPT library's benchmark (Kumar et al., 2016) and the derivative-free algorithms benchmark (Rios & Sahinidis, 2013), we can find that usually one optimization algorithm cannot solve all test problems. Also there is an optimization solution distribution of multiple optimization solver runs due to the algorithm itself or initial solution. Obviously, there is no single ultimate algorithm for all optimization problems, so the substitutable algorithm in the shell of the XANES fit is a good choice.

We notice that the correlation of structure parameters and intermediate fraction has not been studied thoroughly in TR-XAS fitting. FITIT uses Kas's scheme of the conditioning number (Rehr et al., 2005) to deal with the correlation of structure parameters. In the study of pyrazolate-bridged platinum dimer which employed FITIT to fit the fraction (Lockard et al., 2010), the correlation between fraction and structure parameters was not reported. MXAN does not give the correlation between structure parameters. Vorobeva et al. adopted MXAN in the study of excited Ti³⁺ in Al₂O (Vorobeva et al., 2009). They fit the fraction, but they could not distinguish between elongation and compression of the octahedron without a detailed study of the correlation between the structure change and fraction. Nowadays, methods have been developed for studying the correlation between parameters such as the maximal information coefficient (Reshef et al., 2011), distance correlation (Delicado & Smrekar, 2009) and mutual information estimator (Kraskov et al., 2004). It is expected that these methods can be applied in helping to deal with this problem.

3. Application

3.1. Sample and experimental description

We applied our scheme to a model spin crossover iron complex, 1,10-phenanthroline iron(II) sulfate, Fe(II)(phen)₃. Fe(II)(phen)₃ is a low-spin (LS) complex in the ground state, and becomes a high-spin (HS) state upon photo-excitation due to spin crossover. It is important in dynamic magnetic research and has been well studied by TR-XAS. The Nozawa group reconstructed the EXAFS of the excited state of Fe(II)(phen)₃ and found a 0.17 Å bond elongation for the first shell compared with the LS ground state (Nozawa et al., 2010). Fe(II)(phen)₃ was purchased from Alfa Aesar and was used without any further purification. The HS analogue complex, Fe(II)(2-CH₃-phen)₃, was synthesized. Its HS state results from steric hindrance of the methyl group. Fe(II)(2-CH₃-phen)₃ not only provides a HS analogue for structure analysis but also helps to determine the excitation fraction.

TR-XAS measurements were performed at beamline 11-ID-D of the Advanced Photon Source. An experimental description can be found elsewhere (Chen & Zhang, 2013). The Fe(II)(phen)₃ laser-on and laser-off spectra are shown in Fig. 2 along with XAS of Fe(II)(2-CH₃-phen)₃.

Figure 2 TR-XAS of the Fe K-edge of the Fe(II)(phen)3: laser-on (red) and laser-off (black) along with XAS of Fe(II)(2-CH3-phen)3 (blue). TR-XAS was measured in fluorescence mode. The laser excitation wavelength is 400 nm with a repetition rate of 10 kHz. Inset: molecular structure of Fe(II)(phen)3; the iron atom is coordinated to six nitro­gen atoms in octahedral configuration, and each phenanthroline ligand provides two nitro­gen atoms for coordination. Molecular structure drawn by VMD (Humphrey et al., 1996).

In the EXAFS fit, we fit the ground-state Fe(II)(phen)₃ and its HS analogue Fe(II)(2-CH₃-phen)₃. The data were processed by our scheme with the EXAFS equation referring to that in Larch (Newville, 2013). The amplitude reduction factor was fixed to 0.9 and the coordination number was fixed to 6. The variables are the disorder factor σ², the energy shift E₀ and the bond length r. The theoretical amplitude and phase are extracted from the 1.95 Å Fe—N single scattering path calculated by FEFF9. The R-space fit is displayed in Fig. 3 and the fit results are listed in Table 1. The disorder factor and energy shift will be used in the following difference fit.

Table 1
Fit results of the ground state and HS analogue

r is the bond length change related to the reference Fe—N scattering path. σ² is the disorder factor. E₀ is the energy shift.

Sample	r (Å)	σ^2	E0 (eV)	R-factor
Fe(phen)3	1.97	0.009	−4.99	0.005
Fe(II)(2-CH3-phen)3	2.17	0.014	1.81	0.004

Figure 3 Experimental and simulated Fourier transforms of k^2-weighted EXAFS of the ground state (gsexp, gsfit) and reference sample Fe(II)(2-CH3-phen)3 (refexp, reffit).

As stated above, the difference between laser-on and laser-off highlights the excited signal. The difference spectrum (black) is displayed in Fig. 4, compared with the difference obtained from the HS analogue minus the ground state. The discrepancy reflects the excitation fraction; the excitation fraction is 37%. The reconstructed difference spectrum and the difference spectrum obtained from the excitation analogue are in good agreement. In the following we will use our new scheme to fit the difference spectrum in both the XANES and EXAFS regions to obtain the structure change information and excitation fraction.

Figure 4 Difference spectra of the HS excitation state (exc) and excitation analogue reference (ref). The reconstructed difference spectrum (reconstructed exc) is obtained from dividing the experimental spectrum by the excitation fraction (37%).

3.2. EXAFS difference spectrum fitting

The reverse Fourier transform of the experimental difference spectrum is shown in Fig. 5. We find that the 1–1.9 Å R-space window reveals the main features of the experimental spectrum. The k-weight of the fit is 2. We used k-space in the 3–10 Å⁻¹ region and R-space in the 1–1.9 Å region for the fit.

Figure 5 Experimental (exp) and reverse Fourier transform of the difference spectrum (rmax1.9). The R-space range is 1–1.9 Å.

We fit the intermediate's fraction and relative bond length change. The energy shift and the disorder factor, given in Table 1, are fixed during the fitting. The R-space fit is shown in Fig. 6(a). We also performed the reverse Fourier transform of the R-space fit result, as seen in Fig. 6(b). The Fe—N bond length is elongated by 0.17 Å under the excitation state. The excitation fraction is 37%. The results are in agreement with those from the HS analogue in Table 2. The difference EXAFS fit can only provide bond length information, so next we perform the difference XANES fitting.

Table 2
HS state fit result with different calculation core and search methods

Core	Search method	Bond length change (Å)	Fraction	R-factor	Structure	Image of fit result
FEFF	MADS	+0.11	39%	0.0058	Table S9	Fig. 7(b)
FEFF	DIRECT_L	+0.12	36%	0.0060	Table S1	Fig. S1
FEFF	ISRES	+0.09	45%	0.0066	Table S4	Fig. S4
FDMNES (MST)	MADS	+0.10	60.5%	0.0038	Table S10	Fig. 7(b)

Figure 6 (a) Experimental (exp) and simulated (fit) Fourier transforms of the k^2-weighted difference spectrum. Real part of the experimental (expreal) and simulated (fitreal) Fourier transforms of the k^2-weighted difference spectrum. (b) Reverse Fourier transform of the R-space fit result of the difference spectrum (fit), with the k-space (exp) and q-space (qexp) experimental difference spectrum.

3.3. XANES difference fit

First we fit the ground-state structure. We used the FEFF9/FDMNES packages for the XANES calculation and took the NOMAD package to search for the optimized variables. Multiple-scattering theory (MST) was selected in the FDMNES package. The initial structure of the LS ground state was taken from Yan et al. (2000). The radius of the cluster for self-consistent multiple scattering is 6.5 Å, containing all of the atoms in the molecule. The real Hedin–Lundqvist exchange-correlation potential was used in the calculation. The phenanthroline ring is taken as a rigid unit during the fit, and we chose the midpoint of the two ligand N atoms to be the representative of this unit. In total, nine structural variables representing three rigid bodies were fitted, three parameters (r, θ, φ) for each ligand. The best fits of FEFF and FDMNES are shown in Fig. 7(a). All of the features were reproduced by FEFF and they were much better than those reproduced by FDMNES, so we used FEFF's data here. The average Fe—N bond length from FEFF was 1.97 Å. The broadening given in FEFF by default is acceptable, as shown in Fig. 7(a), so no additional broadening was considered here. In the FDMNES fit, we used the arctangent style function representing the energy-dependent Lorentzian broadening width,

where E_F is the Fermi energy, Γ_Hole is the core level broadening width, Γ_m is the maximum of the broadening width of the final state (high energy range) and E_cent is the energy whose width is half of Γ_m. Γ_m/E_larg represents the depth of the center of the arctangent function.

Figure 7 (a) Ground-state XANES calculation of FEFF (FEFF), MST in FDMNES (MST) and FDM in FDMNES (FDM). Experimental XANES of ground state (EXP). (b) Difference XANES calculation of FEFF (FEFF), MST in FDMNES (MST) and FDM in FDMNES (FDM). Experimental difference XANES (EXP).

We fix Γ_Hole to be the default value 1.33 eV of iron in FDMNES. The fit result is E_F = 3.5 eV, Γ_m = 19.6 eV, E_cent = 16.4 eV, E_larg = 47.9 eV. When we use the finite difference full potential method (FDM) in FDMNES to calculate the spectrum based on the fitted structure via MST, the spectrum reproduces the feature better, as shown in Fig. 7. The fitted broadening parameters in FDM are: E_F = −3.6 eV, Γ_m = 23.5 eV, E_cent = 36.7 eV, E_larg = 50.0 eV, and the fitted fraction in FDM is 60%.

In the difference XANES fit, we used the FEFF9/FDMNES (MST) packages for the XANES calculation and took MADS optimization algorithms to search for the best optimized variables. Fit results using FEFF and other optimization algorithms can be found in the supporting information. We fixed the normalization factor and the energy shift of the intermediate state to be the same as those derived from the previous ground-state fit (in the FDMNES case, broadening parameters are also fixed). The fit results are listed in Table 2. Although FDMNES (MST)/NOMAD produces the lower R-factor, the spectrum of the ground state derived by FDMNES is not good enough. Moreover, its intermediate state fraction is much larger than that from the intermediate analogue and that from the difference EXAFS fit. The fit given by FEFF/NOMAD is shown in Fig. 7(b). The main features in the difference spectrum are reproduced. Coordinates of the intermediate structures are listed in the supporting information. The average bond length of the first shell extends to 0.11 Å, and the intermediate fraction is 39%, as shown in Table 3. The fitting results are close to those obtained in the R-space EXAFS fit and reported data (Nozawa et al., 2010).

Table 3
Summary of the difference fitting spectrum compared with the HS analogue result

	First shell bond length change (Å)	Intermediate fraction
HS analogue	+0.17	35%
Difference EXAFS (real and imaginary part)	+0.17	37%
Difference XANES	+0.11	39%

The difference EXAFS result could be utilized in the XANES fitting loop, as listed in the supporting information (Table S0). In `Direct_L Bound', we used the [−0.05 Å, 0.3 Å] boundary limitation of the bond length change. The limitation comes from the difference EXAFS result. The lower bound is −0.05 Å because the EXAFS gives the average bond length extension in which some bonds may contract. For comparison, in `Direct_L', we used the [−0.3 Å, 0.3 Å] boundary limitation of the bond length change. The boundary limitation from EXAFS can improve by about 30% the efficiency of the fit.

4. Conclusion and outlook

We have developed a two-step scheme for data analysis of the TR-XAS. First we fit the difference EXAFS in R-space and obtained the bond length change of the first shell coordination; then we fit the difference XANES to obtain the three-dimensional structure. The scheme is flexible in the XANES fit where both the spectrum calculation package and the global optimization algorithm for the variables are changeable. This scheme was applied in the analysis of a photo-induced spin crossover iron complex. The fit results are in agreement with the structure of the HS state analogue and the reported results. Next, some other XAFS calculation packages, such as XSPECTRA, based on the projector augmented wavefunction (Gougoussis et al., 2009), and OCEAN, based on the Bethe–Salpeter equation (Gilmore et al., 2015), will be included in our scheme. Also the correlation between the structure parameters and the intermediate fraction is a desirable research direction.