1
#' Estimation of Log-likelihood for Linear Mixed Model with Lasso Penalty
2
#'
3
#' @description \code{sigma2lasso} estimates the value of the sigma2 for the
4
#'   linear mixed model with lasso penalty
5
#'
6
#' @seealso \code{\link{ggmix}}
7
#' @param ggmix_object A ggmix_object object of class \code{lowrank} or
8
#'   \code{fullrank}
9
#' @inheritParams ggmix
10
#' @param beta current estimate of the beta parameter including the intercept.
11
#'   this should be of length p+1, where p is the number of variables.
12
#' @param eta current estimate of the eta parameter
13
#' @param sigma2 current estimate of the sigma2 parameter
14
#' @param nt total number of observations
15
#' @param ... Extra parameters. Currently ignored.
16
#' @return A decreasing sequence of tuning parameters
17
#' @note This function isn't meant to be called directly by the user.
18
#' @seealso \code{\link{ggmix_data_object}}
19 4
logliklasso <- function(ggmix_object, ...) UseMethod("logliklasso")
20

21
#' @rdname logliklasso
22
logliklasso.default <- function(ggmix_object, ...) {
23 0
  stop(strwrap("This function should be used with a ggmix object of class
24 0
               lowrank or fullrank"))
25
}
26

27
#' @rdname logliklasso
28
logliklasso.fullrank <- function(ggmix_object,
29
                                 ...,
30
                                 eta, sigma2, beta, nt, x, y) {
31

32
  # this returns the log-likelihood
33

34
  # this allows the user to input the x and y if they want else
35
  # the default is from the ggmix_object
36
  # this flexibility is used to calculate the saturated and intercept loglik
37
  # in the lmmlasso function
38 4
  if (missing(x)) {
39 4
    x <- ggmix_object[["x"]]
40
  }
41

42 4
  if (missing(y)) {
43 4
    y <- ggmix_object[["y"]]
44
  }
45

46 4
  eigenvalues <- ggmix_object[["D"]]
47

48 4
  kernel <- 1 + eta * (eigenvalues - 1)
49

50 4
  -1 * (
51 4
    (nt / 2) * log(2 * pi) +
52 4
      (nt / 2) * log(sigma2) +
53 4
      0.5 * sum(log(kernel)) +
54 4
      (1 / (2 * sigma2)) * sum((y - x %*% beta)^2 / kernel)
55
  )
56
}
57

58

59
# still need to create logliklasso.lowrank, loglikgglasso.fullrank, loglikgglasso.lowrank

Read our documentation on viewing source code .

Loading