Can we decompose an array of algebraic expressions into a product of matrices that separate known variables from unknowns?

Question

Given vectors $\sigma, \tau \in \Bbb R^n \ge0$, and constant $\mu \in \Bbb R^1>0$, let matrix-valued function $C : \Bbb R \to \Bbb R^{2 \times n}$ be defined by $$C (\mu) := \begin{bmatrix} \vec\sigma+\mu(\vec\tau-\vec\sigma) \\\vec\tau-\vec\sigma-\mu(\vec\tau-\vec\sigma) \end{bmatrix}=\begin{bmatrix} \sigma_1+\mu(\tau_1-\sigma_1) & \sigma_2+\mu(\tau_2-\sigma_2) & \cdots & \sigma_n+\mu(\tau_n-\sigma_n)\\ \tau_1-\sigma_1-\mu(\tau_1-\sigma_1) & \tau_2-\sigma_2-\mu(\tau_2-\sigma_2) & \cdots & \tau_n-\sigma_n-\mu(\tau_n-\sigma_n) \end{bmatrix}$$

Can we decompose $C$ into the product of two matrices that separate the unknown $μ$ from the knowns $σ$ and $τ$? Something along the lines of $$C (\mu) = A(\mu) B_{\sigma, \tau}$$

More generally, how does one handle arrays of expressions in problems like this?

$\cdots\cdots$

Context

You might be wondering where this bizarre construction of a matrix comes from. Please bear with me as I try to explain this conflagration of analytical chemistry and linear algebra and why we're interested in decomposing this matrix in this way (if that is even the correct terminology). I should note that my background is analytical chemistry with an emphasis on applications of Raman spectroscopy. Years ago I had an undergraduate linear-algebra course, and later on encountered its application to chemometrics quantitative analysis of spectra as a graduate student. I use Matlab for data analysis, which is based around matrix representations of data and matrix math.

Background

In spectroscopy data analysis, we often have a set of sample spectra, $S \in\Bbb R \to \Bbb R^{m \times n}\ge 0$, that were collected on an instrument where $m$ is the spectral pixel index, e.g. wavelength, frequency, etc and $n$ is the sample index.

For example, here are a set of 100 Raman spectra at ten discrete concentrations: Example spectra. For this set, the Raman shift axis (the $m$-dimension) ranges from 0 to 3700$\,cm^{-1}$ with 1024 gradations (pixels) and there are 100 spectra so the dimensions of $S$ are 1024$\times$100. The first ten spectra are colored blue, the last ten red, and everything in between along the 'jet' color scheme from blue to green to yellow to orange and red.

In many cases, Raman spectra are additive, which is to say that the spectrum of a mixture of gases, for example, is equivalent to the sum of its parts. Thus our set of 100 spectra can be represented with linear algebra as the product of a set of pure component spectra, $P: \Bbb R \to \Bbb R^{m \times p}$ and a matrix of concentrations of each of the $p$ pure components in the $n$ samples, $C: \Bbb R \to \Bbb R^{p \times n}\ge0$: $$S^{m\times n} =P^{m\times p}*C^{p\times n}$$

This expression of the superposition principal is like the linear-algebra form of Beer's Law: the counts in a spectrum are proportional to the concentration and in a mixture, the total counts are proportional to the sum of the concentrations of the components.

In Classical Least Squares analysis, we might use this relationship on a set of spectra from mixtures of known composition to get the least-squares estimate of the pure component spectra using the right-pseudo-inverse of $C$:$$S*[C^\top*(C*C^\top)^{-1}]=P*C*[C^\top*(C*C^\top)^{-1}]$$ $$S*C^\top*(C*C^\top)^{-1}=P*\begin{bmatrix} 1&0&\cdots0\\ 0&1&\cdots0\\ \vdots&\vdots&\ddots0\\ 0&0&\cdots1\\ \end{bmatrix}^{p\times p}=P $$ We could then in turn apply this 'calibration' to a spectrum from an unknown mixture, $U^{m\times 1}$ to estimate its composition, $\Gamma^{p\times 1}$ : $$(P^\top P)^{-1}*P^\top*U = \Gamma$$

