MixSemiRob (version 1.1.0)

mixpf: Profile Likelihood Method for Normal Mixture with Unequal Variance

Description

`mixpf' is used to estimate the following \(C\)-component univariate normal mixture model, using the profile likelihood method (Yao, 2010), with the assumption that the ratio of the smallest variance to the largest variance is \(k\): $$f(x;\boldsymbol{\theta}) = \sum_{j=1}^C\pi_j\phi(x;\mu_j,\sigma_j^2),$$ where \(\boldsymbol{\theta}=(\pi_1,\mu_1,\sigma_1,..,\pi_{C},\mu_C,\sigma_C)^{\top}\) is the parameter to estimate, \(\phi(\cdot;\mu,\sigma^2)\) is the normal density with a mean of \(\mu\) and a standard deviation of \(\sigma\), and \(\pi\)'s are mixing proportions that sum up to 1. Once the results are obtained, one can also find the maximum likelihood estimate (MLE) of \(k\) by plotting the likelihood vs. \(k\) for different \(k\) values and finding the maximum interior mode in the likelihood. See examples below.

Usage

mixpf(x, k = 0.5, C = 2, nstart = 20)

Value

A list containing the following elements:

mu

vector of estimated component means.

sigma

vector of estimated component standard deviations.

pi

vector of estimated mixing proportions.

lik

final likelihood.

Arguments

x

a vector of observations.

k

ratio of the smallest variance to the largest variance. Default is 0.5.

C

number of mixture components. Default is 2.

nstart

number of initializations to try. Default is 20.

References

Yao, W. (2010). A profile likelihood method for normal mixture with unequal variance. Journal of Statistical Planning and Inference, 140(7), 2089-2098.

Examples

set.seed(4)
n = 100
u = runif(n, 0, 1)
x2 = (u <= 0.3) * rnorm(n, 0, 0.5) + (u > 0.3) * rnorm(n, 1.5, 1)

# please set ngrid to 200 to get a smooth likelihood curve
ngrid = 5
grid = seq(from = 0.01, to = 1, length = ngrid)
likelihood = numeric()
for(i in 1:ngrid){
  k = grid[i]
  est = mixpf(x2, k)
  lh = est$lik
  likelihood[i] = lh
}

# visualize likelihood to find the best k
plot(grid, likelihood, type = "l", lty = 2, xlab = "k", ylab = "profile log-likelihood")