Skip to content

Commit

Permalink
Merge branch 'develop' into main
Browse files Browse the repository at this point in the history
  • Loading branch information
Nicholas Brazeau committed Dec 4, 2023
2 parents f5f42a9 + 5240b33 commit 7f96210
Show file tree
Hide file tree
Showing 6 changed files with 132 additions and 62 deletions.
2 changes: 1 addition & 1 deletion R/main.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -162,7 +162,7 @@ deme_inbreeding_spcoef <- function(discdat,
# catch accidental bad M start if user is standardizing distances
if (normalize_geodist & (start_params[names(start_params) == "m"] > 500) ) {
warning("You have selected to normalize geographic distances, but your
migration rate parameter is large. Please consider placing it on a
migration rate start parameter is large. Please consider placing it on a
similar scale to your normalized geographic distances for stability.")
}

Expand Down
3 changes: 3 additions & 0 deletions data-raw/IBD_simulation_data.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -131,6 +131,9 @@ IBD_simulation_data <- sim_IBDIBD(demesize = c(3,3,4), distmat = matrix(c(0,500,
100,750,0), nrow = 3),
rate = 1e-3, Ft = 0.3)

# drop within deme
IBD_simulation_data <- IBD_simulation_data[IBD_simulation_data$geodist > 0, ]

#............................................................
# out
#...........................................................
Expand Down
Binary file modified data/IBD_simulation_data.rda
View file
Binary file not shown.
18 changes: 9 additions & 9 deletions src/main.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -184,16 +184,16 @@ Rcpp::List deme_inbreeding_coef_cpp(Rcpp::List args, Rcpp::List args_functions,
fi_gradtraj[step][i] = fgrad[i];
}

// // get M moments for Adam
// m1t_m[step] = b1 * m1t_m[step-1] + (1-b1) * mgrad;
// v2t_m[step] = b2 * v2t_m[step-1] + (1-b2) * pow(mgrad, 2);
// m1t_m_hat = m1t_m[step] / (1-pow(b1, step));
// v2t_m_hat = v2t_m[step] / (1-pow(b2, step));
//
// // calculate and apply the update for M
// m = m - m_learningrate * (m1t_m_hat / (sqrt(v2t_m_hat) + e));
// get M moments for Adam
m1t_m[step] = b1 * m1t_m[step-1] + (1-b1) * mgrad;
v2t_m[step] = b2 * v2t_m[step-1] + (1-b2) * pow(mgrad, 2);
m1t_m_hat = m1t_m[step] / (1-pow(b1, step));
v2t_m_hat = v2t_m[step] / (1-pow(b2, step));

// calculate and apply the update for M
m = m - m_learningrate * (m1t_m_hat / (sqrt(v2t_m_hat) + e));
// vanilla GD
m = m - m_learningrate * mgrad;
// m = m - m_learningrate * mgrad;
// check bounds on m
// will reflect with normal based on magnitude + standard normal sd it is off to proper interval; NB also want to bound m so that it can only explore distance isolation (repulsion versus attraction)
if (m < m_lowerbound) {
Expand Down
104 changes: 52 additions & 52 deletions tests/testthat/test-AdamByHand.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -54,55 +54,55 @@ test_that("Fi adam by hand", {
})


# test_that("M adam by hand", {
# # sim data
# data("IBD_simulation_data", package = "discent")
# dat <- IBD_simulation_data
#
# # start params
# our_start_params <- rep(0.2, 3)
# names(our_start_params) <- 1:3
# our_start_params <- c(our_start_params, "m" = 1e3)
#
# # tidy date and calculate gradient for F1
# input <- dat %>%
# dplyr::filter(deme1 != deme2) %>%
# dplyr::filter(deme1 == 1 | deme2 == 1)
#
#
# input <- input %>%
# dplyr::mutate(gendist = ifelse(gendist > 0.999, 0.999,
# ifelse(gendist < 0.001, 0.001,
# gendist))) %>% # reasonable bounds on logit
# dplyr::mutate(gendist = logit(gendist))
#
# # now run model
# inputdisc <- dat %>%
# dplyr::filter(deme1 != deme2)
# ret <- deme_inbreeding_spcoef(discdat = inputdisc,
# start_params = our_start_params,
# f_learningrate = 1e-3,
# m_learningrate = 1e-5,
# b1 = 0.9,
# b2 = 0.999,
# e = 1e-8,
# steps = 1e3,
# normalize_geodist = F,
# report_progress = T,
# return_verbose = T)
# # back out for m adam
# b1 <- 0.9
# b2 <- 0.999
# e <- 1e-8
# m_learningrate <- 1e-5
# mt_m <- b1 * 0 + (1-b1) * ret$fi_gradtraj[2,1]
# vt_m <- b2 * 0 + (1-b2) * (ret$fi_gradtraj[2,1]^2)
# mt_mhat <- mt_m / (1 - b1^1)
# vt_mhat <- vt_m / (1 - b2^1)
# mnew1 <- ret$m_run[1] - m_learningrate * (mt_mhat/(sqrt(vt_mhat) + e))
# ret$m_run[1]
#
# # test out
# testthat::expect_equal(mnew1, ret$m_run[1])
#
# })
test_that("M adam by hand", {
# sim data
data("IBD_simulation_data", package = "discent")
dat <- IBD_simulation_data

# start params
our_start_params <- rep(0.2, 3)
names(our_start_params) <- 1:3
our_start_params <- c(our_start_params, "m" = 1e3)

# tidy date and calculate gradient for F1
input <- dat %>%
dplyr::filter(deme1 != deme2) %>%
dplyr::filter(deme1 == 1 | deme2 == 1)


input <- input %>%
dplyr::mutate(gendist = ifelse(gendist > 0.999, 0.999,
ifelse(gendist < 0.001, 0.001,
gendist))) %>% # reasonable bounds on logit
dplyr::mutate(gendist = logit(gendist))

# now run model
inputdisc <- dat %>%
dplyr::filter(deme1 != deme2)
ret <- deme_inbreeding_spcoef(discdat = inputdisc,
start_params = our_start_params,
f_learningrate = 1e-3,
m_learningrate = 1e-5,
b1 = 0.9,
b2 = 0.999,
e = 1e-8,
steps = 1e3,
normalize_geodist = F,
report_progress = T,
return_verbose = T)
# back out for m adam
b1 <- 0.9
b2 <- 0.999
e <- 1e-8
m_learningrate <- 1e-5
mt_m <- b1 * 0 + (1-b1) * ret$fi_gradtraj[2,1]
vt_m <- b2 * 0 + (1-b2) * (ret$fi_gradtraj[2,1]^2)
mt_mhat <- mt_m / (1 - b1^1)
vt_mhat <- vt_m / (1 - b2^1)
mnew1 <- ret$m_run[1] - m_learningrate * (mt_mhat/(sqrt(vt_mhat) + e))
ret$m_run[1]

# test out
testthat::expect_equal(mnew1, ret$m_run[1])

})
67 changes: 67 additions & 0 deletions tests/testthat/test-consistentgradient.R
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,67 @@
test_that("model has deterministic results from same start", {

# sim data
data("IBD_simulation_data", package = "discent")
inputdisc <- IBD_simulation_data
# start params
our_start_params <- rep(0.2, 3)
names(our_start_params) <- 1:3
our_start_params <- c(our_start_params, "m" = 1e3)

# given same start parameters and data and iters
# model should always follow the same gradient/trajectory
mod1 <- deme_inbreeding_spcoef(discdat = inputdisc,
start_params = our_start_params,
f_learningrate = 1e-5,
m_learningrate = 1e-1,
b1 = 0.9,
b2 = 0.999,
e = 1e-8,
steps = 1e4,
normalize_geodist = F,
report_progress = F,
return_verbose = F)

mod2 <- deme_inbreeding_spcoef(discdat = inputdisc,
start_params = our_start_params,
f_learningrate = 1e-5,
m_learningrate = 1e-1,
b1 = 0.9,
b2 = 0.999,
e = 1e-8,
steps = 1e4,
normalize_geodist = F,
report_progress = F,
return_verbose = F)

mod3 <- deme_inbreeding_spcoef(discdat = inputdisc,
start_params = our_start_params,
f_learningrate = 1e-5,
m_learningrate = 1e-1,
b1 = 0.9,
b2 = 0.999,
e = 1e-8,
steps = 1e4,
normalize_geodist = F,
report_progress = F,
return_verbose = F)

mod4 <- deme_inbreeding_spcoef(discdat = inputdisc,
start_params = our_start_params,
f_learningrate = 1e-5,
m_learningrate = 1e-1,
b1 = 0.9,
b2 = 0.999,
e = 1e-8,
steps = 5e4, # different number of steps!!
normalize_geodist = F,
report_progress = F,
return_verbose = F)

testthat::expect_equal(mod1,mod2)
testthat::expect_equal(mod1,mod3)
testthat::expect_equal(mod2,mod2) # trans, so know true, but for complete
testthat::expect_false(identical(mod2,mod4)) # more steps, different results
})


0 comments on commit 7f96210

Please sign in to comment.