# Fitting a 1D distribution using Gaussian Mixtures

**URL:** <https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113>\
**Category:** General Usage\
**Tags:** statistics\
**Created:** [August 28, 2019, 3:00pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113 "2019-08-28T15:00:59Z")\
**Posts on this page:** 19\
**Page:** 1

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 28, 2019, 3:01pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/1 "2019-08-28T15:01:00Z")

</div>

I struggle to fit my data using Gaussian mixtures with Julia. I found the package [GitHub - davidavdav/GaussianMixtures.jl: Large scale Gaussian Mixture Models](https://github.com/davidavdav/GaussianMixtures.jl) which seems to work similarly to the scikit-stuff I used before but the actual fit gets stuck and every mixture component has the same parameters. I thought that’s an easy starting exercise but apparently I don’t grasp it.

Here is what I have:

```julia
using GaussianMixtures
using PGFPlotsX
using StatsBase

gmm = GMM(2, 1) # two components, one dimension

# dummy data
d₁ = randn(20000, 1) .- 5
d₂ = randn(100000, 1) .+ 23
d₊ = vcat(d₁, d₂) # 120000×1 Array{Float64,2}

```

This is how the distribution looks like:

```julia
h = fit(Histogram, d₊[:,1], nbins=40)
fig = @pgf Axis(Plot({"ybar interval"}, Coordinates(h.edges[1][2:end], h.weights)))

```

![distribution](https://global.discourse-cdn.com/julialang/original/3X/6/1/6180a18e43bc279b6e3c6e175bb4a7c3f46f447a.png)

Running the fit should be easy:

```julia
result = em!(gmm, d₊)

```

however, this spits out the following log, indicating that the optimisation got stuck after the first iteration:

```julia
┌ Info: Running 10 iterations EM on diag cov GMM with 2 Gaussians in 1 dimensions
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/train.jl:242
┌ Info: iteration 1, average log likelihood -223.922195
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 2, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 3, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 4, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 5, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 6, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 7, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 8, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 9, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: iteration 10, average log likelihood -3.768908
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71
┌ Info: EM with 120000 data points 10 iterations avll -3.768908
│ 24000.0 data points per parameter
└ @ GaussianMixtures /home/tgal/.julia/packages/GaussianMixtures/RGtTJ/src/gmms.jl:71

```

Looking at the \mu parameter, both Gauss components have the same value:

```julia
> gmm.μ
2×1 Array{Float64,2}:
 18.332104281524813
 18.332104281524813

```

The same holds if I try with more components…

So what is going wrong here?

---

<div class="post-metadata">

**Author:** ![uwbanjoman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/uwbanjoman/32/9998_2.png) [@uwbanjoman](https://discourse.julialang.org/u/uwbanjoman)\
**Post date:** [August 28, 2019, 3:36pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/2 "2019-08-28T15:36:16Z")

</div>

must be something in your setup  
my output is like this:  
julia\> result = em!(gmm, d₊)  
[ Info: Running 10 iterations EM on diag cov GMM with 2 Gaussians in 1 dimensions  
[ Info: iteration 1, average log likelihood -223.973397  
[ Info: iteration 2, average log likelihood -3.768754  
[ Info: iteration 3, average log likelihood -3.768754  
[ Info: iteration 4, average log likelihood -3.768754  
[ Info: iteration 5, average log likelihood -3.768754  
[ Info: iteration 6, average log likelihood -3.768754  
[ Info: iteration 7, average log likelihood -3.768754  
[ Info: iteration 8, average log likelihood -3.768754  
[ Info: iteration 9, average log likelihood -3.768754  
[ Info: iteration 10, average log likelihood -3.768754  
┌ Info: EM with 120000 data points 10 iterations avll -3.768754  
└ 24000.0 data points per parameter  
10-element Array{Float64,1}:  
-223.97339691979332  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004  
-3.768753993438004

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 28, 2019, 3:37pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/3 "2019-08-28T15:37:42Z")

</div>

Uhm, what do you mean? The output seems the same as mine… and I guess the `gmm.\mu` is also the same.

---

<div class="post-metadata">

**Author:** ![uwbanjoman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/uwbanjoman/32/9998_2.png) [@uwbanjoman](https://discourse.julialang.org/u/uwbanjoman)\
**Post date:** [August 28, 2019, 3:44pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/4 "2019-08-28T15:44:21Z")

</div>

What I mean is that on my PC i don’t get the logfile spitted out. I thought that that was your question / problem.

and yes  
julia\> gmm.μ  
2×1 Array{Float64,2}:  
18.335820572982005  
18.335820572982005

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 28, 2019, 3:48pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/5 "2019-08-28T15:48:07Z")

</div>

Small deviations are fine 😉

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 28, 2019, 4:16pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/6 "2019-08-28T16:16:14Z")

</div>

I am not sure what else I am doing wrong. I also tried more components with:

```julia
d₁ = randn(20000, 1) .- 5
d₂ = randn(100000, 1) .+ 23
d₃ = randn(10000, 1) .- 42
d₄ = randn(150000, 1) .+ 69
d₊ = vcat(d₁, d₂, d₃, d₄);

```

and get the same behaviour:

```julia
> gmm.μ
5×1 Array{Float64,2}:
 43.32143763742449
 43.32143763742449
 43.32143763742449
 43.32143763742449
 43.32143764486415

```

Anyways, it seems that I simply don’t understand the API.

This one works (if I do it in one-shot):

```julia
gmm = GMM(4, d₊)

```

Which yields

```julia
4×1 Array{Float64,2}:
  68.9987907346628 
  22.99872660258784
  -4.98850081738524
 -41.99187152150204

```

The problem is still that the sigmas and weights do not match…

---

<div class="post-metadata">

**Author:** ![pistacliffcho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pistacliffcho/32/8518_2.png) [@pistacliffcho](https://discourse.julialang.org/u/pistacliffcho)\
**Post date:** [August 28, 2019, 6:55pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/7 "2019-08-28T18:55:19Z")

</div>

This looks like a pretty significant bug in the Gaussian Mixtures package, but you should have an easy workaround.

From the Git page on how you initialized your model:

> Initialize a GMM with `n` multivariate Gaussians of dimension `d` . The means are all set to **0** (the origin) and the variances to **I** , which is silly by itself.

This is problematic: during the E-step, since each component has the exact same mean/variance, every single row will have the exact same probability of coming from each component…which means during the M-step, it will have the exact same estimate for each component parameter…which means during the next E-step, every row will have the same probability of being from each component, etc.

**Solution:**

This has a very simple solution: just randomly jitter the mean parameters before you start fitting the model, and you will get separation of the data by component, and the EM algorithm should actually work.

This should be included in the initialization of the model, but until that’s fixed you can just do it manually. I have no idea how the structure of the package works, so I don’t know how to do this, but I assume it can’t be that hard.

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 29, 2019, 7:01am UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/8 "2019-08-29T07:01:02Z")

</div>

Makes perfect sense, thanks!

Do you happen to know how to get the expectation value for a given `X`? I only managed to draw random samples with `rand(gmm, ...)` but when converting to a `MixtureModel` with `mm = MixtureModel(gmm)` and then `pdf(mm, x)` I get values which are not well weighted. It seems that some binning/normalisation is lost during the conversion.

However, I should probably open an issue 😉

---

<div class="post-metadata">

**Author:** ![pistacliffcho](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pistacliffcho/32/8518_2.png) [@pistacliffcho](https://discourse.julialang.org/u/pistacliffcho)\
**Post date:** [August 29, 2019, 8:49am UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/9 "2019-08-29T08:49:08Z")

</div>

If you mean in general, is there a closed form solution for the mean of a Gaussian mixture model, sure:

\sum\_{i = 1}^k p\_i \mu\_i

If you mean how does one extract that from the GuassianMixture.jl, that I have no idea; I’ve never used it before.

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 29, 2019, 8:51am UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/10 "2019-08-29T08:51:32Z")

</div>

Yeah I meant the second one. Alright, time for an issue. Thanks!

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 29, 2019, 8:57am UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/11 "2019-08-29T08:57:16Z")

</div>

Btw. do you know any other Julia package which can do GMM fitting?

---

<div class="post-metadata">

**Author:** ![mrVeng](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mrveng/32/8836_2.png) [@mrVeng](https://discourse.julialang.org/u/mrVeng)\
**Post date:** [August 29, 2019, 9:04am UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/12 "2019-08-29T09:04:32Z")

</div>

You would need to work in a Bayesian framework but `Turing` could do that easily for you, see here a tutorial for that:  
[Turing GMM](https://turing.ml/dev/tutorials/1-gaussianmixturemodel/)

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 29, 2019, 9:08am UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/13 "2019-08-29T09:08:19Z")

</div>

Wow that looks quite nice. I’ll try the tutorial and translate it to my 1D problem…  
In the end all I need is nothing else than the parameters of multiple Gaussians which describe the distribution of the data.

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [August 29, 2019, 12:09pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/14 "2019-08-29T12:09:30Z")

</div>

For the sake of completeness, here is the fit and the comparison of the fit to the data using `GaussianMixtures.jl`:

```julia
using GaussianMixtures
using PGFPlotsX
using StatsBase
using Random
Random.seed!(23);

d₁ = randn(20000, 1) .- 5
d₂ = randn(100000, 1) .+ 23
d₊ = vcat(d₁, d₂);

n = length(d₊)

gmm = GMM(2, d₊)

n_bins = 100

h = fit(Histogram, d₊[:,1], nbins=n_bins)
h₀ = fit(Histogram, rand(gmm, length(d₊))[:], nbins=n_bins)

fig = @pgf Axis(
    Plot({"ybar interval"},
        Coordinates(h.edges[1][2:end], h.weights)
    ),
    Plot({color="red"},
        Coordinates(h₀.edges[1][2:end].+(h₀.edges[1].step.hi/2), h₀.weights)
    )
)

```

![gmm_fit](https://global.discourse-cdn.com/julialang/original/3X/0/c/0c8cc1951694d28d00ab1f675f6fd1e2cbec6d36.png)

---

<div class="post-metadata">

**Author:** ![william](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/william/32/14117_2.png) [@william](https://discourse.julialang.org/u/william)\
**Post date:** [April 26, 2020, 2:09pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/15 "2020-04-26T14:09:14Z")

</div>

@tamasgal were you able to try @pistacliffcho’s solution for fitting the GMM? Did it eventually work?

Agreed that it’s bizarre that all parameters would be initialized to be identical. Random initialization works; ive also seen means initialized to random samples from the training data, or initialized with k-means / Lloyd’s algorithm.

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [April 26, 2020, 2:21pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/16 "2020-04-26T14:21:36Z")

</div>

Yes, you can see in my example just above your post how I randomly set the initial values, it works fine!

Edit: ah, my mistake, I did not include that in my example, but I remember that random starting values worked 😉

---

<div class="post-metadata">

**Author:** ![sairus7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sairus7/32/10816_2.png) [@sairus7](https://discourse.julialang.org/u/sairus7)\
**Post date:** [April 12, 2024, 3:10pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/17 "2024-04-12T15:10:51Z")

</div>

Sorry for reviving old thread, but looks like as of 2024 none of this is working.  
I’ve changed the examples to the following one, notice manual setting of different means, without this `em!` fits identical gausses:

```julia
using StatsBase
using StatsPlots
using GaussianMixtures

d1 = randn(20000, 1) .- 3
d2 = randn(100000, 1) .+ 3
dsum = vcat(d1, d2);
edges = -10:0.2:10

h = fit(Histogram, dsum[:,1], edges)
plot(h, alpha = 0.5)

gmm = GMM(2,1,kind=:full)
gmm.μ .= [-1, 1] # ! without this both gausses become identical
em!(gmm, dsum)

n = length(dsum)

h0 = fit(Histogram, rand(gmm, length(dsum))[:], edges)
plot!(h0, alpha = 0.5)

```

![image](https://global.discourse-cdn.com/julialang/original/3X/e/d/edad07c3ffb444b932b2dfab8370b02398b868a1.png)

---

<div class="post-metadata">

**Author:** ![tamasgal](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamasgal/32/27946_2.png) [@tamasgal](https://discourse.julialang.org/u/tamasgal)\
**Post date:** [April 12, 2024, 3:29pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/18 "2024-04-12T15:29:10Z")

</div>

At least I can contradict that the example above works, just like 4 years ago 😉

I was using these packages and the plot looks the same as in [Fitting a 1D distribution using Gaussian Mixtures - #14 by tamasgal](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/14)

```julia
Status `~/tmp/gmm/Project.toml`
  [cc18c42c] GaussianMixtures v0.3.9
  [8314cec4] PGFPlotsX v1.6.1
⌅ [2913bbd2] StatsBase v0.33.21

```

---

<div class="post-metadata">

**Author:** ![sairus7](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/sairus7/32/10816_2.png) [@sairus7](https://discourse.julialang.org/u/sairus7)\
**Post date:** [April 12, 2024, 3:51pm UTC](https://discourse.julialang.org/t/fitting-a-1d-distribution-using-gaussian-mixtures/28113/19 "2024-04-12T15:51:21Z")

</div>

Indeed, re-run this in a fresh env and this worked. looks like I’ve had polluted namespace from other packages, since

```julia
gmm = GMM(2, d₊)

```

gave me  
`ERROR: MethodError: no method matching weights(::Vector{Float64})`