# Constrained non-negative matrix factorization enabling real-time insights of
_in situ_ and high-throughput experiments

Phillip M. Maffettone  National Synchrotron Light Source II, Brookhaven
National Laboratory, Upton, NY 11973, USA  Aidan C. Daly  Center for
Computational Biology, Flatiron Institute, U.S.A  Daniel Olds<EMAIL_ADDRESS>National Synchrotron Light Source II, Brookhaven National Laboratory, Upton,
NY 11973, USA

###### Abstract

Non-negative Matrix Factorization (NMF) methods offer an appealing
unsupervised learning method for real-time analysis of streaming spectral data
in time-sensitive data collection, such as in situ characterization of
materials. However, canonical NMF methods are optimized to reconstruct a full
dataset as closely as possible, with no underlying requirement that the
reconstruction produces components or weights representative of the true
physical processes. In this work, we demonstrate how constraining NMF weights
or components, provided as known or assumed priors, can provide significant
improvement in revealing true underlying phenomena. We present a PyTorch based
method for efficiently applying constrained NMF and demonstrate this on
several synthetic examples. When applied to streaming experimentally measured
spectral data, an expert researcher-in-the-loop can provide and dynamically
adjust the constraints. This set of interactive priors to the NMF model can,
for example, contain known or identified independent components, as well as
functional expectations about the mixing of components. We demonstrate this
application on measured X-ray diffraction and pair distribution function data
from in situ beamline experiments. Details of the method are described, and
general guidance provided to employ constrained NMF in extraction of critical
information and insights during in situ and high-throughput experiments.

##  I Introduction

High-throughputKusne  _et al._ (2014); Gómez-Bombarelli  _et al._ (2016);
Greenaway  _et al._ (2018); Langner  _et al._ (2020); Ludwig (2019); Decker
_et al._ (2017) and _in situ_ Yeung  _et al._ (2016); MacLeod  _et al._ (2020)
analyses have become ubiquitous in modern materials discovery. Automation and
engineering has increased capabilities for data streaming and _in situ_
experimentation in many settings, including university laboratories,Li  _et
al._ (2020) industrial research,Fleischer  _et al._ (2018) and central
facilities.Campbell  _et al._ (2020) This technological advancement opens up
new opportunities for interrogating the dynamic processes of a system
including phase transitionsOlds  _et al._ (2017) and reactive processes.Yeung
_et al._ (2016) Such opportunities bear the burden of an increased rate of
data production, as well as a commensurate increase in the volume of data to
process. For many experiments, these high-volume datasets require physically
meaningful analysis to be performed in real-time. Diffraction experiments are
a prototypical materials characterization technique that produce such data.
X-ray diffraction (XRD) and pair distribution function (PDF) analysis are
necessary for comprehending energy materials across batteries,Hua  _et al._
(2021); Shadike  _et al._ (2021) catalysis,Wang, Li, and Domen (2019); Löffler
_et al._ (2019) photovoltaics,Olds (2020) and multi-ferroics.Clarke  _et al._
(2021) This challenge is most pressing at central facilities, which have the
technological capacity to produce vast amounts of data across variable sample
environments.

With the advances in accelerator technology achieved in recent decades,
available photon flux at synchrotron X-ray lightsources is increasing at a
rate even faster than Moore’s law.Olds (2020) With new beamlines capable of
utilizing these record-breaking levels of brightness and resolution,
scientists can now perform many of the critical high-throughput and _in situ_
studies required to realize future energy technologies. Largely driven by
these new capabilities, facilities are generating more data now than ever
before. For example, at the National Synchrotron Light Source II (NSLS-II)
produced over 1 PB of raw data in the last year, and that rate is expected to
increase as the facility matures.Campbell  _et al._ (2020) Despite these huge
data generation rates, however, approaches to data analysis have barely
changed over time, with users often relying on methods and code that are
decades old. Consequently, data collected in seconds to minutes may take weeks
or months of analysis to understand. Due to such limitations, conventional
knowledge extraction is often divorced from the beamtime measurement process.
This lack of real-time feedback can force users into ‘flying blind’ during
their limited time at the beamline, leading to missed opportunities and
mistakes.

This information bottleneck is most pronounced during diffraction experiments,
where the atomic structure of a sample can be measured in fractions of a
second, in both real (PDF) and reciprocal (XRD) space. Historically, such data
are modeled using iterative refinement procedures such as Rietveld refinements
to extract the relevant physical insights such as crystallographic phase,
lattice constants, atomic positions, site occupancies, and atomic displacement
parameters.Rietveld (1967) As these refinements may take minutes to hours for
a skilled expert to perform, they are daunting for technique novices and often
do not scale well to large datastreams.Giacovazzo (2011) While such full
information extraction from the data can be performed _a posteriori_ , users
distinctly benefit from the extraction of critical information at-the-
beamline, such as if the sample has undergone a sought transition or phase
transformation.

Machine learning (ML) tools offer great promise for monitoring and analyzing
large volumes of streaming data.Campbell  _et al._ (2020) These tools have
been demonstrated for classifying diffraction patterns,Iwasaki, Kusne, and
Takeuchi (2017); Xiong  _et al._ (2017); Long  _et al._ (2009); Oviedo  _et
al._ (2019); Lee  _et al._ (2020); Ziletti  _et al._ (2018); Aguiar  _et al._
(2019); Di Chen  _et al._ (2020) or decomposing large datasets.Takeuchi  _et
al._ (2005); Olds  _et al._ (2017); Stanev  _et al._ (2018); Bai  _et al._
(2017, 2018); Suram  _et al._ (2017) Particularly, unsupervised learning is
useful for offering rapid diagnostics and analysis without prior training.Hua
_et al._ (2021); Geddes  _et al._ (2019); Stanev  _et al._ (2018) Because
unsupervised decomposition—or data reduction—algorithms do not require any
prior data labeling or expected phases, they are exceptionally useful for
interpreting datasets where there are potentially new or unexpected phases.
Two such algorithms, principle component analysis (PCA)Olds  _et al._ (2017)
and non-neagative matrix fatorization (NMF),Geddes  _et al._ (2019); Li  _et
al._ (2020) have been employed to segregate components from diffraction
datasets. The algorithms aim to describe a dataset of $m$ experimental
patterns (XRD or PDF) as a weighted linear composition of $k$ fundamental
components. Of these, NMF is particularly well suited to scattering data
because its decomposition aligns with physics: the component signals are
enforced to remain greater than or equal to zero, and are strictly additive to
one another in the dataset reconstruction.

