Gamma creation

[TTitscher]

As discussed by @aradermacher and @ajafarihub, a reasonable definition of the (prior for the) gamma noise distribution is based on two quantiles. This is the main content of this PR, e.g.

gamma = bayes.vb.Gamma.FromQuantiles(100, 400) # the quantiles itself default to 0.05 and 0.95

In the context of VB, the gamma distribution is used to model the noise precision. In practice, it can be convenient to generate this distribution from the standard deviation (SD). So the line above is identical to

gamma = bayes.vb.Gamma.FromSDQuantiles(0.1, 0.2) # as 100 = 1/0.1**2 and 400 = 1/0.2**2

IMO, this is a bit sketchy as those not actually the quantiles of an inverted squared gamma distribution...

Also there is VBResult.summary:

# result.summary(tablefmt="fancy_grid"):
╒════════╤══════════╤═════════╤════════════╤══════════╤══════════╤═══════════╤═══════════╕
│ name   │   median │    mean │         sd │       5% │      25% │       75% │       95% │
╞════════╪══════════╪═════════╪════════════╪══════════╪══════════╪═══════════╪═══════════╡
│ p0     │   6.9813 │  6.9813 │ 0.0109935  │  6.96322 │  6.97389 │   6.98872 │   6.99939 │
├────────┼──────────┼─────────┼────────────┼──────────┼──────────┼───────────┼───────────┤
│ p1     │  42.0053 │ 42.0053 │ 0.00634869 │ 41.9949  │ 42.0011  │  42.0096  │  42.0158  │
├────────┼──────────┼─────────┼────────────┼──────────┼──────────┼───────────┼───────────┤
│ noise  │  99.0231 │ 99.0891 │ 4.43096    │ 91.9152  │ 96.0656  │ 102.041   │ 106.488   │
╘════════╧══════════╧═════════╧════════════╧══════════╧══════════╧═══════════╧═══════════╛
# or result.summary(gamma_as_sd=True):
name    median                  mean  sd                           5%        25%        75%         95%
------  -----------------  ---------  --------------------  ---------  ---------  ---------  ----------
p0      6.981304431696486   6.9813    0.010993495726152509   6.96322    6.97389    6.98872    6.99939
p1      42.00534689378616  42.0053    0.006348692991699628  41.9949    42.0011    42.0096    42.0158
noise   ?                   0.100459  ?                      0.104305   0.102027   0.098995   0.0969057

[aradermacher]

Thanks for integrating this method here, Thomas. I would prefer to use shape instead of alpha avoiding to mixing these things up again ;) In the resukt summary sd is based on the mean precision or is that the standard deviation of the gamma distribution??

[TTitscher]

Thanks for your comments! Indeed, "sd" in the context of "gamma for the precision" can be confusing. In the summary, it is really 'just' the standard deviation of the gamma distribution, so the sd of the precision.

If you provide summary(gamma_as_sd=True), the quantiles are transformed by 1/x**0.5, such that the mean then refers to the mean of the standard deviation distribution. But mathematically, this may only be correct for the mean, not for the other quantiles.

Do we find a better way do deal with this separation of:

• mathematically well-defined Gamma distribution
• use case of using Gamma to describe the precision of a zero-mean noise
  • where the standard deviation instead of precision is almost always more practical

?