In your model driver truncatedmodel() the fit function looks like:
z...@fit <- function(x, y, w) {
para <- list(mean = mean(x), sd = sd(x), lower = lower, upper= upper)
para$df <- 4
with(para, eval(z...@definecomponent))
}
w are the a-posteriori probabilities and denote the weights with which
observations are currently assigned to each component. These weights
are not used in your fit function and hence, the parameters of each
component are estimated identically using the whole sample. Please
modify the fit function to take the weights into account for example by
z...@fit <- function(x, y, w) {
para <- cov.wt(y, wt = w)[c("center", "cov")]
para$df <- (3 * ncol(y) + ncol(y)^2)/2
if (diagonal) {
para$cov <- diag(diag(para$cov))
para$df <- 2 * ncol(y)
}
with(para, eval(z...@definecomponent))
}
Please also note that your fit function is not appropriate because you
also have to take the truncation into account in the M-step. See for
example for grouped and truncated data:
G. McLachlan and P. Jones (1988) Fitting Mixture Models to Grouped and
Truncated Data via the EM Algorithm. Biometrics, 44(2): 571-578.
Best,
Bettina
Giovanni Luca Ciampaglia wrote:
Hello all,
I am trying to fit a truncated mixture model and I wrote a driver for
flexmix following the example in the vignette, but it doesn't work for
me: it assigns all data points to one component only, e.g.:
> source('bugged.R')
Call:
flexmix(formula = x ~ 1, k = 2, model = truncatedmodel(lower = -4,
upper = 4))
prior size post>0 ratio
Comp.1 0.494 0 1000 0
Comp.2 0.506 1000 1000 1
'log Lik.' -707703.3 (df=9)
AIC: 1415425 BIC: 1415469
What am I doing wrong? Please find my code attached.
cheers
------------------------------------------------------------------------
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