While NMF has the advantages of retaining the additive superposition and non-
negativity of XRD and PDF, it comes with some significant drawbacks. Firstly,
NMF is incapable of interpreting translational variance. Since the NMF model
uses only a linear combination of components, a given component cannot
effectively reproduce data that is identical except for translation.
Unfortunately this is a common occurrence for XRD and PDF experiments due to
thermal expansion, changes in nanostructure, and lattice strain due to
Vergad’s law.Stanev  _et al._ (2018) Efforts to overcome this translation have
included discarding translationally redundant components, although this
requires additional hyperparameters for detecting the redundancy that are not
effective across different experimental challenges.Iwasaki, Kusne, and
Takeuchi (2017); Stanev  _et al._ (2018) Furthermore, there is no unique set
of components and weights that can satisfy the NP-hard optimization for the
NMF decomposition. This means that the NMF result will vary based on random
initialization, and is not guaranteed to converge to a local minimum for
reconstruction error. Lastly, while non-negativity is a sensible constraint,
it doesn’t require the resulting components to correspond to interpretable
patterns nor the weights to carry any meaning with respect to volume fraction.
It can be argued that there is a unique set of components in the analysis of
diffraction data that is physically meaningful. Metropolis matrix
factorization attempts to identify this unique set by constraining the
constraining select components and finding the global optimum for the
constrained decomposition.Geddes  _et al._ (2019); Hua  _et al._ (2021)
Although this approach is powerful, it remains too costly to run real-time
during an measurement that is streaming data.

To successfully apply decomposition methods during a diffraction experiment,
they must operate efficiently and resolve components and weights that are
physically meaningful. Instead of simply finding the decomposition model that
best fits the data, we need the _most likely physical_ model that best fits
the data. By increasing the likelihood that found components and weights are
physically relevant to real chemical process occurring within the sample, it
would present results in an immediately convincing representation to the
researcher. For example, a sudden shift in weights between two components may
suggest to a researcher that a second order phase transition has occurred
during an experiment; however, if the researcher is presented with nonphysical
components to inspect, they are far less likely to understand or trust the
information even if it later proves true. Fortunately, the nature of these _in
situ_ experiments typically follows a standard trajectory across a set of
state variables, under which the sample begins readily known phase to the
researcher that may undergo a transition at some future point into known or
unknown phases.Li  _et al._ (2021) Various combinations of known and unknown
components or weights can be occur in a data series, such as knowledge of two
components at the edges of a dataset with curiosity about a third, or no
knowledge of the components, but expectations for their concentrations.Geddes
_et al._ (2019) Thus, through use of constraints we can employ NMF in highly
productive and insightful ways for diffraction experiments.

In this contribution, we present our advances developing constrained NMF
methods for use on _in situ_ studies at total scattering and powder
diffraction beamlines. The algorithms are developed using the PyTorchPaszke
_et al._ (2019) framework for speed and scalability. We first demonstrate the
capabilities of constrained NMF on several synthetic examples, in which the
limitations of canonical NMF are explored. We then use constrained NMF to
decompose temperature dependent XRD and PDF datasets of a ferroelectric
perovskite and a molten salt, wherein there is no _a priori_ knowledge of
components in the latter. Lastly, we provide best practices on how to apply
these methods during _in situ_ experiments and interpret the results as they
relate to commonly encountered scenarios in the analysis of energy materials.
These developments are extensible to applications in spectroscopy and
microscopy, and can also be deployed at university and industrial
laboratories. Within each of these results, constrained NMF behaves as a
companion agent in concert with the researcher, leveraging the expertise of
the scientist while providing real-time insights.

##  II Results and Discussion

The goal of NMF is to find an $m\times k$ non-negative matrix $\mathbf{W}$ of
weights and a $k\times n$ non-negative matrix $\mathbf{H}$ of components, to
approximate a dataset $\mathbf{X}\approx\mathbf{WH}$, where $m$ is the number
of patterns in a dataset, $n$ is the number of points in a pattern, and $k$ is
the chosen number of components. The details of this optimization are
presented in Section IV. We augmented this approach to allow the user — or an
algorithm — to initialize any column of $\mathbf{W}$ or row of $\mathbf{H}$,
and to subsequently constrain those vectors from being updated during the
optimization. It is worth noting that the algorithm is often presented with
the contents of $\mathbf{W}$ and $\mathbf{H}$ reversed; however, we find more
intuitive to have $\mathbf{W}$ refer to weights.

###  II.1  Synthetic and experimental datasets

To test the capabilities of canonical and constrained NMF, we assembled two
synthetic datasets (Fig. 1). The first approach highlights a common case where
two unique, overlapping, components will mix linearly with coefficients that
sum to one. This occurs in binary mixtures when the measured response function
is a molar weighted sum of the pure response functions, as is the case for
diffraction. We represent this with two Gaussians mixed linearly, and include
random noise in coefficients of mixing and the final functions. Under these
circumstances, it is likely to know the weights of the individual components
_a priori_ , but not necessarily the components themselves.

The second synthetic dataset contains three unique functions that are mixed
non-monotonically as the dataset progresses. The three functions are a non-
overlapping Lorentzian, and a Gaussian and box function that overlap to create
non-differentiable result. This dataset explores a limiting and challenging
use case for the underlying components, that can occur with a known mask or
background augmentation. The weights represent commonly encountered ways that
physical and chemical processes — and therefore components — can evolve
naturally. This family of challenges highlights the scenario wherein _a
priori_ knowledge of one or all components exists, but not necessarily
knowledge the weights.

Figure 1: Example datasets used to demonstrate the capabilities of canonical
and constrained NMF. (a) Two noisy Gaussian curves centered at 2.5 and -2.5
are mixed in linear combinations with coefficients that vary approximately,
but not exactly, linearly across the dataset. (b) A Lorentzian centered at
3.0, a box function of width two centered at 5.0, and a Gaussian centered at
6.0 are mixed as linear combination with functionally varying weights across
the dataset.

We next examined to two distinct experimental challenges, analyzing the
exceedingly subtle temperature driven phase transitions of barium titanate,
and understanding the structural evolution of a salt, NaCl:CrCl3, heated
beyond it’s melting temperature. BaTiO3 is a well studied material that has
been used to probe the mechanism of ferroelecticity because its phase
transitions exist at relatively low temperature.Page  _et al._ (2010); Wegner
_et al._ (2020) At high temperature, BaTiO3 exists as a paraelectric cubic
phase, where on average the Ti atoms are octahedrally coordinated by six O
atoms. At room temperature the average displacement of the Ti elongates along
[001], resulting in a ferroelectric tetragonal phase. Further cooling results
in an orthorhombic phase with polarization along the [011]. At low temperature
the rhombohedral ground state is reached, with all Ti atoms polarized along
the [111]. This symmetry breaking results in subtle peak splitting, which is
further obfuscated by thermal expansion driven peak shifting (Fig. 2). The
differences between patterns is so subtle that classification across these
temperature-induced phase transitions is not possible considering only
Rietveld refinements without expert intervention.Page  _et al._ (2010)

Contrary to changes in BaTiO3, an crystalline–amorphous transition is abrupt
and distinct. We used constrained NMF to understand the behavior of NaCl:CrCl3
across a temperature range of 300–963 K, which included its melting transition
(Fig. 2). These molten salts are of interest for the development of nuclear
reactors.Li  _et al._ (2021) In this experiment there was a residual solid
phase after the melting transition. Without any additional analysis, it
remained unclear to the research team whether that solid was a distinct high-
temperature phase or an experimental artefact. Together, these four datasets
present the clear utility of constrained NMF for on-the-fly and _in situ_
analysis of diffraction data.

