The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

What-if: calibrating exposure

library(ambre)
set.seed(2024)
library(dplyr)
library(purrr)

Two numbers decide how much of the water each person actually meets:

Together with the pathogen concentration, these drive the dose, and the dose drives the risk. This vignette is about the exposure side of the model: where those two numbers come from by default, how to read them, and how to override them to ask “what if?” — including the two most common management measures, personal protective equipment and reduced-contact irrigation.

If you are new to the pipeline, read vignette("a-get-started", package = "ambre") first. The barrier-side levers (treatment plants, on-field practices, their log-reductions) are a different story, told in vignette("b-initial-vs-new-scenario", package = "ambre").

How exposure is set by default

You never set volume, frequency or concentration by hand in a normal run. When you call inflow_concentration() (which run_qmra_intial_situation() and run_qmra_supplementary_process() call for you), it runs three updates in order:

inflow_concentration()
  ├─ update_pathogen()    # pathogen concentration
  ├─ update_frequency()   # number of events / year, per path
  └─ update_volume()      # litres per event, per path
scenario <- create_scenario(
  system.file("input_1culture_2pop.xlsx", package = "ambre")
)
scenario_conc <- inflow_concentration(
  scenario = scenario,
  pathogenName = "Campylobacter jejuni"
)

update_pathogen() look up inflow database. update_frequency() and update_volume() look up the exposure path of each scenario row and copy the path-specific numbers out of the shipped database config_ambre$path into that row’s embedded config$exposure table. So the model does not use one generic volume for every path — each crop × population × path combination gets its own exposure profile.

You can see the result. The config$exposure table has three rows — number_of_repeatings (the Monte-Carlo count), number_of_exposures (events per year) and volume_perEvent:

# Row 1 of this two-row scenario: irrigation-staff droplet ingestion
scenario_conc$config[[1]]$exposure
#> # A tibble: 3 × 10
#>   name                 type    value   min   max mode  mean  sd    meanlog sdlog
#>   <chr>                <chr>   <dbl> <dbl> <dbl> <lgl> <lgl> <lgl> <lgl>   <lgl>
#> 1 number_of_repeatings value    1000 NA    NA    NA    NA    NA    NA      NA   
#> 2 number_of_exposures  value      60 NA    NA    NA    NA    NA    NA      NA   
#> 3 volume_perEvent      triang…    NA  1e-3  1e-3 NA    NA    NA    NA      NA

Two quirks are worth knowing up front, and we return to them at the end:

Inspect the base numbers

Before changing anything, look at what the model is assuming.

Start from the PathogenName, you can see the default value of concentration for the simulated pathogen :

# Row 1 of this two-row scenario: Campylobacter jejuni concentration
dplyr::filter(scenario_conc$config[[1]]$inflow, PathogenName == "Campylobacter jejuni")
#> # A tibble: 1 × 13
#>   PathogenID PathogenName   PathogenGroup simulate type  value   min   max mode 
#>        <dbl> <chr>          <chr>            <dbl> <chr> <lgl> <dbl> <dbl> <lgl>
#> 1          1 Campylobacter… Bacteria             1 unif… NA      100  5000 NA   
#> # ℹ 4 more variables: mean <lgl>, sd <lgl>, meanlog <lgl>, sdlog <lgl>

Then start from a path description and get its PathID with query_exp_path():

pid <- query_exp_path(
  pathName = "Ingestion of water droplets during maintenance of the irrigation system"
)
pid
#> [1] 1

Then read the per-path volume and frequency straight from the database. Both helpers accept a vector of IDs, so you can inspect a whole scenario at once:

scenario$PathID          # the two paths in this input file
#> [1] 1 4

query_volume(scenario$PathID)     # litres per event: min / max
#> # A tibble: 2 × 2
#>     min   max
#>   <dbl> <dbl>
#> 1 0.001 0.001
#> 2 0.001 0.001
query_frequency(scenario$PathID)  # events per year: min / max
#> # A tibble: 2 × 2
#>     min   max
#>   <dbl> <dbl>
#> 1    40    60
#> 2    48    48

These are exactly the numbers update_volume() and update_frequency() inject. The volumes here are tiny — 0.001 litre is one millilitre of incidentally swallowed droplets — the frequency 4060 becomes 60 because only the maximum is used and the concentration is between 100 and 500 log/L.

Override for a what-if

To explore sensitivity you replace those numbers with your own. Two helpers do this, and both expect one value per scenario row.

update_volume_with_desired_value() takes a volume argument that is a data frame (or list) with a min, a max and a type column. Use min == max for a fixed volume, or min < max for a distribution spread according to the type specified:

scenario_ppe <- update_volume_with_desired_value(
  scenario = scenario,
  volume = data.frame(
    min = c(0.0005, 0.0005),   # 0.5 mL per event, one entry per row
    max = c(0.0005, 0.0005),
    type = c("triangle", "triangle")
  )
)

scenario_ppe$config[[1]]$exposure
#> # A tibble: 3 × 10
#>   name                 type    value   min   max mode  mean  sd    meanlog sdlog
#>   <chr>                <chr>   <dbl> <dbl> <dbl> <lgl> <lgl> <lgl> <lgl>   <lgl>
#> 1 number_of_repeatings value    1000 NA    NA    NA    NA    NA    NA      NA   
#> 2 number_of_exposures  value     365 NA    NA    NA    NA    NA    NA      NA   
#> 3 volume_perEvent      triang…    NA  5e-4  5e-4 NA    NA    NA    NA      NA

