From: Benjamin Auder Date: Mon, 28 Oct 2019 21:38:47 +0000 (+0100) Subject: Fix default weights, fix reporting scripts X-Git-Url: https://git.auder.net/%7B%7B%20asset%28%27mixstore/css/img/pieces/scripts/%7B%7B%20pkg.url%20%7D%7D?a=commitdiff_plain;h=0ad4c8de650e9f27ec3754c9cb9b2a03db5aff24;p=morpheus.git Fix default weights, fix reporting scripts --- diff --git a/pkg/R/optimParams.R b/pkg/R/optimParams.R index 85c21e7..06d1684 100644 --- a/pkg/R/optimParams.R +++ b/pkg/R/optimParams.R @@ -59,7 +59,7 @@ optimParams = function(K, link=c("logit","probit"), optargs=list()) weights <- optargs$weights if (is.null(weights)) - weights <- rep(1, K) + weights <- rep(1, 3) # Build and return optimization algorithm object methods::new("OptimParams", "li"=link, "M1"=as.double(M[[1]]), diff --git a/reports/accuracy.R b/reports/accuracy.R index e9c9d1b..91d9c61 100644 --- a/reports/accuracy.R +++ b/reports/accuracy.R @@ -11,11 +11,10 @@ optimBeta <- function(N, n, K, p, beta, b, link, weights, ncores) M <- computeMoments(fargs$X, fargs$Y) fargs$optargs$M <- M mu <- computeMu(fargs$X, fargs$Y, fargs$optargs) - res2 <- NULL + op <- optimParams(K,fargs$optargs$link,fargs$optargs) + x_init <- list(p=rep(1/K,K-1), beta=mu, b=rep(0,K)) tryCatch({ - op <- optimParams(K,fargs$optargs$link,fargs$optargs) - x_init <- list(p=rep(1/K,K-1), beta=mu, b=rep(0,K)) - res2 <- do.call(rbind, op$run(x_init)) + res2 <- do.call(rbind, op$run(x_init)) }, error = function(e) { res2 <- NA }) @@ -50,7 +49,6 @@ optimBeta <- function(N, n, K, p, beta, b, link, weights, ncores) fargs$X = io$X fargs$Y = io$Y fargs$optargs$K = ncol(fargs$beta) - fargs$optargs$M = computeMoments(io$X,io$Y) fargs }, N=N, ncores=ncores, verbose=TRUE) p <- c(p, 1-sum(p)) @@ -119,7 +117,7 @@ if (d == 2) { } mr <- optimBeta(N, n, K, p, beta, b, link, weights, ncores) -mr_params <- list("N"=N, "n"=n, "K"=K, "d"=d, "link"=link, +mr_params <- list("N"=N, "nc"=ncores, "n"=n, "K"=K, "d"=d, "link"=link, "p"=c(p,1-sum(p)), "beta"=beta, "b"=b, "weights"=weights) save("mr", "mr_params", file=paste("res_",n,"_",d,"_",link,"_",strw,".RData",sep="")) diff --git a/reports/multistart.R b/reports/multistart.R index 9971cd2..5fa80ad 100644 --- a/reports/multistart.R +++ b/reports/multistart.R @@ -1,12 +1,76 @@ library(morpheus) +testMultistart <- function(N, n, K, p, beta, b, link, nstart, ncores) +{ + ms <- multiRun( + list(n=n,p=p,beta=beta,b=b,optargs=list(K=K,link=link,nstart=nstart)), + list( + function(fargs) { + # 1 start + library(morpheus) + K <- fargs$optargs$K + op <- optimParams(K, fargs$optargs$link, fargs$optargs) + x_init <- list(p=rep(1/K,K-1), beta=fargs$mu, b=rep(0,K)) + res <- NULL + tryCatch({ + res <- do.call(rbind, op$run(x_init)) + }, error = function(e) { + res <- NA + }) + res + }, + function(fargs) { + # B starts + library(morpheus) + K <- fargs$optargs$K + op <- optimParams(K, fargs$optargs$link, fargs$optargs) + best_val <- Inf + best_par <- list() + for (i in 1:fargs$optargs$nstart) + { + #x_init <- list(p=rep(1/K,K-1), beta=i*fargs$mu, b=rep(0,K)) + M <- matrix(rnorm(d*K), nrow=d, ncol=K) + M <- t(t(M) / sqrt(colSums(M^2))) + x_init <- list(p=rep(1/K,K-1), beta=M, b=rep(0,K)) + tryCatch({ + par <- op$run(x_init) + }, error = function(e) { + par <- NA + }) + if (!is.na(par[0])) + { + val <- op$f( op$linArgs(par) ) + if (val < best_val) + { + best_par <- par + best_val <- val + } + } + } + # Bet that at least one run succeded: + do.call(rbind,best_par) + } + ), + prepareArgs = function(fargs, index) { + library(morpheus) + io = generateSampleIO(fargs$n, fargs$p, fargs$beta, fargs$b, fargs$optargs$link) + fargs$optargs$M <- computeMoments(io$X, io$Y) + mu <- computeMu(io$X, io$Y, fargs$optargs) + fargs$mu <- mu + fargs + }, N=N, ncores=ncores, verbose=TRUE) + for (i in 1:2) + ms[[i]] <- alignMatrices(ms[[i]], ref=rbind(p,beta,b), ls_mode="exact") + ms +} + #model = binomial K <- 2 p <- .5 b <- c(-.2, .5) # Default values: link = "logit" -N <- 100 +N <- 10 d <- 2 n <- 1e4 ncores <- 1 @@ -39,47 +103,8 @@ betas <- list( matrix( c(1,2,-1,0,3,4,-1,-3,0,2, 2,-3,0,1,0,-1,-4,3,2,0), ncol=K ) ) #d=10 beta <- betas[[ ifelse( d==2, 1, ifelse(d==5,2,3) ) ]] -ms <- multiRun( - list(n=n,p=p,beta=beta,b=b,optargs=list(K=K,link=link,nstart=nstart)), list( - function(fargs) { - # 1 start - library(morpheus) - K <- fargs$optargs$K - op <- optimParams(K, fargs$optargs$link, fargs$optargs) - x_init <- list(p=rep(1/K,K-1), beta=fargs$mu, b=rep(0,K)) - do.call(rbind,op$run(x_init)) - }, - function(fargs) { - # B starts - library(morpheus) - K <- fargs$optargs$K - op <- optimParams(K, fargs$optargs$link, fargs$optargs) - best_val <- Inf - best_par <- list() - for (i in 1:fargs$optargs$nstart) - { - x_init <- list(p=rep(1/K,K-1), beta=i*fargs$mu, b=rep(0,K)) - par <- op$run(x_init) - val <- op$f( op$linArgs(par) ) - if (val < best_val) - { - best_par <- par - best_val <- val - } - } - do.call(rbind,best_par) - }), - prepareArgs = function(fargs) { - library(morpheus) - io = generateSampleIO(fargs$n, fargs$p, fargs$beta, fargs$b, fargs$optargs$link) - fargs$optargs$M <- computeMoments(io$X, io$Y) - mu <- computeMu(io$X, io$Y, fargs$optargs) - fargs$mu <- mu - }, N=N, ncores=ncores, verbose=TRUE) -for (i in 1:2) - ms[[i]] <- alignMatrices(ms[[i]], ref=rbind(p,beta,b), ls_mode="exact") - -ms_params <- list("N"=N, "nc"=ncores, "n"=n, "K"=K, "link"=link, - "p"=p, "beta"=beta, "b"=b, "nstart"=nstart) +ms <- testMultistart(N, n, K, p, beta, b, link, nstart, ncores) +ms_params <- list("N"=N, "nc"=ncores, "n"=n, "K"=K, "d"=d, "link"=link, + "p"=c(p,1-sum(p)), "beta"=beta, "b"=b, "nstart"=nstart) -save(ms, ms_params, file="multistart.RData") +save("ms", "ms_params", file="multistart.RData") diff --git a/reports/test.R b/reports/test.R deleted file mode 100644 index 2ce9a44..0000000 --- a/reports/test.R +++ /dev/null @@ -1,58 +0,0 @@ -library(morpheus) -morph <- function(fargs) { - K <- fargs$optargs$K - M <- computeMoments(fargs$X, fargs$Y) - fargs$optargs$M <- M - mu <- computeMu(fargs$X, fargs$Y, fargs$optargs) - res2 <- NULL - tryCatch({ - op <- optimParams(K,link,fargs$optargs) - x_init <- list(p=rep(1/K,K-1), beta=mu, b=rep(0,K)) - res2 <- do.call(rbind, op$run(x_init)) - }, error = function(e) { - res2 <- NA - }) - res2 -} - -#model = binomial; default values: -link = "probit" -N <- 10 -d <- 2 -n <- 1e4 -ncores <- 1 - -if (d == 2) { - K <- 2 - p <- .5 - b <- c(-.2, .5) - beta <- matrix( c(1,-2, 3,1), ncol=K ) -} else if (d == 5) { - K <- 2 - p <- .5 - b <- c(-.2, .5) - beta <- matrix( c(1,2,-1,0,3, 2,-3,0,1,0), ncol=K ) -} else if (d == 10) { - K <- 3 - p <- c(.3, .3) - b <- c(-.2, 0, .5) - beta <- matrix( c(1,2,-1,0,3,4,-1,-3,0,2, 2,-3,0,1,0,-1,-4,3,2,0, -1,1,3,-1,0,0,2,0,1,-2), ncol=K ) -} else if (d == 20) { - K <- 3 - p <- c(.3, .3) - b <- c(-.2, 0, .5) - beta <- matrix( c(1,2,-1,0,3,4,-1,-3,0,2,2,-3,0,1,0,-1,-4,3,2,0, -1,1,3,-1,0,0,2,0,1,-2,1,2,-1,0,3,4,-1,-3,0,2, 2,-3,0,1,0,-1,-4,3,2,0,1,1,2,2,-2,-2,3,1,0,0), ncol=K ) -} - -fargs = list(n=n, p=p, beta=beta, b=b) -fargs$optargs = list(link=link) - -io = generateSampleIO(fargs$n, fargs$p, fargs$beta, fargs$b, fargs$optargs$link) -fargs$X = io$X -fargs$Y = io$Y -fargs$optargs$K = ncol(fargs$beta) -fargs$optargs$M = computeMoments(io$X,io$Y) - -res2 <- morph(fargs) - -save("res2", file="test.RData") diff --git a/reports/test.sh b/reports/test.sh deleted file mode 100644 index b617a09..0000000 --- a/reports/test.sh +++ /dev/null @@ -1,16 +0,0 @@ -#!/bin/bash - -# arg --vanilla maybe possible on cluster -for d in 2 5; do - for link in "logit" "probit"; do - R --slave --args N=10 n=1e3 nc=3 d=$d link=$link out$d$link 2>&1 - done -done - -#for d in 2 5; do -# for n in 5000 10000 100000 500000 1000000; do -# for link in "logit" "probit"; do -# R --slave --args N=1000 n=$n nc=64 d=$d link=$link out_$n$link$d 2>&1 -# done -# done -#done