Figure 2:  Constrained NMF was used to analyze two variable temperature
diffraction datasets. (a) BaTiO3 undergoes three phase transitions across the
temperature range 150–450 K, that are imperceptible even to traditional
refinement techniques. (b) NaCl:CrCl3 diffraction was measured a temperature
range of 300–963 K that included a melting transition and an anomolous high
temperature solid phase (shown by the remaining sharp peaks above 650 K).

###  II.2 Constraining weights with synthetic datasets

With the first synthetic dataset (Fig.1a) we demonstrated the challenge of
retaining physically relevant components from NMF when the weights are readily
available. This is a common situation in diffraction studies of variable
composition or phase mixing, when other experimental information (preparation
conditions or spectroscopy) will indicate the concentration of individual
phases without specifying the diffraction patterns of the phases. We
considered a continuous mixture of two Gaussian functions, where noise is
added to the composition and to the final function. From this we can see the
power of canonical NMF for dimensionality reduction and the benefit offered by
constraints.

In decomposing the Gaussian dataset (Fig 1(a)) using canonical NMF, the
reconstruction loss as defined by squared error (SE) (eq. 3) achieved a mean
value (MSE) of 5.788 over the dataset. The common crystallographic metric of
R-factor (eq. 6) was 0.066 for the reconstruction on average. The
decomposition into two components created a model of the mixtures that fits
the data well. Even though the learned components are similar to the true
functions, their respective magnitudes are not obviously linked to the
magnitudes of those functions (Fig 3d). Furthermore, the weights must
compensate for this increased predicted magnitude, and inevitably over-fit the
noise in the dataset and (Fig 3e). This is the expected behavior of NMF, but
is naive to any prior knowledge or physical constraints we can impose on the
dataset.

Figure 3:  The reconstruction of the (a) first, (b) median, and (c) last
pattern of a dataset of mixed Gaussian functions using canonical NMF. The full
dataset is shown in Figure 1(a). (d) Both the ordering and amplitude of the
learned components fails to match the ground truth functions. (e) The learned
weights compensate for this false magnitude and fit accurately to the non-
linearity in the true weights.

However, when we informed the NMF model with prior knowledge of weights, and
constrained those weights, we constructed the most likely model that best fits
the data. This resultant model successfully retains the physical significance
on the weights and components. In this use case, the weights were limited to
the coordinates of a simplex (that is, they sum to unity), and expected to
vary linearly. By initializing the NMF with these exact weights, and
constraining their update in the optimization algorithm, we created a
decomposition with a MSE of 6.903 and mean R-factor 0.0722. By design these
metrics are slightly worse than those for the canonical NMF because we do not
allow the decomposition to over-fit the noise in the weights, which not
exactly linear. This forces the weights to follow the prior and keeps the
reconstruction robust to outliers. Importantly, the resulting components from
the constrained NMF are physically meaningful with respect to the true
function. In a crystallographic or spectroscopic context, this is important
because peak intensity carries physical significance.

Figure 4:  The reconstruction of the (a) first, (b) median, and (c) last
pattern of a dataset of mixed Gaussian functions using constrained NMF. The
weights of the decomposition were constrained to be linear with respect to the
dataset index. The full dataset is shown in Figure 1(a). (d) Both the ordering
and amplitude of the learned components successfully matches the ground truth
functions. (e) By fixing the weights, they cannot match the non-linearity in
the ground truth but serve as a physically meaningful approximation.

The results from this constrained NMF are obtained rapidly, and can be
produced on-the-fly during measurement. The output components offer immediate
insight into the underlying dataset, while depending on an physical
approximate prior. These components can be used as starting conditions for
exact refinement, such as a Rietveld refinement,Rietveld (1967) of each
pattern or spectrum in a dataset. When analysis is less time constrained,
other methods for further tuning the optimization of constrained NMF can be
explored with continual human input. Some such methods include using non-rigid
constraints and Markov chain Monte Carlo optimization,Geddes  _et al._ (2019)
which would enable a flexible, non-exact, linear prior.

###  II.3 Constraining components with synthetic datasets

As a second controlled example, we created a synthetic dataset from three
diverse functions with smooth, non-monotonically varying weights in a linear
combination (Fig 1b). Those functions included an overlapping box and Gaussian
function, and a Lorentz function which does not overlap either. In this
scenario, knowledge of any of these components significantly constrains the
NP-hard optimization problem of NMF.

Unlike the previous example, this dataset was not effectively decomposed by
canonical NMF. The resulting reconstruction MSE was 25.441 with a mean
R-factor of 0.034. One reason for this failure was the non-differentiable box
function. This property makes decomposition particularly challenging for the
gradient descent methods used in the optimization. As shown in Figure 5d, the
model found a local minimum in the optimization where the learned components
were partial mixtures of the true components. This is both physically
inaccurate and results in a poor reconstruction. In a case such as this, the
decomposition from NMF results in a weak, non-physical, dimensionality
reduction, which offers no new knowledge to the scientist.

Figure 5:  The reconstruction of the (a) first, (b) median, and (c) last
pattern of a dataset of mixed Gaussian, Lorenzian, and box functions using
canonical NMF. The full dataset is shown in Figure 1(b). (d) Canonical NMF
fails to successfully decompose this dataset and produce an accurate
reconstruction. The learned components are mixtures of the underlying
functions. (e) The learned weights compensate for the inaccurate components,
albeit the decomposition is trapped in a local optimum.

Corollary to the previous example where prior knowledge of the weights enabled
a physically meaningful extraction of components, when the components are fed
to NMF as a prior to decompose these composite functions, we were able to
extract the ground truth weights. As shown in Figure 5e, the weights are not
trivial, and prior knowledge of the component functions enabled the model to
learn those weights (Fig. 6). Again we see the rigidity of the priors ignoring
the inherent noise, guiding the model away from overfitting. Yet, the
constrained NMF outperformed canonical NMF even by the standard metrics
(MSE=1.336, $<\mathrm{R-factor}>=$0.013). This approach is also applicable if
only some of the components are known (Fig. S1–6), albeit less robust for this
use case. This suggest that a procedure which allows for iterative inclusion
of more priors through human or automatic intervention would allow for the
most complex datasets to be successfully resolved using constrained NMF.

Figure 6: The reconstruction of the (a) first, (b) median, and (c) last
pattern of a dataset of mixed Gaussian, Lorenzian, and box functions using
constrained NMF. The full dataset is shown in Figure 1(b). (d) Prior knowledge
of the functions fixes the underlying components. (e) The learned weights are
meaningful and create an accurate reconstruction of the dataset.

###  II.4 Constraints enable physically meaningful phase segregation of
BaTiO3.

