ourMELONS/R/computePopulationLogml.R

34 lines
1 KiB
R
Raw Normal View History

2021-02-01 09:22:58 +01:00
computePopulationLogml <- function(pops, adjprior, priorTerm = 0) {
2020-10-19 14:08:25 +02:00
# Palauttaa length(pops)*1 taulukon, jossa on laskettu korikohtaiset
2021-02-01 09:22:58 +01:00
# ======================================================== #
# Limiting COUNTS size #
# ======================================================== #
COUNTS <- COUNTS[seq_len(nrow(adjprior)), seq_len(ncol(adjprior)), pops, drop=FALSE]
2020-10-19 14:08:25 +02:00
x <- size(COUNTS, 1)
y <- size(COUNTS, 2)
z <- length(pops)
2021-02-01 09:22:58 +01:00
# ======================================================== #
# Computation #
# ======================================================== #
term1 <- squeeze(
2021-01-15 11:47:20 +01:00
# FIXME: assumes COUNTS has 3 dims. Where does this come from?
2020-10-19 14:08:25 +02:00
sum(
sum(
reshape(
lgamma(
repmat(adjprior, c(1, 1, length(pops))) +
2021-02-01 09:22:58 +01:00
COUNTS[seq_len(nrow(adjprior)), seq_len(ncol(adjprior)), pops, drop=FALSE]
2020-10-19 14:08:25 +02:00
),
c(x, y, z)
),
1
),
2
)
2021-02-01 09:22:58 +01:00
)
popLogml <- term1 - sum(lgamma(1 + SUMCOUNTS[pops, ]), 2) - priorTerm
2020-10-19 14:08:25 +02:00
return(popLogml)
}