Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

how to pass parameters to control.fixed? #46

Open
ThierryO opened this issue Aug 27, 2018 · 6 comments
Open

how to pass parameters to control.fixed? #46

ThierryO opened this issue Aug 27, 2018 · 6 comments

Comments

@ThierryO
Copy link

Create some example data

set.seed(20180827)
n <- 1000
dataset <- data.frame(
  Area = runif(n, min = 1, max = 4),
  X = rnorm(n, mean = 2, sd = 1)
)
dataset$logArea <- log(dataset$Area)
dataset$eta <- 2 + 0.7 * dataset$X
dataset$lambda <- exp(dataset$eta) * dataset$Area
dataset$Y <- rpois(n, dataset$lambda)

Model it with INLA. First using the defaults, then fixed the parameter of logArea to 1 by setting the prior to mean = 1 and var = 0.

library(INLA)
inla(
  Y ~ X + logArea,
  family = "poisson",
  data = dataset
)$summary.fixed
inla(
  Y ~ X + logArea,
  family = "poisson",
  control.fixed = list(
    mean = list(logArea = 1, default = 0),
    prec = list(logArea = 1e20, default = 0.001)
  ),
  data = dataset
)$summary.fixed

bru() seems to ignore the options when passed as a list. Passing the output of bru.options() yields an error.

library(inlabru)
m <- bru(
  Y ~ X + logArea,
  family = "poisson",
  options = list(
    control.fixed = list(
      mean = list(logArea = 1, default = 0),
      prec = list(logArea = 1e20, default = 0.001)
    )
  ),
  data = dataset
)
summary(m)
bru(
  Y ~ X + logArea,
  family = "poisson",
  options = bru.options(
    control.fixed = list(
      mean = list(logArea = 1, default = 0),
      prec = list(logArea = 1e20, default = 0.001)
    )
  ),
  data = dataset
)

