\journalcode

GuRakitaLan Thatcher Kamm O’Nolan Mcbride Wustrow Neilson Chapman \cauthor[b][email protected] \cauthor[b]Simon J. [email protected] \aff[a]School of Statistics and Data Science, Nankai University, \cityTianjin 300071, \countryChina \aff[b]Department of Applied Physics and Applied Mathematics, Fu Foundation School of Engineering & Applied Sciences, Columbia University, \cityNew York, NY 10025 \countryUSA \aff[c]Department of Chemistry, Stony Brook University, \cityStony Brook, NY 11794, \countryUSA \aff[d]Department of Chemistry, Colorado State University, \city Fort Collins, CO 80523 , \countryUSA

Stretched Non-negative Matrix Factorization

Ran Yevgeny Ling Zach Gabrielle E Daniel Brennan Allison James R Karena W

Abstract

An algorithm is described and tested that carries out a non negative matrix factorization (NMF) ignoring any stretching of the signal along the axis of the independent variable. This extended NMF model is called stretchedNMF. Variability in a set of signals due to this stretching is then ignored in the decomposition. This can be used, for example, to study sets of powder diffraction data collected at different temperatures where the materials are undergoing thermal expansion. It gives a more meaningful decomposition in this case where the component signals resemble signals from chemical components in the sample. The stretchedNMF model introduces a new variable, the stretching factor, to describe any expansion of the signal. To solve stretchedNMF, we discretize it and employ Block Coordinate Descent framework algorithms. The initial experimental results indicate that stretchedNMF model outperforms the conventional NMF for sets of data with such an expansion. A further enhancement to stretchedNMF for the case of powder diffraction data from crystalline materials called sparse-stretchedNMF, which makes use of the sparsity of the powder diffraction signals, allows correct extractions even for very small stretches where stretchedNMF struggles. As well as demonstrating the model performance on simulated PXRD patterns and atomic pair distribution functions (PDFs), it also proved successful when applied to real data taken from an in situ chemical reaction experiment.

keywords:

non-negative matrix factorization

keywords:

data expansion

keywords:

functional optimization

keywords:

pair distribution function

1 Introduction

Nonnegative matrix factorization (NMF) is an unsupervised machine learning method used for decomposing compressed data. NMF extracts distinct components from related signal sets in various research fields, including signal processing [buciuNonnegativeMatrixFactorization2008], biomedical engineering [sraNonnegativeMatrixApproximation2006], pattern recognition [cichockiFastLocalAlgorithms2009], image engineering [buciuNonnegativeMatrixFactorization2008a] and so on. NMF differs from principle component analysis (PCA) [jolliffePrincipalComponentAnalysis2002] by applying positivity constraints on the extracted components and their weights. It is then attractive for attempting to find components that resemble physical signals in the case where the positivity constraints are expected to hold. In crystallography, NMF has demonstrated significant potential in finding physically plausible structural signals from diffraction data collected from in situ chemical reactions [longRapidIdentificationStructural2009a, kusneHighthroughputDeterminationStructural2015a, huaNonequilibriumMetalOxides2021]. Recently, NMF has also been used for in situ time-dependent diffraction measurements [liuValidationNonnegativeMatrix2021a, thatc;aca22] and spatially resolved electron diffraction maps [rakit;am23], single-layer nanosheets [beauvaisResolvingSinglelayerNanosheets2021], integrated multimodal analysis [onolanMultimodalAnalyticalToolkit2021], and metal–organic frameworks [chenNodeDistortionTunable2023, rayderUnveilingUnexpectedModulatorCO22023].

The conventional NMF method assumes that the components remain fixed with respect to time, and therefore can hardly capture changes in the components over time. For example, in temperature series experiments, increased temperature can expand the inter-atomic distances, resulting in the stretching of peak positions in the measured powder diffraction pattern or atomic pair distribution function (PDF) data.

To address this limitation, extended NMF models have been proposed. One such model is the Shifted NMF, which accounts for shifts in the onset of a frequency profile, which can be induced by the Doppler effect for spectrometry data [morupShiftedNonNegativeMatrix2007]. However, Shifted NMF is not able to solve temperature series data problem because the change in the component is a stretch, not a shift. Another approach is to incorporate stretching regression steps into the analysis workflow [rakitaActiveReactionControl2020c].

In this paper, we propose a new extended NMF model called stretchedNMF, to explore a more fundamental aspect of the algorithm itself. We introduce a stretching factor matrix to describe the stretching scales of each component and each component is allowed to have different entire stretching factors at different moments. stretchedNMF can be developed to account for a simple stretching of the measured signal and returns only components that explain variability beyond this stretching.

In this paper, we first develop the mathematical formulas of stretchedNMF in the form of functional optimization. We present the method of discretization and the optimization algorithm. Then using both simulated and real data, we show that stretchedNMF significantly outperforms conventional NMF in the case of diffraction data with thermal expansion. Furthermore, we show that the algorithm may be used to extract different chemical components from the data if there are multiple components that have differential thermal expansivities. This gives an interesting possibility for extracting the components in a multi-phase sample from a temperature dependent measurement of that sample, even when those components are not changing chemically during the measurement. Although we focused on diffraction signals from temperature series data, the algorithm may be used for any case where part of the changes to the signal are exactly, or approximately, a stretch of its dependent variable.

2 Stretched Nonnegative matrix factorization

Nonnegative matrix factorization (NMF) is a mathematical tool to approximate a given matrix $Z\in\mathbb{R}^{N\times M}$ by the product of two low-rank nonnegative matrices,

Z\approx XY,

(1)

where $X\in\mathbb{R}^{N\times K}$ and $Y\in\mathbb{R}^{K\times M}$ , and $K\ll N,M$ [leeLearningPartsObjects1999d]. Its description and use is described in detail in multiple places [berryAlgorithmsApplicationsApproximate2007b, wangNonnegativeMatrixFactorization2013a]. The common NMF model uses the square of Euclidian distance (SED) as the objective function, and the corresponding optimization problem is written as

	$\displaystyle\min\limits_{X\in\mathbb{R}^{N\times K},Y\in\mathbb{R}^{K\times M}}$	$\displaystyle\quad\frac{1}{2}\left\\|XY-Z\right\\|_{F}^{2},$
	s.t.	$\displaystyle\quad X\geq 0\text{ and }Y\geq 0.$		(2)

Similar to principal component analysis [abdiPrincipalComponentAnalysis2010a] , the NMF decomposition will find components that explain variability in the signals in the set of data. Unlike PCA, a constraint of positivity is applied to both the components and the weights. Since many real physical signals, and their weights, obey positivity, NMF is more likely to find components that resemble signals from different physical components contributing to a compound signal coming from multiple sources. As such, it is finding extensive use in scientific applications [renNonnegativeMatrixFactorization2018, liuValidationNonnegativeMatrix2021a, gobinetApplicationNonnegativeMatrix2004].

Here we address a situation where one aspect of the variability, a stretching of the signal on the axis of its independent variable, is not of scientific interest, for example, due to the thermal expansion of a material affecting its diffraction pattern. We formulate an approach named stretchedNMF which extends the conventional NMF decomposition whilst accounting for the stretching in the algorithm.

Suppose the experimental signals, which are columns of $Z$ , $z^{m}$ for $m=1\dots M$ , and the components, which are columns in $X$ , $x_{k}$ for $k=1\dots K$ , are continuous functions of an independent variable $r$ . Then the conventional NMF optimization problem may be written as

\min\limits_{y_{k}^{m}\geq 0,x_{k}\geq 0}\quad\sum\limits_{m=1}^{M}\left\|\sum\limits_{k=1}^{K}y_{k}^{m}x_{k}(r)-z^{m}(r)\right\|_{L_{2}}^{2},

(3)

where $y_{k}^{m}$ is the weight of the $k$ th component at the $m$ th position in the dataset. Now, we assume that there is an $m$ -dependent stretching of the signal along the $r$ axis. The component signals stretch with component-dependent rates that we capture in a stretching factor, $\{a_{k}^{m}\}_{m=1,\dots,M}$ . We add the stretching factors $a_{k}^{m}$ into Eq. 3 and the optimization problem becomes

\min\limits_{a_{k}^{m},y_{k}^{m}\geq 0,x_{k}\geq 0}\quad\sum\limits_{m=1}^{M}\left\|\sum\limits_{k=1}^{K}y_{k}^{m}x_{k}(r/a_{k}^{m})-z^{m}(r)\right\|_{L_{2}}^{2}.

(4)

