FindDistributionParameters of a sum of a mixture distribution

Question

I have this Mixture Distribution of two Normal Distributions:

 mix[p_, μ1_, μ2_, σ_] := 
 MixtureDistribution[{p, 
 1 - p}, {NormalDistribution[μ1, σ], 
 NormalDistribution[μ2, σ]}];

Say $x_1, x_2$ are samples from this mixture. Set $y=x_1+x_2$, so that $y$ is the sum of two samples from the mixture. I now generate 100 $y's$ (with arbitrary p, mu and sigma, each is a sum of two $x$:

mixdat = RandomVariate[mix[0.75, 0.5, -1.5, 0.2], 1000];
mixdatSum = 
Total[RandomVariate[mix[0.75, 0.5, -1.5, 0.2], 2]] & /@ Range[100];

mixdatSum contains these 100 values. Now given the data in mixdatSum I want to estimate the Distribution Parameters. So a natural way to solve this would be:

FindDistributionParameters[mixdatSum, 
TransformedDistribution[x1 + x2, 
Thread[Distributed[{x1, x2}, 
MixtureDistribution[{p, 
  1 - p}, {NormalDistribution[μ1, σ], 
  NormalDistribution[μ2, σ]}]]]]]

Problem is, the sum of a mixture does not follow the mixture distribution, so this does not work and I get lots of errors.

However I want to treat the sum of a mixture AS IF IT WERE following the same mixture, for approximation purposes.

Is there a way to do this?

Answer 1

Since the PDF can be computed in closed form you might have some luck with ProbabilityDistribution and some half-way reasonable starting values.

Generate the data...

mix[p_, m1_, m2_, s_] := MixtureDistribution[{p, 1 - p}, {NormalDistribution[m1, s], 
    NormalDistribution[m2, s]}];

mixdatSum = Plus @@ RandomVariate[mix[0.75, 0.5, -1.5, 0.2], {2, 100}];

Use your TransformedDistribution to get the PDF.

sumdist = TransformedDistribution[x1 + x2, Thread[Distributed[{x1, x2}, 
     MixtureDistribution[{p,  1 - p}, {NormalDistribution[m1, s], 
       NormalDistribution[m2, s]}]]]];

pdf[z_] = FullSimplify[PDF[sumdist, z], DistributionParameterAssumptions[sumdist]]

Create a ProbabilityDistribution from your PDF with the proper assumptions set.

pdist = ProbabilityDistribution[pdf[z], {z, -Infinity, Infinity}, 
  Assumptions -> DistributionParameterAssumptions[sumdist]]

With some starting values in the right ballpark you can get good estimates of the parameters.

FindDistributionParameters[mixdatSum, pdist, {{m1, 1}, {m2, -1}, {s, 1}, {p, .5}}]

(*{m1 -> 0.514195, m2 -> -1.53428, s -> 0.226609, p -> 0.719953}