Introduction

Quantitative assessment of tissue properties from high-dimensional MRI data requires signal decomposition or appropriate parametrization^1,2. Imposing physiological models in data analysis carries the risk of estimating biased parameters. On the other hand, using standard decomposition strategies without presumptive physiological models may result in physically unrealistic components difficult to interpret. This work proposes an extension of the non-negative matrix factorization (NMF), introducing additional constraints on the decomposed signals tailored to MR signals, expected to be of multi-exponential or similar character: the latent signals must be monotonous (inspired by Bhatt and Ayyar³) and have a monotonous first derivative (slope). We present a demonstration of the novel method, referred to as monotonous slope non-negative matrix factorization (msNMF).

Methods

The standard NMF approximates a data matrix, $$$\textbf{X}\in\mathbb{R}^{m\times n}$$$, with a low rank matrix consisting of two factors: a set of non-negative basis functions, $$$\textbf{W}\in\mathbb{R}_{+}^{m\times k}$$$ (signal components), and the associated non-negative weights, $$$\textbf{H} \in \mathbb{R}_{+}^{k \times n}$$$ (mixture maps):
$$\textbf{X} = \textbf{W}\textbf{H}+\textbf{E}$$
The residual matrix ($$$\textbf{E}\in\mathbb{R}^{m\times n}$$$) accounts for noise and non-factorizable signals. The problem was split into two sub-problems to alternately optimize for one of the two factors using Frobenius norm. In the following equations, the sub-problems have been written via standard quadratic programming³ with the new monotonicity constraints formulated as inequality terms:
$$\text{Given }\textbf{W}:\quad\min_\textbf{H}\text{vec}(\textbf{H})^\mathrm{T}(\textbf{I}_n\otimes \textbf{W}^\mathrm{T}\textbf{W})\text{vec}(\textbf{H})-2 \text{vec}(\textbf{X}^\mathrm{T}\textbf{W})^\mathrm{T}\text{vec}(\textbf{H})\quad\text{subject to: }-\textbf{I}_{nk}\text{vec}(\textbf{H})\leq\textbf{0}_{nk}$$
$$\text{Given } \textbf{H}:\quad\min_\textbf{W}\text{vec}(\textbf{W})^\mathrm{T}(\textbf{HH}^\mathrm{T}\otimes \textbf{I}_m)\text{vec}(\textbf{W})-2\text{vec}(\textbf{XH}^\mathrm{T})^\mathrm{T}\text{vec}(\textbf{W}) \quad\text{subject to: }\begin{bmatrix}-\textbf{I}_{mk}\\ \textbf{A}\\\textbf{B}\end{bmatrix}\text{vec}(\textbf{W})\leq\textbf{0}_{mk}$$
Here, $$$\text{vec()}$$$ is the column vectorization of a matrix and $$$\otimes$$$ is the Kronecker tensor product. $$$\textbf{A}$$$ and $$$\textbf{B}$$$ are designed to ensure monotonous behavior of the signal components and their first derivative, respectively.

Preclinical evaluation
The DWI data (published by Lundell et al⁴) consists of multi-dimensional diffusion encoding measurements of a post mortem brain of a Vervet monkey in buffer solution. Data was recorded with a 4.7 T pre-clinical MRI-system using TE=68 ms and 12 b-values (range 242-4832 s/mm²). Two acquisitions with a difference in the implemented gradient waveform generated data with an average diffusion time of ~10 ms (short diffusion time) and ~20 ms (long diffusion time), respectively. Comparisons of these should capture the effect of time-dependent diffusion, and hence be reflected in the decomposition.

The relaxometry data (Harkins et al⁵) was originally used to demonstrate intra-myelinic edema (IME) quantification in rat spinal white matter. It consists of multi-echo T₂-weighted measurements (40 echoes with 9 ms separation) recorded at 9.4 T for three groups of eight rats, each receiving a different amount of IME-inducing toxin in the diet: 0 ppm, 300 ppm and 600 ppm hexachlorophene (HCP). The data was decomposed into three water components to allow comparison to the reference analysis by Harkins et al and the accompanying histology values.

Results and Discussion

Three characteristic signal components were extracted from the DWI measurements (Figure 1b). Signals from the brain tissue has primarily been decomposed into two components (green and purple), where a different composition of these distinguishes gray from white matter (Figure 1c). The third short-lived component (yellow) is distributed as the surrounding watery buffer solution and CSF leftovers. The difference between the two sets of maps (Figure 1c, third column), provided to capture the effect of time-dependent diffusion, uncovers areas of high cell density. Tissues of limited diffusion appear strongly positive for the long-lasting signal component (green) and negative for the intermediate signal decay (purple). This involves white matter, visual cortex (thin arrow), and in particular the cerebellar cortex (thick arrow), which in accordance with the pronounced contrast consists of a very restricted diffusion environment^4,6,7. Such clear contrasts across cerebral and cerebellar cortex have to our knowledge not been found using other MR contrast mechanisms before. This makes the msNMF applied to DWI data a promising and potentially powerful tool for tumor tissue characterization.

The decomposition of the relaxometry data is presented in Figure 2. For all three components, the contrast between gray and white matter increases for the IME-induced rats. Consistent with expectations of relaxation properties (IME causes T₂-hyperintensity^5,8), the fraction of the green component (given by the normalized $$$\textbf{H}$$$) is considered an estimate of the EWF. The msNMF-based EWF estimates in four white matter tracks (Figure 3a) confirm a clear increase of edema with increasing HCP intoxication (Figure 3b). There is a high degree of similarity to the equivalent barplot presented by Harkins et al (Figure 3c), although here EWF estimates are slightly higher considering the HCP-groups. In a comparison of the EWF estimates with the histology-based edema water content fraction, W_ed, (Figure 4a), a high agreement (R² of 0.88, P<0.001) confirms that the decomposition is related to the underlying compartmentalization of the tissue. Similar to the results reported by Harkins et al (Figure 4b), an expected small overestimation of EWF compared to W_ed was observed⁵.

Conclusion

The novel decomposition method (msNMF) separated MRI signals into identifiable components specific to the underlying structure, including imitated pathological conditions, and produced results consistent with those published. The msNMF is insensitive to partial volume effects and potentially applicable to a range of multi-dimensional MRI data. The demonstrated contrast across areas with different cell body density, for example, qualifies it as a candidate for characterization of malignant lesions for better staging, prognostics and management.

Acknowledgements

The authors would like to thank Mark D. Does and Kevin D. Harkins for their cooperativeness and sharing of data, and Nirav Bhatt and Arun Ayyar for providing their software for monotonous NMF. The study was funded by the Danish Cancer Society (grant R167-A10637-17-S2) and Henrik Lundell additionally received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No 804746).