Notice that if $a_{k}^{m}>1$ , $x_{k}$ is stretched, and if $a_{k}^{m}<1$ , $x_{k}$ is compressed. In practice, we consider a finite $r$ range $[0,r_{\max}]$ . Therefore, without loss of generality, we define $x_{k}(r)=0$ for $r\geq r_{\max}$ . Thus, when $a_{k}^{m}>1$ , $x_{k}(r/a_{k}^{m})=0$ for $r\geq r_{\max}/a_{k}^{m}$ . Now we are able to expand the $L_{2}$ norm in equation (4) as an integral over the $r$ range as

	$\displaystyle\min\limits_{a_{k}^{m},y_{k}^{m}\geq 0,x_{k}\geq 0}$	$\displaystyle\quad\sum\limits_{m=1}^{M}\int_{0}^{r_{\max}}\left(\sum\limits_{k=1}^{K}y_{k}^{m}x_{k}(r/a_{k}^{m})-z^{m}(r)\right)^{2}dr,$
	s.t.	$\displaystyle\quad x_{k}(r)=0\text{, if }r\geq r_{\max}.$		(5)

For fixed component $k$ , $\{a_{k}^{m}\}_{m=1,\cdots,M}$ is a series of stretching factors, which usually change smoothly with time $m$ . However, the optimization problem in equation (2) is non-convex, and hence the smoothness of $\{a_{k}^{m}\}_{m=1,\cdots,M}$ may be violated when we solve it numerically. Therefore, we add a regularization term to the objective function to make it favor smooth $a_{k}$ , i.e.,

$\displaystyle\min\limits_{a_{k}^{m},y_{k}^{m}\geq 0,x_{k}\geq 0}$	$\displaystyle\quad\sum\limits_{m=1}^{M}\int_{0}^{r_{\max}}\left(\sum\limits_{k=1}^{K}y_{k}^{m}x_{k}(r/a_{k}^{m})-z^{m}(r)\right)^{2}dr$
	$\displaystyle\hskip 85.35826pt+\rho\sum\limits_{k=1}^{K}\sum\limits_{m=1}^{M-2}(a_{k}^{m}-2a_{k}^{m+1}+a_{k}^{m+2})^{2},$
s.t.	$\displaystyle\quad x_{k}(r)=0\text{, if }r\geq r_{\max},$	(6)

where $\sum\limits_{k=1}^{K}\sum\limits_{m=1}^{M-2}(a_{k}^{m}-2a_{k}^{m+1}+a_{k}^{m+2})^{2}$ is the smoothness regularization and $\rho$ is the parameter to control the effect of regularization. In our numerical testing section, we initiate a large $\rho$ and gradually decrease it in subsequent iterations.

3 Numerical Solution of stretchedNMF

In this section, we describe the numerical implementation of the stretchedNMF.

In order to numerically solve the functional optimization problem (2), we discretize the functionals and solve the corresponding vector optimization problem. Unlike Shifted NMF [morupShiftedNonNegativeMatrix2007], we cannot get benefits from discretizing the frequency domain of the components after applying the Fourier transform. So we choose to discretize the problem in the real $r$ space, without loss of generality, using a uniform grid on $[0,r_{\max}]$ . Since we have introduced the stretching factors, when we discretize the functionals $x_{k}(r/a_{k}^{m})$ , on this uniform grid the arguments $r/a_{k}^{m}$ are actually not on the grid nodes. Therefore, we apply a spline interpolation, that is we approximate $x_{k}(r/a_{k}^{m})$ from $x_{k}(r)$ , where the interpolant is a piecewise polynomial. In terms of the order of the spline, we need at least a quadratic order, i.e., a piecewise quadratic polynomial with continuous derivatives on the grid points. The smoothness of the spline will help the convergence of the discretized optimization solution. In this paper, we use a quadratic spline interpolation to approximate $x_{k}(r/a_{k}^{m})$ in the optimization problem (2). Explicitly, let $0=r_{0}<r_{1}<\dots<r_{n}=r_{max}$ be the uniform grid nodes, resulting in an interval of $h=r_{max}/n$ . The quadratic piecewise polynomial approximation, $S_{i}(r)$ , of $x(r)$ for $r\in[r_{i},r_{i+1}]$ is

\displaystyle S_{i}(r)=q_{i}(r-r_{i})(r-r_{i+1})+\left[x(r_{i+1})-x(r_{i})\right](r-r_{i})/h+x(r_{i}),

(7)

where $q_{i}$ is the quadratic coefficient to be determined. The derivatives of the polynomials $S_{i}(r)$ and $S_{i+1}(r)$ are

	$\displaystyle S_{i}^{\prime}(r)=q_{i}(2r-r_{i}-r_{i+1})+\left[x(r_{i+1})-x(r_{i})\right]/h,$		(8)
	$\displaystyle S_{i+1}^{\prime}(r)=q_{i+1}(2r-r_{i+1}-r_{i+2})+\left[x(r_{i+2})-x(r_{i+1})\right]/h.$		(9)

Notice the fact that the second-order spline should have continuous derivatives over the entire domain, which means that $S_{i}^{\prime}(r_{i+1})=S_{i+1}^{\prime}(r_{i+1})$ at positions $r_{i+1}$ for $i=0,\dots,n-2$ , using (8) and (9), we get

\displaystyle q_{i}+q_{i+1}=[x(r_{i})-2x(r_{i+1})+x(r_{i+2})]/h^{2}.

(10)

Since we have $x(r)=0$ , for $r\geq r_{max}$ , we set $S_{n-1}(r_{n})=0$ and $S^{\prime}_{n-1}(r_{n})=0$ . Then we can write $q$ as

\left(\begin{array}[]{c}q_{0}\\ q_{1}\\ \vdots\\ q_{n-1}\\ \end{array}\right)=\frac{1}{h^{2}}\left(\begin{array}[]{lllll}1&1&&&\\ &1&1&&\\ &&\ddots&\ddots&\\ &&&1&1\\ &&&&1\\ \end{array}\right)^{-1}\left(\begin{array}[]{rrrrr}1&-2&1&&\\ &1&-2&1&\\ &&\ddots&\ddots&\ddots\\ &&&1&-2\\ &&&&1\\ \end{array}\right)\left(\begin{array}[]{c}x(r_{0})\\ x(r_{1})\\ \vdots\\ x(r_{n-1})\\ \end{array}\right).

(11)

Now we can write $x_{k}(r/a_{k}^{m})$ in terms of $x_{k}(r_{i})$ as a linear transformation

x_{k}(r/a_{k}^{m})=q_{i}(r/a_{k}^{m}-r_{i})(r/a_{k}^{m}-r_{i+1})+[x(r_{i+1})-x(r_{i})](r/a_{k}^{m}-r_{i})/h+x(r_{i}),

(12)

if $r/a_{k}^{m}\in[r_{i},r_{i+1}]$ and $x_{k}(r/a_{k}^{m})$ is set to zero if $r/a_{k}^{m}\geq r_{max}$ . Since the leading coefficient $q$ is also linearly dependent on $x$ as shown in equation (11), we can conclude the linear transformation $x_{k}(r_{i}/a_{k}^{m})=b_{i,a_{k}^{m}}^{T}x_{k}$ for $i=0,1,\dots,n$ in a matrix form

\left(\begin{array}[]{c}x_{k}(r_{0}/a_{k}^{m})\\ x_{k}(r_{1}/a_{k}^{m})\\ \vdots\\ x_{k}(r_{n}/a_{k}^{m})\\ \end{array}\right)=\left(\begin{array}[]{c}\makebox[20.00003pt]{$\smash{-}\mkern-7.0mu\cleaders\hbox{$\mkern-2.0mu\smash{-}\mkern-2.0mu$}\hfill\mkern-7.0mu\smash{-}$}\hskip 5.69054ptb_{0,a_{k}^{m}}^{T}\hskip 5.69054pt\makebox[20.00003pt]{$\smash{-}\mkern-7.0mu\cleaders\hbox{$\mkern-2.0mu\smash{-}\mkern-2.0mu$}\hfill\mkern-7.0mu\smash{-}$}\\ \makebox[20.00003pt]{$\smash{-}\mkern-7.0mu\cleaders\hbox{$\mkern-2.0mu\smash{-}\mkern-2.0mu$}\hfill\mkern-7.0mu\smash{-}$}\hskip 5.69054ptb_{1,a_{k}^{m}}^{T}\hskip 5.69054pt\makebox[20.00003pt]{$\smash{-}\mkern-7.0mu\cleaders\hbox{$\mkern-2.0mu\smash{-}\mkern-2.0mu$}\hfill\mkern-7.0mu\smash{-}$}\\ \vdots\\ \makebox[20.00003pt]{$\smash{-}\mkern-7.0mu\cleaders\hbox{$\mkern-2.0mu\smash{-}\mkern-2.0mu$}\hfill\mkern-7.0mu\smash{-}$}\hskip 5.69054ptb_{n,a_{k}^{m}}^{T}\hskip 5.69054pt\makebox[20.00003pt]{$\smash{-}\mkern-7.0mu\cleaders\hbox{$\mkern-2.0mu\smash{-}\mkern-2.0mu$}\hfill\mkern-7.0mu\smash{-}$}\\ \end{array}\right)\left(\begin{array}[]{c}x_{k}(r_{0})\\ x_{k}(r_{1})\\ \vdots\\ x_{k}(r_{n})\\ \end{array}\right),

(13)

and denote the coefficient matrix as $B_{a_{k}^{m}}$ . Now we are ready to write the discretization of the optimization problem in (2) as

