Nothing
#selection
select_drop = function(x, m, K, id0, s.id, BIC.all, EM.iter, nk.min, tol, silent, orders, Parallel, n.cores)
{
p = ncol(x)
n = nrow(x)
BIC = s.id
for (k in 1:K) {
BIC[[k]][which(!BIC[[k]])] = NA
}
#cat("begin dropping, current BIC are ", BIC.all,"\n")
id.true = NULL
flag = NULL
id.true1 = NULL
if((!Parallel)||K ==1)
{for (k in 1:K) {
Id.true = which(s.id[[k]])
id.true = append(id.true,list(Id.true))
for (i in 1:length(Id.true)){
j = Id.true[i]%/%(p^2)+1
j1 = ifelse((Id.true[i]%%(p^2))%%p == 0, p, (Id.true[i]%%(p^2))%%p)
j2 = ifelse((Id.true[i]%%(p^2))%%p == 0,(Id.true[i]%%(p^2))%/%p, (Id.true[i]%%(p^2))%/%p+1 )
s.id.new = s.id
s.id.new[[k]][j1,j2,j] = FALSE
obj.new = EM_c(x = x, m, K, s.id = s.id.new, id0 = id0, EM.iter, tol)
BIC[[k]][j1,j2,j] = obj.new$BIC
if(min(table( obj.new$id)) < nk.min || length(table( obj.new$id))!=K) BIC[[k]][j1,j2,j] = BIC.all+10
else if(is.na(obj.new$llk)||is.infinite(obj.new$llk)) BIC[[k]][j1,j2,j]= BIC.all+10
else if(min(obj.new$s2, na.rm = TRUE)< 1e-07 || any(is.na(obj.new$s2))) BIC[[k]][j1,j2,j]= BIC.all+10
#cat("BIC are ", BIC[[k]][j1,j2,j] ,"\n")
}
flag = append(flag, all(is.na(BIC[[k]])))
}
}else{
chk = Sys.getenv("_R_CHECK_LIMIT_CORES_", "")
if(nzchar(chk) && chk == "TRUE"){
n.cores = 2L
}else{
n.cores = ifelse(n.cores > detectCores(), detectCores(), n.cores)
}
cl = makeCluster(n.cores)
clusterSetRNGStream(cl, 123)
clusterExport(cl, varlist=c("x", "m","K", "p", "s.id","nk.min","id0","BIC","BIC.all","id.true1", "EM.iter", "tol","EM_c"), envir=environment())
result = parLapply(cl, 1:K, function(k)
{
Id.true = which(s.id[[k]])
id.true1 = append(id.true1,list(Id.true))
BIC1 = BIC[[k]]
for (i in 1:length(Id.true)){
j = Id.true[i]%/%(p^2)+1
j1 = ifelse((Id.true[i]%%(p^2))%%p == 0, p, (Id.true[i]%%(p^2))%%p)
j2 = ifelse((Id.true[i]%%(p^2))%%p == 0,(Id.true[i]%%(p^2))%/%p, (Id.true[i]%%(p^2))%/%p+1 )
s.id.new = s.id
s.id.new[[k]][j1,j2,j] = FALSE
obj.new = EM_c(x = x, m, K, s.id = s.id.new, id0 = id0, EM.iter, tol)
BIC1[j1,j2,j] = obj.new$BIC
if(min(table( obj.new$id)) < nk.min || length(table( obj.new$id))!=K) BIC1[j1,j2,j] = BIC.all+10
else if(is.na(obj.new$llk)||is.infinite(obj.new$llk)) BIC1[j1,j2,j]= BIC.all+10
else if(min(obj.new$s2, na.rm = TRUE)< 1e-07 || any(is.na(obj.new$s2))) BIC1[j1,j2,j]= BIC.all+10
#cat("BIC are ", BIC[[k]][j1,j2,j] ,"\n")
}
return(list(BIC1 = BIC1, id.true1 = id.true1))
})
stopCluster(cl)
for (k in 1:K){
BIC[[k]] = result[[k]]$BIC1
id.true = append(id.true, result[[k]]$id.true1)
flag = append(flag, all(is.na(BIC[[k]])))
}
}
if(!all(flag)){
id = numeric(K)
bic = numeric(K)
for (k in 1:K) {
id[k] = which.min(c(BIC.all, BIC[[k]][which(!is.na(BIC[[k]]))]))
bic[k] = min(c(BIC.all, BIC[[k]][which(!is.na(BIC[[k]]))]), na.rm = TRUE)
}
if(!max(id) == 1){
k = which.min(bic)
BIC.all = min(bic)
i = id.true[[k]][id[k]-1]
j = i%/%(p^2)+1
j1 = ifelse((i%%(p^2))%%p == 0, p, (i%%(p^2))%%p)
j2 = ifelse((i%%(p^2))%%p == 0, (i%%(p^2))%/%p, (i%%(p^2))%/%p+1 )
s.id[[k]][j1,j2,j] = FALSE
obj = EM_c(x = x, m, K, s.id = s.id, id0 = id0, EM.iter, tol)
id0 = obj$id
if(!silent) cat(paste0("Current BIC is ", BIC[[k]][j1,j2,j], "."), paste0("In cluster ", k, ", variable"), paste0("x", orders[j2],"^",j), "for response", paste0("x", orders[j1]), "is excluded", "\n")
}else{
if(!silent) cat(paste0("Current BIC is ", BIC.all,". No variable is excluded"), "\n")
obj = EM_c(x = x, m, K, s.id = s.id, id0 = id0, EM.iter, tol)
id0 = obj$id
}
}else{
if(!silent) cat(paste0("Current BIC is ", BIC.all,". No variable is excluded"), "\n")
obj = EM_c(x = x, m, K, s.id = s.id, id0 = id0, EM.iter, tol)
id0 = obj$id
}
return(list(BIC = BIC.all, llk = obj$llk, id = id0, Pi = obj$Pi, s.id = s.id, beta = obj$beta, s2 = obj$s2, tau = obj$tau, n_pars = obj$n_pars ))
}
select_add = function(x, m, K, id0, s.id, BIC.all, EM.iter, nk.min, tol, silent, orders, Parallel, n.cores){
n = nrow(x)
p = ncol(x)
n = nrow(x)
BIC = s.id
for (k in 1:K) {
BIC[[k]][which(!BIC[[k]])] = NA
BIC[[k]][which(BIC[[k]])] = NA
}
#cat("begin adding BIC are ", BIC.all, "\n")
s.id.1 = array(FALSE, c(p, p, m))
for (i in 1:m) {
s.id.1[lower.tri(s.id.1[,,i])] = TRUE
}
id.true = NULL
id.true1 = NULL
flag = NULL
if((!Parallel)||K ==1)
{for (k in 1:K)
{
Id.true = which(s.id.1)[!which(s.id.1)%in% which(s.id[[k]])]
id.true = append(id.true,list(Id.true))
for (i in 1:length(Id.true)){
j = Id.true[i]%/%(p^2)+1
j1 = ifelse((Id.true[i]%%(p^2))%%p == 0, p, (Id.true[i]%%(p^2))%%p)
j2 = ifelse((Id.true[i]%%(p^2))%%p == 0,(Id.true[i]%%(p^2))%/%p, (Id.true[i]%%(p^2))%/%p+1 )
s.id.new = s.id
s.id.new[[k]][j1,j2,j] = TRUE
obj.new = EM_c(x = x, m, K, s.id = s.id.new, id0 = id0, EM.iter, tol)
BIC[[k]][j1,j2,j] = obj.new$BIC
if(min(table( obj.new$id)) < nk.min || length(table( obj.new$id))!=K) BIC[[k]][j1,j2,j] = BIC.all+10
else if(is.na(obj.new$llk)||is.infinite(obj.new$llk)) BIC[[k]][j1,j2,j]= BIC.all+10
else if(min(obj.new$s2, na.rm = TRUE)< 1e-07 || any(is.na(obj.new$s2))) BIC[[k]][j1,j2,j]= BIC.all+10
# cat(" BIC are ", BIC[[k]][j1,j2,j], "\n")
}
flag = append(flag, all(is.na(BIC[[k]])))
}
}else{
chk = Sys.getenv("_R_CHECK_LIMIT_CORES_", "")
if(nzchar(chk) && chk == "TRUE"){
n.cores = 2L
}else{
n.cores = ifelse(n.cores > detectCores(), detectCores(), n.cores)
}
cl = makeCluster(n.cores)
clusterSetRNGStream(cl, 123)
clusterExport(cl, varlist=c("x", "m","n","K", "p", "s.id","s.id.1","nk.min","id0","BIC","BIC.all","id.true1", "EM.iter", "tol","EM_c"), envir=environment())
result = parLapply(cl, 1:K, function(k){
BIC1 = BIC[[k]]
Id.true = which(s.id.1)[!which(s.id.1)%in% which(s.id[[k]])]
id.true1 = append(id.true1,list(Id.true))
for (i in 1:length(Id.true)){
j = Id.true[i]%/%(p^2)+1
j1 = ifelse((Id.true[i]%%(p^2))%%p == 0, p, (Id.true[i]%%(p^2))%%p)
j2 = ifelse((Id.true[i]%%(p^2))%%p == 0,(Id.true[i]%%(p^2))%/%p, (Id.true[i]%%(p^2))%/%p+1 )
s.id.new = s.id
s.id.new[[k]][j1,j2,j] = TRUE
obj.new = EM_c(x = x, m, K, s.id = s.id.new, id0 = id0, EM.iter, tol)
BIC1[j1,j2,j] = obj.new$BIC
if(min(table( obj.new$id)) < nk.min || length(table( obj.new$id))!=K) BIC1[j1,j2,j]= BIC.all+10
else if(is.na(obj.new$llk)||is.infinite(obj.new$llk)) BIC1[j1,j2,j]= BIC.all+10
else if(min(obj.new$s2, na.rm = TRUE)< 1e-07 || any(is.na(obj.new$s2))) BIC1[j1,j2,j]= BIC.all+10
# cat(" BIC are ", BIC[[k]][j1,j2,j], "\n")
}
return(list(BIC1 = BIC1, id.true1 = id.true1))
})
stopCluster(cl)
for (k in 1:K){
BIC[[k]] = result[[k]]$BIC1
id.true = append(id.true, result[[k]]$id.true1)
flag = append(flag, all(is.na(BIC[[k]])))
}
}
if(!all(flag)){
id = numeric(K)
bic = numeric(K)
for (k in 1:K) {
id[k] = which.min(c(BIC.all, BIC[[k]][which(!is.na(BIC[[k]]))]))
bic[k] = min(c(BIC.all, BIC[[k]][which(!is.na(BIC[[k]]))]), na.rm = TRUE)
}
if(!max(id) == 1){
k = which.min(bic)
BIC.all = min(bic)
i = id.true[[k]][id[k]-1]
j = i%/%(p^2)+1
j1 = ifelse((i%%(p^2))%%p == 0, p, (i%%(p^2))%%p)
j2 = ifelse((i%%(p^2))%%p == 0, (i%%(p^2))%/%p, (i%%(p^2))%/%p+1 )
s.id[[k]][j1,j2,j] = TRUE
obj = EM_c(x = x, m, K, s.id = s.id, id0 = id0, EM.iter, tol)
id0 = obj$id
if(!silent) cat(paste0("Current BIC is ", BIC[[k]][j1,j2,j], "."), paste0("In cluster ", k, ", variable"), paste0("x", orders[j2],"^",j), "for response", paste0("x", orders[j1]), "is included", "\n")
}else{
if(!silent) cat(paste0("Current BIC is ", BIC.all,". No variable is included"), "\n")
obj = EM_c(x = x, m, K, s.id = s.id, id0 = id0, EM.iter, tol)
id0 = obj$id
}
}else{
if(!silent) cat(paste0("Current BIC is ", BIC.all,". No variable is included"), "\n")
obj = EM_c(x = x, m, K, s.id = s.id, id0 = id0, EM.iter, tol)
id0 = obj$id
}
return(list(BIC = BIC.all, llk = obj$llk, id = id0, Pi = obj$Pi, s.id = s.id, beta = obj$beta, s2 = obj$s2, tau = obj$tau, n_pars = obj$n_pars ))
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.