Last updated on 2020-02-19 10:49:00 CET.
| Flavor | Version | Tinstall | Tcheck | Ttotal | Status | Flags |
|---|---|---|---|---|---|---|
| r-devel-linux-x86_64-debian-clang | 0.3 | 5.97 | 35.74 | 41.71 | ERROR | |
| r-devel-linux-x86_64-debian-gcc | 0.3 | 5.13 | 28.04 | 33.17 | ERROR | |
| r-devel-linux-x86_64-fedora-clang | 0.3 | 50.58 | ERROR | |||
| r-devel-linux-x86_64-fedora-gcc | 0.3 | 49.29 | ERROR | |||
| r-devel-windows-ix86+x86_64 | 0.3 | 9.00 | 54.00 | 63.00 | OK | |
| r-devel-windows-ix86+x86_64-gcc8 | 0.3 | 13.00 | 75.00 | 88.00 | OK | |
| r-patched-linux-x86_64 | 0.3 | 4.56 | 36.90 | 41.46 | OK | |
| r-patched-solaris-x86 | 0.3 | 77.80 | OK | |||
| r-release-linux-x86_64 | 0.3 | 5.22 | 37.57 | 42.79 | OK | |
| r-release-windows-ix86+x86_64 | 0.3 | 9.00 | 49.00 | 58.00 | OK | |
| r-release-osx-x86_64 | 0.3 | OK | ||||
| r-oldrel-windows-ix86+x86_64 | 0.3 | 4.00 | 47.00 | 51.00 | OK | |
| r-oldrel-osx-x86_64 | 0.3 | OK |
Version: 0.3
Check: examples
Result: ERROR
Running examples in 'mmtfa-Ex.R' failed
The error most likely occurred in:
> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: mmtfa
> ### Title: mmtfa: Function for Model-Based Clustering and Classification
> ### with Mixtures of Multivariate t Factor Analyzers
> ### Aliases: mmtfa
>
> ### ** Examples
>
> ###Note that only one model is run for each example
> ###in order to reduce computation time
>
> #Clustering iris data with hard random start
> tirisr <- mmtfa(iris[,-5], models="UUUU", Gs=1:3, Qs=1, init="hard")
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
mmtfa
--- call from context ---
mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
--- call from argument ---
if (class(duptest) == "try-error") {
break
}
--- R stacktrace ---
where 1: mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
where 2: mmtfa(iris[, -5], models = "UUUU", Gs = 1:3, Qs = 1, init = "hard")
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
{
modgen <- modelgen()
modold <- modgen$modold
p <- ncol(x)
n <- nrow(x)
zlist3 <- ll <- dff <- it <- store <- meanlist <- siglist <- list()
if (clas > 0) {
testindex <- sample(1:n, ceiling(n * (clas/100)))
kno <- vector(mode = "numeric", length = n)
kno[testindex] <- 1
unkno <- (kno - 1) * (-1)
Gs <- length(unique(known))
}
gvec <- 1:max(Gs)
qvec <- 1:max(Qs)
gstuff <- paste("G=", gvec, sep = "")
qstuff <- paste("Q=", qvec, sep = "")
bic <- rands <- icl <- logls <- array(-Inf, dim = c(length(models),
max(Qs), max(Gs)))
meansave <- meansave2 <- sigsave <- zmatsave <- sigsave2 <- zmatsave2 <- NA
oldmodvec <- modold[match(models, modgen$allmodels)]
if (class(init) != "list" && !(init %in% c("kmeans", "hard",
"disc", "cont", "soft", "uniform"))) {
stop("'init' must be one of 'kmeans', 'hard', 'soft', 'uniform' or a list. See ?mmtfa.")
}
zmatin <- list()
for (G in Gs) {
if (G == 1) {
zmatin[[G]] <- matrix(1, n, 1)
}
else {
if (is.character(init)) {
if (init == "hard") {
zmatin[[G]] <- discrandz(n, G)
}
if (init == "soft") {
zmatin[[G]] <- contrandz(n, G)
}
if (init == "uniform") {
if (clas > 0) {
zmatin[[G]] <- uniformz(n, G, clas, kno,
known)
}
else {
stop("Uniform initialization not available for clustering.")
return(NULL)
}
}
if (init == "kmeans") {
zmatin[[G]] <- kmeansz(x, n, G)
}
}
else {
zmatin[[G]] <- givenz(n, G, init[[G]])
}
}
}
for (modnum in 1:length(models)) {
modnew <- models[modnum]
mod <- modold[which(modgen$allmodels == modnew)]
for (G in Gs) {
delta <- matrix(0, n, G)
mug <- matrix(0, G, p)
om <- rep(0, G)
yg <- sigma <- sigmainv <- sg <- array(0, dim = c(p,
p, G))
tri <- array(0, dim = c(p, p, G))
w <- matrix(0, n, G)
for (q in Qs) {
lg <- array(0, dim = c(p, q, G))
betag <- array(0, dim = c(q, p, G))
thetag <- array(0, dim = c(q, q, G))
singular <- 0
breakit <- 0
if (G == 1) {
CCCCgroup <- c("UUCU", "UUCC", "UCCU", "UCCC",
"CUCU", "CUCC", "CCCU", "CCCC")
if (any(mod == CCCCgroup)) {
cccdum <- oldmodvec[oldmodvec %in% CCCCgroup]
if (length(cccdum) > 0) {
if (mod != cccdum[1]) {
breakit <- 1
}
}
}
CCUCgroup <- c("UUUU", "UUUC", "UCUU", "UCUC",
"CUUU", "CUUC", "CCUU", "CCUC", "Mt1U", "Mt1C",
"Mt2U", "Mt2C", "Mt3U", "Mt3C", "Mt4U", "Mt4C")
if (any(mod == CCUCgroup)) {
ccudum <- oldmodvec[oldmodvec %in% CCUCgroup]
if (length(ccudum) > 0) {
if (mod != ccudum[1]) {
breakit <- 1
}
}
}
}
if (breakit == 0) {
zmat <- zmatin[[G]]
vg <- vginit(dfstart, G)
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- muginit(G, p, x, zmat, ng)
sg <- sginit(p, G, x, mug, zmat, n, ng)
sgc <- sginitc(G, sg, pig, p, n, x)
if (substring(mod, 1, 3) == "CCC") {
for (g in 1:G) {
sg[, , g] <- sgc
}
}
if (substring(mod, 1, 1) == "U" | substring(mod,
3, 3) == "2" | substring(mod, 3, 3) == "4") {
lg <- lginitu(p, q, G, sg)
}
if (substring(mod, 1, 1) == "C" | substring(mod,
3, 3) == "1" | substring(mod, 3, 3) == "3") {
lg <- lginitc(p, q, G, sgc)
dumg <- lginitu(p, q, G, sg)
}
if (substring(mod, 3, 3) == "U") {
yg <- yginitu(p, G, sg, lg, mod, pig, dumg)
}
if (substring(mod, 3, 3) == "C") {
yg <- yginitc(p, G, sg, lg, sgc, mod, pig)
}
if (substring(mod, 1, 1) == "M") {
yg <- array(0, dim = c(p, p, G))
om <- rep(0, G)
tri <- array(0, dim = c(p, p, G))
ygst <- yginitu(p, G, sg, lg, mod, pig, dumg)
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
for (g in 1:G) {
om[g] <- det(ygst[, , g])^(1/p)
}
}
if (substring(mod, 3, 3) == "3" | substring(mod,
3, 3) == "4") {
dom <- 0
for (g in 1:G) {
dom <- dom + pig[g] * det(ygst[, , g])^(1/p)
}
om[] <- dom
}
for (g in 1:G) {
tri[, , g] <- ygst[, , g]/(det(ygst[, ,
g])^(1/p))
}
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
av <- diag(p) - diag(p)
for (g in 1:G) {
av <- av + pig[g] * tri[, , g]
}
tri[, , ] <- av
}
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
w <- winit(x, n, G, mug, sigmainv, vg, p, sg,
zmat)
cycle <- 0
dfnewg <- vg
}
conv <- 0
num <- matrix(0, n, G)
ft <- matrix(0, n, G)
logl <- NaN
while (conv != 1) {
if (breakit == 1) {
break
}
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- mugupdate(G, zmat, w, x, p, mug, n)
if (dfupdate == "approx") {
testing <- try(dfnewg <- dfupdatefun2(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
if (dfupdate == "numeric") {
testing <- try(dfnewg <- dfupdatefun(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
ng <- ngupdate(zmat)
sg <- sgupdate(p, G, n, x, mug, zmat, w, ng,
mod, pig, sg)
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
testing <- try(lg <- lgupdate(mod, p, q, G,
ng, yginv, sg, betag, thetag, om, tri, lg),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
if (substring(mod, 1, 1) == "M") {
om <- omupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri)
tri <- triupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri, ng)
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
else {
yg <- ygupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig)
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
duptest <- try(delta <- deltaup(x, mug, sigma,
sigmainv, G, n, delta), silent = TRUE)
if (class(duptest) == "try-error") {
break
}
suppressWarnings(zup <- zupdate(x, G, pig,
dfnewg, p, yg, q, betag, lg, mug, sigmainv,
n, clas, kno, known, unkno, delta))
zmat <- zup$zmat
if (any(is.nan(zmat))) {
break
}
w <- wupdate(x, n, G, mug, sigmainv, dfnewg,
p, delta, w)
cycle <- cycle + 1
logl[cycle] <- sum(log(rowSums(zup$num))) -
sum(zup$kcon)
if (is.na(logl[cycle])) {
break
}
if (cycle > 3) {
if (is.finite(logl[cycle - 2])) {
ak <- (logl[cycle] - logl[cycle - 1])/(logl[cycle -
1] - logl[cycle - 2])
linf <- logl[cycle - 1] + (logl[cycle] -
logl[cycle - 1])/(1 - ak)
if (abs(linf - logl[cycle - 1]) < eps) {
conv <- 1
}
if ((logl[cycle] - logl[cycle - 1]) < 0) {
break
}
}
else {
break
}
}
}
if (conv == 1) {
bic[modnum, q, G] <- bicdum <- BICcalc(conv,
G, p, mod, q, logl, n, gauss)
icl[modnum, q, G] <- icldum <- ICLcalc(conv,
n, zmat, bic, modnum, q, G)
if (bicdum == max(bic)) {
meansave <- mug
sigsave <- sigma
zmatsave <- zmat
dfsave <- dfnewg
itsave <- cycle
llsave <- logl[cycle]
}
if (icldum == max(icl)) {
meansave2 <- mug
sigsave2 <- sigma
zmatsave2 <- zmat
dfsave2 <- dfnewg
itsave2 <- cycle
llsave2 <- logl[cycle]
}
}
}
}
}
dimnames(bic) <- list(models, qstuff, gstuff)
dimnames(icl) <- list(models, qstuff, gstuff)
maxes <- which(bic == max(bic), arr.ind = TRUE)
maxicl <- which(icl == max(icl), arr.ind = TRUE)
if (nrow(maxes) > 1) {
message("WARNING: Maximum BIC tie between two or more models")
bestmodnum <- maxes[1:nrow(maxes), 1]
bestmod <- models[bestmodnum]
bestq <- maxes[1:nrow(maxes), 2]
bestg <- maxes[1:nrow(maxes), 3]
itf <- "MULTIPLE"
dff1 <- "MULTIPLE"
bestz <- "MULTIPLE"
bestzmap <- "MULTIPLE"
tab <- "MULTIPLE"
blogl <- "Multiple"
}
if (nrow(maxes) == 1) {
bestmodnum <- maxes[1]
bestmod <- models[bestmodnum]
bestq <- maxes[2]
bestg <- maxes[3]
bestz <- zmatsave
dff1 <- dfsave
itf <- itsave
blogl <- llsave
bestzmap <- apply(bestz, 1, which.max)
if (clas > 0) {
newmap <- bestzmap
newmap[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tab <- table(known, newmap)
}
else {
if (!is.null(known)) {
tab <- table(known, bestzmap)
}
else {
tab <- NULL
}
}
}
if (nrow(maxicl) > 1) {
message("WARNING: Maximum ICL tie between two or more models")
bestmodnumicl <- maxicl[1:nrow(maxicl), 1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[1:nrow(maxicl), 2]
bestgicl <- maxicl[1:nrow(maxicl), 3]
dff1icl <- "MULTIPLE"
bestzicl <- "MULTIPLE"
bestzmapicl <- "MULTIPLE"
itficl <- "MULTIPLE"
tabicl <- "MULTIPLE"
bloglicl <- "MULTIPLE"
}
if (nrow(maxicl) == 1) {
bestmodnumicl <- maxicl[1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[2]
bestgicl <- maxicl[3]
bestzicl <- zmatsave2
dff1icl <- dfsave2
itficl <- itsave2
bloglicl <- llsave2
bestzmapicl <- apply(bestzicl, 1, which.max)
if (clas > 0) {
newmapicl <- bestzmapicl
newmapicl[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tabicl <- table(known, newmapicl)
}
else {
if (!is.null(known)) {
tabicl <- table(known, bestzmapicl)
}
else {
tabicl <- NULL
}
}
}
iclresults <- list()
par <- list()
paricl <- list()
par[["mean"]] <- meansave
par[["sigma"]] <- sigsave
par[["df"]] <- dff1
paricl[["sigma"]] <- sigsave2
paricl[["mean"]] <- meansave2
paricl[["df"]] <- dff1icl
store[["parameters"]] <- par
store[["allbic"]] <- bic[, Qs, Gs]
iclresults[["allicl"]] <- icl[, Qs, Gs]
store[["bic"]] <- max(bic)
iclresults[["icl"]] <- max(icl)
store[["modelname"]] <- bestmod
store[["bestmodel"]] <- paste("The best model (BIC of ",
round(max(bic), 2), ") is ", bestmod, " with G=", bestg,
sep = "")
store[["Q"]] <- bestq
store[["G"]] <- bestg
store[["classification"]] <- bestzmap
iclresults[["bestmodel"]] <- paste("The best model (ICL of ",
round(max(icl), 2), ") is ", bestmodicl, " with G=",
bestgicl, sep = "")
iclresults[["modelname"]] <- bestmodicl
iclresults[["Q"]] <- bestqicl
iclresults[["G"]] <- bestgicl
iclresults[["fuzzy"]] <- bestzicl
iclresults[["logl"]] <- bloglicl
iclresults[["classification"]] <- bestzmapicl
iclresults[["parameters"]] <- paricl
store[["tab"]] <- tab
iclresults[["tab"]] <- tabicl
store[["iter"]] <- itf
iclresults[["iter"]] <- itficl
store[["x"]] <- x
store[["fuzzy"]] <- bestz
store[["logl"]] <- blogl
store[["iclresults"]] <- iclresults
store
}
<bytecode: 0x33fbc70>
<environment: namespace:mmtfa>
--- function search by body ---
Function mmtfaEM in namespace mmtfa has this body.
----------- END OF FAILURE REPORT --------------
Error in if (class(duptest) == "try-error") { :
the condition has length > 1
Calls: mmtfa -> mmtfaEM
Execution halted
Flavor: r-devel-linux-x86_64-debian-clang
Version: 0.3
Check: examples
Result: ERROR
Running examples in ‘mmtfa-Ex.R’ failed
The error most likely occurred in:
> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: mmtfa
> ### Title: mmtfa: Function for Model-Based Clustering and Classification
> ### with Mixtures of Multivariate t Factor Analyzers
> ### Aliases: mmtfa
>
> ### ** Examples
>
> ###Note that only one model is run for each example
> ###in order to reduce computation time
>
> #Clustering iris data with hard random start
> tirisr <- mmtfa(iris[,-5], models="UUUU", Gs=1:3, Qs=1, init="hard")
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
mmtfa
--- call from context ---
mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
--- call from argument ---
if (class(duptest) == "try-error") {
break
}
--- R stacktrace ---
where 1: mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
where 2: mmtfa(iris[, -5], models = "UUUU", Gs = 1:3, Qs = 1, init = "hard")
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
{
modgen <- modelgen()
modold <- modgen$modold
p <- ncol(x)
n <- nrow(x)
zlist3 <- ll <- dff <- it <- store <- meanlist <- siglist <- list()
if (clas > 0) {
testindex <- sample(1:n, ceiling(n * (clas/100)))
kno <- vector(mode = "numeric", length = n)
kno[testindex] <- 1
unkno <- (kno - 1) * (-1)
Gs <- length(unique(known))
}
gvec <- 1:max(Gs)
qvec <- 1:max(Qs)
gstuff <- paste("G=", gvec, sep = "")
qstuff <- paste("Q=", qvec, sep = "")
bic <- rands <- icl <- logls <- array(-Inf, dim = c(length(models),
max(Qs), max(Gs)))
meansave <- meansave2 <- sigsave <- zmatsave <- sigsave2 <- zmatsave2 <- NA
oldmodvec <- modold[match(models, modgen$allmodels)]
if (class(init) != "list" && !(init %in% c("kmeans", "hard",
"disc", "cont", "soft", "uniform"))) {
stop("'init' must be one of 'kmeans', 'hard', 'soft', 'uniform' or a list. See ?mmtfa.")
}
zmatin <- list()
for (G in Gs) {
if (G == 1) {
zmatin[[G]] <- matrix(1, n, 1)
}
else {
if (is.character(init)) {
if (init == "hard") {
zmatin[[G]] <- discrandz(n, G)
}
if (init == "soft") {
zmatin[[G]] <- contrandz(n, G)
}
if (init == "uniform") {
if (clas > 0) {
zmatin[[G]] <- uniformz(n, G, clas, kno,
known)
}
else {
stop("Uniform initialization not available for clustering.")
return(NULL)
}
}
if (init == "kmeans") {
zmatin[[G]] <- kmeansz(x, n, G)
}
}
else {
zmatin[[G]] <- givenz(n, G, init[[G]])
}
}
}
for (modnum in 1:length(models)) {
modnew <- models[modnum]
mod <- modold[which(modgen$allmodels == modnew)]
for (G in Gs) {
delta <- matrix(0, n, G)
mug <- matrix(0, G, p)
om <- rep(0, G)
yg <- sigma <- sigmainv <- sg <- array(0, dim = c(p,
p, G))
tri <- array(0, dim = c(p, p, G))
w <- matrix(0, n, G)
for (q in Qs) {
lg <- array(0, dim = c(p, q, G))
betag <- array(0, dim = c(q, p, G))
thetag <- array(0, dim = c(q, q, G))
singular <- 0
breakit <- 0
if (G == 1) {
CCCCgroup <- c("UUCU", "UUCC", "UCCU", "UCCC",
"CUCU", "CUCC", "CCCU", "CCCC")
if (any(mod == CCCCgroup)) {
cccdum <- oldmodvec[oldmodvec %in% CCCCgroup]
if (length(cccdum) > 0) {
if (mod != cccdum[1]) {
breakit <- 1
}
}
}
CCUCgroup <- c("UUUU", "UUUC", "UCUU", "UCUC",
"CUUU", "CUUC", "CCUU", "CCUC", "Mt1U", "Mt1C",
"Mt2U", "Mt2C", "Mt3U", "Mt3C", "Mt4U", "Mt4C")
if (any(mod == CCUCgroup)) {
ccudum <- oldmodvec[oldmodvec %in% CCUCgroup]
if (length(ccudum) > 0) {
if (mod != ccudum[1]) {
breakit <- 1
}
}
}
}
if (breakit == 0) {
zmat <- zmatin[[G]]
vg <- vginit(dfstart, G)
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- muginit(G, p, x, zmat, ng)
sg <- sginit(p, G, x, mug, zmat, n, ng)
sgc <- sginitc(G, sg, pig, p, n, x)
if (substring(mod, 1, 3) == "CCC") {
for (g in 1:G) {
sg[, , g] <- sgc
}
}
if (substring(mod, 1, 1) == "U" | substring(mod,
3, 3) == "2" | substring(mod, 3, 3) == "4") {
lg <- lginitu(p, q, G, sg)
}
if (substring(mod, 1, 1) == "C" | substring(mod,
3, 3) == "1" | substring(mod, 3, 3) == "3") {
lg <- lginitc(p, q, G, sgc)
dumg <- lginitu(p, q, G, sg)
}
if (substring(mod, 3, 3) == "U") {
yg <- yginitu(p, G, sg, lg, mod, pig, dumg)
}
if (substring(mod, 3, 3) == "C") {
yg <- yginitc(p, G, sg, lg, sgc, mod, pig)
}
if (substring(mod, 1, 1) == "M") {
yg <- array(0, dim = c(p, p, G))
om <- rep(0, G)
tri <- array(0, dim = c(p, p, G))
ygst <- yginitu(p, G, sg, lg, mod, pig, dumg)
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
for (g in 1:G) {
om[g] <- det(ygst[, , g])^(1/p)
}
}
if (substring(mod, 3, 3) == "3" | substring(mod,
3, 3) == "4") {
dom <- 0
for (g in 1:G) {
dom <- dom + pig[g] * det(ygst[, , g])^(1/p)
}
om[] <- dom
}
for (g in 1:G) {
tri[, , g] <- ygst[, , g]/(det(ygst[, ,
g])^(1/p))
}
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
av <- diag(p) - diag(p)
for (g in 1:G) {
av <- av + pig[g] * tri[, , g]
}
tri[, , ] <- av
}
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
w <- winit(x, n, G, mug, sigmainv, vg, p, sg,
zmat)
cycle <- 0
dfnewg <- vg
}
conv <- 0
num <- matrix(0, n, G)
ft <- matrix(0, n, G)
logl <- NaN
while (conv != 1) {
if (breakit == 1) {
break
}
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- mugupdate(G, zmat, w, x, p, mug, n)
if (dfupdate == "approx") {
testing <- try(dfnewg <- dfupdatefun2(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
if (dfupdate == "numeric") {
testing <- try(dfnewg <- dfupdatefun(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
ng <- ngupdate(zmat)
sg <- sgupdate(p, G, n, x, mug, zmat, w, ng,
mod, pig, sg)
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
testing <- try(lg <- lgupdate(mod, p, q, G,
ng, yginv, sg, betag, thetag, om, tri, lg),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
if (substring(mod, 1, 1) == "M") {
om <- omupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri)
tri <- triupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri, ng)
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
else {
yg <- ygupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig)
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
duptest <- try(delta <- deltaup(x, mug, sigma,
sigmainv, G, n, delta), silent = TRUE)
if (class(duptest) == "try-error") {
break
}
suppressWarnings(zup <- zupdate(x, G, pig,
dfnewg, p, yg, q, betag, lg, mug, sigmainv,
n, clas, kno, known, unkno, delta))
zmat <- zup$zmat
if (any(is.nan(zmat))) {
break
}
w <- wupdate(x, n, G, mug, sigmainv, dfnewg,
p, delta, w)
cycle <- cycle + 1
logl[cycle] <- sum(log(rowSums(zup$num))) -
sum(zup$kcon)
if (is.na(logl[cycle])) {
break
}
if (cycle > 3) {
if (is.finite(logl[cycle - 2])) {
ak <- (logl[cycle] - logl[cycle - 1])/(logl[cycle -
1] - logl[cycle - 2])
linf <- logl[cycle - 1] + (logl[cycle] -
logl[cycle - 1])/(1 - ak)
if (abs(linf - logl[cycle - 1]) < eps) {
conv <- 1
}
if ((logl[cycle] - logl[cycle - 1]) < 0) {
break
}
}
else {
break
}
}
}
if (conv == 1) {
bic[modnum, q, G] <- bicdum <- BICcalc(conv,
G, p, mod, q, logl, n, gauss)
icl[modnum, q, G] <- icldum <- ICLcalc(conv,
n, zmat, bic, modnum, q, G)
if (bicdum == max(bic)) {
meansave <- mug
sigsave <- sigma
zmatsave <- zmat
dfsave <- dfnewg
itsave <- cycle
llsave <- logl[cycle]
}
if (icldum == max(icl)) {
meansave2 <- mug
sigsave2 <- sigma
zmatsave2 <- zmat
dfsave2 <- dfnewg
itsave2 <- cycle
llsave2 <- logl[cycle]
}
}
}
}
}
dimnames(bic) <- list(models, qstuff, gstuff)
dimnames(icl) <- list(models, qstuff, gstuff)
maxes <- which(bic == max(bic), arr.ind = TRUE)
maxicl <- which(icl == max(icl), arr.ind = TRUE)
if (nrow(maxes) > 1) {
message("WARNING: Maximum BIC tie between two or more models")
bestmodnum <- maxes[1:nrow(maxes), 1]
bestmod <- models[bestmodnum]
bestq <- maxes[1:nrow(maxes), 2]
bestg <- maxes[1:nrow(maxes), 3]
itf <- "MULTIPLE"
dff1 <- "MULTIPLE"
bestz <- "MULTIPLE"
bestzmap <- "MULTIPLE"
tab <- "MULTIPLE"
blogl <- "Multiple"
}
if (nrow(maxes) == 1) {
bestmodnum <- maxes[1]
bestmod <- models[bestmodnum]
bestq <- maxes[2]
bestg <- maxes[3]
bestz <- zmatsave
dff1 <- dfsave
itf <- itsave
blogl <- llsave
bestzmap <- apply(bestz, 1, which.max)
if (clas > 0) {
newmap <- bestzmap
newmap[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tab <- table(known, newmap)
}
else {
if (!is.null(known)) {
tab <- table(known, bestzmap)
}
else {
tab <- NULL
}
}
}
if (nrow(maxicl) > 1) {
message("WARNING: Maximum ICL tie between two or more models")
bestmodnumicl <- maxicl[1:nrow(maxicl), 1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[1:nrow(maxicl), 2]
bestgicl <- maxicl[1:nrow(maxicl), 3]
dff1icl <- "MULTIPLE"
bestzicl <- "MULTIPLE"
bestzmapicl <- "MULTIPLE"
itficl <- "MULTIPLE"
tabicl <- "MULTIPLE"
bloglicl <- "MULTIPLE"
}
if (nrow(maxicl) == 1) {
bestmodnumicl <- maxicl[1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[2]
bestgicl <- maxicl[3]
bestzicl <- zmatsave2
dff1icl <- dfsave2
itficl <- itsave2
bloglicl <- llsave2
bestzmapicl <- apply(bestzicl, 1, which.max)
if (clas > 0) {
newmapicl <- bestzmapicl
newmapicl[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tabicl <- table(known, newmapicl)
}
else {
if (!is.null(known)) {
tabicl <- table(known, bestzmapicl)
}
else {
tabicl <- NULL
}
}
}
iclresults <- list()
par <- list()
paricl <- list()
par[["mean"]] <- meansave
par[["sigma"]] <- sigsave
par[["df"]] <- dff1
paricl[["sigma"]] <- sigsave2
paricl[["mean"]] <- meansave2
paricl[["df"]] <- dff1icl
store[["parameters"]] <- par
store[["allbic"]] <- bic[, Qs, Gs]
iclresults[["allicl"]] <- icl[, Qs, Gs]
store[["bic"]] <- max(bic)
iclresults[["icl"]] <- max(icl)
store[["modelname"]] <- bestmod
store[["bestmodel"]] <- paste("The best model (BIC of ",
round(max(bic), 2), ") is ", bestmod, " with G=", bestg,
sep = "")
store[["Q"]] <- bestq
store[["G"]] <- bestg
store[["classification"]] <- bestzmap
iclresults[["bestmodel"]] <- paste("The best model (ICL of ",
round(max(icl), 2), ") is ", bestmodicl, " with G=",
bestgicl, sep = "")
iclresults[["modelname"]] <- bestmodicl
iclresults[["Q"]] <- bestqicl
iclresults[["G"]] <- bestgicl
iclresults[["fuzzy"]] <- bestzicl
iclresults[["logl"]] <- bloglicl
iclresults[["classification"]] <- bestzmapicl
iclresults[["parameters"]] <- paricl
store[["tab"]] <- tab
iclresults[["tab"]] <- tabicl
store[["iter"]] <- itf
iclresults[["iter"]] <- itficl
store[["x"]] <- x
store[["fuzzy"]] <- bestz
store[["logl"]] <- blogl
store[["iclresults"]] <- iclresults
store
}
<bytecode: 0x558621a76b48>
<environment: namespace:mmtfa>
--- function search by body ---
Function mmtfaEM in namespace mmtfa has this body.
----------- END OF FAILURE REPORT --------------
Error in if (class(duptest) == "try-error") { :
the condition has length > 1
Calls: mmtfa -> mmtfaEM
Execution halted
Flavor: r-devel-linux-x86_64-debian-gcc
Version: 0.3
Check: examples
Result: ERROR
Running examples in ‘mmtfa-Ex.R’ failed
The error most likely occurred in:
> ### Name: mmtfa
> ### Title: mmtfa: Function for Model-Based Clustering and Classification
> ### with Mixtures of Multivariate t Factor Analyzers
> ### Aliases: mmtfa
>
> ### ** Examples
>
> ###Note that only one model is run for each example
> ###in order to reduce computation time
>
> #Clustering iris data with hard random start
> tirisr <- mmtfa(iris[,-5], models="UUUU", Gs=1:3, Qs=1, init="hard")
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
mmtfa
--- call from context ---
mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
--- call from argument ---
if (class(duptest) == "try-error") {
break
}
--- R stacktrace ---
where 1: mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
where 2: mmtfa(iris[, -5], models = "UUUU", Gs = 1:3, Qs = 1, init = "hard")
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
{
modgen <- modelgen()
modold <- modgen$modold
p <- ncol(x)
n <- nrow(x)
zlist3 <- ll <- dff <- it <- store <- meanlist <- siglist <- list()
if (clas > 0) {
testindex <- sample(1:n, ceiling(n * (clas/100)))
kno <- vector(mode = "numeric", length = n)
kno[testindex] <- 1
unkno <- (kno - 1) * (-1)
Gs <- length(unique(known))
}
gvec <- 1:max(Gs)
qvec <- 1:max(Qs)
gstuff <- paste("G=", gvec, sep = "")
qstuff <- paste("Q=", qvec, sep = "")
bic <- rands <- icl <- logls <- array(-Inf, dim = c(length(models),
max(Qs), max(Gs)))
meansave <- meansave2 <- sigsave <- zmatsave <- sigsave2 <- zmatsave2 <- NA
oldmodvec <- modold[match(models, modgen$allmodels)]
if (class(init) != "list" && !(init %in% c("kmeans", "hard",
"disc", "cont", "soft", "uniform"))) {
stop("'init' must be one of 'kmeans', 'hard', 'soft', 'uniform' or a list. See ?mmtfa.")
}
zmatin <- list()
for (G in Gs) {
if (G == 1) {
zmatin[[G]] <- matrix(1, n, 1)
}
else {
if (is.character(init)) {
if (init == "hard") {
zmatin[[G]] <- discrandz(n, G)
}
if (init == "soft") {
zmatin[[G]] <- contrandz(n, G)
}
if (init == "uniform") {
if (clas > 0) {
zmatin[[G]] <- uniformz(n, G, clas, kno,
known)
}
else {
stop("Uniform initialization not available for clustering.")
return(NULL)
}
}
if (init == "kmeans") {
zmatin[[G]] <- kmeansz(x, n, G)
}
}
else {
zmatin[[G]] <- givenz(n, G, init[[G]])
}
}
}
for (modnum in 1:length(models)) {
modnew <- models[modnum]
mod <- modold[which(modgen$allmodels == modnew)]
for (G in Gs) {
delta <- matrix(0, n, G)
mug <- matrix(0, G, p)
om <- rep(0, G)
yg <- sigma <- sigmainv <- sg <- array(0, dim = c(p,
p, G))
tri <- array(0, dim = c(p, p, G))
w <- matrix(0, n, G)
for (q in Qs) {
lg <- array(0, dim = c(p, q, G))
betag <- array(0, dim = c(q, p, G))
thetag <- array(0, dim = c(q, q, G))
singular <- 0
breakit <- 0
if (G == 1) {
CCCCgroup <- c("UUCU", "UUCC", "UCCU", "UCCC",
"CUCU", "CUCC", "CCCU", "CCCC")
if (any(mod == CCCCgroup)) {
cccdum <- oldmodvec[oldmodvec %in% CCCCgroup]
if (length(cccdum) > 0) {
if (mod != cccdum[1]) {
breakit <- 1
}
}
}
CCUCgroup <- c("UUUU", "UUUC", "UCUU", "UCUC",
"CUUU", "CUUC", "CCUU", "CCUC", "Mt1U", "Mt1C",
"Mt2U", "Mt2C", "Mt3U", "Mt3C", "Mt4U", "Mt4C")
if (any(mod == CCUCgroup)) {
ccudum <- oldmodvec[oldmodvec %in% CCUCgroup]
if (length(ccudum) > 0) {
if (mod != ccudum[1]) {
breakit <- 1
}
}
}
}
if (breakit == 0) {
zmat <- zmatin[[G]]
vg <- vginit(dfstart, G)
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- muginit(G, p, x, zmat, ng)
sg <- sginit(p, G, x, mug, zmat, n, ng)
sgc <- sginitc(G, sg, pig, p, n, x)
if (substring(mod, 1, 3) == "CCC") {
for (g in 1:G) {
sg[, , g] <- sgc
}
}
if (substring(mod, 1, 1) == "U" | substring(mod,
3, 3) == "2" | substring(mod, 3, 3) == "4") {
lg <- lginitu(p, q, G, sg)
}
if (substring(mod, 1, 1) == "C" | substring(mod,
3, 3) == "1" | substring(mod, 3, 3) == "3") {
lg <- lginitc(p, q, G, sgc)
dumg <- lginitu(p, q, G, sg)
}
if (substring(mod, 3, 3) == "U") {
yg <- yginitu(p, G, sg, lg, mod, pig, dumg)
}
if (substring(mod, 3, 3) == "C") {
yg <- yginitc(p, G, sg, lg, sgc, mod, pig)
}
if (substring(mod, 1, 1) == "M") {
yg <- array(0, dim = c(p, p, G))
om <- rep(0, G)
tri <- array(0, dim = c(p, p, G))
ygst <- yginitu(p, G, sg, lg, mod, pig, dumg)
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
for (g in 1:G) {
om[g] <- det(ygst[, , g])^(1/p)
}
}
if (substring(mod, 3, 3) == "3" | substring(mod,
3, 3) == "4") {
dom <- 0
for (g in 1:G) {
dom <- dom + pig[g] * det(ygst[, , g])^(1/p)
}
om[] <- dom
}
for (g in 1:G) {
tri[, , g] <- ygst[, , g]/(det(ygst[, ,
g])^(1/p))
}
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
av <- diag(p) - diag(p)
for (g in 1:G) {
av <- av + pig[g] * tri[, , g]
}
tri[, , ] <- av
}
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
w <- winit(x, n, G, mug, sigmainv, vg, p, sg,
zmat)
cycle <- 0
dfnewg <- vg
}
conv <- 0
num <- matrix(0, n, G)
ft <- matrix(0, n, G)
logl <- NaN
while (conv != 1) {
if (breakit == 1) {
break
}
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- mugupdate(G, zmat, w, x, p, mug, n)
if (dfupdate == "approx") {
testing <- try(dfnewg <- dfupdatefun2(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
if (dfupdate == "numeric") {
testing <- try(dfnewg <- dfupdatefun(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
ng <- ngupdate(zmat)
sg <- sgupdate(p, G, n, x, mug, zmat, w, ng,
mod, pig, sg)
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
testing <- try(lg <- lgupdate(mod, p, q, G,
ng, yginv, sg, betag, thetag, om, tri, lg),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
if (substring(mod, 1, 1) == "M") {
om <- omupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri)
tri <- triupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri, ng)
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
else {
yg <- ygupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig)
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
duptest <- try(delta <- deltaup(x, mug, sigma,
sigmainv, G, n, delta), silent = TRUE)
if (class(duptest) == "try-error") {
break
}
suppressWarnings(zup <- zupdate(x, G, pig,
dfnewg, p, yg, q, betag, lg, mug, sigmainv,
n, clas, kno, known, unkno, delta))
zmat <- zup$zmat
if (any(is.nan(zmat))) {
break
}
w <- wupdate(x, n, G, mug, sigmainv, dfnewg,
p, delta, w)
cycle <- cycle + 1
logl[cycle] <- sum(log(rowSums(zup$num))) -
sum(zup$kcon)
if (is.na(logl[cycle])) {
break
}
if (cycle > 3) {
if (is.finite(logl[cycle - 2])) {
ak <- (logl[cycle] - logl[cycle - 1])/(logl[cycle -
1] - logl[cycle - 2])
linf <- logl[cycle - 1] + (logl[cycle] -
logl[cycle - 1])/(1 - ak)
if (abs(linf - logl[cycle - 1]) < eps) {
conv <- 1
}
if ((logl[cycle] - logl[cycle - 1]) < 0) {
break
}
}
else {
break
}
}
}
if (conv == 1) {
bic[modnum, q, G] <- bicdum <- BICcalc(conv,
G, p, mod, q, logl, n, gauss)
icl[modnum, q, G] <- icldum <- ICLcalc(conv,
n, zmat, bic, modnum, q, G)
if (bicdum == max(bic)) {
meansave <- mug
sigsave <- sigma
zmatsave <- zmat
dfsave <- dfnewg
itsave <- cycle
llsave <- logl[cycle]
}
if (icldum == max(icl)) {
meansave2 <- mug
sigsave2 <- sigma
zmatsave2 <- zmat
dfsave2 <- dfnewg
itsave2 <- cycle
llsave2 <- logl[cycle]
}
}
}
}
}
dimnames(bic) <- list(models, qstuff, gstuff)
dimnames(icl) <- list(models, qstuff, gstuff)
maxes <- which(bic == max(bic), arr.ind = TRUE)
maxicl <- which(icl == max(icl), arr.ind = TRUE)
if (nrow(maxes) > 1) {
message("WARNING: Maximum BIC tie between two or more models")
bestmodnum <- maxes[1:nrow(maxes), 1]
bestmod <- models[bestmodnum]
bestq <- maxes[1:nrow(maxes), 2]
bestg <- maxes[1:nrow(maxes), 3]
itf <- "MULTIPLE"
dff1 <- "MULTIPLE"
bestz <- "MULTIPLE"
bestzmap <- "MULTIPLE"
tab <- "MULTIPLE"
blogl <- "Multiple"
}
if (nrow(maxes) == 1) {
bestmodnum <- maxes[1]
bestmod <- models[bestmodnum]
bestq <- maxes[2]
bestg <- maxes[3]
bestz <- zmatsave
dff1 <- dfsave
itf <- itsave
blogl <- llsave
bestzmap <- apply(bestz, 1, which.max)
if (clas > 0) {
newmap <- bestzmap
newmap[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tab <- table(known, newmap)
}
else {
if (!is.null(known)) {
tab <- table(known, bestzmap)
}
else {
tab <- NULL
}
}
}
if (nrow(maxicl) > 1) {
message("WARNING: Maximum ICL tie between two or more models")
bestmodnumicl <- maxicl[1:nrow(maxicl), 1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[1:nrow(maxicl), 2]
bestgicl <- maxicl[1:nrow(maxicl), 3]
dff1icl <- "MULTIPLE"
bestzicl <- "MULTIPLE"
bestzmapicl <- "MULTIPLE"
itficl <- "MULTIPLE"
tabicl <- "MULTIPLE"
bloglicl <- "MULTIPLE"
}
if (nrow(maxicl) == 1) {
bestmodnumicl <- maxicl[1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[2]
bestgicl <- maxicl[3]
bestzicl <- zmatsave2
dff1icl <- dfsave2
itficl <- itsave2
bloglicl <- llsave2
bestzmapicl <- apply(bestzicl, 1, which.max)
if (clas > 0) {
newmapicl <- bestzmapicl
newmapicl[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tabicl <- table(known, newmapicl)
}
else {
if (!is.null(known)) {
tabicl <- table(known, bestzmapicl)
}
else {
tabicl <- NULL
}
}
}
iclresults <- list()
par <- list()
paricl <- list()
par[["mean"]] <- meansave
par[["sigma"]] <- sigsave
par[["df"]] <- dff1
paricl[["sigma"]] <- sigsave2
paricl[["mean"]] <- meansave2
paricl[["df"]] <- dff1icl
store[["parameters"]] <- par
store[["allbic"]] <- bic[, Qs, Gs]
iclresults[["allicl"]] <- icl[, Qs, Gs]
store[["bic"]] <- max(bic)
iclresults[["icl"]] <- max(icl)
store[["modelname"]] <- bestmod
store[["bestmodel"]] <- paste("The best model (BIC of ",
round(max(bic), 2), ") is ", bestmod, " with G=", bestg,
sep = "")
store[["Q"]] <- bestq
store[["G"]] <- bestg
store[["classification"]] <- bestzmap
iclresults[["bestmodel"]] <- paste("The best model (ICL of ",
round(max(icl), 2), ") is ", bestmodicl, " with G=",
bestgicl, sep = "")
iclresults[["modelname"]] <- bestmodicl
iclresults[["Q"]] <- bestqicl
iclresults[["G"]] <- bestgicl
iclresults[["fuzzy"]] <- bestzicl
iclresults[["logl"]] <- bloglicl
iclresults[["classification"]] <- bestzmapicl
iclresults[["parameters"]] <- paricl
store[["tab"]] <- tab
iclresults[["tab"]] <- tabicl
store[["iter"]] <- itf
iclresults[["iter"]] <- itficl
store[["x"]] <- x
store[["fuzzy"]] <- bestz
store[["logl"]] <- blogl
store[["iclresults"]] <- iclresults
store
}
<bytecode: 0x3d36bd0>
<environment: namespace:mmtfa>
--- function search by body ---
Function mmtfaEM in namespace mmtfa has this body.
----------- END OF FAILURE REPORT --------------
Error in if (class(duptest) == "try-error") { :
the condition has length > 1
Calls: mmtfa -> mmtfaEM
Execution halted
Flavor: r-devel-linux-x86_64-fedora-clang
Version: 0.3
Check: examples
Result: ERROR
Running examples in ‘mmtfa-Ex.R’ failed
The error most likely occurred in:
> ### Name: mmtfa
> ### Title: mmtfa: Function for Model-Based Clustering and Classification
> ### with Mixtures of Multivariate t Factor Analyzers
> ### Aliases: mmtfa
>
> ### ** Examples
>
> ###Note that only one model is run for each example
> ###in order to reduce computation time
>
> #Clustering iris data with hard random start
> tirisr <- mmtfa(iris[,-5], models="UUUU", Gs=1:3, Qs=1, init="hard")
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
mmtfa
--- call from context ---
mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
--- call from argument ---
if (class(duptest) == "try-error") {
break
}
--- R stacktrace ---
where 1: mmtfaEM(x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
where 2: mmtfa(iris[, -5], models = "UUUU", Gs = 1:3, Qs = 1, init = "hard")
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (x, Gs, Qs, clas, init, scale, models, dfstart, dfupdate,
gauss, eps, known)
{
modgen <- modelgen()
modold <- modgen$modold
p <- ncol(x)
n <- nrow(x)
zlist3 <- ll <- dff <- it <- store <- meanlist <- siglist <- list()
if (clas > 0) {
testindex <- sample(1:n, ceiling(n * (clas/100)))
kno <- vector(mode = "numeric", length = n)
kno[testindex] <- 1
unkno <- (kno - 1) * (-1)
Gs <- length(unique(known))
}
gvec <- 1:max(Gs)
qvec <- 1:max(Qs)
gstuff <- paste("G=", gvec, sep = "")
qstuff <- paste("Q=", qvec, sep = "")
bic <- rands <- icl <- logls <- array(-Inf, dim = c(length(models),
max(Qs), max(Gs)))
meansave <- meansave2 <- sigsave <- zmatsave <- sigsave2 <- zmatsave2 <- NA
oldmodvec <- modold[match(models, modgen$allmodels)]
if (class(init) != "list" && !(init %in% c("kmeans", "hard",
"disc", "cont", "soft", "uniform"))) {
stop("'init' must be one of 'kmeans', 'hard', 'soft', 'uniform' or a list. See ?mmtfa.")
}
zmatin <- list()
for (G in Gs) {
if (G == 1) {
zmatin[[G]] <- matrix(1, n, 1)
}
else {
if (is.character(init)) {
if (init == "hard") {
zmatin[[G]] <- discrandz(n, G)
}
if (init == "soft") {
zmatin[[G]] <- contrandz(n, G)
}
if (init == "uniform") {
if (clas > 0) {
zmatin[[G]] <- uniformz(n, G, clas, kno,
known)
}
else {
stop("Uniform initialization not available for clustering.")
return(NULL)
}
}
if (init == "kmeans") {
zmatin[[G]] <- kmeansz(x, n, G)
}
}
else {
zmatin[[G]] <- givenz(n, G, init[[G]])
}
}
}
for (modnum in 1:length(models)) {
modnew <- models[modnum]
mod <- modold[which(modgen$allmodels == modnew)]
for (G in Gs) {
delta <- matrix(0, n, G)
mug <- matrix(0, G, p)
om <- rep(0, G)
yg <- sigma <- sigmainv <- sg <- array(0, dim = c(p,
p, G))
tri <- array(0, dim = c(p, p, G))
w <- matrix(0, n, G)
for (q in Qs) {
lg <- array(0, dim = c(p, q, G))
betag <- array(0, dim = c(q, p, G))
thetag <- array(0, dim = c(q, q, G))
singular <- 0
breakit <- 0
if (G == 1) {
CCCCgroup <- c("UUCU", "UUCC", "UCCU", "UCCC",
"CUCU", "CUCC", "CCCU", "CCCC")
if (any(mod == CCCCgroup)) {
cccdum <- oldmodvec[oldmodvec %in% CCCCgroup]
if (length(cccdum) > 0) {
if (mod != cccdum[1]) {
breakit <- 1
}
}
}
CCUCgroup <- c("UUUU", "UUUC", "UCUU", "UCUC",
"CUUU", "CUUC", "CCUU", "CCUC", "Mt1U", "Mt1C",
"Mt2U", "Mt2C", "Mt3U", "Mt3C", "Mt4U", "Mt4C")
if (any(mod == CCUCgroup)) {
ccudum <- oldmodvec[oldmodvec %in% CCUCgroup]
if (length(ccudum) > 0) {
if (mod != ccudum[1]) {
breakit <- 1
}
}
}
}
if (breakit == 0) {
zmat <- zmatin[[G]]
vg <- vginit(dfstart, G)
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- muginit(G, p, x, zmat, ng)
sg <- sginit(p, G, x, mug, zmat, n, ng)
sgc <- sginitc(G, sg, pig, p, n, x)
if (substring(mod, 1, 3) == "CCC") {
for (g in 1:G) {
sg[, , g] <- sgc
}
}
if (substring(mod, 1, 1) == "U" | substring(mod,
3, 3) == "2" | substring(mod, 3, 3) == "4") {
lg <- lginitu(p, q, G, sg)
}
if (substring(mod, 1, 1) == "C" | substring(mod,
3, 3) == "1" | substring(mod, 3, 3) == "3") {
lg <- lginitc(p, q, G, sgc)
dumg <- lginitu(p, q, G, sg)
}
if (substring(mod, 3, 3) == "U") {
yg <- yginitu(p, G, sg, lg, mod, pig, dumg)
}
if (substring(mod, 3, 3) == "C") {
yg <- yginitc(p, G, sg, lg, sgc, mod, pig)
}
if (substring(mod, 1, 1) == "M") {
yg <- array(0, dim = c(p, p, G))
om <- rep(0, G)
tri <- array(0, dim = c(p, p, G))
ygst <- yginitu(p, G, sg, lg, mod, pig, dumg)
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
for (g in 1:G) {
om[g] <- det(ygst[, , g])^(1/p)
}
}
if (substring(mod, 3, 3) == "3" | substring(mod,
3, 3) == "4") {
dom <- 0
for (g in 1:G) {
dom <- dom + pig[g] * det(ygst[, , g])^(1/p)
}
om[] <- dom
}
for (g in 1:G) {
tri[, , g] <- ygst[, , g]/(det(ygst[, ,
g])^(1/p))
}
if (substring(mod, 3, 3) == "1" | substring(mod,
3, 3) == "2") {
av <- diag(p) - diag(p)
for (g in 1:G) {
av <- av + pig[g] * tri[, , g]
}
tri[, , ] <- av
}
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
w <- winit(x, n, G, mug, sigmainv, vg, p, sg,
zmat)
cycle <- 0
dfnewg <- vg
}
conv <- 0
num <- matrix(0, n, G)
ft <- matrix(0, n, G)
logl <- NaN
while (conv != 1) {
if (breakit == 1) {
break
}
ng <- ngupdate(zmat)
pig <- pigupdate(ng, n)
mug <- mugupdate(G, zmat, w, x, p, mug, n)
if (dfupdate == "approx") {
testing <- try(dfnewg <- dfupdatefun2(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
if (dfupdate == "numeric") {
testing <- try(dfnewg <- dfupdatefun(mod,
dfnewg, ng, zmat, w, G, p, n, x, mug, sigmainv),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
}
ng <- ngupdate(zmat)
sg <- sgupdate(p, G, n, x, mug, zmat, w, ng,
mod, pig, sg)
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
testing <- try(lg <- lgupdate(mod, p, q, G,
ng, yginv, sg, betag, thetag, om, tri, lg),
silent = TRUE)
if (!all(is.finite(testing))) {
break
}
if (substring(mod, 1, 1) == "M") {
om <- omupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri)
tri <- triupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig, om, tri, ng)
for (g in 1:G) {
yg[, , g] <- om[g] * tri[, , g]
}
}
else {
yg <- ygupdate(mod, q, G, yg, p, sg, lg,
betag, thetag, pig)
}
yginv <- yginvup(p, G, yg)
sigma <- sigmaup(p, G, lg, yg, sigma)
testing <- try(sigmainv <- sigmainvup(p, G,
yginv, lg, q, sigmainv), silent = TRUE)
if (!all(is.finite(testing))) {
break
}
betag <- betagup(q, p, G, lg, sigmainv, betag)
thetag <- thetagup(q, G, betag, lg, sg, thetag)
duptest <- try(delta <- deltaup(x, mug, sigma,
sigmainv, G, n, delta), silent = TRUE)
if (class(duptest) == "try-error") {
break
}
suppressWarnings(zup <- zupdate(x, G, pig,
dfnewg, p, yg, q, betag, lg, mug, sigmainv,
n, clas, kno, known, unkno, delta))
zmat <- zup$zmat
if (any(is.nan(zmat))) {
break
}
w <- wupdate(x, n, G, mug, sigmainv, dfnewg,
p, delta, w)
cycle <- cycle + 1
logl[cycle] <- sum(log(rowSums(zup$num))) -
sum(zup$kcon)
if (is.na(logl[cycle])) {
break
}
if (cycle > 3) {
if (is.finite(logl[cycle - 2])) {
ak <- (logl[cycle] - logl[cycle - 1])/(logl[cycle -
1] - logl[cycle - 2])
linf <- logl[cycle - 1] + (logl[cycle] -
logl[cycle - 1])/(1 - ak)
if (abs(linf - logl[cycle - 1]) < eps) {
conv <- 1
}
if ((logl[cycle] - logl[cycle - 1]) < 0) {
break
}
}
else {
break
}
}
}
if (conv == 1) {
bic[modnum, q, G] <- bicdum <- BICcalc(conv,
G, p, mod, q, logl, n, gauss)
icl[modnum, q, G] <- icldum <- ICLcalc(conv,
n, zmat, bic, modnum, q, G)
if (bicdum == max(bic)) {
meansave <- mug
sigsave <- sigma
zmatsave <- zmat
dfsave <- dfnewg
itsave <- cycle
llsave <- logl[cycle]
}
if (icldum == max(icl)) {
meansave2 <- mug
sigsave2 <- sigma
zmatsave2 <- zmat
dfsave2 <- dfnewg
itsave2 <- cycle
llsave2 <- logl[cycle]
}
}
}
}
}
dimnames(bic) <- list(models, qstuff, gstuff)
dimnames(icl) <- list(models, qstuff, gstuff)
maxes <- which(bic == max(bic), arr.ind = TRUE)
maxicl <- which(icl == max(icl), arr.ind = TRUE)
if (nrow(maxes) > 1) {
message("WARNING: Maximum BIC tie between two or more models")
bestmodnum <- maxes[1:nrow(maxes), 1]
bestmod <- models[bestmodnum]
bestq <- maxes[1:nrow(maxes), 2]
bestg <- maxes[1:nrow(maxes), 3]
itf <- "MULTIPLE"
dff1 <- "MULTIPLE"
bestz <- "MULTIPLE"
bestzmap <- "MULTIPLE"
tab <- "MULTIPLE"
blogl <- "Multiple"
}
if (nrow(maxes) == 1) {
bestmodnum <- maxes[1]
bestmod <- models[bestmodnum]
bestq <- maxes[2]
bestg <- maxes[3]
bestz <- zmatsave
dff1 <- dfsave
itf <- itsave
blogl <- llsave
bestzmap <- apply(bestz, 1, which.max)
if (clas > 0) {
newmap <- bestzmap
newmap[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tab <- table(known, newmap)
}
else {
if (!is.null(known)) {
tab <- table(known, bestzmap)
}
else {
tab <- NULL
}
}
}
if (nrow(maxicl) > 1) {
message("WARNING: Maximum ICL tie between two or more models")
bestmodnumicl <- maxicl[1:nrow(maxicl), 1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[1:nrow(maxicl), 2]
bestgicl <- maxicl[1:nrow(maxicl), 3]
dff1icl <- "MULTIPLE"
bestzicl <- "MULTIPLE"
bestzmapicl <- "MULTIPLE"
itficl <- "MULTIPLE"
tabicl <- "MULTIPLE"
bloglicl <- "MULTIPLE"
}
if (nrow(maxicl) == 1) {
bestmodnumicl <- maxicl[1]
bestmodicl <- models[bestmodnumicl]
bestqicl <- maxicl[2]
bestgicl <- maxicl[3]
bestzicl <- zmatsave2
dff1icl <- dfsave2
itficl <- itsave2
bloglicl <- llsave2
bestzmapicl <- apply(bestzicl, 1, which.max)
if (clas > 0) {
newmapicl <- bestzmapicl
newmapicl[testindex] <- NA
newknown <- known
newknown[testindex] <- NA
tabicl <- table(known, newmapicl)
}
else {
if (!is.null(known)) {
tabicl <- table(known, bestzmapicl)
}
else {
tabicl <- NULL
}
}
}
iclresults <- list()
par <- list()
paricl <- list()
par[["mean"]] <- meansave
par[["sigma"]] <- sigsave
par[["df"]] <- dff1
paricl[["sigma"]] <- sigsave2
paricl[["mean"]] <- meansave2
paricl[["df"]] <- dff1icl
store[["parameters"]] <- par
store[["allbic"]] <- bic[, Qs, Gs]
iclresults[["allicl"]] <- icl[, Qs, Gs]
store[["bic"]] <- max(bic)
iclresults[["icl"]] <- max(icl)
store[["modelname"]] <- bestmod
store[["bestmodel"]] <- paste("The best model (BIC of ",
round(max(bic), 2), ") is ", bestmod, " with G=", bestg,
sep = "")
store[["Q"]] <- bestq
store[["G"]] <- bestg
store[["classification"]] <- bestzmap
iclresults[["bestmodel"]] <- paste("The best model (ICL of ",
round(max(icl), 2), ") is ", bestmodicl, " with G=",
bestgicl, sep = "")
iclresults[["modelname"]] <- bestmodicl
iclresults[["Q"]] <- bestqicl
iclresults[["G"]] <- bestgicl
iclresults[["fuzzy"]] <- bestzicl
iclresults[["logl"]] <- bloglicl
iclresults[["classification"]] <- bestzmapicl
iclresults[["parameters"]] <- paricl
store[["tab"]] <- tab
iclresults[["tab"]] <- tabicl
store[["iter"]] <- itf
iclresults[["iter"]] <- itficl
store[["x"]] <- x
store[["fuzzy"]] <- bestz
store[["logl"]] <- blogl
store[["iclresults"]] <- iclresults
store
}
<bytecode: 0x26194f8>
<environment: namespace:mmtfa>
--- function search by body ---
Function mmtfaEM in namespace mmtfa has this body.
----------- END OF FAILURE REPORT --------------
Error in if (class(duptest) == "try-error") { :
the condition has length > 1
Calls: mmtfa -> mmtfaEM
Execution halted
Flavor: r-devel-linux-x86_64-fedora-gcc