greta-dev / greta.dynamics

a greta extension for modelling dynamical systems
https://greta-dev.github.io/greta.dynamics/
Other
6 stars 2 forks source link

add optional max and min state values in iterate_dynamic_function() to help prevent numerical over/underflow #29

Closed goldingn closed 8 months ago

goldingn commented 8 months ago

only on greta_2 branch for now

goldingn commented 8 months ago

e.g. ideally this wouldn't overflow to Infinity

library(greta.dynamics)
#> Loading required package: greta
#> 
#> Attaching package: 'greta'
#> The following objects are masked from 'package:stats':
#> 
#>     binomial, cov2cor, poisson
#> The following objects are masked from 'package:base':
#> 
#>     %*%, apply, backsolve, beta, chol2inv, colMeans, colSums, diag,
#>     eigen, forwardsolve, gamma, identity, rowMeans, rowSums, sweep,
#>     tapply
fun <- function(state, iter, r) {
  state * r
}
n_times <- 200
r <- normal(0, 100, truncation = c(0, Inf))
#> ℹ Initialising python and checking dependencies, this may take a moment.
#> ✔ Initialising python and checking dependencies ... done!
#> 
results <- iterate_dynamic_function(fun,
                                    initial_state = 1,
                                    tol = 0,
                                    niter = n_times,
                                    r = r)
set.seed(1)
vals <- calculate(results$all_states,
                  nsim = 1)
tail(vals[[1]][1, , ], 30)
#>  [1] 6.717551e+296 3.656307e+298 1.990097e+300 1.083193e+302 5.895729e+303
#>  [6] 3.208996e+305 1.746630e+307           Inf           Inf           Inf
#> [11]           Inf           Inf           Inf           Inf           Inf
#> [16]           Inf           Inf           Inf           Inf           Inf
#> [21]           Inf           Inf           Inf           Inf           Inf
#> [26]           Inf           Inf           Inf           Inf           Inf

Created on 2024-03-13 with reprex v2.0.2

goldingn commented 8 months ago

30 looks like this now:

library(greta.dynamics)
#> Loading required package: greta
#> 
#> Attaching package: 'greta'
#> The following objects are masked from 'package:stats':
#> 
#>     binomial, cov2cor, poisson
#> The following objects are masked from 'package:base':
#> 
#>     %*%, apply, backsolve, beta, chol2inv, colMeans, colSums, diag,
#>     eigen, forwardsolve, gamma, identity, rowMeans, rowSums, sweep,
#>     tapply
fun <- function(state, iter, r) {
  state * r
}
n_times <- 200
r <- normal(0, 100, truncation = c(0, Inf))
#> ℹ Initialising python and checking dependencies, this may take a moment.
#> ✔ Initialising python and checking dependencies ... done!
#> 
results <- iterate_dynamic_function(fun,
                                    initial_state = 1,
                                    tol = 0,
                                    niter = n_times,
                                    r = r)
set.seed(1)
vals <- calculate(results$all_states,
                  nsim = 1)
tail(vals[[1]][1, , ], 30)
#>  [1] 6.717551e+296 3.656307e+298 1.990097e+300 1.083193e+302 5.895729e+303
#>  [6] 3.208996e+305 1.746630e+307           Inf           Inf           Inf
#> [11]           Inf           Inf           Inf           Inf           Inf
#> [16]           Inf           Inf           Inf           Inf           Inf
#> [21]           Inf           Inf           Inf           Inf           Inf
#> [26]           Inf           Inf           Inf           Inf           Inf

r <- normal(0, 100, truncation = c(0, Inf))
results <- iterate_dynamic_function(fun,
                                    initial_state = 1,
                                    tol = 0,
                                    niter = n_times,
                                    r = r,
                                    state_limits = c(0, 1e3))
set.seed(1)
vals <- calculate(results$all_states,
                  nsim = 1)
tail(vals[[1]][1, , ], 30)
#>  [1] 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000
#> [16] 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000 1000

Created on 2024-03-13 with reprex v2.0.2