\min\limits_{a_{k}^{m},\,y_{k}^{m},\,x_{k}\geq 0}\sum\limits_{m=1}^{M}\left\|\sum\limits_{k=1}^{K}y_{k}^{m}B_{a_{k}^{m}}x_{k}-z^{m}\right\|^{2}+\rho\sum\limits_{k=1}^{K}\sum\limits_{m=1}^{M-2}(a_{k}^{m}-2a_{k}^{m+1}+a_{k}^{m+2})^{2},

(14)

where $y_{k}^{m}$ , $x_{k}$ and $z^{m}$ are discretized functionals on the uniform grid $0=r_{0}<r_{1}<\cdots<r_{n}=r_{\max}$ .

If the theoretical convergence is neglected, linear interpolation may be used as an approximation. In this case, we set $q_{i}=0$ in (12). The final form of the optimization problem is still (14), but with a different $B_{a_{k}^{m}}$ with higher sparsity.

Among the existing methods, a popular approach to solve the conventional NMF is alternating non-negative least squares (ANLS) [paateroPositiveMatrixFactorization1994, linProjectedGradientMethods2007, kimNonnegativeMatrixFactorization2008, guanNeNMFOptimalGradient2012, huangQuadraticRegularizationProjected2015]. ANLS alternatively adjusts $X$ and $Y$ to minimize the objective function and each subproblem can be solved by the non-negative linear least square method. In fact, this framework is also called the block coordinate descent (BCD) method with two blocks. In our problem Eq. 14, which can be simplified as

\min\limits_{A,Y,X\geq 0}f(A,Y,X),

(15)

there are three blocks $A$ , $Y$ and $X$ . Applying the BCD method with three blocks, we can solve the problem Eq. 15 using algorithm 2.

Algorithm 2 Block Coordinate Descent (BCD) Method

1:for

t=1,2,\cdots

A:=\arg\min_{A\geq 0}\langle\nabla_{A}f(\hat{A},Y,X),A\rangle

Y:=\arg\min_{Y\geq 0}\langle\nabla_{Y}f(A,\hat{Y},X),Y\rangle

X:=\arg\min_{X\geq 0}\langle\nabla_{X}f(A,Y,\hat{X}),X\rangle

5:end for

Similar to conventional NMF, the subproblems of $Y$ and $X$ in Lines 3 and 4 are convex quadratic programming problems that can be easily solved by existing solvers. But the subproblem of $A$ in Line 2 is highly non-convex and therefore we cannot efficiently solve it for its global minimum. In practice, we use a subspace trust-region method [colemanInteriorTrustRegion1996] to find a local minimum.

The convergence of the BCD method for 3 blocks is not guaranteed [grippoConvergenceBlockNonlinear2000]. Here we use an algorithm that can guarantee its convergence for a quadratic spline approximation that is called the linearized block coordinate descent method [xuGloballyConvergentAlgorithm2017]. The outline of the framework is presented in Algorithm 4, where $\alpha_{t}$ is the step size and $\hat{A}/\hat{X}/\hat{Y}$ are the extrapolations of the current $A/X/Y$ and previous $A/X/Y$ . In each iteration, the algorithm randomly chooses one block and minimizes the corresponding linear approximation and a proximal term. One can refer to [xuGloballyConvergentAlgorithm2017] for more information about parameter selections.

Algorithm 4 linearized block coordinate descent method

1:for

t=1,2,\cdots

2: pick one of the following to implement in a deterministic or random manner;

A:=\arg\min_{A\geq 0}\langle\nabla_{A}f(\hat{A},Y,X),A\rangle+\frac{1}{\alpha_{t}}\|A-\hat{A}\|^{2}

Y:=\arg\min_{Y\geq 0}\langle\nabla_{Y}f(A,\hat{Y},X),Y\rangle+\frac{1}{\alpha_{t}}\|Y-\hat{Y}\|^{2}

X:=\arg\min_{X\geq 0}\langle\nabla_{X}f(A,Y,\hat{X}),X\rangle+\frac{1}{\alpha_{t}}\|X-\hat{X}\|^{2}

6:end for

4 Diffraction Use Case

4.1 Introduction

Here we test the approach using simulated and also real x-ray powder diffraction (PXRD) data [pecharskyFundamentalsPowderDiffraction2008, dinnebierPowderDiffractionTheory2008a], and atomic pair distribution function (PDF) [egamiBraggPeaksStructural2012d] data. PXRD and PDF patterns are continuous 1D signals that encode the 3D arrangement of atoms in a material. We assume a situation where the PDF and PXRD patterns have been measured for samples as a function of temperature and are undergoing thermal expansion, where the thermal expansion coefficient of each sample is different. The thermal expansion causes Bragg peaks in the PXRD, and peaks in the PDF, to change their positions. In principle, thermal expansion can be different along different directions of the crystal, but often it is quite isotropic and appears as a stretching of the pattern where the peak shifts increase with increasing distance along the independent variable axis as required for this algorithm to work. This makes it an interesting use-case for stretchedNMF, though we note that the stretchedNMF may be applied to any series of signals where one aspect of the variability is a continuous stretching on the axis of the independent variable.

The goal of our testing use-case is to see if we can use NMF in general, but stretchedNMF in particular, to separate the chemical components in a binary chemical mixture where the two components have different thermal expansion coefficients. For example, this could be used by a material scientist to discover the chemical components in a synthesis product by measuring the mixture as a function of temperature and running stretchedNMF on the mixture, where the algorithm returns mathematical components that resemble the PXRD or PDF signals of the actual chemical components. These mathematical components could then be given to algorithms such as the structureMining [yangStructureminingScreeningStructure2020g] or spacegroupMining [liuUsingMachineLearning2019] algorithms that are implemented as a service on the PDFitc.com website [yangCloudPlatformAtomic2021c]. These algorithms, given an uploaded PDF, will return a rank ordered list of candidate structures consistent with that PDF. Our test of the algorithm will therefore consist of taking either simulated or actual measured data over a wide temperature range from binary mixtures where the components have different thermal expansion coefficients. These signals will be fed to NMF and stretchedNMF to extract two components which will then be analyzed to see if they resemble the signals from the actual chemical components. In the case of the PDF, an interesting test of this is to take the extracted mathematical components and giving them to the structureMining algorithm to see if it correctly identifies the chemical component from the stretchedNMF and conventional NMF extraction.

4.2 Data

To evaluate the performance of the stretchedNMF and sparse-stretchedNMF algorithms, we tested them on the following test-case datasets:

4.2.1 Simulated PXRD and PDF data with increasing lattice parameters

