R/target-u.R

In BrandonEdwards/spmm: Generalized spatial Poisson mixture models

#' Target u
#'
#' @importFrom pracma sqrtm
#'

target_u <- function(D = NULL,
                     W = NULL,
                     v = NULL,
                     C = NULL,
                     R = NULL,
                     kappa = NULL,
                     k = NULL,
                     B = NULL,
                     u = NULL,
                     constant_B = FALSE)
{
  zeta <- kappa
  xi <- eigen(pracma::sqrtm(D)$Binv %*% W %*% pracma::sqrtm(D)$Binv)$values
  #
  # eigen_check <- (1/min(xi) < zeta) && (zeta < 1)
  # if (!all(eigen_check))
  # {
  #   return(-Inf)
  # }

  prob_u_given_B_1 <- 0

  if (isFALSE(constant_B))
  {
    for (i in 1:C)
    {
      for (j in 1:R)
      {
        prob_u_given_B_1 <- prob_u_given_B_1 + log(1 - (zeta[i] * xi[j]))
      }
    }
  }


  prob_u_given_B_2 <- 0

  for (i in 1:C)
  {
    for (j in 1:R)
    {
      for (l in 1:R)
      {
        prob_u_given_B_2 <- prob_u_given_B_2 +
          B[i,i] * u[i,j] * W[j,l] * u[i,l]
      }
    }
  }

  prob_u_given_B_3 <- 0
  for (i in 1:(C - 1))
  {
    for (j in (i+1):C)
    {
      for (k in 1:R)
      {
        for (l in 1:R)
        {
          prob_u_given_B_3 <- prob_u_given_B_3 +
            B[i,i] * u[i,j] * W[j,l] * u[i,l]
        }
      }
    }
  }

  prob_u_given_B <- 0.5*prob_u_given_B_1 + 0.5*prob_u_given_B_2 + prob_u_given_B_3

  return(prob_u_given_B)
}