update_frequency_with_desired_value() takes a plain numeric vector of integer, one number of events per year per row:

scenario_night <- update_frequency_with_desired_value(
  scenario = scenario,
  frequency = c(30L, 30L)   # cap both paths at 30 events / year
)

scenario_night$config[[1]]$exposure
#> # A tibble: 3 × 10
#>   name                 type    value   min   max mode  mean  sd    meanlog sdlog
#>   <chr>                <chr>   <dbl> <dbl> <dbl> <lgl> <lgl> <lgl> <lgl>   <lgl>
#> 1 number_of_repeatings value    1000  NA      NA NA    NA    NA    NA      NA   
#> 2 number_of_exposures  value      30  NA      NA NA    NA    NA    NA      NA   
#> 3 volume_perEvent      triang…    NA   0.5     3 NA    NA    NA    NA      NA

The volume_perEvent row now reads 0.0005, and number_of_exposures reads 30: your what-if values have replaced the ones the database would have injected.

update_concentration() takes a data.frame with PathogenName, min, max and type of ditribution law values.

concentration_custom <- data.frame(PathogenName = c("Campylobacter jejuni"),
                                   min = c(1),
                                   max = c(2),
                                   type = c("uniform"))

scenario_pathogen <- update_pathogen(scenario = scenario, pathoName = concentration_custom$PathogenName)
scenario_low_pathogen <- update_concentration(scenario = scenario_pathogen ,
                     concentration = concentration_custom)

dplyr::filter(scenario_low_pathogen$config[[1]]$inflow, PathogenName == "Campylobacter jejuni")
#> # A tibble: 1 × 13
#>   PathogenID PathogenName PathogenGroup simulate value mode  mean  sd    meanlog
#>        <dbl> <chr>        <chr>            <dbl> <lgl> <lgl> <lgl> <lgl> <lgl>  
#> 1          1 Campylobact… Bacteria             1 NA    NA    NA    NA    NA     
#> # ℹ 4 more variables: sdlog <lgl>, min <dbl>, max <dbl>, type <chr>

Modelling management measures

Read those two overrides as interventions. This is the exposure-side complement to the barrier log-reductions in vignette("b-treatment-vs-multibarrier", package = "ambre"): instead of removing pathogens from the water, you change how much of the water each person meets.

Both are assumptions you choose and should defend — they are not credited from a barrier database. (The barriere_path / barriere_specific tables in config_ambre that would encode route-specific exposure reductions are shipped data, not yet wired into the engine.)

A mini sensitivity study

Let us quantify how much the volume assumption matters. One catch first: run_qmra_initial_situation() and run_qmra_supplementary_process() re-derives the per-path volumes from the database on every call (via inflow_concentration()), so it would overwrite any override you set. To inject a what-if value you should use run_qmra_custom function:

Run it twice on the same file and pathogen — once with a pessimistic bare-hand volume, once with a PPE volume ten times smaller:

sc <- create_scenario(system.file("input_1culture_2pop.xlsx", package = "ambre"))

regulation_reduction <- config_ambre$regulation$regulation_value |> 
  dplyr::filter(Country == "France") |>
  dplyr::select(-c(Concentration, Country, RegulationID))
regulation_concentration <- config_ambre$regulation$regulation_value |> 
  dplyr::filter(Country == "France") |>
  dplyr::select(-c(Country, RegulationID, Reduction))

concentration_custom <- data.frame(PathogenName = c("Rotavirus"),
                                   min = c(1000),
                                   max = c(2000),
                                   type = c("uniform"))

set.seed(2024)
bare <- run_qmra_custom(scenario = sc,
                  concentration = concentration_custom,
                  volume = data.frame(min=c(0.005,0.005), max = c(0.01, 0.01), type = c("triangle", "triangle")),
                  frequency = c(48L, 60L),
                  regulationLog = regulation_reduction,
                  regulationConcentration = regulation_concentration,
                  initialSituation = TRUE) # min 5 mL, max 10 mL

set.seed(2024)
ppe  <- run_qmra_custom(scenario = sc,
                  concentratio = concentration_custom,
                  volume = data.frame(min=c(0.0005,0.0005), max = c(0.001, 0.001), type = c("triangle", "triangle")),
                  frequency = c(48L, 60L),
                  regulationLog = regulation_reduction,
                  regulationConcentration = regulation_concentration,
                  initialSituation = TRUE) # min 0.5 mL min 0.1 mL

As the the other functions run_qmra_* this custom function return 2 plots and to formattable. See d-interpreting-riskfor more detail on the output.

library(ggplot2)
cowplot::plot_grid(
      bare$dalys$Rotavirus +
        labs(subtitle = "min 5 mL, max 10 mL") +
        theme(plot.subtitle = element_text(hjust = 0.5)),
      ppe$dalys$Rotavirus +
        labs(subtitle = "min 0.5 mL min 0.1 mL") +
        theme(plot.subtitle = element_text(hjust = 0.5)),
      ncol = 2,
      align = "h"
    )

The risk tracks the volume almost proportionally: a tenfold cut in swallowed volume gives roughly a tenfold cut in DALYs. More tellingly, the 95th percentile crosses the line — the bare-hand upper tail sits above 1e-6 while the PPE tail drops below it. Deciding whether a scenario meets the target can come down to this single exposure assumption. For how to read these ranges against the regulatory line, see vignette("d-interpreting-risk", package = "ambre"); for why every number is a distribution rather than a point, see vignette("g-monte-carlo-engine", package = "ambre").

Technical notes and caveats

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.