From 88cce42e8849a033615bdd3ab291a2bc10fbe8e0 Mon Sep 17 00:00:00 2001 From: Sean Hackett Date: Tue, 28 May 2024 16:46:43 -0700 Subject: [PATCH 1/3] exposed the number of iterations to run optimizer since some current timecourses were not converging --- R/data.R | 11 ++++ R/impulse_fitting.R | 14 +++- data/example_timecourse.rda | Bin 0 -> 733 bytes man/estimate_one_timecourse_params_tf.Rd | 3 + man/estimate_timecourse_params_tf.Rd | 3 + man/example_timecourse.Rd | 21 ++++++ ...st-likelihood.R => test-impulse_fitting.R} | 60 +++++++++++++++++- 7 files changed, 108 insertions(+), 4 deletions(-) create mode 100644 R/data.R create mode 100644 data/example_timecourse.rda create mode 100644 man/example_timecourse.Rd rename tests/testthat/{test-likelihood.R => test-impulse_fitting.R} (72%) diff --git a/R/data.R b/R/data.R new file mode 100644 index 0000000..47030b2 --- /dev/null +++ b/R/data.R @@ -0,0 +1,11 @@ +#' Example Timecourse +#' +#' An example timecourse of fluorescent reporter dynamics in yeast. +#' +#' @format A tibble with 24 rows and 3 columns: +#' \describe{ +#' \item{tc_id}{A number representing a distinct timecourse} +#' \item{time}{Time (in divisions)} +#' \item{abundance}{Relative expression of the reporter} +#' } +"example_timecourse" diff --git a/R/impulse_fitting.R b/R/impulse_fitting.R index 682eef4..afdcecb 100644 --- a/R/impulse_fitting.R +++ b/R/impulse_fitting.R @@ -28,6 +28,7 @@ #' } #' @param fit_intercept If TRUE, the intercept will be fit, if FALSE, the intercept will be constrainted to zero #' @param learning_rate learning rate for the Adams optimizer +#' @param n_iterations the number of iterations to run the optimizer #' @param verbose if TRUE then print additional information #' #' @return a timecourse list: @@ -77,6 +78,7 @@ estimate_timecourse_params_tf <- function( "time_scale" = 15), fit_intercept = FALSE, learning_rate = 0.1, + n_iterations = 100, verbose = FALSE ) { @@ -105,6 +107,7 @@ estimate_timecourse_params_tf <- function( checkmate::assertLogical(use_prior, len = 1) checkmate::assertLogical(fit_intercept, len = 1) checkmate::assertNumber(learning_rate, lower = 0) + checkmate::assertNumber(n_iterations, lower = 1) checkmate::assertLogical(verbose, len = 1) # test parameters supplied for parameter initialization / priors @@ -153,6 +156,7 @@ estimate_timecourse_params_tf <- function( prior_pars = prior_pars, fit_intercept = fit_intercept, learning_rate = learning_rate, + n_iterations = n_iterations, verbose = verbose ) @@ -197,6 +201,7 @@ estimate_one_timecourse_params_tf <- function ( prior_pars = NULL, fit_intercept = FALSE, learning_rate = 0.1, + n_iterations = 100, verbose = FALSE ) { @@ -207,6 +212,7 @@ estimate_one_timecourse_params_tf <- function ( checkmate::assertNumber(learning_rate, lower = 0) checkmate::assertLogical(verbose, len = 1) checkmate::assertLogical(fit_intercept, len = 1) + checkmate::assertNumber(n_iterations, lower = 1) tfp <- reticulate::import("tensorflow_probability") @@ -368,6 +374,12 @@ estimate_one_timecourse_params_tf <- function ( ) } else { # constant noise + + noise_value <- ifelse(use_prior, 0.2, 1) + if (verbose) { + cli::cli_alert("No noise variable is present so setting noise to {.field {noise_value}}") + } + noise_entry <- matrix( ifelse(use_prior, 0.2, 1), nrow = nrow(one_timecourse), @@ -407,7 +419,7 @@ estimate_one_timecourse_params_tf <- function ( continue <- TRUE while (continue) { # train - for (i in 1:100) { + for (i in 1:n_iterations) { sess$run( train, feed_dict = timecourse_dict diff --git a/data/example_timecourse.rda b/data/example_timecourse.rda new file mode 100644 index 0000000000000000000000000000000000000000..16730877fc1b45403dd8d2181a68d5e9e94b9f76 GIT binary patch literal 733 zcmV<30wVoFT4*^jL0KkKS!ztf$N&Ml|NsC0+xY9h{eJDV64n2&-r>ON#X&d{?&voL zPl0#M7m>gM90Y`jh*N0>hJ!%R4FRTr0BH3z0iYhCp`+3O05lC54^Yr)w9;tysp@)& z88Ajn2-6VL2*fnN1Ykx244Py#WHbgQ0%8~mh|^3=34sG92+4sOVj5u>hL`}12*5#; zOoohxfW*K|Ljf@wX^DX_2_i_+dV!-uLTwP!KpJQUhJb0H217xkAO?WY(@h#_qiSe0 z00vAVa9LNjk56wtM9iwBAQf&Rwv|MsHA9IU-&7%;o=RY5fIW(uk}xeVhc8-?u9}2X zFeZYRE-JC1_m!6wL^% + dplyr::select(variable, value) %>% + tidyr::spread(variable, value), + timepts = unique(example_timecourse$time), + model = "sigmoid" + ) + + combined_data <- example_timecourse %>% + dplyr::left_join(fitted_timecourses, by = "time") + + #ggplot(combined_data, aes(x = time, y = abundance)) + + # geom_point() + + # geom_line(aes(y = fit), color = "red") + + # ggtitle("Example timecourse with sigmoid fit") + + # theme_minimal() + + timecourse_logLik <- combined_data %>% + dplyr::mutate( + normal_logLik = dnorm(abundance, fit, sd = 0.2, log = TRUE) + ) %>% + dplyr::summarize(logLik = sum(normal_logLik)) + + testthat::expect_equal( + object = timecourse_logLik$logLik, + expected = consensus_kinetics$loss$logLik[1], + tolerance = 0.01 + ) }) From e583b0fe8d036073b7d75296cde1e4babecda290 Mon Sep 17 00:00:00 2001 From: Sean Hackett Date: Tue, 28 May 2024 16:48:05 -0700 Subject: [PATCH 2/3] minor --- tests/testthat/test-impulse_fitting.R | 8 -------- 1 file changed, 8 deletions(-) diff --git a/tests/testthat/test-impulse_fitting.R b/tests/testthat/test-impulse_fitting.R index 1dd9dbc..31add82 100644 --- a/tests/testthat/test-impulse_fitting.R +++ b/tests/testthat/test-impulse_fitting.R @@ -112,14 +112,6 @@ test_that("Bayesian impulse minimizes MAP estimate", { test_that("Sigmoid Likelihood is correct", { - #CONDA_ENV <- "base397" - #CONDA_PATH <- "/scratch4/alex/python39-for-sean/miniforge3/bin/conda" - #reticulate::use_condaenv( - # CONDA_ENV, - # conda = CONDA_PATH, - # required=TRUE - #) - fitted_kinetics <- impulse::estimate_timecourse_params_tf( measurements = example_timecourse, prior_pars = c(v_sd = 1.2, rate_shape = 2, rate_scale = 0.25, time_shape = 2, time_scale = 15), From 2b71cf85ae918bb107d55ea56483e815f7558600 Mon Sep 17 00:00:00 2001 From: Sean Hackett Date: Tue, 28 May 2024 16:51:44 -0700 Subject: [PATCH 3/3] updating r cmd check action --- .github/workflows/R-CMD-check.yaml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/workflows/R-CMD-check.yaml b/.github/workflows/R-CMD-check.yaml index 8aff012..a837faf 100644 --- a/.github/workflows/R-CMD-check.yaml +++ b/.github/workflows/R-CMD-check.yaml @@ -27,9 +27,9 @@ jobs: GITHUB_PAT: ${{ secrets.GITHUB_TOKEN }} steps: - - uses: actions/checkout@v2 + - uses: actions/checkout@v4 - - uses: r-lib/actions/setup-r@v1 + - uses: r-lib/actions/setup-r@v2 with: r-version: ${{ matrix.config.r }}