# Bayesian inference on a sum of iid random variables with known distribution

Let $X_1$, $X_2$, ..., $X_n$ be iid RV's following a mixture distribution of two lognormals such that the pdf of each $X_i$ is $f_{mix}(x)=pf_1(x) + (1-p)f_2(x)$ where $f_1(x)$ and $f_2(x)$ are lognormal pdfs with parameters $\mu_1,\sigma$ and $\mu_2,\sigma$, respectively.

Define $S_i$ as a sum of 10 $X$s, e.g. $S_1 = X_1 +..+ X_{10}, \ \ S_2 = X_{11} +..+X_{20},...$

I am given only $S_1,S_2,\dots,S_{n/10}$.

How can I infer the mixing proportion parameter $p$ here? That means, I want to know what the proportion of the two lognormals in my mixture is among the 10 samples is (but I only have the sum measurement).

Note that the sum of lognormals is not lognormally distributed.

My problem could be related to this one, but I am not sure: Bayesian inference on a sum of iid real-valued random variables

update: assume that I have many samples of $S_i$. If that helps we can assume that the $\mu$s are well separated, making the mixture nicely bimodal. Generally in mixtures it is that the mixing proportions sum to 1. This is also why I thought about the Dirichlet process. I have no restrictions on the values of $p$ otherwise.

• Do you have many $S_{n}$'s? Do you know more about the ranges of $p$, $\mu_{1}$, $\mu_{2}$ or $\sigma$? – Ytsen de Boer Apr 13 '17 at 7:23
• @YtsendeBoer I added more information in the question. There are no restrictions on the ranges of $p, \mu, \sigma$, except all of them are positive of course because of the lognormal. – spore234 Apr 13 '17 at 7:36

Each $X_i$ comes from either one of the two lognormals with probabilities $p$ and $1-p$. Let $Z_j$ be the number of $X_i$'s in the sum $S_j=\sum_{i=10j-9}^{10j}X_i$ that comes from the first lognormal. Clearly $Z_j \sim \mbox{bin}(10,p)$. Conditional on $Z_j$, each $S_j$ is a sum of $Z_j$ lognormals with parameters $\mu_1,\sigma^2$ and $10-Z_j$ lognormals with parameters $\mu_2,\sigma^2$, with pdf $$f_{S|Z=z}(s;\mu_1,\mu_2,\sigma^2), \tag{1}$$ given by a complicated convolution integral. The unconditional distribution of each $S_j$ is the 11-component mixture $$f_S(s;\mu_1,\mu_2,\sigma^2,p)=\sum_{z=0}^{10}{10 \choose z}p^z(1-p)^{10-z}f_{S|Z=z}(s;\mu_1,\mu_2,\sigma^2). \tag{2}$$ If $\sigma^2$ is sufficiently small or $\mu_1$ and $\mu_2$ sufficiently different this is going to be a 11-modal distribution with modes located near $10e^{\mu_1+\sigma^2},9e^{\mu_1+\sigma^2}+e^{\mu_2+\sigma^2},\dots,10e^{\mu_2+\sigma^2}$. This suggest that all five parameters in principle are identifiable given enough data.

Perhaps you can approximate the convolution in (1) by a single moment-matched lognormal as discussed here, use (2) to compute the likelihood and then compute approximate maximum likelihood estimates by maximising the resulting log likelihood numerically. Or you could do approximate Bayesian inference using this approximate likelihood function. This option would allow using informative priors on some of the parameters which might be necessary in practice if there is too much overlap between each component of (2).

• thanks, this is very impressive. Actually, I already implemented the exact same solution you pointed out as your first suggestion. I did numerically approximate the likelihood you gave (using optim). Of course I could also estimate this likelihood using MCMC, but actually I was hoping for a more inherently Bayesian solution. So could you elaborate on the approximate Bayesian thing a bit? – spore234 Apr 13 '17 at 16:32
• @spore234 well, you could give MCMC a try and report back with the observation that estimates of $\mu_1$, $\mu_2$ and $p$ oscillate and do not converge. Perhaps that would convince OP that the model is not identifiable? – matus Apr 13 '17 at 17:22
• @spore234 I was thinking you could just run some generic MCMC sampler like MCMCmetrop1R in MCMCpack, possibly building any prior knowledge you might have into the prior for $p, \mu_1,\mu_2$ and $\sigma$. If the approximation based on moment-matching gives bad result, you could perhaps use ABC using the above MLEs as summary statistics. You would then obtain samples from the exact posterior distribution, not conditional on all the data, but conditional on these summary statistics. – Jarle Tufto Apr 13 '17 at 17:46

A simulation based answer is to treat the $X_i$'s as latent variables and include these in a global MCMC sampler. At iteration $t$, it could proceed as follows

1. Simulate the $X_i^t$'s given the $S_j$'s and the current value of the parameters ${\theta}^t$
2. Given the $X_i^t$'s and their log-Normal mixture distribution, simulate the next value of the parameters, $\theta^{t+1}$

Step 2 is straightforward in that this is equivalent to simulate the posterior of a Normal sample. Step 1 can be decomposed in a sequence of Gibbs steps where each $X_i$ is generated conditional on the $8$ other $X_k$'s and the corresponding $S_j$. Meaning any MCMC move targeting the distribution $$\{pf_1(x_i;\mu_1,\sigma_1)+(1-p)f_2(x_i;\mu_2,\sigma_2)\}\times \{pf_1(s_j-x_{10(j-1)+1}-\ldots-x_{10j};\mu_1,\sigma_1)+(1-p)f_2(s_j-x_{10(j-1)+1}-\ldots-x_{10j};\mu_2,\sigma_2)\}$$

The parameter $p$ is not identifiable. This can be seen by writing

$f_{mix}=p\cdot \mathcal{LN}(\mu_1,\sigma) + (1-p)\cdot \mathcal{LN}(\mu_2,\sigma)$

as

$f_{mix}=\mathcal{LN}(\nu_1,\sigma) +\mathcal{LN}(\nu_2,\sigma)$

with $\nu_1=\log p + \mu_1$ and $\nu_2=\log (1-p) + \mu_2$

Then we can always choose $\mu_1=\nu_1-\log p$ and $\mu_2=\nu_2-\log (1-p)$. Then $$p(f_{mix}|\mu_1,\mu_2,\sigma,p)=p(f_{mix}|\nu_1,\nu_2,\sigma,p)=p(f_{mix}|\nu_1,\nu_2,\sigma)$$ and hence $$p(S_n|\nu_1,\nu_2,\sigma,p)=p(S_n|\nu_1,\nu_2,\sigma)$$ In other words the data likelihood is independent of $p$ given $\mu_1$ and $\mu_2$.

• -1: $f_{mix}$ is not a random variable but the pdf of a random variable which is lognormal with parameters $\mu_1,\sigma^2$ with probability $p$ and lognormal with parameters $\mu_2,\sigma^2$ with probability $1-p$. Yes, this is perhaps not the best notation. – Jarle Tufto Apr 13 '17 at 17:37
• The question asks explicitly for bayesian solution. In bayesian framework $p$ and any other parameter is treated as random variable. As a consequence $f_{mix}$ is a random variable. – matus Apr 14 '17 at 10:29