We used a set of 20 simulated PXRD and PDF patterns. The PXRD and PDF were from a weighted sum of a simulated cubic perovskite $\mathrm{BaTiO_{3}}$ phase and a cubic wurtzite $\mathrm{ZnSe}$ . The $\mathrm{BaTiO_{3}}$ : $\mathrm{ZnSe}$ phase-fraction was set to 1:1, which corresponded to an atomic concentration ratio 0.61:1.00 when initial lattice parameters of $\mathrm{BaTiO_{3}}$ =4.18Å and $\mathrm{ZnSe}$ =5.62Å were used. The CIF files used for $\mathrm{BaTiO_{3}}$ and $\mathrm{ZnSe}$ were from structures reported in [keler_reaction_1960] and [andreev_synthesis_1995], respectively, and downloaded from the Springer Materials database (https://materials.springer.com/isp/crystallographic/docs/sd\_0304044, https://materials.springer.com/isp/crystallographic/docs/sd\_1929775).

To simulate lattice expansion we assumed a constant expansion coefficient, $\alpha_{BTO}$ and $\alpha_{ZS}$ , for each component. For a more complete exploration of the performance of conventional NMF with that of stretchedNMF and sparse-stretchedNMF, we did not restrict ourselves to expansion coefficients that resemble actual thermal expansivities. However, since the linear thermal expansion coefficient of $\mathrm{BaTiO_{3}}$ is approximately twice of that of $\mathrm{ZnSe}$ [bland_thermal_1959, su_thermal_2009], we always set the expansion coefficient of $\mathrm{BaTiO_{3}}$ to be twice that of $\mathrm{ZnSe}$ .

We note that to test purely the effects of stretching, which is the basis of the current NMF modification, we fixed and did not vary the atomic displacement parameters (ADPs) that would result in changes in the attenuation of the PXRD Bragg peaks and broadening of peaks in the simulated PDF. Such effects are likely to be present in real data and may require a further modification to the NMF algorithm in the future but this is beyond the scope of the current paper. This set of simulations assumed no phase transition or chemical reaction to be occurring and the relative weights of the components were not varied in the computed dataset.

The PXRD patterns were simulated using Dans-Diffraction [porter_danporterdans_diffraction_2020]. Pseudo-Voight lineshapes were used. The PDFs were simulated using Diffpy-CMI [juhas_complex_2015]. The code used to generate the PDFs is reproduced in the supplementary information and can be found at https://github.com/yevgenyr/diffpysim. The static set of parameters used for the simulations is reproduced in Table 1.

Table 1: The static set of parameters that were used for PXRD and PDF simulations.

	Parameter	Value	Units
PDF	steps	20
	qmin	0.1	$\mathrm{\AA^{-1}}$
	qmax	30	$\mathrm{\AA^{-1}}$
	qdamp	0.03
	Uiso	0.007	$\mathrm{\AA^{2}}$
	rmin	0
	rmax	120	$\mathrm{\AA}$
	rstep	0.01	$\mathrm{\AA}$
PXRD	voigt_profile.sigma	1.5
	voigt_profile.gamma	1.5
	qmax	30	$\mathrm{\AA^{-1}}$
	Uiso	0.007	$\mathrm{\AA^{2}}$

Representative PXRD and PDF patterns are shown in Fig. 1.

Refer to caption — Figure 1: Example simulated signals used in the tests. These were from a linear combination of the PDFs of \chBaTiO3 and \chZnSe where each was computed with a linearly expanding cubic lattice parameter. The PDFs are plotted offset from one another expanding from the bottom to the top of the figure. In the curves shown an overall expansion values of 20% and 10% was used from the first to the last curve. Top panel shows the simulated XRD and the bottom panel shows the PDF.

4.2.2 Experimental PXRD data - thermal expansion

To test the stretchedNMF and sparse-stretchedNMF algorithms on real data, we use part of an in situ solid-state synthesis reaction dataset where no phase-transition or chemical reaction occurred but which spanned a rather broad temperature range. This allows us to evaluate how the algorithms perform for the effect of thermal expansion of a phase mixture from real data. The PXRD experiment was done at the 28-ID-2 beamline (XPD instrument) at the NSLS-II facility at Brookhaven National Laboratory. A large area 2D Perkin Elmer detector was used to acquire the data. To gain high spectral resolution in the PXRD, the distance between the sample and the detector was set to 144 cm. The beam wavelength was 0.1949 Å.

A stoichiometric mixture of 2:1 $\mathrm{YOCl}$ ( $>$ 98% tetragonal phase) and $\mathrm{MgMn_{2}O_{4}}$ (spinel phase) was uniformly mixed and sealed in a quartz capillary. It was then heated in a gradient furnance, meaning that each location on the quartz tube had a different temperature.[onolan_thermal-gradient_2020] The absolute temperatures at each point along the sample were calibrated from the lattice expansion of a known calibration material, Ni. The data went from a low temperature of 368 $\mathrm{{}^{o}C}$ to a highest temperature of 668 $\mathrm{{}^{o}C}$ with a total of 20 individual temperature points. Using ‘pyFAI’ [ashiotis_fast_2015], the collected 2D diffraction patterns were then cleaned by masking the beam-stop and over-bright/dead pixels, followed by an azimuthal integration to gain 1D PXRD patterns. The 1D PXRD data was then used as inputs to the different NMF algorithms.

4.2.3 Experimental PXRD data - thermal expansion and reaction

We also tested the NMF algorithms on another PXRD dataset, but this time, where a solid-state chemical reaction happened together with the thermal expansion so that the weights of the components as well as the thermal expansion were varying during the experiment. The data were measured as the temperature changed from 28 C to 370 C in 215 steps during the reaction of

\ch{CuCl2}+\ch{Na2Se2}\longrightarrow\ch{CuSe2}+2\ch{NaCl}.

Here the components involved in the reaction are NaCl, CuSe, \chCu2Se, Se, pyrite, and marcasite, as determined by a multi-phase Rietveld refinement on the full dataset carried out previously. The full details of experiment are published in [martinolichPolymorphSelectivitySuperconducting2015].

4.3 Algorithm Developments

In the case of PDF data, we apply stretchedNMF to time-series data according to the workflow shown in the chart in Fig. 2.

A common experimental function (for example, the output of xPDFsuite [yang;arxiv14] and PDFgetX3 [juhas2013pdfgetx3], is the $G(r)$ function [farrowRelationshipAtomicPair2009a]. This function oscillates above and below zero. NMF works on the basis that signals are positive and in order to avoid the loss of signal where the function goes negative, we need to modify the signal into a non-negative form. Here we use an offset method, by taking the smallest of all data values and adding its absolute value to all data. This approach has the advantage of being simple and has been successfully applied to the deep learning method of PDF analysis [liuValidationNonnegativeMatrix2021a].

After running the NMF solvers, we must restore the components to valid $G(r)$ functions (oscillating around zero). To do this we utilize the solved weights and stretching factors to recover the components according to

\min\limits_{x_{k}}\sum\limits_{m=1}^{M}\left\|\sum\limits_{k=1}^{K}y_{k}^{m}B_{a_{k}^{m}}x_{k}-z^{m}\right\|^{2},

(16)

where, $z^{m}$ is the original data rather than the data after the offset pre-processing and the other symbols are described alongsided Eq. 14. The weight, $y$ , and stretching factors, $a$ , are fixed to be those obtained from the NMF solution, and we remove the constraint that the components must be non-negative. Functions resembling $G(r)$ are then recovered from the NMF components and may be fit using standard PDF modeling protocols. This is reasonable because it is based on our trust in the weights and stretching factors of the NMF solver’s solution of the preprocessed data. This approach is highly automated and can be applied to both conventional NMF and stretchedNMF, because the stretching factor of the conventional NMF is always $1$ .

For the case of PXRD data from highly crystalline samples, we have the additional observation that the spectrum consists of a sparse set of sharp peaks. That is, the function value is zero in between the Bragg peaks (neglecting backgrounds and any diffuse scattering). We can make use of this property to enhance our ability to decompose signals by adding a sparse regularization term to the optimization problem. For the case where there are smooth backgrounds in experimental PXRD data, the background can be easily and automatically eliminated to make the data sparse. In this case we make two modifications to the optimization problem in Eq. (14). The first is adding the $l_{1/2}$ sparse regularization term to $x$ [xuL1Regularization2010]. The second is adding an upper bound on $y$ , in order to prevent $x$ from collapsing to zero as a whole, resulting in

	$\displaystyle\min\limits_{a_{k}^{m}>0,\,0\leq y_{k}^{m}\leq 1,\,x_{k}\geq 0}$	$\displaystyle\sum\limits_{m=1}^{M}\left\\|\sum\limits_{k=1}^{K}y_{k}^{m}B_{a_{k}^{m}}x_{k}-z^{m}\right\\|^{2}+$		(17)
		$\displaystyle\rho\sum\limits_{k=1}^{K}\sum\limits_{m=1}^{M-2}(a_{k}^{m}-2a_{k}^{m+1}+a_{k}^{m+2})^{2}+\eta\sum\limits_{k=1}^{K}\sum\limits_{i=1}^{n}(x_{k,i})^{\frac{1}{2}}.$

We refer to this as sparse-stretchedNMF.

In this optimization model, there are two regularization parameters, $\rho$ and $\eta$ . From experience, the smoothness parameter $\rho$ is not sensitive and is usually adjusted by multiplying by 10. The sparsity parameter $\eta$ can be adjusted by doubling.

5 Numerical Results

5.1 Results on Simulated PDF

First, we compare the performance of the conventional NMF and the stretchedNMF on simulated PDF data. The PDFs were generated by a combination of two components, namely simulated \chBaTiO3 and ZnSe. The weight coefficients for each component were set as constants. We also assigned different linearly increasing rates for the thermal expansion of \chBaTiO3 and ZnSe. Specifically, we used artificially generated rates such that \chBaTiO3 and \chZnSe linearly expands from the first PDF to the last with 20% and 10% expansions, respectively.

We then applied the conventional NMF and stretchedNMF methods to extract two components from the simulated PDF data. These could then be compared with the ground-truth PDFs. In principal, any of the ground-truth PDFs could be picked as we apply a stretching factor to the NMF component signal before the comparison. In this study, we selected the ground-truth PDF that resulted in the minimal residual when using only the scale factor variable. We further optimized the agreement between the NMF component and the selected ground-truth PDF by varying both the scale-factor and stretch factor variables.

We first evaluated the outcomes of the conventional NMF approach. These findings are illustrated in Fig. 3(a-d) and Tab. 2.

Table 2: Results of the comparison between the NMF extracted components and the ground-truth PDFs on simulated PDF data test.

R_{w}

and PC are the residual and the Pearson Correlation values between the numerical solution and the ground truth.

Method	Component	Ground Truth	$R_{w}$ (PC)
Conventional NMF	1	\chBaTiO3	0.9488 (0.3157)
	2	\chBaTiO3	0.9050 (0.4255)
	1	ZnSe	0.9186 (0.3953)
	2	ZnSe	0.8462 (0.5328)
stretchedNMF	1	\chBaTiO3	0.1357 (0.9911)
	2	\chBaTiO3	0.9750 (0.2225)
	1	ZnSe	0.8600 (0.5104)
	2	ZnSe	0.1162 (0.9937)

Fig. 3(a-d) depicts the resulting PDFs in a matrix layout, with the NMF extracted components being represented as rows (in red) and the ground-truth PDFs as columns (in blue). The difference curves (ground-truth - NMF component) are plotted below in green. Large residuals and large $R_{w}$ factors are evident between all the NMF components and the ground-truth curves and the NMF extraction has failed to produce components that resemble the actual signals. This is not surprising since the weights of the two components are not varying in the test.

The same test was applied using the stretchedNMF algorithm and the results are shown in Fig. 3(e-h) and Tab. 2. In this case we can see that the stretchedNMF extracted signal I is closely related to the ground truth component I and likewise for the component II. This is evident as a very flat difference curve in Fig. 3(e) and (h) and small $R_{w}$ for these pairings in Tab. 2.

This shows that even in the absence of changes in component weights the stretchedNMF algorithm can extract components just from a differential stretching of the structure signal.

5.2 Results on Simulated PXRD

We carry out the same comparison of NMF vs stretchedNMF for the case of powder diffraction signals. Similar to the simulated PDF case, the data comprise of a combination of simulated \chBaTiO3 and ZnSe, where \chBaTiO3 and ZnSe have 20% and 10% linearly varying expansions, respectively.

The results of the comparison are presented in Fig. 4 and Tab. 3.

Table 3: Comparison between the NMF extracted components and the ground-truth PXRDs on the simulated PXRD data set.

R_{w}

and PC are the residual and the Pearson Correlation values, respectively, between the numerical solutions and the ground truth.

Method	Component	Ground Truth	$R_{w}$ (PC)
Conventional NMF	1	\chBaTiO3	0.8506 (0.5040)
	2	\chBaTiO3	0.7196 (0.6830)
	1	ZnSe	0.7976 (0.5979)
	2	ZnSe	0.5655 (0.8225)
stretchedNMF	1	\chBaTiO3	0.0437 (0.9990)
	2	\chBaTiO3	0.8427 (0.5290)
	1	ZnSe	0.8400 (0.5326)
	2	ZnSe	0.0272 (0.9996)

As is evident in Fig. 4(a-d), none of the extracted conventional NMF components resemble ground-truth curves. Again, this is not a surprise because the weights of the components are not changing. However, for the stretchedNMF extraction we see that the first extracted component (Comp I) corresponds well the \chBaTiO3 pattern (Truth I), and the second extracted component (Comp II) corresponds well to the \chZnSe diffraction pattern (Truth II) (Fig. 4(e) and (h)).

As with the simulated PDF data, the stretchedNMF algorithm can extract components resembling the physical signals from a phase mixture where the weights are not changing but there is a variable thermal expansion.

5.3 Results on simulated PDF and PXRD data with small expansion coefficients

The tests above show that even in the presence of large stretches of signals stretchedNMF can automatically extract signals that resemble real physical signals whereas conventional NMF cannot, at least in the case where the component weights are not changing.

We now would like to see how well stretchedNMF can perform for smaller stretching factors, for example, for magnitudes that might occur in physical systems due to thermal expansion. The simulated data is still taken as the combination of \chBaTiO3 and \chZnSe. The weights are set to constants as before. However, in this example we set the thermal expansion rates of \chBaTiO3 and \chZnSe to 4% and 2%, respectively. Both simulated PDF and PXRD are tested.

First, we compare the performance of the conventional NMF and the stretchedNMF on simulated PDF data. The results are presented in Tab. 4 and Fig. 5.

Table 4: Results of the comparison between the NMF extracted components and the ground-truth PDFs on simulated PDF data sets with 2% expansions on ZnSe.

R_{w}

and PC are the residual and the Pearson Correlation values between the numerical solution and the ground truth.

Method	Componenpt	Ground Truth	$R_{w}$ (PC)
Conventional NMF	1	\chBaTiO3	0.7922 (0.6103)
	2	\chBaTiO3	0.8268 (0.5625)
	1	ZnSe	0.7607 (0.6491)
	2	ZnSe	0.7861 (0.6180)
stretchedNMF	1	\chBaTiO3	0.0960 (0.9961)
	2	\chBaTiO3	0.9732 (0.2299)
	1	ZnSe	0.9725 (0.2333)
	2	ZnSe	0.1237 (0.9935)

Unlike the previous figures we just plot the agreement of the extracted component and the ground-truth curve that shows the best agreement. The poor performance of the conventional NMF is evident in Fig. 5(a) and (b), whereas again, even for this much smaller stretch, the stretchedNMF algorithm still gives a good extraction of the physical components (Fig. 5(c) and (d)).

We get the same overall result for the test on simulated PXRD data as for the PDF data. The results are shown in Fig. 6 and Tab. 5. Again, stretchedNMF gives a very good extraction of the physical components even for this small relative expansion coefficient (Fig. 6(c) and (d)) whereas conventional NMF does not (Fig. 6(a) and (b))

Table 5: Results of the comparison between the NMF extracted components and the ground-truth PXRDs on simulated PXRD data test with 2% expansions on ZnSe.

R_{w}

and PC are the residual and the Pearson Correlation values between the numerical solution and the ground truth.

Method	Component	Ground Truth	$R_{w}$ (PC)
Conventional NMF	1	\chBaTiO3	0.8764 (0.4645)
	2	\chBaTiO3	0.8869 (0.4449)
	1	ZnSe	0.5165 (0.8550)
	2	ZnSe	0.5003 (0.8648)
stretchedNMF	1	\chBaTiO3	0.0376 (0.9993)
	2	\chBaTiO3	0.9896 (0.1271)
	1	ZnSe	0.9991 (0.0217)
	2	ZnSe	0.0305 (0.9995)

The results are less ideal when the expansion rates are reduced further to \chBaTiO3 and ZnSe changing linearly from 1 to 1.02 and 1.01, respectively. The results are summarized in Tab. 6, and Fig. 7.

Table 6: Results of the comparison between the NMF extracted components and the ground-truth PXRDs on simulated PXRD data test with a 1% differential expansion between the components.

R_{w}

and PC are the residual and the Pearson Correlation values between the numerical solution and the ground truth.

Method	Component	Ground Truth	$R_{w}$ (PC)
Conventional NMF	1	\chBaTiO3	0.8430 (0.5238)
	2	\chBaTiO3	0.8491 (0.5142)
	1	ZnSe	0.5248 (0.8493)
	2	ZnSe	0.5059 (0.8609)
stretchedNMF	1	\chBaTiO3	0.7584 (0.6413)
	2	\chBaTiO3	0.8815 (0.4578)
	1	ZnSe	0.6482 (0.7574)
	2	ZnSe	0.4490 (0.8921)
sparse-stretchedNMF	1	\chBaTiO3	0.0765 (0.9971)
	2	\chBaTiO3	0.9902 (0.1243)
	1	ZnSe	0.9951 (0.0823)
	2	ZnSe	0.0556 (0.9985)

At this level of expansion, even the stretchedNMF is not correctly extracting the physical components. For example, it incorrectly assigns peaks in the spectrum of its extracted components in red at around $Q=1.5$ , 2 and 2.5 Å^-1 (Fig. 7(c) and (d)). These same peaks are partially misassigned by the conventional NMF algorithm.

However, the sparse-stretchedNMF algorithm does a good job of extracting physical components from the powder PXRD simulations (Fig. 7(e) and (f)) even in this challenging case with a relatively small (1%) differential expansion. The components of sparse-stretchedNMF are close to ground truths. This indicates that sparse-stretchedNMF can enhance the performance of stretchedNMF.

These tests show that the stretchedNMF algorithm is able to extract physically meaningful PDF and PXRD signals from sets of data where the signals are unchanged except for a different relative stretch between the two curves. If there is a large differential change in lattice parameter across the dataset stretchedNMF can still extract ground-truth PDF and PXRD signals. For relative stretches of a few percent, comparable to what might be expected for a mixture of compounds with a differential thermal expansion, this is also true for both PDF and PXRD data. When the differential thermal expansion gets to around 1%, stretchedNMF starts to struggle to extract physical components. However, for PXRD data the sparse-stretchedNMF algorithm still performs well. We note that the PDF data is not sparse, and therefore sparse-stretchedNMF algorithm is applied only on PXRD data.

We should note that in these ground-truth tests on simulated data we wanted to test how well stretchedNMF can handle datasets that contain stretches, for example, as might come from thermal expansion. We therefore did not include in the simulation other effects of temperature changes such as increases in atomic displacement factors (ADPs). In principle, we would like to develop a new algorithm that can eliminate changes in ADP in the same way as stretchedNMF eliminates stretches. This problem will be left for a future paper. Preliminary tests on simulated data with combined stretching and increased-ADP effects indicate that stretchedNMF and sparse-stretchedNMF still perform reasonably well and clearly outperform the conventional NMF algorithm, but with larger errors than in the constant-ADP tests reported here. Despite this known shortcoming, we would still like to see whether stretchedNMF and sparse-stretchedNMF can perform well on experimental data from a variable temperature experiment, and this is discussed in the following section.

5.4 Results on measured PXRD data I

Here, we test the NMF algorithms on measured PXRD data. The data are from the in situ chemical reaction experiment described in Section 4.2.2.

Multi-component Rietveld refinements were carried out and indicate that the chemical components in this reaction are \chMgMn2O4, orthorhombic \chYMnO3, and rhombohedral and tetragonal \chYOCl (\chrYOCl and \chtYOCl, respectively) where \chMgMn2O4 and \chtYOCl are the dominant phases. The results of the Rietveld refinements for the two majority phases were used as ground truth against which to compare the performance of the NMF algorithms.

The results are shown in Fig. 8 and the resulting $R_{w}$ and PC are listed in Tab. 7.

Table 7: Results of the comparison between the NMF extracted components and the ground truth from Rietveld refinement on real PXRD data test.

R_{w}

and PC are the residual and the Pearson Correlation values between the numerical solution and the ground truth.

Method	Component	Ground Truth	$R_{w}$ (PC)
Conventional NMF	1	\chMgMn2O4	0.9214(0.3402)
	2	\chMgMn2O4	0.9175(0.3519)
	1	\chtYOCl	0.4498(0.8863)
	2	\chtYOCl	0.4575(0.8819)
stretchedNMF	1	\chMgMn2O4	0.8020(0.5655)
	2	\chMgMn2O4	0.9831(0.1246)
	1	\chtYOCl	0.6433(0.7492)
	2	\chtYOCl	0.3493(0.9330)
sparse-stretchedNMF	1	\chMgMn2O4	0.4851(0.8682)
	2	\chMgMn2O4	0.9977(0.0137)
	1	\chtYOCl	0.9646(0.2254)
	2	\chtYOCl	0.3273(0.9422)

In Fig. 8 the blue curves in the top row (a, b, e) are from the diffraction pattern of \chMgMn2O4 and the blue curves in the bottom row (c, d, f) are from t\chYOCl. The red curves in each panel show the relevant extracted component from the NMF algorithm used. The columns are sorted by the NMF algorithm used. The first column (a, c) used regular NMF, the second (b,d) used the stretchedNMF algorithm, and the third column (e, f) used the sparse-stretchedNMF algorithm.

All NMF solvers give reasonable results for the \chtYOCl chemical component. The peak positions are consistent with the ground truth, and the inconsistency of intensity is acceptable. But for \chMgMn2O4, the NMF and stretchedNMF derived components are poor. They are much better using the sparse-stretchedNMF algorithm, which gives better agreement both visually and in terms of the $R_{w}$ between the ground-truth and the extracted components. For this case, from the perspective of separation ability, sparse-stretchedNMF is superior to stretchedNMF which is superior to the conventional NMF in this test.

The scaled weights from all NMF solvers are compared to the weights from Rietveld refinement which can be considered as ground-truth. The results are shown in Fig. 9.

The weights of the chemical components are not changing during the experiment and so we would expect the weights to be largely independent of temperature. The conventional NMF clearly does not return constant weights and is getting confused by the thermal expansion in the data. The stretchedNMF and sparse-stretchedNMF methods do yield almost constant weights. Rietveld refinements were carried out on these data-sets and can be treated as a ground-truth. The results of the Rietveld refinement are shown as the dashed curve. stretchedNMF is doing quite well, but sparse-stretchedNMF is doing very well in reproducing the results of the Rietveld refinement.

5.5 Results on real PXRD data II

We also tested the NMF algorithms on a real PXRD dataset from an in situ chemical reaction experiment, which was published in [martinolichPolymorphSelectivitySuperconducting2015]. The data were measured as the temperature changed from 28 K to 370 K in 215 steps during the reaction

\ch{CuCl2}+\ch{Na2Se2}\longrightarrow\ch{CuSe2}+2\ch{NaCl},

where chemical components that found to appear during reaction are NaCl, CuSe, \chCu2Se, Se, pyrite, and marcasite.

The top panel in Fig. 10 shows the measured PXRD data during the in situ reaction experiment.

The curves obtained by a multi-phase Rietveld refinement fit [martinolichPolymorphSelectivitySuperconducting2015] are shown in blue in the panels below. The Rietveld refined phase weights are shown in blue in the right hand column below [martinolichPolymorphSelectivitySuperconducting2015]. The components extracted from a stretchedNMF decomposition are shown in red, plotted on top of the ground-truth components, and the extracted weights are shown in red on top of the Rietveld extracted weights in the right hand column. The results are very good and indicate that, except for \chSe, the components obtained from sparse-stretchedNMF matched well with the ground truth, as do the extracted weights.

This shows that the stretchedNMF algorithm can be used as a rapid way to extract reliable components and weights from data collected at different temperatures. This approach can be very helpful looking at large amounts of data very rapidly as it is being collected to look for known phases and unknown phases without having to carry out a complex multicomponent Rietveld campaign in real time.

6 Discussion and Conclusion

This paper presents a novel functional optimization model called stretchedNMF, which is an extension to the traditional NMF model. The proposed model introduces a new variable, the stretching factor, that enables the components to undergo stretching transformations. Furthermore, a regularization term is incorporated to ensure the stretching factors are smooth over time or temperature. To solve the optimization problem, we discretize it and employ Block Coordinate Descent (BCD) framework algorithms. The initial experimental results indicate that for data where stretches in the signal are observed, such as diffraction data where thermal expansion has taken place, the proposed stretchedNMF model outperforms the conventional NMF. This is true even for PXRD and PDF data with small stretching degrees corresponding to realistic thermal expansivities. However, a further enhancement to stretchedNMF, which makes use of the sparsity of powder diffraction patterns, called sparse-stretchedNMF allows correct extractions even for very small stretches where stretchedNMF struggled.

Despite the utility described here, there are some limitations to the stretchedNMF model. One is where the stretching is anisotropic in a material. This would require a model-dependent correction to account for different stretches in different crystallographic directions.

We also note that experimental noise can affect the outcome. This has not been studied in detail in this paper, but we note that we obtained good results from real data that included noise. To further address the noise issue, different regularization techniques can be utilized.

We note that the current model only considers stretching, adding shift transforms it into a first-order polynomial transformation. In this case, only a new block is added to the computation, but a better approximation can be obtained. Incorporating higher-order polynomial transformations could further balance the computational and approximative accuracy of the model. Further research is needed to investigate and optimize the stretchedNMF model’s potential in overcoming these challenges.

Finally, we note that although the motivation for the development, and all the tests, were on diffraction data where underlying structures have undergone thermal expansion, the stretchedNMF algorithm will work on any signal decomposition that smooth continuous variations in a stretching fact as a characteristic of the signal and it is not limited to use on diffraction data.

\ack

Acknowledgment

We would like to thank Dr. Daniel Olds, for assistance during the measurements of the experimental PDF data. The work described here was funded by the Next Generation Synthesis Center (GENESIS), an Energy Frontier Research Center funded by the U.S. Department of Energy, Office of Science, Basic Energy Sciences under Award Number DE-SC0019212. X-ray PDF measurements were conducted on beamline 28-ID-2 of the National Synchrotron Light Source II, a US DOE Office of Science User Facility operated for the DOE Office of Science by Brookhaven National Laboratory under contract No. DESC0012704. Qiang Du is also partially supported by DOE-ASCRDE-SC0022317. GEK received training and support as a part of QuADS: Quantitative Analysis of Dynamic Structures National Science Foundation Research Traineeship Program, grant number NSF DGE 1922639.

References

[1] \harvarditemAbdi \harvardand Williams2010abdiPrincipalComponentAnalysis2010a Abdi, H. \harvardand Williams, L. J. \harvardyearleft2010\harvardyearright. WIREs Computational Statistics, \volbf2(4), 433–459.
[2] \harvarditem[Andreev et al.]Andreev, Bulanyi, Hayward \harvardand Mozharovsikii1995andreev_synthesis_1995 Andreev, A., Bulanyi, M., Hayward, S. \harvardand Mozharovsikii, L. \harvardyearleft1995\harvardyearright. (Russian) Journal of Inorganic Chemistry (translated from Zhurnal Neorganicheskoi Khimii), \volbf40, 1079–1082.
[3] \harvarditem[Ashiotis et al.]Ashiotis, Deschildre, Nawaz, Wright, Karkoulis, Picca \harvardand Kieffer2015ashiotis_fast_2015 Ashiotis, G., Deschildre, A., Nawaz, Z., Wright, J. P., Karkoulis, D., Picca, F. E. \harvardand Kieffer, J. \harvardyearleft2015\harvardyearright. Journal of Applied Crystallography, \volbf48(2), 510–519. Publisher: International Union of Crystallography.
\harvardurl//scripts.iucr.org/cgi-bin/paper?fv5028
[4] \harvarditem[Beauvais et al.]Beauvais, Chupas, O’Nolan, Parise \harvardand Chapman2021beauvaisResolvingSinglelayerNanosheets2021 Beauvais, M. L., Chupas, P. J., O’Nolan, D., Parise, J. B. \harvardand Chapman, K. W. \harvardyearleft2021\harvardyearright. ACS MATERIALS LETTERS, \volbf3(6), 698–703.
[5] \harvarditem[Berry et al.]Berry, Browne, Langville, Pauca \harvardand Plemmons2007berryAlgorithmsApplicationsApproximate2007b Berry, M. W., Browne, M., Langville, A. N., Pauca, V. P. \harvardand Plemmons, R. J. \harvardyearleft2007\harvardyearright. Computational Statistics & Data Analysis, \volbf52(1), 155–173.
[6] \harvarditemBland1959bland_thermal_1959 Bland, J. A. \harvardyearleft1959\harvardyearright. Canadian Journal of Physics, \volbf37(4), 417–421. Publisher: NRC Research Press.
\harvardurlhttps://cdnsciencepub.com/doi/10.1139/p59-046
[7] \harvarditemBuciu2008buciuNonnegativeMatrixFactorization2008a Buciu, I. \harvardyearleft2008\harvardyearright. International Journal of Computers, Communications & Control (IJCCC), \volbf3(3), 67–74.
[8] \harvarditem[Buciu et al.]Buciu, Nikolaidis \harvardand Pitas2008buciuNonnegativeMatrixFactorization2008 Buciu, I., Nikolaidis, N. \harvardand Pitas, I. \harvardyearleft2008\harvardyearright. IEEE Transactions on Neural Networks, \volbf19(6), 1090–1100.
[9] \harvarditem[Chen et al.]Chen, Stroscio, Liu, Lu, Hupp, Gagliardi \harvardand Chapman2023chenNodeDistortionTunable2023 Chen, Z., Stroscio, G. D., Liu, J., Lu, Z., Hupp, J. T., Gagliardi, L. \harvardand Chapman, K. W. \harvardyearleft2023\harvardyearright. Journal of the American Chemical Society, \volbf145(1), 268–276.
[10] \harvarditemCichocki \harvardand Phan2009cichockiFastLocalAlgorithms2009 Cichocki, A. \harvardand Phan, A.-H. \harvardyearleft2009\harvardyearright. IEICE Transactions on Fundamentals of Electronics, Communications and Computer Sciences, \volbfE92-A(3), 708–721.
[11] \harvarditemColeman \harvardand Li1996colemanInteriorTrustRegion1996 Coleman, T. F. \harvardand Li, Y. \harvardyearleft1996\harvardyearright. SIAM Journal on Optimization, \volbf6(2), 28.
[12] \harvarditemDinnebier \harvardand Billinge2008dinnebierPowderDiffractionTheory2008a Dinnebier, R. E. \harvardand Billinge, S. J. L. \harvardyearleft2008\harvardyearright. Powder Diffraction: Theory and Practice. Royal Society of Chemistry.
[13] \harvarditemEgami \harvardand Billinge2012egamiBraggPeaksStructural2012d Egami, T. \harvardand Billinge, S. J. L. \harvardyearleft2012\harvardyearright. Underneath the Bragg Peaks: Structural Analysis of Complex Materials. No. 16 in Pergamon Materials Series. Amsterdam: Elsevier, 2nd ed.
[14] \harvarditemFarrow \harvardand Billinge2009farrowRelationshipAtomicPair2009a Farrow, C. L. \harvardand Billinge, S. J. L. \harvardyearleft2009\harvardyearright. Acta Crystallographica Section A: Foundations of Crystallography, \volbf65(3), 232–239.
[15] \harvarditem[Gobinet et al.]Gobinet, Perrin \harvardand Huez2004gobinetApplicationNonnegativeMatrix2004 Gobinet, C., Perrin, E. \harvardand Huez, R. \harvardyearleft2004\harvardyearright. In 2004 12th European Signal Processing Conference, pp. 1095–1098.
[16] \harvarditemGrippo \harvardand Sciandrone2000grippoConvergenceBlockNonlinear2000 Grippo, L. \harvardand Sciandrone, M. \harvardyearleft2000\harvardyearright. Operations Research Letters, \volbf26(3), 127–136.
[17] \harvarditem[Guan et al.]Guan, Tao, Luo \harvardand Yuan2012guanNeNMFOptimalGradient2012 Guan, N., Tao, D., Luo, Z. \harvardand Yuan, B. \harvardyearleft2012\harvardyearright. IEEE Transactions on Signal Processing, \volbf60(6), 2882–2898.
[18] \harvarditem[Hua et al.]Hua, Allan, Gong, Chater, Schmidt, Geddes, Robertson, Bruce \harvardand Goodwin2021huaNonequilibriumMetalOxides2021 Hua, X., Allan, P. K., Gong, C., Chater, P. A., Schmidt, E. M., Geddes, H. S., Robertson, A. W., Bruce, P. G. \harvardand Goodwin, A. L. \harvardyearleft2021\harvardyearright. Nature Communications, \volbf12(1), 561.
[19] \harvarditem[Huang et al.]Huang, Liu \harvardand Zhou2015huangQuadraticRegularizationProjected2015 Huang, Y., Liu, H. \harvardand Zhou, S. \harvardyearleft2015\harvardyearright. Data Mining and Knowledge Discovery, \volbf29(6), 1665–1684.
[20] \harvarditemJolliffe2002jolliffePrincipalComponentAnalysis2002 Jolliffe, I. T. \harvardyearleft2002\harvardyearright. Principal Component Analysis. Springer Series in Statistics. New York: Springer, 2nd ed.
[21] \harvarditem[Juhás et al.]Juhás, Davis, Farrow \harvardand Billinge2013juhas2013pdfgetx3 Juhás, P., Davis, T., Farrow, C. L. \harvardand Billinge, S. J. \harvardyearleft2013\harvardyearright. Journal of Applied Crystallography, \volbf46(2), 560–566.
[22] \harvarditem[Juhás et al.]Juhás, Farrow, Yang, Knox \harvardand Billinge2015juhas_complex_2015 Juhás, P., Farrow, C., Yang, X., Knox, K. \harvardand Billinge, S. \harvardyearleft2015\harvardyearright. Acta Crystallographica Section A: Foundations and Advances, \volbf71(6), 562–568. Number: 6 Publisher: International Union of Crystallography.
\harvardurlhttp://scripts.iucr.org/cgi-bin/paper?ae5008
[23] \harvarditemKeler \harvardand Karpenko1960keler_reaction_1960 Keler, E. \harvardand Karpenko, N. \harvardyearleft1960\harvardyearright. (Russian) Journal of Inorganic Chemistry (translated from Zhurnal Neorganicheskoi Khimii), (5), 322–325.
[24] \harvarditemKim \harvardand Park2008kimNonnegativeMatrixFactorization2008 Kim, H. \harvardand Park, H. \harvardyearleft2008\harvardyearright. SIAM Journal on Matrix Analysis and Applications, \volbf30(2), 18.
[25] \harvarditem[Kusne et al.]Kusne, Keller, Anderson, Zaban \harvardand Takeuchi2015kusneHighthroughputDeterminationStructural2015a Kusne, A. G., Keller, D., Anderson, A., Zaban, A. \harvardand Takeuchi, I. \harvardyearleft2015\harvardyearright. Nanotechnology, \volbf26(44), 444002.
[26] \harvarditemLee \harvardand Seung1999leeLearningPartsObjects1999d Lee, D. D. \harvardand Seung, H. S. \harvardyearleft1999\harvardyearright. Nature, \volbf401(6755), 788–791.
[27] \harvarditemLin2007linProjectedGradientMethods2007 Lin, C.-J. \harvardyearleft2007\harvardyearright. Neural Computation, \volbf19(10), 2756–2779.
[28] \harvarditem[Liu et al.]Liu, Tao, Hsu, Du \harvardand Billinge2019liuUsingMachineLearning2019 Liu, C.-H., Tao, Y., Hsu, D. J., Du, Q. \harvardand Billinge, S. J. L. \harvardyearleft2019\harvardyearright. Acta Crystallogr. A, \volbf75, 633–643.
[29] \harvarditem[Liu et al.]Liu, Wright, Gu, Bandi, Wustrow, Todd, O’Nolan, Beauvais, Neilson, Chupas, Chapman \harvardand Billinge2021liuValidationNonnegativeMatrix2021a Liu, C.-H., Wright, C. J., Gu, R., Bandi, S., Wustrow, A., Todd, P. K., O’Nolan, D., Beauvais, M. L., Neilson, J. R., Chupas, P. J., Chapman, K. W. \harvardand Billinge, S. J. L. \harvardyearleft2021\harvardyearright. Journal of Applied Crystallography, \volbf54(3).
[30] \harvarditem[Long et al.]Long, Bunker, Li, Karen \harvardand Takeuchi2009longRapidIdentificationStructural2009a Long, C. J., Bunker, D., Li, X., Karen, V. L. \harvardand Takeuchi, I. \harvardyearleft2009\harvardyearright. Review of Scientific Instruments, \volbf80(10), 103902.
[31] \harvarditem[Martinolich et al.]Martinolich, Kurzman \harvardand Neilson2015martinolichPolymorphSelectivitySuperconducting2015 Martinolich, A. J., Kurzman, J. A. \harvardand Neilson, J. R. \harvardyearleft2015\harvardyearright. Journal of the American Chemical Society, \volbf137(11), 3827–3833.
[32] \harvarditem[Morup et al.]Morup, Madsen \harvardand Hansen2007morupShiftedNonNegativeMatrix2007 Morup, M., Madsen, K. H. \harvardand Hansen, L. K. \harvardyearleft2007\harvardyearright. In 2007 IEEE Workshop on Machine Learning for Signal Processing, pp. 139–144.
[33] \harvarditem[O’Nolan et al.]O’Nolan, Huang, Kamm, Grenier, Liu, Todd, Wustrow, Thinh Tran, Montiel, Neilson, Billinge, Chupas, Thornton \harvardand Chapman2020onolan_thermal-gradient_2020 O’Nolan, D., Huang, G., Kamm, G. E., Grenier, A., Liu, C.-H., Todd, P. K., Wustrow, A., Thinh Tran, G., Montiel, D., Neilson, J. R., Billinge, S. J. L., Chupas, P. J., Thornton, K. S. \harvardand Chapman, K. W. \harvardyearleft2020\harvardyearright. Journal of Applied Crystallography, \volbf53(3), 662–670. _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1107/S160057672000415X.
\harvardurlhttps://onlinelibrary.wiley.com/doi/abs/10.1107/S160057672000415X
[34] \harvarditem[O’Nolan et al.]O’Nolan, Zhao, Chen, Grenier, Beauvais, Newton, Nenoff, Chupas \harvardand Chapman2021onolanMultimodalAnalyticalToolkit2021 O’Nolan, D., Zhao, H., Chen, Z., Grenier, A., Beauvais, M. L., Newton, M. A., Nenoff, T. M., Chupas, P. J. \harvardand Chapman, K. W. \harvardyearleft2021\harvardyearright. Chemical Science, \volbf12(41), 13836–13847.
[35] \harvarditemPaatero \harvardand Tapper1994paateroPositiveMatrixFactorization1994 Paatero, P. \harvardand Tapper, U. \harvardyearleft1994\harvardyearright. Environmetrics, \volbf5(2), 111–126.
[36] \harvarditemPecharsky \harvardand Zavalij2008pecharskyFundamentalsPowderDiffraction2008 Pecharsky, V. \harvardand Zavalij, P. \harvardyearleft2008\harvardyearright. Fundamentals of Powder Diffraction and Structural Characterization of Materials, Second Edition. Springer Science & Business Media.
[37] \harvarditemPorter2020porter_danporterdans_diffraction_2020 Porter, D., \harvardyearleft2020\harvardyearright. DanPorter/dans_diffraction. Language: eng.
\harvardurlhttps://zenodo.org/record/3859501
[38] \harvarditem[Rakita et al.]Rakita, Hart, Das, Shahrezaei, Foley, Mathaudhu, Nicolopoulos, Taheri \harvardand Billinge2023rakit;am23 Rakita, Y., Hart, J. L., Das, P. P., Shahrezaei, S., Foley, D. L., Mathaudhu, S. N., Nicolopoulos, S., Taheri, M. L. \harvardand Billinge, S. J. L. \harvardyearleft2023\harvardyearright. Acta Materialia, \volbf242, 118426.
[39] \harvarditem[Rakita et al.]Rakita, O’Nolan, McAuliffe, Veith, Chupas, Billinge \harvardand Chapman2020rakitaActiveReactionControl2020c Rakita, Y., O’Nolan, D., McAuliffe, R. D., Veith, G. M., Chupas, P. J., Billinge, S. J. L. \harvardand Chapman, K. W. \harvardyearleft2020\harvardyearright. Journal of the American Chemical Society, \volbf142(44), 18758–18762.
[40] \harvarditem[Rayder et al.]Rayder, Formalik, Vornholt, Frank, Lee, Alzayer, Chen, Sengupta, Islamoglu, Paesani, Chapman, Snurr \harvardand Farha2023rayderUnveilingUnexpectedModulatorCO22023 Rayder, T. M., Formalik, F., Vornholt, S. M., Frank, H., Lee, S., Alzayer, M., Chen, Z., Sengupta, D., Islamoglu, T., Paesani, F., Chapman, K. W., Snurr, R. Q. \harvardand Farha, O. K. \harvardyearleft2023\harvardyearright. Journal of the American Chemical Society, \volbf145(20), 11195–11205.
[41] \harvarditem[Ren et al.]Ren, Pueyo, Zhu, Debes \harvardand Duchêne2018renNonnegativeMatrixFactorization2018 Ren, B., Pueyo, L., Zhu, G. B., Debes, J. \harvardand Duchêne, G. \harvardyearleft2018\harvardyearright. The Astrophysical Journal, \volbf852(2), 104.
[42] \harvarditemSra \harvardand Dhillon2006sraNonnegativeMatrixApproximation2006 Sra, S. \harvardand Dhillon, I. \harvardyearleft2006\harvardyearright. Nonnegative Matrix Approximation: Algorithms and Applications. Tech. rep. Computer Science Department, University of Texas at Austin.
[43] \harvarditem[Su et al.]Su, Feth \harvardand Lehoczky2009su_thermal_2009 Su, C.-H., Feth, S. \harvardand Lehoczky, S. L. \harvardyearleft2009\harvardyearright. Materials Letters, \volbf63(17), 1475–1477.
\harvardurlhttps://www.sciencedirect.com/science/article/pii/S0167577X09002638
[44] \harvarditem[Thatcher et al.]Thatcher, Liu, Yang, McBride, Thinh Tran, Wustrow, Karlsen, Neilson, Ravnsbæk \harvardand Billinge2022thatc;aca22 Thatcher, Z., Liu, C.-H., Yang, L., McBride, B. C., Thinh Tran, G., Wustrow, A., Karlsen, M. A., Neilson, J. R., Ravnsbæk, D. B. \harvardand Billinge, S. J. L. \harvardyearleft2022\harvardyearright. Acta Crystallographica Section A: Foundations and Advances, \volbf78(3).
[45] \harvarditemWang \harvardand Zhang2013wangNonnegativeMatrixFactorization2013a Wang, Y.-X. \harvardand Zhang, Y.-J. \harvardyearleft2013\harvardyearright. IEEE Transactions on Knowledge and Data Engineering, \volbf25(6), 1336–1353.
[46] \harvarditemXu \harvardand Yin2017xuGloballyConvergentAlgorithm2017 Xu, Y. \harvardand Yin, W. \harvardyearleft2017\harvardyearright. Journal of Scientific Computing, \volbf72(2), 700–734.
[47] \harvarditem[Xu et al.]Xu, Zhang, Wang, Chang \harvardand Liang2010xuL1Regularization2010 Xu, Z., Zhang, H., Wang, Y., Chang, X. \harvardand Liang, Y. \harvardyearleft2010\harvardyearright. Science China Information Sciences, \volbf53(6), 1159–1169.
[48] \harvarditem[Yang et al.]Yang, Culbertson, Thomas, Vuong, Kjær, Jensen, Tucker \harvardand Billinge2021yangCloudPlatformAtomic2021c Yang, L., Culbertson, E. A., Thomas, N. K., Vuong, H. T., Kjær, E. T. S., Jensen, K. M. Ø., Tucker, M. G. \harvardand Billinge, S. J. L. \harvardyearleft2021\harvardyearright. Acta Crystallographica Section A: Foundations and Advances, \volbf77(1), 2–6.
[49] \harvarditem[Yang et al.]Yang, Juhás, Terban, Tucker \harvardand Billinge2020yangStructureminingScreeningStructure2020g Yang, L., Juhás, P., Terban, M. W., Tucker, M. G. \harvardand Billinge, S. J. L. \harvardyearleft2020\harvardyearright. Acta Crystallographica Section A: Foundations and Advances, \volbf76(3), 395–409.
[50] \harvarditem[Yang et al.]Yang, Juhas, Farrow \harvardand Billinge2015yang;arxiv14 Yang, X., Juhas, P., Farrow, C. L. \harvardand Billinge, S. J. L. \harvardyearleft2015\harvardyearright. arXiv:1402.3163 [cond-mat].
[51]

Stretched Non-negative Matrix Factorization

Abstract

keywords:

keywords:

keywords:

keywords:

1 Introduction

2 Stretched Nonnegative matrix factorization

3 Numerical Solution of stretchedNMF

4 Diffraction Use Case

4.1 Introduction

4.2 Data

4.2.1 Simulated PXRD and PDF data with increasing lattice parameters

4.2.2 Experimental PXRD data - thermal expansion

4.2.3 Experimental PXRD data - thermal expansion and reaction

4.3 Algorithm Developments

5 Numerical Results

5.1 Results on Simulated PDF

5.2 Results on Simulated PXRD

5.3 Results on simulated PDF and PXRD data with small expansion coefficients

5.4 Results on measured PXRD data I

5.5 Results on real PXRD data II

6 Discussion and Conclusion

Acknowledgment

References

7 Supplementary Information

7.1 Simulated PDF

7.2 Simulated PXRD

7.3 Simulated PDF and PXRD data with small expansion coefficients

7.4 Real PXRD I