Skip to content

Commit

Permalink
gamma site models have either a gamma cat count less than two and no …
Browse files Browse the repository at this point in the history 
…distribution, or a gamma cat count of two or more and a distribution. Thanks @bjoelle. Fix #50
  • Loading branch information
richelbilderbeek committed Oct 25, 2018
1 parent 1b8715c commit 911f560
Show file tree
Hide file tree
Showing 12 changed files with 198 additions and 71 deletions.
35 changes: 25 additions & 10 deletions R/create_gamma_site_model.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,14 @@
#' @param gamma_shape gamma curve shape parameter
#' @param prop_invariant the proportion invariant, must be a value
#' from 0.0 to 1.0
#' @param gamma_shape_prior_distr the distribution of the gamma shape prior
#' @param gamma_shape_prior_distr the distribution of the gamma shape prior.
#' \code{gamma_shape_prior_distr} must be \code{NA} for
#' a \code{gamma_cat_count} of zero or one.
#' For a \code{gamma_cat_count} of two or more,
#' leaving \code{gamma_shape_prior_distr} equal to its default
#' value of \code{NA}, a default distribution is used.
#' Else \code{gamma_shape_prior_distr} must be a
#' distribution, as can be created by \code{\link{create_distr}}
#' @param freq_equilibrium the frequency in which the rates are at equilibrium
#' are either \code{estimated}, \code{empirical} or \code{all_equal}.
#' \code{get_freq_equilibrium_names} returns the possible values
Expand All @@ -29,13 +36,7 @@ create_gamma_site_model <- function(
gamma_cat_count = "0",
gamma_shape = "1.0",
prop_invariant = "0.0",
gamma_shape_prior_distr = create_exp_distr(
id = NA,
mean = create_mean_param(
id = 0, # TODO: id must be NA by default
value = "1.0" # string to match XML
)
),
gamma_shape_prior_distr = NA,
freq_equilibrium = "estimated"
) {
if (gamma_cat_count < 0) {
Expand All @@ -47,10 +48,24 @@ create_gamma_site_model <- function(
if (prop_invariant < 0.0 || prop_invariant > 1.0) {
stop("'prop_invariant' must be in range [0.0, 1.0]")
}
if (!is_distr(gamma_shape_prior_distr)) {
if (gamma_cat_count >= 2 && is.na(gamma_shape_prior_distr)) {
gamma_shape_prior_distr <- create_exp_distr(
id = NA,
mean = create_mean_param(
id = 0, # TODO: id must be NA by default
value = "1.0" # string to match XML
)
)
}
if (!is.na(gamma_shape_prior_distr) && !is_distr(gamma_shape_prior_distr)) {
stop("'gamma_shape_prior_distr' must be a distribution")
}

if (gamma_cat_count < 2 && !is.na(gamma_shape_prior_distr)) {
stop(
"'gamma_shape_prior_distr' must be NA ",
"for a 'gamma_cat_count' of less than two"
)
}

gamma_site_model <- list(
gamma_cat_count = gamma_cat_count,
Expand Down
3 changes: 2 additions & 1 deletion R/create_random.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -193,7 +193,8 @@ create_rnd_gamma_site_model <- function() {
whitelist <- c(
"'gamma_cat_count' must be positive",
"'gamma_shape' must be positive",
"'prop_invariant' must be in range \\[0\\.0, 1\\.0\\]"
"'prop_invariant' must be in range \\[0\\.0, 1\\.0\\]",
"'gamma_shape_prior_distr' must be NA for a 'gamma_cat_count' of less than two" # nolint indeed long error message, preferred this over using paste0
)
testit::assert(
is_in_patterns(line = error$message, patterns = whitelist)
Expand Down
107 changes: 59 additions & 48 deletions R/init_site_models.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -73,17 +73,19 @@ init_gtr_site_model <- function(
) {
testit::assert(is_gtr_site_model(gtr_site_model))

# Gamma site model
if (!is_init_distr(gtr_site_model$gamma_site_model$gamma_shape_prior_distr)) { # nolint internal function
gtr_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
gtr_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
gtr_site_model$gamma_site_model$gamma_shape_prior_distr
)
# Initialize gamma site model, if any
if (!is_one_na(gtr_site_model$gamma_site_model$gamma_shape_prior_distr)) {
if (!is_init_distr(gtr_site_model$gamma_site_model$gamma_shape_prior_distr)) { # nolint internal function
gtr_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
gtr_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
gtr_site_model$gamma_site_model$gamma_shape_prior_distr
)
}
}

if (!is_init_distr(gtr_site_model$rate_ac_prior_distr)) { # nolint internal function
Expand Down Expand Up @@ -204,17 +206,19 @@ init_hky_site_model <- function(
) {
testit::assert(is_hky_site_model(hky_site_model))

# Initialize gamma site model
if (!is_init_distr(hky_site_model$gamma_site_model$gamma_shape_prior_distr)) {
hky_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
hky_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
hky_site_model$gamma_site_model$gamma_shape_prior_distr
)
# Initialize gamma site model, if any
if (!is_one_na(hky_site_model$gamma_site_model$gamma_shape_prior_distr)) {
if (!is_init_distr(hky_site_model$gamma_site_model$gamma_shape_prior_distr)) {
hky_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
hky_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
hky_site_model$gamma_site_model$gamma_shape_prior_distr
)
}
}

# kappa_prior_distr
Expand Down Expand Up @@ -251,20 +255,24 @@ init_jc69_site_model <- function(
param_id = 0
) {
testit::assert(is_jc69_site_model(jc69_site_model))
# Initialize gamma site model
if (
!is_init_distr(jc69_site_model$gamma_site_model$gamma_shape_prior_distr)
) {
jc69_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
jc69_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
jc69_site_model$gamma_site_model$gamma_shape_prior_distr
)

# Initialize gamma site model, if any
if (!is_one_na(jc69_site_model$gamma_site_model$gamma_shape_prior_distr)) {
if (
!is_init_distr(jc69_site_model$gamma_site_model$gamma_shape_prior_distr)
) {
jc69_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
jc69_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
jc69_site_model$gamma_site_model$gamma_shape_prior_distr
)
}
}

testit::assert(is_init_gamma_site_model(jc69_site_model$gamma_site_model)) # nolint internal function
testit::assert(is_init_jc69_site_model(jc69_site_model)) # nolint internal function
jc69_site_model
Expand All @@ -287,19 +295,22 @@ init_tn93_site_model <- function(
param_id = 0
) {
testit::assert(is_tn93_site_model(tn93_site_model))
# Initialize gamma site model
if (
!is_init_distr(tn93_site_model$gamma_site_model$gamma_shape_prior_distr)
) {
tn93_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
tn93_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
tn93_site_model$gamma_site_model$gamma_shape_prior_distr
)

# Initialize gamma site model, if any
if (!is_one_na(tn93_site_model$gamma_site_model$gamma_shape_prior_distr)) {
if (
!is_init_distr(tn93_site_model$gamma_site_model$gamma_shape_prior_distr)
) {
tn93_site_model$gamma_site_model$gamma_shape_prior_distr <- init_distr(
tn93_site_model$gamma_site_model$gamma_shape_prior_distr,
distr_id = distr_id,
param_id = param_id
)
distr_id <- distr_id + 1
param_id <- param_id + get_distr_n_params(
tn93_site_model$gamma_site_model$gamma_shape_prior_distr
)
}
}

# kappa_1_prior_distr
Expand Down
3 changes: 2 additions & 1 deletion R/is_gamma_site_model.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,7 @@ is_gamma_site_model <- function(x) {
if (x$prop_invariant < 0.0) return(FALSE)
if (x$prop_invariant > 1.0) return(FALSE)
if (!"gamma_shape_prior_distr" %in% names(x)) return(FALSE)
if (!is_distr(x$gamma_shape_prior_distr)) return(FALSE)
if (!is.na(x$gamma_shape_prior_distr) &&
!is_distr(x$gamma_shape_prior_distr)) return(FALSE)
TRUE
}
1 change: 1 addition & 0 deletions R/is_init_gamma_site_model.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,7 @@ is_init_gamma_site_model <- function(
x
) {
if (!is_gamma_site_model(x)) return(FALSE)
if (is_one_na(x$gamma_shape_prior_distr)) return(TRUE) # nolint internal function
if (!is_init_distr(x$gamma_shape_prior_distr)) return(FALSE) # nolint internal function
TRUE
}
12 changes: 9 additions & 3 deletions man/create_gamma_site_model.Rd
View file

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

4 changes: 2 additions & 2 deletions tests/testthat/test-create_beast2_input_by_reproducing_files_clock.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -200,12 +200,12 @@ test_that("aco_nd2_strict_rln_2_4.xml, strict RLN, example 10", {
site_models = list(
create_jc69_site_model(
gamma_site_model = create_gamma_site_model(
gamma_shape_prior_distr = create_exp_distr(id = 0)
gamma_cat_count = 0
)
),
create_jc69_site_model(
gamma_site_model = create_gamma_site_model(
gamma_shape_prior_distr = create_exp_distr(id = 0)
gamma_cat_count = 0
)
)
),
Expand Down
45 changes: 45 additions & 0 deletions tests/testthat/test-create_gamma_site_model.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,6 +17,30 @@ test_that("Can specify HKY proportion invariant", {

})

test_that("zero gamma count categories has no distribution", {

gamma_site_model <- create_gamma_site_model(gamma_cat_count = 0)
expect_false(
is_distr(gamma_site_model$gamma_shape_prior_distr)
)
})

test_that("one gamma count categories has no distribution", {

gamma_site_model <- create_gamma_site_model(gamma_cat_count = 1)
expect_false(
is_distr(gamma_site_model$gamma_shape_prior_distr)
)
})

test_that("two gamma count categories has a distribution", {

gamma_site_model <- create_gamma_site_model(gamma_cat_count = 2)
expect_true(
is_distr(gamma_site_model$gamma_shape_prior_distr)
)
})


test_that("abuse", {

Expand Down Expand Up @@ -51,4 +75,25 @@ test_that("abuse", {
"'gamma_shape_prior_distr' must be a distribution"
)

testthat::expect_error(
create_gamma_site_model(
gamma_cat_count = 0,
gamma_shape_prior_distr = create_exp_distr()
),
paste0(
"'gamma_shape_prior_distr' must be NA ",
"for a 'gamma_cat_count' of less than two"
)
)
testthat::expect_error(
create_gamma_site_model(
gamma_cat_count = 1,
gamma_shape_prior_distr = create_exp_distr()
),
paste0(
"'gamma_shape_prior_distr' must be NA ",
"for a 'gamma_cat_count' of less than two"
)
)

})
2 changes: 1 addition & 1 deletion tests/testthat/test-init_gamma_site_model.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3,11 +3,11 @@ context("init_gamma_site_model")
test_that("use", {

gamma_site_model <- create_gamma_site_model(
gamma_cat_count = 2,
gamma_shape_prior_distr = create_one_div_x_distr(id = NA)
)
testit::assert(!beautier:::is_init_gamma_site_model(gamma_site_model))
gamma_site_model <- beautier:::init_gamma_site_model(gamma_site_model)
testthat::expect_true(beautier:::is_init_gamma_site_model(gamma_site_model))


})
28 changes: 24 additions & 4 deletions tests/testthat/test-init_site_models.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,13 @@ context("init_site_models")

test_that("initialize JC69 site model", {

gamma_site_model <- create_gamma_site_model(
gamma_cat_count = 2,
gamma_shape_prior_distr = create_exp_distr()
)

id <- "a"
before <- list(create_jc69_site_model())
before <- list(create_jc69_site_model(gamma_site_model = gamma_site_model))
testit::assert(beautier:::is_jc69_site_model(before[[1]]))
testit::assert(!beautier:::are_init_site_models(before))
after <- beautier:::init_site_models(before, ids = id)
Expand All @@ -14,8 +19,13 @@ test_that("initialize JC69 site model", {

test_that("initialize HKY site model", {

gamma_site_model <- create_gamma_site_model(
gamma_cat_count = 2,
gamma_shape_prior_distr = create_exp_distr()
)

id <- "a"
before <- list(create_hky_site_model())
before <- list(create_hky_site_model(gamma_site_model = gamma_site_model))
testit::assert(beautier:::is_hky_site_model(before[[1]]))
testit::assert(!beautier:::are_init_site_models(before))
after <- beautier:::init_site_models(before, ids = id)
Expand All @@ -26,8 +36,13 @@ test_that("initialize HKY site model", {

test_that("initialize TN93 site model", {

gamma_site_model <- create_gamma_site_model(
gamma_cat_count = 2,
gamma_shape_prior_distr = create_exp_distr()
)

id <- "a"
before <- list(create_tn93_site_model())
before <- list(create_tn93_site_model(gamma_site_model = gamma_site_model))
testit::assert(beautier:::is_tn93_site_model(before[[1]]))
testit::assert(!beautier:::are_init_site_models(before))
after <- beautier:::init_site_models(before, ids = id)
Expand All @@ -38,8 +53,13 @@ test_that("initialize TN93 site model", {

test_that("initialize GTR site model", {

gamma_site_model <- create_gamma_site_model(
gamma_cat_count = 2,
gamma_shape_prior_distr = create_exp_distr()
)

id <- "a"
before <- list(create_gtr_site_model())
before <- list(create_gtr_site_model(gamma_site_model = gamma_site_model))
testit::assert(beautier:::is_gtr_site_model(before[[1]]))
testit::assert(!beautier:::are_init_site_models(before))
after <- beautier:::init_site_models(before, ids = id)
Expand Down
Loading

1 comment on commit 911f560

@lintr-bot
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

R/create_site_model.R:76:3: warning: no visible global function definition for ‘check_gamma_site_model’

​  check_gamma_site_model(gamma_site_model)
  ^~~~~~~~~~~~~~~~~~~~~~

R/get_site_model_n_distrs.R:33:32: warning: no visible global function definition for ‘get_gamma_site_model_n_distrs’

​  gamma_site_model_n_distrs <- get_gamma_site_model_n_distrs(gamma_site_model)
                               ^~~~~~~~~~~~~~~~~~~~~~~~~~~~~

R/init_site_models.R:211:1: style: Lines should not be more than 80 characters.

​    if (!is_init_distr(hky_site_model$gamma_site_model$gamma_shape_prior_distr)) {
^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

tests/testthat/test-get_site_models_n_distrs.R:8:1: style: Lines should not be more than 80 characters.

​  site_model_0_distrs <- create_jc69_site_model(gamma_site_model = gamma_site_model)
^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

tests/testthat/test-get_site_models_n_distrs.R:9:1: style: Lines should not be more than 80 characters.

​  site_model_1_distrs <- create_hky_site_model(gamma_site_model = gamma_site_model)
^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

tests/testthat/test-get_site_models_n_distrs.R:10:1: style: Lines should not be more than 80 characters.

​  site_model_2_distrs <- create_tn93_site_model(gamma_site_model = gamma_site_model)
^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Please sign in to comment.