Parametric terms with binning: fitting


Birthweight model including covariates at municipality level with linear effects.

Load packages, read data and source custom scripts

Paths are defined relative to the git repository location.

rm(list = ls())
library(bamlss)
#> Loading required package: coda
#> Loading required package: colorspace
#> Loading required package: mgcv
#> Loading required package: nlme
#> This is mgcv 1.8-31. For overview type 'help("mgcv-package")'.
#> 
#> Attaching package: 'bamlss'
#> The following object is masked from 'package:mgcv':
#> 
#>     smooth.construct
library(gamlss.dist)
#> Loading required package: MASS
path_proj <- day2day::git_path()
path_data <- file.path(path_proj, "data")
path_processed <- file.path(path_data, "processed")
path_modelled <- file.path(path_data, "modelled")

path_modelled_data <- file.path(path_modelled, "bw-muni-03-covs-param-bin.rds")
path_modelled_sink <- gsub("\\.rds$", "\\.txt", path_modelled_data)
path_modelled_form <- gsub("(\\.rds)$", "-form\\1", path_modelled_data)

bwdata_model <- fst::read_fst(file.path(path_processed, "bwdata_41_model.fst"))

Define formula for our model

Now we define the same models as in the previous study.

form_sigma <- sigma ~ 1

form_mu <- born_weight ~ poly(remoteness, 2) + poly(prop_tap_toilet, 2)

form <- list(form_mu, form_sigma)

Run the model of interest and save results