Following this validation of constrained NMF in controlled systems, we sought
to test its ability to construct a meaningful decomposition in a subtle and
challenging crystallographic problem with measured data. The BaTiO3 dataset
presented in Figure 2 contains four distinct phases across the temperature
range measured. These transitions are subtle and not detectable considering
traditional methods without expert intervention. Between these second order
phase transitions, there is continuous peak shifting due to thermal expansion.

Canonical NMF is poorly suited to address translational variance (peak
shifting) in decomposing a dataset into composite patterns. Stanev  _et al._
(2018); Bai  _et al._ (2018); Suram  _et al._ (2017) However, it was also
incapable of identifying the clear second order phase transitions within the
data. Instead it learned mixtures of the true components as its pure
components. From Figure 7b, one might falsely conclude that the high
temperature phase (blue) increases approximately linearly in composition with
temperature, or that there are three components which make up the low
temperature phase. Despite this, the reconstruction loss from NMF was quite
low with an MSE of 0.150 and a mean R-factor of 0.039. This demonstrates that
canonical NMF is not well suited for crystallographic datasets with subtle
transitions, a common case for diffraction studies.

Figure 7:  Canonical NMF fails to produce a physically meaningful
decomposition of the variable temperature BaTiO3 dataset (Fig. 2(a)). (a) The
learned components are randomized mixtures of real patterns. (b) The resulting
weights carry no significance regarding the composition or any other property
of the physical system.

We next considered constraining the components of NMF to extract a better,
more physically realistic, decomposition of the BaTiO3 dataset. We sought to
approach this in a way that made no presumption of a specific phase, and only
operated under the knowledge that the dataset was measured along a single
state variable (temperature). Our synthetic experiment mixing diverse
functions demonstrated that constraining all of the components led to the
optimal solution, albeit partial constraints enabled the extraction of
physically sensible intermediates. This concept has been demonstrated
previously using different approaches on diffraction data.Gu, Du, and Billinge
(2021); Geddes  _et al._ (2019) To construct a physically meaningful solution
without explicit knowledge of the underlying phases, we designed an iterative
algorithm that wraps our constrained NMF to handle datasets over state
variables such as temperature or composition.

The iterative application of constrained NMF automatically selects components
to constrain, thus depending only on the user selected number of components.
It begins by assuming the end members of a dataset (in this case lowest and
highest temperature patterns) are unique, fixing those members as components,
and optimizing a partially constrained NMF. In the limiting case of a rank 2
NMF, this amounts to combinatorial decomposition, which has been successful in
interpreting PDF data.Olds  _et al._ (2017) When more components are desired,
the algorithm chooses the most populated variable component, then finds the
location where this component is most prevalent in the dataset. The next
component is then fixed to the member at this location, and the NMF
optimization proceeds with one additional constraint. This process continues
until all members of the NMF are constrained. Critically for total scattering
datasets, this results in weights that are bounded between 0 and 1, thus
possessing a meaningful mapping to composition.

The automated, iterative, constrained NMF approached a similar reconstruction
accuracy as the canonical NMF (MSE=0.349, $<\mathrm{R-factor}>=$0.051).
However, the resulting components and weights are immediately interpretable by
the scientist. As the fixed components were chosen directly from the measured
dataset, these can be used directly as reference diffraction patterns for pure
phases (Fig. 8a). Furthermore, the weights can be interpreted as molar
fractions, showing phase mixtures during transitions (Fig. 8b). The
applications of constraints — even without prior knowledge of what those
constraints should be — enabled a substantial improvement over canonical NMF,
providing both physically meaningful components and weights.

Figure 8:  Automated constrained NMF produces a physically meaningful
decomposition of the variable temperature BaTiO3 dataset (Fig. 2(a)). (a) The
learned components correspond to the pure phases at the median temperatures
where they are observed. (b) The resulting weights approximate mixing between
phases and thermal expansion causing translational variance in the XRD
patterns.

There is still some room for improvement of this method. The transitions
provided by the constrained NMF are still too gradual for second order phase
transitions. While they occur in approximately the correct location, they
indicate more mixing than is actually observed. This is a direct result of
translation dependence of NMF: peak shifts in pure phases are likely cast as
mixtures between minimally and maximally shifted components that are otherwise
identical. Supervised or semi-supervised convolutional deep learning methods
have exceeded this performance.Oviedo  _et al._ (2019); Lee  _et al._ (2020);
Maffettone  _et al._ (2021) Nonetheless, to the authors’ knowledge this is the
best performance of any fully unsupervised ML approach at segregating the
phases in such a dataset.

###  II.5 Constraints enable physically meaningful _in situ_ analysis of
melting salts.

We then used the constrained NMF to monitor an experimental transition into an
amorphous state with an unknown number of component phases. Canonical NMF
behaves as expected with the same shortcomings we noted in the study of
BaTiO3. It elucidated the solid to amorphous phase boundary with a clear
second order phase transition shown at 676 K (Fig. 9). The reconstruction MSE
was 3.355 with a mean R-factor of 0.020. By integrating the NMF into data
acquisition pipeline,Allan  _et al._ (2019) this enables the researcher to
note the transition temperature on-the-fly; however, the canonical approach is
unable to extract meaningful components. Without additional processing, the
solid material that existed past the transition temperature was unable to be
characterized. During such an experiment—especially at a central facility—it
is crucial to know whether the high temperature solid phase is unique and
mandates further investigation at the time of measurement. In the case of any
experiment that is irreversible, destructive, or expensive to reproduce,
immediate access to this information is paramount.

Figure 9:  Canonical NMF produces a confusing decomposition, offering limited
clarity to the scientist. (a) The learned components are non-physical and
cannot be interpreted as XRD patterns. (b) The learned weights provide little
information beyond an abrupt change around 400 K.

The constrained NMF approach offered this immediate analysis capability. Using
the elbow methodThorndike (1953) to examine the MSE and KL-divergence loss, we
noted diminishing returns in the NMF reconstruction at around 4 components
(Fig. 10). The 4-component MSE was 9.193 with a mean R-factor of 0.0334.
Considering the region around the inflection point in these learning curves,
the uncertain effects of peak shifting from thermal expansion, and potential
secondary phases, we considered a constrained NMF that allowed for 3, 4, and 5
components. Using the same automatic iterative constraining procedure from
above, the constrained NMF elucidated a significantly more meaningful
decomposition than canonical NMF.

Figure 10:  Following the loss curve of NMF with respect to the number of
components shows diminishing returns around 4 components. Since this heuristic
identification is subjective in nature, multiple instances of NMF can be run
with a range of components.

Using multiple constrained NMFs of different rank elucidated the difference
between a unique and prominent component, and a shifting component due to
thermal expansion. In the rank-3 constrained NMF, two solid components and one
amorphous component arose in the decomposition (Fig. 11a–b). The solid
components were nearly identical and showed a linear variation in weight
(composition), which was 0.5/0.5 at the midpoint in the temperature range
prior to the second order phase transition. A small and decaying amount of the
high temperature solid phase was present early after the transition. This
linear exchange between the two similar components below 400 K arises from the
translational dependence of NMF: the algorithm compensates by considering two
identical [besides shift] patterns at both ends of the relevant temperature
range. From the rank-3 constrained NMF alone, we can see the gradual impact of
thermal expansion followed by a second order transition where some of the
material remained solid.

