Geoscience Reference
In-Depth Information
of the distribution. h e maximization step computes the parameters that
maximize the expected logarithmic likelihood computed in the expectation
step. h e function
mgdistribution.fit
creates an object of the
gmdistribution
class (see Section 2.5 and MATLAB Help on object-oriented programming
for details on objects and classes). h e function
gmdistribution.fit
treats
NaN values as missing data: rows of
data
with NaN values are excluded from
the i t. We can now determine the Gaussian mixture distribution with two
components in a single dimension.
gmfit = gmdistribution.fit(data,2)
Gaussian mixture distribution with 2 components in 1 dimensions
Component 1:
Mixing proportion: 0.513162
Mean: 13.0942
Component 2:
Mixing proportion: 0.486838
Mean: 6.4730
h us we obtain the means and relative mixing proportions of both
distributions. In our example both normal distributions with means of
6.4730 and 13.0942 contribute ~50% (~0.49 and ~0.51, respectively) to the
mixture distribution. h e object
gmfit
contains several layers of information,
including the mean
gmfit.mu
and the standard deviation
gmfit.Sigma
that we
use to calculate the probability density function
y
of the mixed distribution.
x = 0 : 1/30 : 20;
y1 = normpdf(x,gmfit.mu(1,1),gmfit.Sigma(:,:,1));
y2 = normpdf(x,gmfit.mu(2,1),gmfit.Sigma(:,:,2));
h eobject
gmfit
also contains information on the relative mixing proportions
of the two distributions in the layer
gmfit.PComponents
. We can use this
information to scale
y1
and
y2
to the correction proportions relative to each
other.
y1 = gmfit.PComponents(1,1) * y1/trapz(x,y1);
y2 = gmfit.PComponents(1,2) * y2/trapz(x,y2);
We can now superimpose the two scaled probability density functions
y1
and
y2
, and scale the result
y
to the same integral of 200 as the original data.
h e integral of the original data is determined using the function
trapz
to
perform a trapezoidal numerical integration.
y = y1 + y2;
y = trapz(v,n) * y/trapz(x(1:10:end),y(1:10:end));