{
    sink(path_modelled_sink)
    bamlss_model <- bamlss(
        form, data = bwdata_model, binning = TRUE,
        n.iter = 1000, burnin = 0, cores = 4, combine = FALSE, light = TRUE
    )
    sink()
}
readLines(path_modelled_sink)
#>  [1] "AICc 4458216. logPost -2229216 logLik -2229102 edf 6.0000 eps 0.0299 iteration   1"
#>  [2] "AICc 4458216. logPost -2229216 logLik -2229102 edf 6.0000 eps 0.0000 iteration   2"
#>  [3] "AICc 4458216. logPost -2229216 logLik -2229102 edf 6.0000 eps 0.0000 iteration   2"
#>  [4] "elapsed time:  0.39sec"                                                            
#>  [5] "Starting the sampler..."                                                           
#>  [6] "Starting the sampler..."                                                           
#>  [7] "Starting the sampler..."                                                           
#>  [8] "Starting the sampler..."                                                           
#>  [9] ""                                                                                  
#> [10] "|                    |   0%  5.61min"                                              
#> [11] "|                    |   0%  5.58min"                                              
#> [12] "|                    |   0%  5.66min"                                              
#> [13] "|                    |   0%  5.68min"                                              
#> [14] "|*                   |   5%  4.83min 15.25sec"                                     
#> [15] "|*                   |   5%  4.99min 15.77sec"                                     
#> [16] "|*                   |   5%  4.99min 15.76sec"                                     
#> [17] "|*                   |   5%  5.11min 16.13sec"                                     
#> [18] "|**                  |  10%  4.47min 29.78sec"                                     
#> [19] "|**                  |  10%  4.54min 30.28sec"                                     
#> [20] "|**                  |  10%  4.68min 31.17sec"                                     
#> [21] "|**                  |  10%  4.75min 31.67sec"                                     
#> [22] "|***                 |  15%  4.14min 43.86sec"                                     
#> [23] "|***                 |  15%  4.23min 44.78sec"                                     
#> [24] "|***                 |  15%  4.32min 45.73sec"                                     
#> [25] "|***                 |  15%  4.48min 47.45sec"                                     
#> [26] "|****                |  20%  3.89min 58.40sec"                                     
#> [27] "|****                |  20%  3.97min 59.48sec"                                     
#> [28] "|****                |  20%  4.01min  1.00min"                                     
#> [29] "|****                |  20%  4.17min  1.04min"                                     
#> [30] "|*****               |  25%  3.63min  1.21min"                                     
#> [31] "|*****               |  25%  3.70min  1.23min"                                     
#> [32] "|*****               |  25%  3.73min  1.24min"                                     
#> [33] "|*****               |  25%  3.85min  1.28min"                                     
#> [34] "|******              |  30%  3.39min  1.45min"                                     
#> [35] "|******              |  30%  3.46min  1.48min"                                     
#> [36] "|******              |  30%  3.49min  1.49min"                                     
#> [37] "|******              |  30%  3.56min  1.53min"                                     
#> [38] "|*******             |  35%  3.15min  1.70min"                                     
#> [39] "|*******             |  35%  3.20min  1.72min"                                     
#> [40] "|*******             |  35%  3.23min  1.74min"                                     
#> [41] "|*******             |  35%  3.28min  1.77min"                                     
#> [42] "|********            |  40%  2.90min  1.93min"                                     
#> [43] "|********            |  40%  2.95min  1.97min"                                     
#> [44] "|********            |  40%  2.96min  1.97min"                                     
#> [45] "|********            |  40%  3.01min  2.00min"                                     
#> [46] "|*********           |  45%  2.66min  2.17min"                                     
#> [47] "|*********           |  45%  2.69min  2.20min"                                     
#> [48] "|*********           |  45%  2.71min  2.22min"                                     
#> [49] "|*********           |  45%  2.75min  2.25min"                                     
#> [50] "|**********          |  50%  2.41min  2.41min"                                     
#> [51] "|**********          |  50%  2.45min  2.45min"                                     
#> [52] "|**********          |  50%  2.46min  2.46min"                                     
#> [53] "|**********          |  50%  2.49min  2.49min"                                     
#> [54] "|***********         |  55%  2.17min  2.65min"                                     
#> [55] "|***********         |  55%  2.20min  2.69min"                                     
#> [56] "|***********         |  55%  2.21min  2.70min"                                     
#> [57] "|***********         |  55%  2.23min  2.73min"                                     
#> [58] "|************        |  60%  1.92min  2.89min"                                     
#> [59] "|************        |  60%  1.95min  2.93min"                                     
#> [60] "|************        |  60%  1.96min  2.94min"                                     
#> [61] "|************        |  60%  1.98min  2.97min"                                     
#> [62] "|*************       |  65%  1.69min  3.13min"                                     
#> [63] "|*************       |  65%  1.71min  3.18min"                                     
#> [64] "|*************       |  65%  1.72min  3.19min"                                     
#> [65] "|*************       |  65%  1.73min  3.22min"                                     
#> [66] "|**************      |  70%  1.45min  3.38min"                                     
#> [67] "|**************      |  70%  1.47min  3.42min"                                     
#> [68] "|**************      |  70%  1.47min  3.43min"                                     
#> [69] "|**************      |  70%  1.48min  3.46min"                                     
#> [70] "|***************     |  75%  1.21min  3.62min"                                     
#> [71] "|***************     |  75%  1.22min  3.67min"                                     
#> [72] "|***************     |  75%  1.22min  3.67min"                                     
#> [73] "|***************     |  75%  1.23min  3.70min"                                     
#> [74] "|****************    |  80% 57.92sec  3.86min"                                     
#> [75] "|****************    |  80% 58.68sec  3.91min"                                     
#> [76] "|****************    |  80% 58.72sec  3.91min"                                     
#> [77] "|****************    |  80% 59.10sec  3.94min"                                     
#> [78] "|*****************   |  85% 43.42sec  4.10min"                                     
#> [79] "|*****************   |  85% 43.94sec  4.15min"                                     
#> [80] "|*****************   |  85% 43.98sec  4.15min"                                     
#> [81] "|*****************   |  85% 44.27sec  4.18min"                                     
#> [82] "|******************  |  90% 28.95sec  4.34min"                                     
#> [83] "|******************  |  90% 29.27sec  4.39min"                                     
#> [84] "|******************  |  90% 29.31sec  4.40min"                                     
#> [85] "|******************  |  90% 29.48sec  4.42min"                                     
#> [86] "|******************* |  95% 14.47sec  4.58min"                                     
#> [87] "|******************* |  95% 14.63sec  4.63min"                                     
#> [88] "|******************* |  95% 14.63sec  4.63min"                                     
#> [89] "|******************* |  95% 14.72sec  4.66min"                                     
#> [90] "|********************| 100%  0.00sec  4.83min"                                     
#> [91] ""                                                                                  
#> [92] "|********************| 100%  0.00sec  4.87min"                                     
#> [93] ""                                                                                  
#> [94] "|********************| 100%  0.00sec  4.87min"                                     
#> [95] ""                                                                                  
#> [96] "|********************| 100%  0.00sec  4.88min"
system.time(saveRDS(bamlss_model, file = path_modelled_data))
#>    user  system elapsed 
#>   0.571   0.000   0.571
saveRDS(form, file = path_modelled_form)

Time to execute the task

Only useful when executed with Rscript.

proc.time()
#>     user   system  elapsed 
#> 1199.991    5.477  310.706