The Big Question

Consider a scenario where we don't know $P$ and we don't have all of the information about $C$: We have a stock mixture with unknown amounts of two components. We were able to spike that stock mixture with known amounts of one of the components to generate new mixtures from which we collected spectra (i.e. the method of standard additions). Alas, we neglected to measure the spectrum from the pure component and ran out of it! (Hmmm. This scenario is embarrassingly specific...) Nevertheless, from this data set Can we estimate the composition of the stock mixture,$\mu$, and obtain the two pure-component spectra $P^{m \times 2}$?

We might consider that the concentration matrix, $[C]$, can be represented as an array of algebraic expressions containing what we know about the compositions of the samples we generated and measured spectra from. We introduce a variable $\mu$ which is the unknown fraction of the spiked component in the stock mixture. Also we introduce $\sigma_n$, the concentration of the pure component added to the stock mixture to generate a given sample, and finally, $\tau_n$, the total concentration of both components in a given sample. The concentration of the spiked component in a given sample is the sum of the spike $\sigma_n$ and the contribution from the stock, $\mu*(\tau_n-\sigma_n)$. The concentration of the second component is simply $\tau_n$ minus the concentration of the first component, which leads to the unconventional array of algebraic expressions above.

Where we're going with this

Having decomposed $C$ into the product of matrices of unknowns, $A(\mu)^{2\times q}$, separated from knowns, $B_{\sigma, \tau}\,^{q \times n}$, we will use again the right-pseudo inverse to move what we know to the left side of the matrix equation: $$S*B^+=P*A$$ I suspect that doing so will enable us to obtain the least-squares estimates of the pure component spectra and the composition of the stock mixture.

Many thanks to @Rodrigodeazevedo, who first answered this question and patiently worked with me to clear up the notation

Answer 1

Let ${\rm a} := \begin{bmatrix} 1 \\ -1\end{bmatrix}$ and ${\rm b} := \begin{bmatrix} 0 \\ 1\end{bmatrix}$. Hence,

$$\rm C (\mu) := a \sigma^\top + b \tau^\top + \mu a (\tau - \sigma)^\top = \begin{bmatrix} | & | & |\\ \mathrm a & \mathrm b & \mathrm a\\ | & | & |\end{bmatrix} \begin{bmatrix} 1 & & \\ & 1 & \\ & & \mu\end{bmatrix} \begin{bmatrix} | & | & |\\ \sigma & \tau & \tau - \sigma\\ | & | & |\end{bmatrix}^\top$$

---

>>> from sympy import *
>>> sigma1, sigma2, sigma3 = symbols('sigma1, sigma2, sigma3')
>>> tau1, tau2, tau3 = symbols('tau1, tau2, tau3')
>>> mu = symbols('mu')
>>> a = Matrix([1,-1])
>>> b = Matrix([0,1])
>>> latex( a * Matrix([sigma1, sigma2, sigma3]).T + b * Matrix([tau1, tau2, tau3]).T + mu * a * Matrix([tau1 - sigma1, tau2 - sigma2, tau3 - sigma3]).T )

produces

$$\left[\begin{matrix}\mu \left(- \sigma_{1} + \tau_{1}\right) + \sigma_{1} & \mu \left(- \sigma_{2} + \tau_{2}\right) + \sigma_{2} & \mu \left(- \sigma_{3} + \tau_{3}\right) + \sigma_{3}\\- \mu \left(- \sigma_{1} + \tau_{1}\right) - \sigma_{1} + \tau_{1} & - \mu \left(- \sigma_{2} + \tau_{2}\right) - \sigma_{2} + \tau_{2} & - \mu \left(- \sigma_{3} + \tau_{3}\right) - \sigma_{3} + \tau_{3}\end{matrix}\right]$$