I have a large dataset coming from absorption spectroscopy experiment and i'm facing a problem to extract relevant informations. My matrix look like this :
$ S_n(x) = \left( \frac{a_n}{\sigma_1\sqrt{2\pi}} e^{\frac{-(x-\mu_1)^2}{2\sigma_1^2}} -\frac{a_n}{\sigma_1\sqrt{2\pi}} e^{\frac{-(x-\mu'_1)^2}{2\sigma_1^2}}\right) + \left( \frac{b_n}{\sigma_2\sqrt{2\pi}} e^{\frac{-(x-\mu_2)^2}{2\sigma_2^2}} -\frac{b_n}{\sigma_2\sqrt{2\pi}} e^{\frac{-(x-\mu'_2)^2}{2\sigma_2^2}}\right) + \left( \frac{c_n}{\sigma_3\sqrt{2\pi}} e^{\frac{-(x-\mu_3)^2}{2\sigma_3^2}} -\frac{c_n}{\sigma_3\sqrt{2\pi}} e^{\frac{-(x-\mu'_3)^2}{2\sigma_3^2}}\right) $
Where $S_n(x)$ is the one-dimensionnal spectra (and it's a know value).
$n$ value is at least 10 but it can reach 50 for example and I know the value of $a_n, b_n$ and $c_n$. I'm looking for the value of $\sigma_1, \sigma_2$ and $\sigma_3$ and for all the $\mu$ value (1, 2, 3 and ' value).
I have a very basic level in analysis or algebra but is there any way to find those parameters ?
Thank you in advance for your help !
Wow, that's a very interesting model! I can think of two approaches to your problem. One is simply to use the Solver routine in Excel. The other is more theoretical (at first), but might yield a good solution.
First we can simplify your model expression a tad: \begin{align*} S_n(x) =& \frac{a_n}{\sigma_1\sqrt{2\pi}}\left[ \exp\left(\frac{-(x-\mu_1)^2}{2\sigma_1^2}\right) - \exp\left(\frac{-(x-\mu'_1)^2}{2\sigma_1^2}\right)\right]\\ + &\frac{b_n}{\sigma_2\sqrt{2\pi}}\left[ \exp\left(\frac{-(x-\mu_2)^2}{2\sigma_2^2}\right) - \exp\left(\frac{-(x-\mu'_2)^2}{2\sigma_2^2}\right)\right]\\ + &\frac{c_n}{\sigma_3\sqrt{2\pi}}\left[ \exp\left(\frac{-(x-\mu_3)^2}{2\sigma_3^2}\right) - \exp\left(\frac{-(x-\mu'_3)^2}{2\sigma_3^2}\right)\right]. \end{align*} (With exponentials that complicated, it's usually better typesetting practice to use the
expfunction.)Now, for each $n,$ we have a series of data points $\{x_j, S_{nj}\}$ to which we wish to fit the model. We'll try least-squares regression (this is not a linear model, because the coefficients we're trying to find don't show up linearly in the model) with a cost function $$C_n=\sum_{j=1}^N(S_n(x_j)-S_{nj})^2. $$ Here $S_n(x_j)$ refers to the model evaluated at the $x_j$ data point, and $S_{nj}$ refers to the known (measured) $S$ value of the $n$th spectrum at $x_j$.
To solve for $\{\sigma_i,\mu_i, \mu_i'\},\; 1\le i\le 3,$ the usual procedure is to take the partial derivatives of $C_n$ and set them equal to zero. This will be a system of nine equations. Now, we will need to assemble these complicated ingredients one step at a time. Here are some intermediate calculations: \begin{align*} \frac{\partial S_n(x_j)}{\partial \mu_1}&=\frac{a_n}{\sigma_1\sqrt{2\pi}}\,\exp\left(\frac{-(x-\mu_1)^2}{2\sigma_1^2}\right)\left(\frac{2(x-\mu_1)}{2\sigma_1^2}\right)=\frac{a_n(x-\mu_1)}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x-\mu_1)^2}{2\sigma_1^2}\right) \\ \frac{\partial S_n(x_j)}{\partial \mu_1'}&=-\frac{a_n(x-\mu_1')}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x-\mu_1')^2}{2\sigma_1^2}\right). \end{align*} The other partials w.r.t. the means $\mu_2, \mu_2', \mu_3,$ and $\mu_3'$ will be analogous. The standard deviation partials are much more complicated. First, let's absorb non-$\sigma$ constants in one term together, like this: $$_{1}S_{n}(x)=\frac{_1A_n}{\sigma_1}\left[ \exp\left(-\frac{B_1}{\sigma_1^2}\right) - \exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right],$$ where $_1A_n=a_n/\sqrt{2\pi}, \; B_1=(x-\mu_1)^2/2,$ and $B_1'=(x-\mu_1')^2/2.$ I'm using the notation $_1S_n(x)$ to denote the first term in $S_n(x)$. Now we calculate the partial derivative thus: \begin{align*} \frac{\partial (_1S_n(x))}{\partial \sigma_1}&=-\frac{_1A_n}{\sigma_1^2}\left[ \exp\left(-\frac{B_1}{\sigma_1^2}\right) - \exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] +\frac{_1A_n}{\sigma_1}\left[ \left(\frac{2B_1}{\sigma_1^3}\right)\exp\left(-\frac{B_1}{\sigma_1^2}\right) - \left(\frac{2B_1'}{\sigma_1^3}\right)\exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] \\ &=-\frac{\sigma_1^2\,_1A_n}{\sigma_1^4}\left[ \exp\left(-\frac{B_1}{\sigma_1^2}\right) - \exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] +\frac{2\,_1A_n}{\sigma_1^4}\left[ B_1\exp\left(-\frac{B_1}{\sigma_1^2}\right) - B_1'\exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] \\ &=\frac{_1A_n}{\sigma_1^4} \left[(2B_1-\sigma_1^2)\exp\left(-\frac{B_1}{\sigma_1^2}\right)-(2B_1'-\sigma_1^2)\exp\left(-\frac{B_1'}{\sigma_1^2}\right)\right] \\ &=\frac{a_n}{\sigma_1^4\sqrt{2\pi}} \left[((x-\mu_1)^2-\sigma_1^2)\exp\left(-\frac{(x-\mu_1)^2}{2\sigma_1^2}\right)-((x-\mu_1')^2-\sigma_1^2)\exp\left(-\frac{(x-\mu_1')^2}{2\sigma_1^2}\right)\right]. \end{align*} A more standard (pun intended, of course) way to write this would probably (pun intended) be: $$\frac{\partial (_1S_n(x))}{\partial \sigma_1}= \frac{a_n}{\sigma_1^2\sqrt{2\pi}} \left[\left(\frac{(x-\mu_1)^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x-\mu_1)^2}{2\sigma_1^2}\right)-\left(\frac{(x-\mu_1')^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x-\mu_1')^2}{2\sigma_1^2}\right)\right] $$ Whew! And those were only the intermediate calculations! What we really need to calculate is $\partial C_n/\partial\mu_1,$ and so on. Thankfully, we've already done most of the hard work. We have \begin{align*} \frac{\partial C_n}{\partial\mu_1}&=2\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})\frac{a_n(x_j-\mu_1)}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x_j-\mu_1)^2}{2\sigma_1^2}\right)\right] \\ \frac{\partial C_n}{\partial\mu_1'}&=-2\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})\frac{a_n(x_j-\mu_1')}{\sigma_1^3\sqrt{2\pi}}\,\exp\left(\frac{-(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right] \\ \frac{\partial C_n}{\partial\sigma_1}&=2\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})\frac{a_n}{\sigma_1^2\sqrt{2\pi}} \left[\left(\frac{(x_j-\mu_1)^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1)^2}{2\sigma_1^2}\right)-\left(\frac{(x_j-\mu_1')^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right]\right]. \end{align*} The formatting there is off, because MathJax on M.SE only has so much width it can use. Thankfully, these expressions simplify a bit when we set them equal to zero: \begin{align*} 0&=\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})(x_j-\mu_1)\,\exp\left(\frac{-(x_j-\mu_1)^2}{2\sigma_1^2}\right)\right] \\ 0&=\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj})(x_j-\mu_1')\,\exp\left(\frac{-(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right] \\ 0&=\sum_{j=1}^{N}\left[(S_n(x_j)-S_{nj}) \left[\left(\frac{(x_j-\mu_1)^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1)^2}{2\sigma_1^2}\right)-\left(\frac{(x_j-\mu_1')^2}{\sigma_1^2}-1\right)\exp\left(-\frac{(x_j-\mu_1')^2}{2\sigma_1^2}\right)\right]\right]. \end{align*} Plus, of course, we have the six other analogous equations.
These are certainly complicated formulae, but they're no insurmountable. It would be possible to code it all up in Python, for example, and use, say, a gradient descent algorithm to solve the system. See Andrew Ng's Machine Learning course for more information.