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
--
Giovanni L. Ciampaglia
PhD Student
University of Lugano, MACS Lab
library(flexmix)
logf <- .Primitive("log")
dtnorm <- function(
x,
lower = -Inf,
upper = Inf,
mean = 0,
sd = 1,
log = FALSE)
{
cmin = pnorm(lower,mean,sd)
cmax = pnorm(upper,mean,sd)
d = array(0, dim=length(x))
d[(x >= lower) & (x <= upper)] = dnorm(x,mean,sd) / (cmax - cmin)
if (log)
return(logf(d))
return(d)
}
truncatedmodel <- function (formula = . ~ ., lower = -Inf, upper = Inf)
{
require(flexmix)
z <- new("FLXMC", weighted = TRUE, formula = formula, dist = "tnorm",
name = "Truncated bi-modal fitting")
z...@definecomponent <- expression({
logLik <- function(x, y)
dtnorm(y, lower=lower, upper=upper, mean=mean, sd=sd,
log=TRUE)
predict <- function(x, ...) array(mean, dim = length(x))
new("FLXcomponent", parameters = list(mean = mean, sd = sd, lower =
lower, upper = upper), df = df, logLik = logLik,
predict = predict)
})
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))
}
z
}
n <- 1000
x <- c()
while ( length(x) < n) {
x <- c(x, rnorm( n - length(x) ))
x <- x[ (x >= -2) & (x <= 2) ]
}
x[1:500] <- x[1:500] - 2
x[501:1000] <- x[501:1000] + 2
m1 <- flexmix( x ~ 1, k = 2, model = truncatedmodel(lower=-4,upper=4))
print(summary(m1))
______________________________________________
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.