# How to model compositional data?

What is the best way to model compositional data problems?

Compositional data is when each example or sample is a vector that sums to 1 (or 100%). In my case, I am interested in the composition of minerals in a rock and I have sensors that tell me the sum of the minerals but not the components that make up the sum.

For example, lets say I have two minerals, $m_1$ and $m_2$, that are made up of 3 elements (like copper and other elements from the periodic table) which form a vector of length 3:

m1 = [0.1, 0.3, 0.6]
m2 = [0.6, 0.2, 0.2]


If a rock has 25% of $m_1$ and 75% of $m_2$, the sensor reading produces the sum of the two minerals (shown in bottom-left subplot below):

\begin{align} &0.25*m_1 + 0.75*m_2 \\ =&0.25*[0.1, 0.3, 0.6] + 0.75*[0.6, 0.2, 0.2] \\ =&[0.475, 0.225, 0.3] \end{align}

I would like to know how to model and solve the problem of unmixing a composition into its underlying components, where the sum of the elements is normalized to 100% (e.g. $0.25m_1 + 0.75m_2$ has the same composition as $0.50m_1 + 1.50m_2$).

Furthermore, my example is simplistic; in reality a composition can have more than just 2 minerals (up to 3000) and each mineral is made up of 118 elements, not just 3 (all the elements of periodic table - though many elements will be zero). The elemental composition of a mineral is assumed to be known (definition of $m_1$ and $m_2$ in the example). Also, the sensor reading is noisy - each element of the observed composition is assumed to have Gaussian noise.

First normalize the result vector. E.g. [.95, .45, .6] by dividing by 2 (sum of the members); giving [.475, .225, .3].

Let $x$ be the share of the first mineral, than $(x-1)$ is the share of the second mineral.

Solve the three linear equation which must give the same result.

$.1 * x + .6 * (1-x) = .475$

$.3 * x + .2 * (1-x) = .225$

$.6 * x + .2 * (1-x) = .3$

Result is $x = 1/4$ as expected.

UPDATE

The above proposed solution works of course only if the number of minerals is less or equal the number of elements. In the update of the question is clearly stated, that this is not the case (3000 minerals and 118 elements).

Lets simulate this oposite case on a small example with 3 minarals and 2 elements.

m1 <- c(0.2, 0.8)
m2 <- c(0.4, 0.6)
m3 <- c(0.9, 0.1)


and with the mix of minerals

 x <- c(.25, .65, .1)


which produce a measurement of

 t(matrix(c(m1,m2,m3),3,2, byrow= T))  %*% x

[,1]
[1,]  0.4
[2,]  0.6


This gives following linear equations

$m_1 + m_2 + m_3 = 1$

$.2 m_1 + .4 m_2 + .9 m_3 = .4$

$.8 m_1 + .6 m_2 + .1 m_3 = .6$

The solution of the equations is not unique

$m_3 \in <0, 1 / 3.5>$

$m_1 = 2 - 2 * m_2 - 4.5 * m_3$

$m_2 = 1 - 3.5 * m_3$

Some alternative solutions are provided below

  0 1 0
0.125 0.825 0.050
0.375 0.475 0.150
0.5 0.3 0.2
0.625 0.125 0.250
0.71428571 0.00000000 0.28571429


This is of course oversimplified example, but shows that you should carefully select the optimization goal, as there could be more "equaly good" solutions. For example with promoting sparsity you will find the solution [0 1 0] which is far from our used mix.

• Where did [.95, .45, .6] come from? The way I defined the problem, the result is already normalized. Your solution works for the example I gave - thanks for that - but it does not generalize to when there are more than 2 minerals (I revised my question to make this more clear) Nov 10 '15 at 23:32
• This is the result of 0.50m1+1.50m2. And yes, I know it can't be what you are looking for; I post it only to sharp your question. Good luck! Nov 11 '15 at 10:44

If I understand correctly, in your problem the mineral-element vector $m1$ and $m2$ is known.

If that is the case, I would suggest LASSO regression. Suppose you have $n$ minerals $m_1$, $m_2$, ..., $m_n\in\mathbb{R}^p$ arranged into a matrix $M=[m_1;m_2,...,m_n]\in\mathbb{R}^{n\times p}$, and an unknown mixing vector $x\in\mathbb{R}^n$ (in your case, $x=[0.25; 0.75]$). Your sensor measures $y=M^\top x\in\mathbb{R}^n$ and you want to recover $x$ from $y$.

The LASSO regression aims to solve $\hat{x}=\arg\min_x ||y-M^\top x||_2 + \lambda||x||_1$. The $\ell_1$ penalty promotes sparsity in $x$ and could be beneficial when you have a lot of candidate minerals.

I don't think you need any special handling due to the compositional nature of the data, other than normalizing the final $\hat{x}$ if needed.

• This seems reasonable, I will give this a try. However my feeling is that it may find solutions that are impossible, such as a negative amount of a mineral, and also that I may get better performance from a model that is built to constrain the sum to 1 because the search space will be smaller. Nov 10 '15 at 23:41