R/likelihood_expanded.R

In BriceonWiley/IntegratedLikelihood.R:

likelihood_expanded <- function(p1, gam, data1, data2) {

  p2 <- p1 / ((1 - p1) * exp(gam) + p1)
  likelihood_1 <- Vectorize(likelihood_p, "p")(p1, data1)
  likelihood_2 <- Vectorize(likelihood_p, "p")(p2, data2)
  matrix(likelihood_1 * likelihood_2, ncol = length(p1))

  # FOR INTEGRATING OVER CLOSED FORMS
  # p2 <- p1 / ((1 - p1) * exp(gam) + p1)
  #
  # likelihood_1 <- Vectorize(likelihood_p_closed, "p")(p1, data1)
  # likelihood_2 <- Vectorize(likelihood_p_closed, "p")(p2, data2)
  #
  # matrix(likelihood_1 * likelihood_2, ncol = length(p1))

  # USE BELOW IF CAN FIND A WAY TO INTEGRATE FASTER
  # p1 <- vec[1, ]
  # phi1 <- vec[2, ]
  # phi2 <- vec[3, ]
  # theta1 <- vec[4, ]
  # theta2 <- vec[5, ]
  #
  # n100 <- data1[1]; n101 <- data1[2]; n110 <- data1[3]
  # n111 <- data1[4]; X1 <- data1[5]; Y1 <- data1[6]
  # n200 <- data2[1]; n201 <- data2[2]; n210 <- data2[3]
  # n211 <- data2[4]; X2 <- data2[5]; Y2 <- data2[6]
  #
  # p2 <- p1 / ((1 - p1) * exp(gam) + p1)
  # pi1 <- p1 * (1 - theta1) + (1 - p1) * phi1
  # pi2 <- p2 * (1 - theta2) + (1 - p2) * phi2
  #
  # # Likelihood 1
  # log_likelihood_11 <- n100 * (log(1 - p1) + log(1 - phi1))
  # # if(any(is.nan(log_likelihood_11))) log_likelihood_11 <- 0
  #
  # log_likelihood_12 <- n101 * (log(1 - p1) + log(phi1))
  # # if(any(is.nan(log_likelihood_12))) log_likelihood_12 <- 0
  #
  # log_likelihood_13 <- n110 * (log(p1) + log(theta1))
  # # if(any(is.nan(log_likelihood_13))) log_likelihood_13 <- 0
  #
  # log_likelihood_14 <- n111 * (log(p1) + log(1 - theta1))
  # # if(any(is.nan(log_likelihood_14))) log_likelihood_14 <- 0
  #
  # log_likelihood_15 <- X1 * log(pi1)
  # # if(any(is.nan(log_likelihood_15))) log_likelihood_15 <- 0
  #
  # log_likelihood_16 <- Y1 * log(1 - pi1)
  # # if(any(is.nan(log_likelihood_16))) log_likelihood_16 <- 0
  #
  # # Likelihood 2
  # log_likelihood_21 <- n200 * (log(1 - p2) + log(1 - phi2))
  # # if(any(is.nan(log_likelihood_21))) log_likelihood_21 <- 0
  #
  # log_likelihood_22 <- n201 * (log(1 - p2) + log(phi2))
  # # if(any(is.nan(log_likelihood_22))) log_likelihood_22 <- 0
  #
  # log_likelihood_23 <- n210 * (log(p2) + log(theta2))
  # # if(any(is.nan(log_likelihood_23))) log_likelihood_23 <- 0
  #
  # log_likelihood_24 <- n211 * (log(p2) + log(1 - theta2))
  # # if(any(is.nan(log_likelihood_24))) log_likelihood_24 <- 0
  #
  # log_likelihood_25 <- X2 * log(pi2)
  # # if(any(is.nan(log_likelihood_25))) log_likelihood_25 <- 0
  #
  # log_likelihood_26 <- Y2 * log(1 - pi2)
  # # if(any(is.nan(log_likelihood_26))) log_likelihood_26 <- 0
  #
  # log_likelihood <-
  #   log_likelihood_11 + log_likelihood_12 +
  #   log_likelihood_13 + log_likelihood_14 +
  #   log_likelihood_15 + log_likelihood_16 +
  #   log_likelihood_21 + log_likelihood_22 +
  #   log_likelihood_23 + log_likelihood_24 +
  #   log_likelihood_25 + log_likelihood_26
  # matrix(exp(log_likelihood), ncol = length(log_likelihood))
}