bayesian-demos

--- title: bayesian-demos tags: statistical-software-demos, statistical-software robots: noindex, nofollow --- Bayesian and Monte Carlo Software - Demonstration Application of Standards ========================================================================== This file demonstrates the application of [rOpenSci](https://ropensci.org) ’s [standards for statistical software](https://ropenscilabs.github.io/statistical-software-review-book/standards.html) to two [Bayesian software](https://ropenscilabs.github.io/statistical-software-review-book/standards.html#bayesian-and-monte-carlo-software) packages. These applications are not intended to represent or reflect evaluations or assessment of the packages, and particularly not of the extent to which they fail to meet standards. Rather, the demonstrations are intended to highlight aspects of the software which could be productively improved by adhering to the standards, and thereby more generally to demonstrate the general usefulness of these standards in advancing and improving software quality. [`fmcmc`](https://github.com/USCbiostats/fmcmc) ----------------------------------------------- ### 1 [General Standards](https://ropenscilabs.github.io/statistical-software-review-book/standards.html#general-standards-for-statistical-software) #### 1.1 Documentation - [x] **G1.0** The primary reference is to a JOSS article, with sufficient secondary references documented throughout. **Statistical Terminology** - [ ] **G1.1** Statistical terminology is generally used without sufficient clarification or definition. **Function-level Documentation** - [x] **G1.2** *Software should use [`roxygen`](https://roxygen2.r-lib.org/) to documentation all functions.* - [x] **G1.2a** internal functions are documented in standard [`roxygen`](https://roxygen2.r-lib.org/) format **Supplementary Documentation** - [ ] **G1.3** No code included to reproduce results which form the basis of performance claims made in associated publications. - [ ] **G1.4** No code included to compare performance claims with alternative implementations in other R packages. #### 1.2 Input Structures **Uni-variate (Vector) Input** - [x] **G2.0** Lengths of inputs generally asserted, and explicit secondary documentation of expected lengths provided. - [x] **G2.1** Types of inputs generally asserted, and explicit secondary documentation of expected types provided. - [ ] **G2.2** Submission of multivariate input to parameters expected to be univariate is possible, and not caught in pre-processing - [x] **G2.3** No functions have single-valued character inputs, so not applicable. - [ ] **G2.4** Provide appropriate mechanisms to convert between different *data types*, potentially including: - [ ] **G2.4a** There is no explicit conversion to `integer` via `as.integer()` - [x] **G2.4b** explicit conversion to continuous uses `as.numeric()` - [x] **G2.4c** explicit conversion to character uses `as.character()` - [x] **G2.4d**–i**G2.4e** No `factor` data used, so not applicable - [x] **G2.5** No `factor` data used, so not applicable **Tabular Input** - [x] **G2.6**–**G2.9** Pre-processing of rectangular input expected to be explicitly controlled through user specification of a log-likelihood function, so not applicable. **Missing or Undefined Values** - [ ] **G2.10** There are no checks for missing data - [ ] **G2.11** There are no options for users to specify how to handle missing (`NA`) data - [ ] **G2.12** Functions assume non-missingness, without explicitly stating such - [ ] **G2.13** There is no processing of other (non-`NA`) undefined values. #### 1.3 Output Structures - [x] **G4.0** Nothing is written to local files, so not applicable #### 1.4 Testing **Test Data Sets** - [x] **G5.0** No standard data sets used in tests (rather, all tests use random data simulated with fixed random seeds) - [x] **G5.1** Data sets are created within software, but not in a way that is appropriate for export, so not applicable. **Responses to Unexpected Input** - [ ] **G5.2** Tests do not cover all error and warning conditions - [x] **G5.2a** Every message produced `stop()`, `warning()`, `message()`, or equivalent is unique - [ ] **G5.2b** Tests do not cover every one of those messages - [ ] **G5.3** Absence of missing or undefined values in return objects not explicitly tested. **Algorithm Tests** - [x] **G5.4** Correctness tests are appropriately implemented - [x] **G5.5** Correctness tests are run with a fixed random seed - [x] **G5.6** Parameter recovery tests appropriately implemented. - [x] **G5.7** Algorithm performance tests appropriately implemented. - [x] **G5.8** Edge condition tests appropriately implemented. - [ ] **G5.9** No noise susceptibility tests implemented. **Extended tests** - [x] **G5.10**–**G5.12** No extended tests included, so not applicable ------------------------------------------------------------------------ ### 2 [Bayesian Standards](https://ropenscilabs.github.io/statistical-software-review-book/standards.html#bayesian-and-monte-carlo-software) #### 2.1 Documentation of Inputs - [ ] **BS1.0** Uses the term “parameter” to refer to “hyperparameters”, although does so entirely consistently - [x] **BS1.1** Extensive descriptions of how to enter data, primarily via vignettes. - [ ] **BS1.2** Descriptions of how to specify prior distributions provided - [ ] **BS1.2a** Not in `README` - [x] **BS1.2b** In vignette - [ ] **BS1.2c** Not in Function-level documentation - [x] **BS1.3** All parameters which control the computational process extensively described. - [x] **BS1.3a** Documents how to use the output of previous simulations as starting points of subsequent simulations. - [x] **BS1.3b** Does not document how to use different sampling algorithms because only implements one. - [x] **BS1.4** Explicitly describes and provides examples of convergence checkers. - [x] **BS1.5** Differences between convergence checkers are explicitly tested. #### 2.2 Input Data Structures and Validation **Input Data** - [x] **BS2.0** Vector inputs are appropriately pre-processed regardless of class structures. - [ ] **BS2.1** Some pre-processing routines are implemented for alternative columns - [x] **BS2.1a** Re-classed columns are processed via `as.data.frame` - [ ] **BS2.1b** List columns fail - [ ] **BS2.2** No *pre-checks* are implemented to ensure input data are commensurate, rather unhelpful errors are issued instead. **Prior Distributions, Model Specifications, and Hyperparameters** - [ ] **BS2.3** Hyperparameters are neither validated nor pre-processed prior to submitting to analytic routines, with errors passed through to multiple parallel computational chains (`burnin`, for example, is only used *after* chains have been calculated, to determine which parts of result are to be discarded). - [ ] **BS2.4** Lengths of hyperparameter vectors are not explicitly checked, rather passed through to attempted matrix multiplication, and then issuing unhelpful error message (“Error in coeffs \* x. : non-conformable arrays”). - [x] **BS2.5** Passing to “coeffs \* x.” described above ensures that lengths of hyperparameter vectors are commensurate with expected model input - [ ] **BS2.6** There are no pre-processing checks to validate appropriateness of numeric values submitted for hyperparameters; for example, hyperparameters defining second-order moments and given negative values are passed through, and lead to unhelpful errors (“Error in if (R\[i\] < klogratio) { : missing value where TRUE/FALSE needed”). **Computational Parameters** - [ ] **BS2.7** There are no checks that values for parameters are positive (see BS2.3 above; negative `burnin` values still lead to chains being calculated). - [ ] **BS2.8** There are no check on lengths of inputs, for example, passing a length = 2 vector to `burnin` generates warnings on various conditional checks expecting single values. - [ ] **BS2.9** There are no checks that arguments are of expected classes or types, rather passing unexpected classes or types generally fails at some stage with uninformative error messages. - [ ] **BS2.10** Parameters of inappropriate type are not automatically rejected, rather usually fail as described above (BS2.9). **Seed Parameters** - [ ] **BS2.11** Seeds are not able to be passed as a parameter, rather random generation is determined by the system `seed` value. - [x] **BS2.12** Results of previous runs can be used as starting points - [x] **BS2.13** Each chain is started with a different seed by default - [ ] **BS2.14** Seeds can not be passed, so no diagnostic messages can be issued when identical seeds are passed to distinct computational chains - [ ] **BS2.15** There is no advice *not* to use `set.seed()`, rather it is the suggested way to control random generation. - [x] **BS2.16** The relevant parameter is called “initial”, which require no pluralization. **Output Verbosity** - [ ] **BS2.17** There is no parameter controlling the verbosity of output, rather just a single `progress` parameter to switch the progress bar on and off, while convergence checkers have no verbosity control at all. - [x] **BS2.18** It is possible to suppress progress indicators while retaining verbosity of warnings and errors, through the `progress` parameter described above (**BS2.17**). - [ ] **BS2.19** It is not possible to suppress warnings, and many uninformative warnings may be triggered as described above. - [ ] **BS2.20** There is no explicit way to catch errors #### 2.3 Pre-processing and Data Transformation **Missing Values** - [ ] **BS3.0** No missing values are permitted, and no checks are performed, rather missing values are passed through unchecked, leading to unhelpful error messages. - [ ] **BS3.1** There are no routines to pre-process missing values prior to passing data through to main computational algorithms. **Perfect Collinearity** - [ ] **BS3.2** There are no pre-processing routines to diagnose perfect collinearity, and there are no diagnostic messages or warnings - [ ] **BS3.3** There are no distinct routines for processing perfectly collinear data. #### 2.4 Analytic Algorithms - [x] **BS4.0** Sampling algorithms are documented, including literary citation - [ ] **BS4.1** There are no comparisons with external samplers - [ ] **BS4.2** There are no methods to validate posterior estimates (other than convergence checkers used during simulation chains). - [x] **BS4.3** Several different kinds of convergence checkers are implemented and documented. - [x] **BS4.4** Computations are able to be stopped on convergence - [x] **BS4.5** Appropriate mechanisms are provided for models which do not converge, by stopping after specified numbers of iterations regardless of convergence. - [ ] **BS4.6** There are no tests to confirm that results with convergence checker are statistically equivalent to results from equivalent fixed number of samples without convergence checking. - [ ] **BS4.7** Effects of parameters of convergence checkers are not tested. #### 2.5 Return Values - [ ] **BS5.0** Return objects do not include information on seed(s) or starting value(s) - [ ] **BS5.1** Return objects do not include metadata on types (or classes) and dimensions of input data - [ ] **BS5.2** Software returns neither their input function of distributional parameters, nor enables subsequent access to such. - [ ] **BS5.3** Return object does not include convergence statistics - [ ] **BS5.4** Return object does not include details of convergence checker used - [ ] **BS5.5** Return object does not include diagnostic statistics to indicate absence of convergence #### 2.6 Additional Functionality - [ ] **BS6.0** Return object does not Implement a default `print` method - [x] **BS6.1** Return object does Implement a default `plot` method - [ ] **BS6.2** Software provides and documents abilities to plot sequences of posterior samples, although burn-in periods are simply excluded from plots. - [x] **BS6.3** There are straightforward abilities to plot posterior distributional estimates - [x] **BS6.4** Provide `summary` methods for return objects are provided - [x] **BS6.5** It is possible to plot both sequences of posterior samples and distributional estimates together in single graphic #### 2.7 Tests **Parameter Recovery Tests** - [x] **BS7.0** Tests for recovery of parametric estimates of a prior distribution - [x] **BS7.1** Tests for recovery of prior distribution in the absence of additional data - [x] **BS7.2** Tests for recovery of expected posterior distribution given a specified prior and input data **Algorithmic Scaling Tests** - [ ] **BS7.3** There are no tests for scaling of algorithmic efficiency with sizes of input data **Scaling of Input to Output Data** - [x] **BS7.4** Tests confirm that predicted or fitted values are on (approximately) the same scale as input values. ------------------------------------------------------------------------ [`rstan`](https://github.com/stan-dev/rstan) -------------------------------------------- ### 3 [General Standards](https://ropenscilabs.github.io/statistical-software-review-book/standards.html#general-standards-for-statistical-software) #### 3.1 Documentation - [x] **G1.0** Primary reference is provided. **Statistical Terminology** - [x] **G1.1** Statistical terminology clarified and unambiguously defined. **Function-level Documentation** - [ ] **G1.2** [`roxygen`](https://roxygen2.r-lib.org/) not used to documentation functions. - [ ] **G1.2a** Internal functions generally not documented in any form **Supplementary Documentation** - [x] **G1.3** Code is provided to reproduce results which form the basis of performance claims made in associated publications (for example, for `sbc()`). - [ ] **G1.4** Code is provided to compare performance claims with alternative implementations in other R packages (for example, for `sbc()`). #### 3.2 Input Structures **Uni-variate (Vector) Input** - [ ] **G2.0** There is no explicit documentation of expected lengths of inputs - [ ] **G2.1** There is no explicit documentation of expected *data types* of vector inputs - [ ] **G2.2** There is no explicit prohibition or restriction on submission of multivariate input to parameters expected to be univariate. - [ ] **G2.3** Some control of single-valued character input provides - [x] **G2.3a** `match.arg` is used - [ ] **G2.3b** `tolower` is not used, nor is their explicit mention of strict case sensitivity - [x] **G2.4** Provide appropriate mechanisms to convert between different *data types*, potentially including: - [x] **G2.4a** explicit conversion to `integer` via `as.integer()` - [x] **G2.4b** explicit conversion to continuous via `as.numeric()` - [x] **G2.4c** explicit conversion to character via `as.character()` (and not `paste` or `paste0`) - [x] **G2.4d**–**G2.4e** There is no explicit `factor` input, so not applicable - [x] **G2.5** There is no explicit `factor` input, so not applicable **Tabular Input** - [x] **G2.6** Restrictions on tabular input very carefully documented - [x] **G2.7** There are conversion routines for multivariate input, for reasons explained in documentation. - [x] **G2.8** There is no type conversion in which information is lost, so not applicable. - [x] **G2.9** Restrictions on tabular input very carefully documented. **Missing or Undefined Values** - [x] **G2.10**–**2.13** Missing data is explicitly prohibited, as extensively documented, and all appropriate checks are also implemented. #### 3.3 Output Structures - [ ] **G4.0** Processing of names of local files explicitly documented #### 3.4 Testing **Test Data Sets** - [x] **G5.0** Tests use standard data sets - [x] **G5.1** Data sets created to test package are not exported, but are made generally available through direct inclusion within package sub-directory. **Responses to Unexpected Input** - [ ] G4.2 Tests do not cover all error and warning conditions - [ ] G4.3 Absence of missing or undefined values in return objects not explicitly tested. **Algorithm Tests** - [x] **G5.3** Correctness tests appropriately implemented. - [x] **G5.5** Correctness tests are run with a fixed random seed - [x] **G5.6** Parameter recovery tests appropriately implemented. - [x] **G5.7** Algorithm performance tests appropriately implemented. - [x] **G5.8** Edge condition tests appropriately implemented. - [ ] **G5.9** No noise susceptibility tests implemented **Extended tests** - [x] G4.10 Extended tests appropriately implemented (notably for underlying C++ code). - [x] G4.11 Extended tests require no large data sets, so not applicable - [ ] G4.12 Conditions necessary to run extended tests not appropriately documented ------------------------------------------------------------------------ ### 4 [Bayesian Standards](https://ropenscilabs.github.io/statistical-software-review-book/standards.html#bayesian-and-monte-carlo-software) #### 4.1 Documentation of Inputs - [ ] **BS1.0** There terms “parameter” and “hyperparameter” are used interchangeably and ambiguously in much documentation (for example, in [*Prior Choice Recommendations*](https://github.com/stan-dev/stan/wiki/Prior-Choice-Recommendations)). - [x] **BS1.1** There is extensive description of how to enter data - [x] **BS1.2** There is extensive Description of how to specify prior distributions - [ ] **BS1.2a** Not in the main package `README`, although that seems okay here, as the stan project is too big for that to serve a simple demonstration purpose. - [x] **BS1.2b** In at least one package vignette - [x] **BS1.2c** In function-level documentation, preferably with code included in examples - [x] **BS1.3** All parameters which control the computational process are extensively described - [x] **BS1.3a** There is sufficient documentation of how to use the output of previous simulations as starting points of subsequent simulations. - [x] **BS1.3b** There is sufficient documentation of how to use different sampling algorithms. - [x] **BS1.4** `rstan` does not directly implement convergence checkers, rather provides convergence statistics throughout, and exports several well-documented helper functions to aid subsequent analyses of such. - [x] **BS1.5** The documentation for the [`optimizing()` function](https://mc-stan.org/rstan/reference/stanmodel-method-optimizing.html) describes several ways of controlling convergence to point estimates. #### 4.2 Input Data Structures and Validation **Input Data** - [ ] **BS2.0** one-dimensional input is not appropriately pre-processed, with submission of `units` vectors erroring (“both operands of the expression should be ‘units’ objects”). - [x] **BS2.1** `rstan` does not convert all plausible forms of two-dimension input to standard forms, but does document very explicitly the expected forms for such, and strongly discourages submitting non-standard formats or classes. - [x] **BS2.1a** `data.frame` or equivalent objects which have columns with non-standard class attributes are simply coerced to numeric, as clearly documented. - [x] **BS2.1b** Submitting list columns is explicitly and strongly discouraged in the [documentation](https://mc-stan.org/rstan/reference/stan.html#passing-data-to-stan), although possible as described there. - [x] **BS2.2** Pre-processing routines ensure all input data is dimensionally commensurate **Prior Distributions, Model Specifications, and Hyperparameters** - [x] **BS2.3** Validation and pre-processing of hyperparameters is expected to be part of user-provide model specifications, as clearly documented. - [x] **BS2.4**–**BS2.6** Hyperparameter vectors are not submitted to `rstan` like most other Bayesian software, and so these standards do not apply. **Computational Parameters** - [ ] **BS2.7** There are some checks that values for parameters are positive where expected, but these are only implemented within each chain, not as single pre-processing step. - [x] **BS2.8** Lengths of inputs are checked and error appropriately - [x] **BS2.9** Arguments of unexpected classes or types are appropriately pre-converted - [ ] **BS2.10** Parameters of inappropriate type are not automatically rejected; rather are submitted to sub-routines, triggering unhelpful errors. **Seed Parameters** - [x] **BS2.11** Seeds can be passed as a parameter - [x] **BS2.12** Results of previous runs can be used as starting points for subsequent runs - [x] **BS2.13** Each chain is always started with a different seed by default - [x] **BS2.14** It is not possible to pass identical seeds to distinct computational chains, so standard not applicable - [x] **BS2.15** Advice *not* to use `set.seed()` is explicitly documented (and has no effect). - [x] **BS2.16** The parameter controlling starting values is called `init`, so appropriately not pluralised. **Output Verbosity** - [x] **BS2.17** There are multiple ways of controlling the verbosity of output - [x] **BS2.18** It is possible to suppress messages and progress indicators, while retaining verbosity of warnings and errors. - [x] **BS2.19** It is possible to suppress warnings - [ ] **BS2.20** No errors are able to be either caught or converted to warnings #### 4.3 Pre-processing and Data Transformation **Missing Values** - [x] **BS3.0** Documentation and error message explicitly state that no missing values are permitted - [x] **BS3.1** No missing values are permitted, no pre-processing standards not applicable **Perfect Collinearity** - [ ] **BS3.2** There are no pre-processing routines to diagnose perfect collinearity, no are diagnostic messages or warnings issued - [ ] **BS3.3** There are no distinct routines for processing perfectly collinear data, rather all data are passed directly to sampling algorithms #### 4.4 Analytic Algorithms - [x] **BS4.0** Sampling algorithms are explicitly documented and referenced. - [ ] **BS4.1** There are no explicit comparisons with external samplers - [x] **BS4.2** A method to validate posterior estimates is implemented (via the [`sbc()` function](https://mc-stan.org/rstan/reference/sbc.html)). - [x] **BS4.3** No convergence checker as such implemented, but extensive convergence statistics are provided, documented, and able to be extracted. - [x] **BS4.4** Computations are able to be stopped on convergence (via appropriate post-processing functions such as [`optimizing()`](https://mc-stan.org/rstan/reference/stanmodel-method-optimizing.html)). - [x] **BS4.5** Appropriate mechanisms are provided for models which do not converge. - [ ] **BS4.6** Convergence statistics are not explicitly tested within the R package (even though the more general `stan` code is checked very comprehensively). - [ ] **BS4.7** Effects of parameters passed to convergence checkers are not tested. #### 4.5 Return Values - [x] **BS5.0** Return objects include both seeds and starting values (able to be extracted via `get_seed()` and `get_inits()` functions). - [ ] **BS5.1** Return objects DO NOT include appropriate metadata on types or dimensions of input data - [x] **BS5.2** `rstan` enables direct access to prior specification through the `get_stancode()` and `get_stanmod()` functions. - [x] **BS5.3** Convergence statistics or equivalent are returned or able to be extracted. - [x] **BS5.4** Where multiple checkers are enabled, return details of convergence checker used - [x] **BS5.5** Appropriate diagnostic statistics are returned (or able to be accessed) to indicate absence of convergence #### 4.6 Additional Functionality - [x] **BS6.0** Default `print` method for return objects is implemented - [x] **BS6.1** Default `plot` method for return objects is implemented - [x] **BS6.2** Sequences of posterior samples are able to be plotted, although burn-in periods are omitted from plot - [x] **BS6.3** Posterior distributional estimates able to be directed plotted - [x] **BS6.4** `summary` methods are provided for return objects - [ ] **BS6.5** It is not possible to plot both sequences of posterior samples and distributional estimates together in single graphic #### 4.7 Tests **Parameter Recovery Tests** - [x] **BS7.0** There are tests for recovery of parametric estimates of prior distribution - [x] **BS7.1** There are tests for recovery of prior distribution in the absence of additional data - [x] **BS7.2** There are tests for recovery of posterior distribution given specified prior and input data **Algorithmic Scaling Tests** - [ ] **BS7.3** There are no tests for scaling of algorithmic efficiency with sizes of input data **Scaling of Input to Output Data** - [x] **BS7.4** There are tests that predicted or fitted values are on the same scale as input values.

Syntax	Example	Reference
# Header	Header	基本排版
- Unordered List	Unordered List
1. Ordered List	Ordered List
- [ ] Todo List	Todo List
> Blockquote	Blockquote
Bold font	Bold font
Italics font	Italics font
~~Strikethrough~~	~~Strikethrough~~
19^th^	19^th
H~2~O	H₂O
++Inserted text++	Inserted text
==Marked text==	Marked text
[link text](https:// "title")	Link
![image alt](https:// "title")	Image
`Code`	`Code`	在筆記中貼入程式碼
```javascript var i = 0; ```	`var i = 0;`
:smile:		Emoji list
{%youtube youtube_id %}	Externals
$L^aT_eX$	L^aT_eX
:::info This is a alert area. :::	This is a alert area.