Figure 11:  Automated constrained NMF was run using 3-6 components
asynchronously on the molten salt dataset. The learned components for the (a)
3, (c) 4, and (e) 6 component decompostions provide insight into thermal
expansion, a phases present during the temperature ramp. The learned weights
for the (b) 3, (d) 4, and (f) 6 component decompositions demonstrate the
second-order phase transition, the characteristic response to thermal
expansion, and persistence of a solid phase after 400 K.

Considering higher orders of constrained NMF showed the uniqueness of the high
temperature solid material, and further validated the assertion that smooth,
linear variations of similar patterns are an indication of peak shifting. The
rank-4 NMF extracts similar components as the rank-3, albeit with a new high
temperature component to describe the fine structure observed after the
transition temperature (Fig. 11c–d). This pattern (orange in Figure 11c) was a
result of an amorphous background (later fully pronounced at high temperature
in dark grey) and persistent, strained, crystallites of the original solid
that had a preferred orientation with respect to the X-ray beam. Lastly, the
rank-6 NMF (Fig. 11e–f) decomposed the dataset into the same components with
an additional subdivision of the low temperature solid phase. This subdivision
occurred at the midpoint of the previous two end members, demonstrating the
behavior of the constrained NMF when compensating for pattern shift. This
response to the imperceptible shift in the data is a useful and predictable
benefit to the iterative and constrained approach: without translational
invariance, this behavior enables the user to confidently note when components
belong to the same, but shifting, phase.

The results of the constrained NMF for decomposing a variable temperature
diffraction dataset are immediately more interpretable than a similar
decomposition using cannonical NMF (Fig 9). As with the BaTiO3 example, the
weights correspond to meaningful compositions. Without solid-solid phase
transitions as seen with BaTiO3, the linear response to pattern shifting is
apparent. While this ameliorates a problem with NMF, a more ideal solution
would be one which maps only one component to a single phase regardless of the
translational variations which occur during an experiment. Moreover, the
recognition of persistent solids by the NMF is useful for adjusting an
experimental protocol on the fly. Through user or algorithm intervention, the
identification of a second order transition could be used to slow the rate of
temperature change, allowing the system time to equilibrate. In such
experiments there is often a lag between the measured temperature and the
material temperature, especially around a transition that involves latent
heat. Indicating this transition on-the-fly enables a more robust and
informative experiment.

##  III Conclusions

Here we developed a constrained NMF approach that interprets streaming
spectral data on-the-fly, promptly providing critical insight and analysis
from time-sensitive experiments such as XRD and PDF. This approach efficiently
handles the drastic increase in data generation rates — critical at modern
beamlines — where conventional analysis approaches cannot keep up with the
experimental pace required for optimized data collection. We demonstrated the
advantages of using constrained NMF to aid in revealing subtle phase
transitions of a ferroelectric perovskite and the identification of
experimental complexities in molten salts. This approach leverages the
combined capabilities of the computer and the researcher. Constrained NMF
processes large datasets in an unsupervised fashion with an interactive prior:
the intuition and decision making of the human-in-the-loop. Although applying
constraints may cause the NMF model to create a worse reconstruction of the
dataset, the additional prior information provided by constraints enables for
a model that is more physically meaningful. These constraints can be chosen
explicitly by the researcher, or by a implicitly from the dataset by an
algorithm. This gives the streaming analysis added value that cannot be
expected from canonical NMF. Some overall limitations of NMF that still apply
to constrained approaches are the inability to retain the partial
translational invariance of XRD and PDF patterns, as well as other common
aberrations that impact a pattern without impacting the underlying
phase.Oviedo  _et al._ (2019); Maffettone  _et al._ (2021) Constrained NMF can
be extended to a broad range of spectral applications such as X-ray
photoelectron spectra, X-ray absorption near edge spectra, photoluminescence
spectra, nuclear magnetic resonance spectra, or mass spectra.Gu, Du, and
Billinge (2021) The outputs can be interpreted asynchronously with supervised
ML, and be used to inform adaptive experiments.Allan  _et al._ (2019)
Constrained NMF enables human-in-the-loop, informed, unsupervised learning. It
maintains the agency of the researcher through leveraging of expertise, while
producing real-time analysis and actionable results during an experiment.

##  IV Methods

###  IV.1 NMF

Given a non-negative data matrix X of dimensions $m\times n$ — where $m$ is
the number of independent measurements and $n$ is the length of the feature
dimension — the goal of NMF is to find a factorization such that:

$\textbf{X}\approx\textbf{W}\textbf{H},$ (1)

where W and H are non-negative matrices of dimension $m\times k$ and $k\times
n$, respectively, and $k$ is the user-specified number of components. As such,
the problem reduces to solving:

$\min_{\textbf{W}\in\mathbb{R}^{m\times c},\textbf{H}\in\mathbb{R}^{c\times
n}}D(\textbf{X},\textbf{W}\textbf{H})\quad\text{s.t.}\quad\textbf{W}\geq
0,\textbf{H}\geq 0,$ (2)

where $D$ is a separable, positive measure of fit and the notation
$\textbf{M}\geq 0$ indicates non-negativity of the entries of matrix $M$.

We employ the $\beta$-divergence $D_{\beta}$ as a loss function during
minimization. The $\beta$-divergence is a popular cost function for NMF
minimization parameterized by a single shape parameter $\beta$ (Févotte and
Idier, 2011). For the purposes of this study, we consider two special cases of
the $\beta$-divergence — squared Euclidean distance, also known as the squared
error (SE) ($\beta=2$):

$D_{2}(\textbf{X},\textbf{X}^{\prime})=\frac{1}{2}\sum_{i=1}^{m}\sum_{j=1}^{n}(\textbf{X}_{i,j}-\textbf{X}^{\prime}_{i,j})^{2},$
(3)

and Kullback-Leibler divergence ($\beta=1$):

$D_{1}(\textbf{X},\textbf{X}^{\prime})=\sum_{i=1}^{m}\sum_{j=1}^{n}\left(\textbf{X}_{i,j}\log\frac{\textbf{X}_{i,j}}{\textbf{X}^{\prime}_{i,j}}-\textbf{X}_{i,j}+\textbf{X}^{\prime}_{i,j}\right),$
(4)

though the methods presented here may generalize to any choice of $\beta$.

We employ the alternating non-negative least squares (ANLS) algorithm in order
to minimize the $\beta$-divergence subject to NMF constraints. Briefly, this
algorithm operates by alternately minimizing
$D_{\beta}(\textbf{X},\textbf{WH})$ with respect to W and H under the
constraints that all matrices remain non-negative. This can be solved by
gradient descent with non-negativity constraints on parameter updates, as
described in Algorithm 1:

Algorithm 1 Alternating non-negative least squares

