paul-buerkner / brms

brms R package for Bayesian generalized multivariate non-linear multilevel models using Stan
https://paul-buerkner.github.io/brms/
GNU General Public License v2.0
1.25k stars 177 forks source link

Fix p limits for some q* distribution functions; fix qgen_extreme_value #1639

Closed venpopov closed 3 months ago

venpopov commented 3 months ago

Closes #1637

This allows quantiles of 0 and 1 to be used for qstudent_t, qskew_normal, qfrechet, qgen_extreme_value and qasym_alplace.

Check if you like it the way I did. In addition, rather than throwing an error when p <0 or >1, the functions produce NaN like all other q* type distribution functions. It probably rarely matters, but I think the consistency is nice. The check itself is still necessary, because some distributions like qfrechet would produce non-meaningful but real values for p > 1.

I added a couple of new tests to check the indended results with p at or beyond the limits

Incorrect output of qgen_extreme_value

While working on the above, I noticed that qgen_extreme_value produces incorrect results when xi != 0. There was no test for qgen_extreme_value, so I initially added the following test (similar to what was there for qfrechet):

  n <- 10
  x <- rgamma(n, 10, 3)
  mu <- rnorm(n)
  res <- pgen_extreme_value(x, mu = mu, sigma = 1:n, xi = 3)
  q <- qgen_extreme_value(res, mu = mu, sigma = 1:n, xi = 3)
  expect_equal(x, q)
#> Error: `x` not equal to `q`.
#> 10/10 mismatches (average diff: 2.06)
#> [1] 3.385 -  0.296 == 3.090
#> [2] 2.194 - -0.797 == 2.990
#> [3] 2.563 -  2.489 == 0.074
#> [4] 2.340 -  1.075 == 1.265
#> [5] 3.846 - -0.274 == 4.119
#> [6] 0.751 -  0.243 == 0.508
#> [7] 2.213 -  1.127 == 1.086
#> [8] 4.624 -  2.251 == 2.373
#> [9] 4.628 -  2.405 == 2.223
#> ...

This line was responsible, so I fixed it and now the new test passes. Let me know if you prefer to have this fix in a separate PR.

paul-buerkner commented 3 months ago

Thank you! Looks good! Merging now.