I am searching for a way to do the following steps below. I am stumped and I have stumped a few profs too. Any help would be appreciated.

  1. I would like to import data from, say, excel or a .txt file
  2. specify a certain number of data points (say 1 to N) to start out with and calculate: a) the coefficients and parameter estimates for a mixture or normal distributions. b) Each mixture component would have the same mean for those data points selected.
  3. This process would be carried out now for data points 1 to N+1 until the end of the data set.
  4. The values of the mixture coefficients as well as all mean and variance parameters are listed as output.

Essentially I am trying to have a mixture of specific normal distributions estimated a subset of my data one time step at a time.

Are there any gurus out there who know how to do this?

I appreciate any assistance!

  • 1
    $\begingroup$ Welcome to Mathematica.SE! I suggest the following: 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Read the faq! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! $\endgroup$
    – user9660
    Nov 22 '14 at 4:12

Here is how I would do it.

First, import the data using Import. This is straightforward from the documentation, using either Excel or txt files, so for the purposes of this exercise, I'm going to generate the data in Mathematica instead.

myD = MixtureDistribution[{2, 1}, {NormalDistribution[2, 2], 
NormalDistribution[2, 1/2]}];

For reference, here is the PDF.

Plot[PDF[myD, x], {x, -1, 5}]

enter image description here

Here is some fake data generated from that distribution, and an example of a subset being estimated. I understood from your question that you were assuming that the means of the component distributions were equal. It's not hard to change that, of course (just make the second mu into something else like othermu).

pretenddata = RandomVariate[myD, 500];
EstimatedDistribution[Take[pretenddata, 100], 
 MixtureDistribution[{w1, w2}, {NormalDistribution[mu, sig1], 
   NormalDistribution[mu, sig2]}]]

(* MixtureDistribution[{0.703728, 0.296272}, 
 {NormalDistribution[1.99722, 1.96856], NormalDistribution[1.99722, 0.460136]}] *)

A recursive estimation of the kind you want is most easily done as a table. It was quite slow, taking a few minutes on my machine with 100 total sample starting from a subsample of 60.

Table[EstimatedDistribution[Take[pretenddata, i], 
  MixtureDistribution[{w1, w2}, {NormalDistribution[mu, sig1], 
    NormalDistribution[mu, sig2]}]], {i, 400, 500}]

Personally I think that might be overkill, unless you have very small sample sizes - which makes me wonder whether you can justify such a sophisticated distribution. Perhaps you could get just as good results by partitioning the data into subsets.

Now to estimate across separate subsets, you could use Map (/@) to map the estimation function onto subsets created using Partition, like this (where pretenddata has 10,000 points, not 500 as above):

results = EstimatedDistribution[#, 
   MixtureDistribution[{w1, w2}, {NormalDistribution[mu, sig1], 
     NormalDistribution[mu, sig2]}]] & /@ Partition[pretenddata, 1000]

This takes rather less time than the recursive approach, which would be infeasible for such a large data set. Which ever approach you take, though, once you are done, though, you can use standard Part (shorthand [[]]) notation to extract different parts of the results. For example here is a bar chart of the estimated weight on the first component. In the estimated distributions, the component weights are normalised to sum to 1.

BarChart[results[[All, 1, 1]]]

enter image description here

The estimated parameters of the component distributions can be accessed as, for example, results[[All, 2, 1, 1]] - in any case you need a subpart of results[[All,2]].

I do not have the statistical smarts to work out how to do this if you do not know how many component distributions there are in the mixture. Presumably there is some likelihood ratio test you can use to check how many components you need.


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.