Error in (function (data, model, stackmaker, n = 10, result = NULL, family, :
INLA returned message: unused argument (inla.options = list(list(TRUE, TRUE, TRUE), list("auto"), list(list(1, 0), list(1e+20, 0.001))))
In addition: Warning message:
In (function (formula, family = "gaussian", contrasts = NULL, data, :

Error in (function (data, model, stackmaker, n = 10, result = NULL, family, :
INLA returned message: unused argument (inla.options = list(list(TRUE, TRUE, TRUE), list("auto"), list(list(1, 0), list(1e+20, 0.001))))

SessionInfo

Session info --------------------------------------------------------------------------------------------------
 setting  value                       
 version  R version 3.5.1 (2018-07-02)
 system   x86_64, linux-gnu           
 ui       RStudio (1.1.453)           
 language nl:en                       
 collate  nl_NL.UTF-8                 
 tz       Europe/Brussels             
 date     2018-08-27                  

Packages ------------------------------------------------------------------------------------------------------
 package      * version date       source                         
 assertthat     0.2.0   2017-04-11 CRAN (R 3.5.1)                 
 cli            1.0.0   2017-11-05 CRAN (R 3.5.1)                 
 colorspace     1.3-2   2016-12-14 CRAN (R 3.5.1)                 
 crayon         1.3.4   2017-09-16 CRAN (R 3.5.1)                 
 dichromat      2.0-0   2013-01-24 CRAN (R 3.5.1)                 
 digest         0.6.15  2018-01-28 CRAN (R 3.5.1)                 
 fansi          0.2.3   2018-05-06 CRAN (R 3.5.1)                 
 ggplot2      * 3.0.0   2018-07-03 CRAN (R 3.5.1)                 
 glue           1.3.0   2018-07-17 CRAN (R 3.5.1)                 
 graphics     * 3.5.1   2018-07-03 local                          
 grDevices    * 3.5.1   2018-07-03 local                          
 grid           3.5.1   2018-07-03 local                          
 gtable         0.2.0   2016-02-26 CRAN (R 3.5.1)                 
 inlabru      * 2.1.9   2018-08-23 Github (fbachl/inlabru@d77d689)
 labeling       0.3     2014-08-23 CRAN (R 3.5.1)                 
 lattice        0.20-35 2017-03-25 CRAN (R 3.5.0)                 
 lazyeval       0.2.1   2017-10-29 CRAN (R 3.5.1)                 
 magrittr       1.5     2014-11-22 CRAN (R 3.5.1)                 
 MASS           7.3-50  2018-04-30 CRAN (R 3.5.0)                 
 Matrix       * 1.2-14  2018-04-09 CRAN (R 3.5.0)                 
 methods      * 3.5.1   2018-07-03 local                          
 mgcv           1.8-24  2018-06-18 CRAN (R 3.5.0)                 
 munsell        0.5.0   2018-06-12 CRAN (R 3.5.1)                 
 nlme           3.1-137 2018-04-07 CRAN (R 3.5.0)                 
 pillar         1.3.0   2018-07-14 CRAN (R 3.5.1)                 
 plyr           1.8.4   2016-06-08 CRAN (R 3.5.1)                 
 R6             2.2.2   2017-06-17 CRAN (R 3.5.1)                 
 RColorBrewer   1.1-2   2014-12-07 CRAN (R 3.5.1)                 
 Rcpp           0.12.17 2018-05-18 CRAN (R 3.5.1)                 
 reshape2       1.4.3   2017-12-11 CRAN (R 3.5.1)                 
 rgdal          1.3-4   2018-08-03 cran (@1.3-4)                  
 rgeos          0.3-28  2018-06-08 CRAN (R 3.5.1)                 
 rlang          0.2.1   2018-05-30 CRAN (R 3.5.1)                 
 scales         0.5.0   2017-08-24 CRAN (R 3.5.1)                 
 sp           * 1.3-1   2018-06-05 CRAN (R 3.5.1)                 
 stats        * 3.5.1   2018-07-03 local                          
 stringi        1.2.3   2018-06-12 CRAN (R 3.5.1)                 
 stringr      * 1.3.1   2018-05-10 CRAN (R 3.5.1)                 
 tibble       * 1.4.2   2018-01-22 CRAN (R 3.5.1)                 
 tools          3.5.1   2018-07-03 local                          
 utf8           1.1.4   2018-05-24 CRAN (R 3.5.1)                 
 utils        * 3.5.1   2018-07-03 local                          
 viridisLite    0.3.0   2018-02-01 CRAN (R 3.5.1)                 
 withr          2.1.2   2018-03-15 CRAN (R 3.5.1)    
@finnlindgren
Copy link
Collaborator

I have an email thread with another user from a few weeks back that looked at this in some detail, including providing some solutions. I’ll try to extract the relevant bits and put here, awaiting a full solution in the package.

@finnlindgren
Copy link
Collaborator

finnlindgren commented Aug 27, 2018

Even if inlabru would recognise control.fixed and pass it on to inla(), this would not have the desired effect. The reason is that internally, fixed effect terms are re-interpreted as f(..., model="linear") terms before the inla() call, and control.fixed doesn't apply to such model components.

Until we come up with a way to automate it by using the information in control.fixed in the re-interpretation step, the following workaround can be used:
Replace

Y ~ X + logArea

by

Y ~ X(map=X, model="linear", mean.linear=0, prec.linear=0.001) +
  logArea(map=logArea, model="linear", mean.linear=1, prec.linear=1e20)

Full example attached:
code46.zip

Note that this also applies to the intercept term, so to control the prior for that you need to include it explicitly, with
+ Intercept(map=Inter, model="linear", prec.linear=...)
and have Inter=1 in the data data.frame.

Dev notes:
The sensible thing seems to be if we can parse control.fixed, but I'm not sure if the point in the processing where the inla formula is constructed is too early for that to work. (If it is, then it needs to be adjusted in the code before calling inla)

@ThierryO
Copy link
Author

The workaround works fine.

I would suggest to update the documentation that mentions control.fixed (?bru and ?bru.options) rather than putting a lot of effort in the parsing.

This issue can be closed.

@finnlindgren
Copy link
Collaborator

Ah, I didn’t realise that’s in the documentation! Will fix it!

@finnlindgren
Copy link
Collaborator

I've updated the documentation in commit 82764e7 and 48cd2cb, which also enables only partially overriding control.* parameters.
I'll keep this open as a wishlist enhancement item.

@finnlindgren
Copy link
Collaborator

Updated code (removing no longer needed control.fixed from the correct bru() call, and more precise comment on the output): code46.zip

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants