diff --git a/R/add_population_age_to.R b/R/add_population_age_to.R index aca63bb1..bfc649df 100644 --- a/R/add_population_age_to.R +++ b/R/add_population_age_to.R @@ -31,7 +31,9 @@ add_population_age_to <- function(contact_data, population = get_polymod_population()) { # get function to interpolate population age distributions to 1y bins - age_population_function <- get_age_population_function(population) + age_population_function <- get_age_population_function( + population + ) # add the population in each 'to' age for the survey context contact_data %>% diff --git a/R/age_ranges.R b/R/age_ranges.R new file mode 100644 index 00000000..7d47064f --- /dev/null +++ b/R/age_ranges.R @@ -0,0 +1,15 @@ +#' @title +#' @param ages vector of ages +#' @noRd +#' @return vector length 1 of min and max ages used in calculating age +#' @author njtierney +age_ranges <- function(ages) { + + min_age <- min(ages) + bin_widths <- diff(ages) + final_bin_width <- bin_widths[length(bin_widths)] + age_max_integration <- max(ages) + final_bin_width + + c(min_age, age_max_integration) + +} \ No newline at end of file diff --git a/R/aggregate_predicted_contacts.R b/R/aggregate_predicted_contacts.R index 2459970d..0ce88102 100644 --- a/R/aggregate_predicted_contacts.R +++ b/R/aggregate_predicted_contacts.R @@ -73,8 +73,8 @@ aggregate_predicted_contacts <- function(predicted_contacts_1y, .groups = "drop" ) %>% # *average* the total contacts within the 'from' contacts, weighted by the - # population distribution (to get contacts for the population-average ember of - # that age group) + # population distribution (to get contacts for the population-average + # member of that age group) dplyr::mutate( pop_age_from = age_population_function(age_from), age_group_from = cut( diff --git a/R/make_bam_predictions.R b/R/make_bam_predictions.R new file mode 100644 index 00000000..2248e42a --- /dev/null +++ b/R/make_bam_predictions.R @@ -0,0 +1,39 @@ +#' .. content for \description{} (no empty lines) .. +#' +#' .. content for \details{} .. +#' +#' @title + +#' @return +#' @author njtierney +#' @export +make_bam_predictions <- function(x) { + + x %>% + # add on prediction features, setting the population to predict to + add_modelling_features( + population = population + ) %>% + dplyr::mutate( + # prediction + contacts = mgcv::predict.bam( + model, + newdata = ., + type = "response" + ), + # uncertainty + se_contacts = mgcv::predict.bam( + model, + newdata = ., + type = "response", + se.fit = TRUE + )$se.fit + ) %>% + dplyr::select( + age_from, + age_to, + contacts, + se_contacts + ) + +} diff --git a/R/predict_contacts.R b/R/predict_contacts.R index 5935c435..5f8a7122 100644 --- a/R/predict_contacts.R +++ b/R/predict_contacts.R @@ -71,12 +71,11 @@ predict_contacts <- function(model, population <- population %>% dplyr::arrange(!!age) - # this could be changed to a function for lower age limit - age_min_integration <- min(population[[age_var]]) - bin_widths <- diff(population[[age_var]]) - final_bin_width <- bin_widths[length(bin_widths)] - age_max_integration <- max(population[[age_var]]) + final_bin_width - + # get the age ranges plus the final age bin (which is the bin width) + the_age_ranges <- age_ranges(population$lower.age.limit) + age_min_integration <- min(the_age_ranges) + age_max_integration <- max(the_age_ranges) + # need to check we are not predicting to 0 populations (interpolator can # predict 0 values, then the aggregated ages get screwed up) pop_fun <- get_age_population_function(population) diff --git a/R/predict_contacts_1y.R b/R/predict_contacts_1y.R index f97949dd..f71f6779 100644 --- a/R/predict_contacts_1y.R +++ b/R/predict_contacts_1y.R @@ -29,9 +29,9 @@ #' contacts and standard error around the prediction. #' @examples #' -#' fairfield <- abs_age_lga("Fairfield (C)") +#' fairfield_abs_data <- abs_age_lga("Fairfield (C)") #' -#' fairfield +#' fairfield_abs_data #' #' # predict the contact rates in 1 year blocks to Fairfield data #' @@ -42,37 +42,19 @@ #' age_max = 2 #' ) #' @export -predict_contacts_1y <- function(model, population, age_min = 0, age_max = 100) { +predict_contacts_1y <- function( + model, + population, + age_min = 0, + age_max = 100 +) { all_ages <- age_min:age_max - # predict contacts to all integer years, adjusting for the population in a given place + # predict contacts to all integer years, adjusting for the population in a + # given place tidyr::expand_grid( age_from = all_ages, age_to = all_ages, ) %>% - # add on prediction features, setting the population to predict to - add_modelling_features( - population = population - ) %>% - dplyr::mutate( - # prediction - contacts = mgcv::predict.bam( - model, - newdata = ., - type = "response" - ), - # uncertainty - se_contacts = mgcv::predict.bam( - model, - newdata = ., - type = "response", - se.fit = TRUE - )$se.fit - ) %>% - dplyr::select( - age_from, - age_to, - contacts, - se_contacts - ) + make_bam_predictions() } diff --git a/R/predict_contacts_flexible.R b/R/predict_contacts_flexible.R new file mode 100644 index 00000000..cbb116c7 --- /dev/null +++ b/R/predict_contacts_flexible.R @@ -0,0 +1,16 @@ +# taken from predict_contacts_1y +#' @export +predict_contacts_flexible <- function( + model, + population, + age_vector +) { + + # predict contacts to all given years + tidyr::expand_grid( + age_from = age_vector, + age_to = age_vector, + ) %>% + make_bam_predictions() + +} diff --git a/R/predict_setting_contacts.R b/R/predict_setting_contacts.R index 0d2cb376..bc54045b 100644 --- a/R/predict_setting_contacts.R +++ b/R/predict_setting_contacts.R @@ -76,7 +76,7 @@ predict_setting_contacts <- function(population, .x = contact_model, .f = predict_contacts, population = population, - age_breaks = age_breaks, + age_breaks = age_breaks .options = furrr::furrr_options(seed = TRUE) )