Merge remote-tracking branch 'origin/tau_nesmr2' into nesmr-logdet

# Conflicts:
#	DESCRIPTION
#	NAMESPACE
#	R/checks.R
#	R/esmr.R
#	R/helpers.R
#	R/likelihood_funcs.R
#	man/esmr.Rd

DESCRIPTION

@@ -1,14 +1,14 @@
 Package: esmr
 Type: Package
 Title: Empirical Shrinkage Mendelian Randomization
-Version: 0.2.1.0176
+Version: 0.2.1.0194
 Author: Jean Morrison
 Maintainer: Jean Morrison <[email protected]>
 Description: Perform univariable or multivariable MR with empirical shrinkage priors.
 License: GPL (>= 3)
 Encoding: UTF-8
 LazyData: true
-RoxygenNote: 7.2.3
+RoxygenNote: 7.3.1
 Suggests:
     knitr,
     rmarkdown
@@ -21,8 +21,6 @@ Imports:
     matrixStats,
     purrr,
     GFA,
-    igraph,
-    numDeriv,
-    rstackdeque
+    igraph
 Remotes:
     jean997/GFA

NAMESPACE

@@ -1,22 +1,13 @@
 # Generated by roxygen2: do not edit by hand
 
+S3method(logLik,esmr)
 export(calc_ell2)
 export(esmr)
 export(estimate_G)
-export(h_det)
-export(h_det_grad)
-export(h_det_grad_vec)
-export(h_det_hessian)
-export(h_det_hessian_vec)
-export(h_det_vec)
-export(nesmr_all_permn)
-export(nesmr_backselect)
-export(nesmr_complete)
-export(nesmr_complete_mvmr)
-export(project_to_DAG)
+export(optimize_lpy2)
 export(snp_beta_rb)
 export(update_beta_joint)
-export(update_l_k)
+export(update_lj)
 import(ashr)
 import(dplyr)
 import(ebnm)

R/best_combinations.R

@@ -0,0 +1,71 @@
+
+get_top_combinations <- function(x, max_logsumexp, nmax = Inf){
+  c <- apply(x, 1, which.max)
+
+  combination_list <- list(c)
+  v <- sapply(1:nrow(x), function(i){x[i,c[i]]}) %>% sum()
+  values <- c(v)
+  l <- Inf
+  i <- 1
+  done <- FALSE
+  #check_ix <- c(1)
+  while(!done){
+    #i <- check_ix[1]
+    #check_ix <- check_ix[-1]
+    #cat(i, " ")
+    v <- values[i]
+    c <- combination_list[[i]]
+    new_comb <- best_valid_one_move(c, x, combination_list)
+    losses <- new_comb[[2]]
+    #cat(losses, "\n")
+    if(all(is.na(losses))){
+      if(i == length(values)) done <- TRUE
+    }else{
+      new_combs <- new_comb[[1]][!is.na(losses)]
+      losses <- losses[!is.na(losses)]
+      combination_list <- c(combination_list, new_combs)
+      new_vals <- v + losses
+      values <- c(values, new_vals)
+      o <- order(values, decreasing = T)
+      values <- values[o]
+      combination_list <- combination_list[o]
+    }
+    i <- i + 1
+    w <- matrixStats::logSumExp(values[1:i])
+    if(i > nmax | w > max_logsumexp){
+      done <- TRUE
+    }
+  }
+  top_combinations <- combination_list[1:i] %>% unlist() %>%
+                         matrix(ncol = nrow(x), byrow = T)
+  return(list(combs = top_combinations, values = values[1:i]))
+}
+
+best_valid_one_move <- function(c, x, combination_list){
+  n <- length(c)
+  best_one_move <- lapply(1:n, function(i){
+    cix <- lapply(combination_list, function(cc){
+      if(all(cc[-i] == c[-i])){
+        return(cc[i])
+      }else{
+        return(NULL)
+      }
+    }) %>% unlist()
+    cix <- c(c[i], cix)
+
+    xx <- x[i,]
+    xx[cix] <- -Inf
+    l <- xx-x[i, c[i]]
+
+    if(!any(is.finite(l))){
+      return(c(NA, NA))
+    }
+    newc <- c
+    newc[i] <- which.max(l)
+    return(list(newc, max(l)))
+  })
+  combinations <- map(best_one_move, 1)
+  losses <- map(best_one_move, 2) %>% unlist()
+  return(list(combinations, losses))
+}
+