Skip to contents

Fit a Model

Source: R/bage_mod-methods.R
fit.bage_mod.Rd

Calculate the posterior distribution for a model.

Usage

# S3 method for class 'bage_mod'
fit(
  object,
  method = c("standard", "inner-outer"),
  vars_inner = NULL,
  optimizer = c("nlminb", "BFGS", "CG"),
  quiet = TRUE,
  start_oldpar = FALSE,
  ...
)

Arguments

object

A bage_mod object, created with mod_pois(), mod_binom(), or mod_norm().

method

Estimation method. Current choices are "standard" (the default) and "inner-outer". See below for details.

vars_inner

Names of variables to use for inner model when method is "inner-outer". If NULL(the default)vars_inner` is the age, sex/gender, and time variables.

optimizer

Which optimizer to use. Current choices are "nlminb" (stats::nlminb()) and "BFGS" (stats::optim() using method "BFGS"), and "GC" (stats::optim() using method "CG") . Default is "nlminb".

quiet

Whether to suppress warnings and progress messages from the optimizer. Default is TRUE.

start_oldpar

Whether the optimizer should start at previous estimates. Used only when fit() is being called on a fitted model. Default is FALSE.

...

Not currently used.

Value

A bage_mod object

Estimation methods

  • "standard" All parameters, other than the lowest-level rates, probabilities, or means are jointly estimated within TMB. The default.

  • "inner-outer". Multiple-stage estimation, which can be faster than "standard" for models with many parameters. In Step 1, the data is aggregated across all dimensions other than those specified in var_inner, and a model for the inner variables is fitted to the data. In Step 2, the data is aggregated across the remaining variables, and a model for the outer variables is fitted to the data. In Step 3, values for dispersion are calculated. Parameter estimtes from steps 1, 2, and 3 are then combined. "inner-outer" methods are still experimental, and may change in future, eg dividing calculations into chunks in Step 2.

See also

Examples

## specify model
mod <- mod_pois(injuries ~ age + sex + year,
                data = nzl_injuries,
                exposure = popn)

## examine unfitted model
mod
#> 
#>     ------ Unfitted Poisson model ------
#> 
#> 
#>    injuries ~ age + sex + year
#> 
#>   exposure = popn
#> 
#> 
#>         term  prior along n_par n_par_free
#>  (Intercept) NFix()     -     1          1
#>          age   RW()   age    12         11
#>          sex NFix()     -     2          2
#>         year   RW()  year    19         18
#> 
#> 
#>  n_draw pr_mean_disp var_time var_age var_sexgender
#>    1000            1     year     age           sex
#> 

## fit model
mod <- fit(mod)

## examine fitted model
mod
#> 
#>     ------ Fitted Poisson model ------
#> 
#> 
#>    injuries ~ age + sex + year
#> 
#>   exposure = popn
#> 
#> 
#>         term  prior along n_par n_par_free std_dev
#>  (Intercept) NFix()     -     1          1       -
#>          age   RW()   age    12         11    0.76
#>          sex NFix()     -     2          2    0.72
#>         year   RW()  year    19         18    0.08
#> 
#> 
#>  n_draw pr_mean_disp var_time var_age var_sexgender optimizer
#>    1000            1     year     age           sex    nlminb
#> 
#> 
#>  time_total time_optim time_draws iter converged                    message
#>        0.18       0.07       0.00   12      TRUE   relative convergence (4)
#> 

## extract rates
aug <- augment(mod)
aug
#> # A tibble: 912 × 9
#>    age   sex    ethnicity  year injuries  popn .observed
#>    <fct> <chr>  <chr>     <int>    <int> <int>     <dbl>
#>  1 0-4   Female Maori      2000       12 35830 0.000335 
#>  2 5-9   Female Maori      2000        6 35120 0.000171 
#>  3 10-14 Female Maori      2000        3 32830 0.0000914
#>  4 15-19 Female Maori      2000        6 27130 0.000221 
#>  5 20-24 Female Maori      2000        6 24380 0.000246 
#>  6 25-29 Female Maori      2000        6 24160 0.000248 
#>  7 30-34 Female Maori      2000       12 22560 0.000532 
#>  8 35-39 Female Maori      2000        3 22230 0.000135 
#>  9 40-44 Female Maori      2000        6 18130 0.000331 
#> 10 45-49 Female Maori      2000        6 13770 0.000436 
#> # ℹ 902 more rows
#> # ℹ 2 more variables: .fitted <rdbl<1000>>, .expected <rdbl<1000>>

## extract hyper-parameters
comp <- components(mod)
comp
#> # A tibble: 37 × 4
#>    term        component level                 .fitted
#>    <chr>       <chr>     <chr>            <rdbl<1000>>
#>  1 (Intercept) effect    (Intercept) -5.8 (-6.9, -4.7)
#>  2 age         effect    0-4                  0 (0, 0)
#>  3 age         effect    5-9         -1.3 (-1.5, -1.2)
#>  4 age         effect    10-14       -0.87 (-1, -0.75)
#>  5 age         effect    15-19       0.86 (0.76, 0.97)
#>  6 age         effect    20-24           1 (0.91, 1.1)
#>  7 age         effect    25-29       0.86 (0.75, 0.96)
#>  8 age         effect    30-34       0.75 (0.65, 0.86)
#>  9 age         effect    35-39       0.76 (0.66, 0.87)
#> 10 age         effect    40-44       0.76 (0.65, 0.86)
#> # ℹ 27 more rows