1:procedure ANLS(X, $\textbf{W}_{0}$, $\textbf{H}_{0}$, $\beta$,
$\lambda_{1}$, $\lambda_{2}$, iter, tol)

2: $i\leftarrow 0$, $L_{\text{prev}}\leftarrow\infty$,
$\Delta\leftarrow\infty$

3: $\textbf{W}\leftarrow\textbf{W}_{0}$, $\textbf{H}\leftarrow\textbf{H}_{0}$

4: while $i\leq\text{iter}$ and $\Delta>\text{tol}$ do

5: $\textbf{X}^{\prime}\leftarrow$ Threshold(WH, 0, $10^{-8}$)

6: $L\leftarrow D_{\beta}(\textbf{X},\textbf{X}^{\prime})$

7:
$\mu_{\textbf{W}}\leftarrow\left\\{\begin{matrix}\sum_{j=1}^{n}\textbf{H}_{.,j}&\text{if
}\beta=1\\\
\textbf{X}^{\prime(\beta-1)}\textbf{W}^{\intercal}&\text{otherwise}\end{matrix}\right.\,$

8:
$\textbf{W}\leftarrow\textbf{W}\cdot\max\left(\mu_{\textbf{W}}-\frac{\partial
L}{\partial\textbf{W}},0\right)/(\mu_{\textbf{W}}+\lambda_{1}+\lambda_{2}\textbf{W})$
$\triangleright$ Non-negative W update

9: $\textbf{X}^{\prime}\leftarrow$ Threshold(WH, 0, $10^{-8}$)

10: $L\leftarrow D_{\beta}(\textbf{X},\textbf{X}^{\prime})$

11:
$\mu_{\textbf{H}}\leftarrow\left\\{\begin{matrix}\sum_{i=1}^{m}\textbf{W}_{i,.}&\text{if
}\beta=1\\\
\textbf{W}^{\intercal}\textbf{X}^{\prime(\beta-1)}&\text{otherwise}\end{matrix}\right.\,$

12:
$\textbf{H}\leftarrow\textbf{H}\cdot\max\left(\mu_{\textbf{H}}-\frac{\partial
L}{\partial\textbf{H}},0\right)/(\mu_{\textbf{H}}+\lambda_{1}+\lambda_{2}\textbf{H})$
$\triangleright$ Non-negative H update

13: $i\leftarrow i+1$, $\Delta\leftarrow L-L_{\text{prev}}$,
$L_{\text{prev}}\leftarrow L$

14: end while

15: return W, H

16:end procedure

17:function Threshold($M$, $t$, $\epsilon$)$\triangleright$ Threshold each
element in matrix $M$

18: for all $M_{i,j}$ in $M$ do

19: $M_{i,j}\leftarrow\left\\{\begin{matrix}M_{i,j}&\text{if }M_{i,j}>t\\\
\epsilon&\text{otherwise}\end{matrix}\right.\,$

20: end for

21:end function

To allow for import of prior knowledge, Algorithm 1 allows for the user to
supply initial values for W and H in the form of inputs $\textbf{W}_{0}$ and
$\textbf{H}_{0}$. To fix certain columns of W or rows of H throughout
optimization in order to reflect ground-truth knowledge, one may fix the
gradients for known components to zero in the following manner:

$\frac{\partial
D_{\beta}}{\partial\textbf{M}_{i,j}}=0\quad\forall\quad\textbf{M}_{i,j}\in
S_{\textbf{M}},$ (5)

where $S_{\textbf{M}}$ is the set of all known components of matrix M.

During the analysis of experimental results, we presented the MSE alongside
the mean R-factor

$R=\frac{\sum_{i=1}^{n}{|\textbf{X}-\textbf{WH}|}}{\sum_{i=1}^{n}{|\textbf{X}|}},$
(6)

as a commonly used crystallographic metric.

Algorithm 1 has been implemented in PyTorch Paszke  _et al._ (2019), and is
publicly available at https://github.com/bnl/pub-Maffettone_2021_04.

###  IV.2 Total scattering data collection

Total scattering measurements were performed at the Pair-Distribution Funcion
(PDF) beamline at the National Synchrotron Light Source II at Brookhaven
National Laboratory. The wavelength of incident X-rays for both measurements
was found to be 0.1671 Å (74.21 keV). Data was collected in transmission mode
with a flat-panel Perkin-Elmer detector. Sample-to-detector distance and
detector orientation was calibrated using Ni powder standard. All raw 2D
images were background subtracted using a dark image collected immediately
prior to each scan. Beamstop and detector aberrations were removed with
masking prior to radial integration that was performed using pyFAI,Ashiotis
_et al._ (2015) with PDFs then being generated using PDFgetX3.Juhás  _et al._
(2013)

For the molten salt measurements, the sample-to-detector distance was found to
be 282.74 mm. The sample was held in 1 mm quartz capillaries resistevely
heated using Nichrome wire, with temperatures measured via a K-type
thermocouple. Measurements were collected with 10 second data collection time
per temperature. For the BaTiO3 measurments, samples were held in 1 mm quartz
capillaries, with temperature controlled through a liquid nitrogen cryostream
directed at the sample. Diffraction (XRD) measurements were collected at 1 m
from the sample with 60 second data collection time per temperature.

###### Acknowledgements.

We acknowledge financial support from the BNL Laboratory Directed Research and
Development (LDRD) projects 20-032 “Accelerating materials discovery with
total scattering via machine learning” (P.M.M., D.O.), and the Simons
Foundation Center for Computational Biology (A.D). This research utilized the
PDF (28-ID-1) Beamline and resources of the National Synchrotron Light Source
II, a U.S. Department of Energy (DOE) Office of Science User Facility operated
for the DOE Office of Science by Brookhaven National Laboratory under Contract
No. DE-SC0012704.

##  V Data Availability

The data that support the findings of this study are available from the
corresponding author upon reasonable request.

## References

  * Kusne  _et al._ (2014) A. G. Kusne, T. Gao, A. Mehta, L. Ke, M. C. Nguyen, K.-M. Ho, V. Antropov, C.-Z. Wang, M. J. Kramer, C. Long, and I. Takeuchi, “On-the-fly machine-learning for high-throughput experiments: search for rare-earth-free permanent magnets,” Scientific Reports 4, 6367 EP – (2014).
  * Gómez-Bombarelli  _et al._ (2016) R. Gómez-Bombarelli, J. Aguilera-Iparraguirre, T. D. Hirzel, D. Duvenaud, D. Maclaurin, M. A. Blood-Forsythe, H. S. Chae, M. Einzinger, D.-G. Ha, T. Wu, G. Markopoulos, S. Jeon, H. Kang, H. Miyazaki, M. Numata, S. Kim, W. Huang, S. I. Hong, M. Baldo, R. P. Adams, and A. Aspuru-Guzik, “Design of efficient molecular organic light-emitting diodes by a high-throughput virtual screening and experimental approach,” Nat. Mater. 15, 1120 EP – (2016).
  * Greenaway  _et al._ (2018) R. L. Greenaway, V. Santolini, M. J. Bennison, B. M. Alston, C. J. Pugh, M. A. Little, M. Miklitz, E. G. B. Eden-Rump, R. Clowes, A. Shakil, H. J. Cuthbertson, H. Armstrong, M. E. Briggs, K. E. Jelfs, and A. I. Cooper, “High-throughput discovery of organic cages and catenanes using computational screening fused with robotic synthesis,” Nat. Commun. 9, 2849 (2018).
  * Langner  _et al._ (2020) S. Langner, F. Häse, J. D. Perea, T. Stubhan, J. Hauch, L. M. Roch, T. Heumueller, A. Aspuru-Guzik, and C. J. Brabec, “Beyond ternary opv: High-throughput experimentation and self-driving laboratories optimize multicomponent systems,” Advanced Materials 32, 1907801 (2020).
  * Ludwig (2019) A. Ludwig, “Discovery of new materials using combinatorial synthesis and high-throughput characterization of thin-film materials libraries combined with computational methods,” npj Comput. Mater. 5, 70 (2019).
  * Decker  _et al._ (2017) P. Decker, D. Naujoks, D. Langenkämper, C. Somsen, and A. Ludwig, “High-throughput structural and functional characterization of the thin film materials system Ni-Co-Al,” ACS Comb. Sci. 19, 618–624 (2017).
  * Yeung  _et al._ (2016) H. H.-M. Yeung, Y. Wu, S. Henke, A. K. Cheetham, D. O’Hare, and R. I. Walton, “In-situ observation of successive crystallizations and metastable intermediates in the formation of metal–organic frameworks,” Angewandte Chemie International Edition 55, 2012–2016 (2016).
  * MacLeod  _et al._ (2020) B. P. MacLeod, F. G. L. Parlane, T. D. Morrissey, F. Häse, L. M. Roch, K. E. Dettelbach, R. Moreira, L. P. E. Yunker, M. B. Rooney, J. R. Deeth, V. Lai, G. J. Ng, H. Situ, R. H. Zhang, M. S. Elliott, T. H. Haley, D. J. Dvorak, A. Aspuru-Guzik, J. E. Hein, and C. P. Berlinguette, “Self-driving laboratory for accelerated discovery of thin-film materials,” Science Advances 6 (2020).
  * Li  _et al._ (2020) Z. Li, M. A. Najeeb, L. Alves, A. Z. Sherman, V. Shekar, P. Cruz Parrilla, I. M. Pendleton, W. Wang, P. W. Nega, M. Zeller, J. Schrier, A. J. Norquist, and E. M. Chan, “Robot-accelerated perovskite investigation and discovery,” Chemistry of Materials 32, 5650–5663 (2020).
  * Fleischer  _et al._ (2018) H. Fleischer, D. Baumann, S. Joshi, X. Chu, T. Roddelkopf, M. Klos, and K. Thurow, “Analytical measurements and efficient process generation using a dual–arm robot equipped with electronic pipettes,” Energies 11 (2018).
  * Campbell  _et al._ (2020) S. Campbell, D. B. Allan, A. Barbour, D. Olds, M. Rakitin, R. Smith, and S. B. Wilkins, “Outlook for artificial intelligence and machine learning at the nsls-ii,” Machine Learning: Science and Technology (2020).
  * Olds  _et al._ (2017) D. Olds, P. F. Peterson, M. K. Crawford, J. R. Neilson, H.-W. Wang, P. S. Whitfield, and K. Page, “Combinatorial appraisal of transition states for in situ pair distribution function analysis,” J. Appl. Cryst. 50, 1744–1753 (2017).
  * Hua  _et al._ (2021) X. Hua, A. S. Eggeman, E. Castillo-Martínez, R. Robert, H. S. Geddes, Z. Lu, C. J. Pickard, W. Meng, K. M. Wiaderek, N. Pereira, G. G. Amatucci, P. A. Midgley, K. W. Chapman, U. Steiner, A. L. Goodwin, and C. P. Grey, “Revisiting metal fluorides as lithium-ion battery cathodes,” Nature Materials (2021).
  * Shadike  _et al._ (2021) Z. Shadike, H. Lee, O. Borodin, X. Cao, X. Fan, X. Wang, R. Lin, S.-M. Bak, S. Ghose, K. Xu, C. Wang, J. Liu, J. Xiao, X.-Q. Yang, and E. Hu, “Identification of lih and nanocrystalline lif in the solid–electrolyte interphase of lithium metal anodes,” Nature Nanotechnology (2021).
  * Wang, Li, and Domen (2019) Z. Wang, C. Li, and K. Domen, “Recent developments in heterogeneous photocatalysts for solar-driven overall water splitting,” Chem. Soc. Rev. 48, 2109–2125 (2019).
  * Löffler  _et al._ (2019) T. Löffler, A. Savan, A. Garzón-Manjón, M. Meischein, C. Scheu, A. Ludwig, and W. Schuhmann, “Toward a paradigm shift in electrocatalysis using complex solid solution nanoparticles,” ACS Energy Letters 4, 1206–1214 (2019).
  * Olds (2020) D. Olds, “Synchrotron x-ray diffraction for energy and environmental materials: The current role and future directions of total scattering beamlines in the functional material scientific ecosystem,” Synchrotron Radiation News 33, 4–10 (2020).
  * Clarke  _et al._ (2021) G. Clarke, C. Ablitt, J. Daniels, S. Checchia, and M. S. Senn, “In situ X-ray diffraction investigation of electric-field-induced switching in a hybrid improper ferroelectric,” Journal of Applied Crystallography 54 (2021).
  * Rietveld (1967) H. M. Rietveld, “Line profiles of neutron powder-diffraction peaks for structure refinement,” Acta Cryst. 22, 151–152 (1967).
  * Giacovazzo (2011) C. Giacovazzo, ed.,  _Fundamentals of Crystallography_ , 3rd ed. (Oxford University Press, 2011).
  * Iwasaki, Kusne, and Takeuchi (2017) Y. Iwasaki, A. G. Kusne, and I. Takeuchi, “Comparison of dissimilarity measures for cluster analysis of x-ray diffraction data from combinatorial libraries,” npj Comput. Mater. 3, 1–9 (2017).
  * Xiong  _et al._ (2017) Z. Xiong, Y. He, J. R. Hattrick-Simpers, and J. Hu, “Automated phase segmentation for large-scale x-ray diffraction data using a graph-based phase segmentation (gphase) algorithm,” ACS Comb. Sci. 19, 137–144 (2017).
  * Long  _et al._ (2009) C. J. Long, D. Bunker, X. Li, V. L. Karen, and I. Takeuchi, “Rapid identification of structural phases in combinatorial thin-film libraries using x-ray diffraction and non-negative matrix factorization,” Rev. Sci. Instrum. 80, 103902 (2009).
  * Oviedo  _et al._ (2019) F. Oviedo, Z. Ren, S. Sun, C. Settens, Z. Liu, N. T. P. Hartono, S. Ramasamy, B. L. DeCost, S. I. P. Tian, and G. Romano, “Fast and interpretable classification of small x-ray diffraction datasets using data augmentation and deep neural networks,” npj Comput. Mater. 5, 1–9 (2019).
  * Lee  _et al._ (2020) J.-W. Lee, W. B. Park, J. H. Lee, S. P. Singh, and K.-S. Sohn, “A deep-learning technique for phase identification in multiphase inorganic compounds using synthetic xrd powder patterns,” Nat. Commun. 11, 86 (2020).
  * Ziletti  _et al._ (2018) A. Ziletti, D. Kumar, M. Scheffler, and L. M. Ghiringhelli, “Insightful classification of crystal structures using deep learning,” Nat. Commun. 9, 2775 (2018).
  * Aguiar  _et al._ (2019) J. A. Aguiar, M. L. Gong, R. R. Unocic, T. Tasdizen, and B. D. Miller, “Decoding crystallography from high-resolution electron imaging and diffraction datasets with deep learning,” Sci. Adv. 5 (2019).
  * Di Chen  _et al._ (2020) Di Chen, Y. Bai, W. Zhao, S. Ament, J. M. Gregoire, and C. P. Gomes, “Deep reasoning networks for unsupervised pattern de-mixing with constraint reasoning,” in  _Proceedings of the 37th International Conference on Machine Learning_ , Proceedings of Machine Learning Research, Vol. 119, edited by F. Bach and D. Blei (PMLR, 2020).
  * Takeuchi  _et al._ (2005) I. Takeuchi, C. J. Long, O. O. Famodu, M. Murakami, J. Hattrick-Simpers, G. W. Rubloff, M. Stukowski, and K. Rajan, “Data management and visualization of x-ray diffraction spectra from thin film ternary composition spreads,” Rev. Sci. Instrum. 76, 062223 (2005).
  * Stanev  _et al._ (2018) V. Stanev, V. V. Vesselinov, A. G. Kusne, G. Antoszewski, I. Takeuchi, and B. S. Alexandrov, “Unsupervised phase mapping of x-ray diffraction data by nonnegative matrix factorization integrated with custom clustering,” npj Computational Materials 4, 43 (2018).
  * Bai  _et al._ (2017) J. Bai, J. Bjorck, Y. Xue, S. K. Suram, J. Gregoire, and C. Gomes, “Relaxation methods for constrained matrix factorization problems: Solving the phase mapping problem in materials discovery,” in  _Integration of AI and OR Techniques in Constraint Programming_ , edited by D. Salvagnin and M. Lombardi (Springer International Publishing, Cham, 2017) pp. 104–112.
  * Bai  _et al._ (2018) J. Bai, Y. Xue, J. Bjorck, R. Le Bras, B. Rappazzo, R. Bernstein, S. K. Suram, R. B. van Dover, J. M. Gregoire, and C. P. Gomes, “Phase mapper: Accelerating materials discovery with ai,” AI Magazine 39, 15–26 (2018).
  * Suram  _et al._ (2017) S. K. Suram, Y. Xue, J. Bai, R. Le Bras, B. Rappazzo, R. Bernstein, J. Bjorck, L. Zhou, R. B. van Dover, C. P. Gomes, and J. M. Gregoire, “Automated phase mapping with agilefd and its application to light absorber discovery in the v–mn–nb oxide system,” ACS Combinatorial Science 19, 37–46 (2017).
  * Geddes  _et al._ (2019) H. S. Geddes, H. Blade, J. F. McCabe, L. P. Hughes, and A. L. Goodwin, “Structural characterisation of amorphous solid dispersions via metropolis matrix factorisation of pair distribution function data,” Chem. Commun. 55, 13346–13349 (2019).
  * Li  _et al._ (2021) Q.-J. Li, D. Sprouster, G. Zheng, J. C. Neuefeind, A. D. Braatz, J. Mcfarlane, D. Olds, S. Lam, J. Li, and B. Khaykovich, “Complex structure of molten nacl–crcl3 salt: Cr–cl octahedral network and intermediate-range order,” ACS Applied Energy Materials (2021).
  * Paszke  _et al._ (2019) A. Paszke, S. Gross, F. Massa, A. Lerer, J. Bradbury, G. Chanan, T. Killeen, Z. Lin, N. Gimelshein, L. Antiga, A. Desmaison, A. Kopf, E. Yang, Z. DeVito, M. Raison, A. Tejani, S. Chilamkurthy, B. Steiner, L. Fang, J. Bai, and S. Chintala, “Pytorch: An imperative style, high-performance deep learning library,” in  _Advances in Neural Information Processing Systems 32_ , edited by H. Wallach, H. Larochelle, A. Beygelzimer, F. d’Alché Buc, E. Fox, and R. Garnett (Curran Associates, Inc., 2019) pp. 8024–8035.
  * Page  _et al._ (2010) K. Page, T. Proffen, M. Niederberger, and R. Seshadri, “Probing local dipoles and ligand structure in batio3 nanoparticles,” Chemistry of Materials 22, 4386–4391 (2010).
  * Wegner  _et al._ (2020) M. Wegner, H. Gu, R. D. James, and E. Quandt, “Correlation between phase compatibility and efficient energy conversion in zr-doped barium titanate,” Scientific Reports 10, 3496 (2020).
  * Gu, Du, and Billinge (2021) R. Gu, Q. Du, and S. J. L. Billinge, “A fast two-stage algorithm for non-negative matrix factorization in streaming data,” (2021).
  * Maffettone  _et al._ (2021) P. M. Maffettone, L. Banko, P. Cui, Y. Lysogorskiy, M. A. Little, D. Olds, A. Ludwig, and A. I. Cooper, “Crystallography companion agent for high-throughput materials discovery,” (2021).
  * Allan  _et al._ (2019) D. Allan, T. Caswell, S. Campbell, and M. Rakitin, “Bluesky’s ahead: A multi-facility collaboration for an a la carte software project for data acquisition and management,” Synchrotron Radiation News 32, 19–22 (2019).
  * Thorndike (1953) R. L. Thorndike, “Who belongs in the family?” Psychometrika 18, 267–276 (1953).
  * Févotte and Idier (2011) C. Févotte and J. Idier, “Algorithms for nonnegative matrix factorization with the beta-divergence,” Neural Computation 23, 2421–2456 (2011).
  * Ashiotis  _et al._ (2015) G. Ashiotis, A. Deschildre, Z. Nawaz, J. P. Wright, D. Karkoulis, F. E. Picca, and J. Kieffer, “The fast azimuthal integration python library: pyfai,” Journal of applied crystallography 48, 510–519 (2015).
  * Juhás  _et al._ (2013) P. Juhás, T. Davis, C. L. Farrow, and S. J. Billinge, “Pdfgetx3: a rapid and highly automatable program for processing powder diffraction data into total scattering pair distribution functions,” Journal of Applied Crystallography 46, 560–566 (2013).