`logDensity` returns `0.0` for a beta distribution with shape parameters `1e200` at `0.5`

Đang mở
#237 0 bình luận 0 reaction 0 người được giao Xem trên GitHub

Chưa có ai nhận issue này.

Đánh giá

Độ khó
4/5
Thời gian dự kiến
3-5 ngày
Mức phù hợp với người mới
48/100
Loại issue
Lỗi
Độ rõ ràng
Đặc tả rõ ràng
Mức độ hoạt động
Sôi nổi
Công nghệ
haskell
Lĩnh vực
data

Hướng nghiên cứu

Start by locating logDensity for the beta distribution and the logBeta implementation used by math-functions; the report identifies cancellation between these calculations for shape parameters 1e200 at 0.5. Check existing beta distribution tests and add coverage for this reproduction. Done means the computed value is close to 230.37929153703981 rather than 0.0.

Do mô hình lập chỉ mục viết ra từ nội dung của issue.

Mô tả

Hello,

Unfortunately, as said in the title: logDensity returns 0.0 for a beta distribution with both shape parameters equal to 1e200, evaluated at 0.5. To reproduce it, save the following as Main.hs and run runghc Main.hs with statistics-0.16.5.0 exposed:

import qualified Statistics.Distribution as D
import Statistics.Distribution.Beta (betaDistr)

main :: IO ()
main = print (D.logDensity (betaDistr 1e200 1e200) 0.5)

It prints:

0.0

I get this on both GHC 9.6.7 and 9.14.1, with statistics-0.16.5.0 and math-functions-0.3.4.4.

The value should be approximately 230.37929153703981. Here is how to check it with mpmath (700 decimal digits of precision, starting from the exact binary floating-point values passed for the shape parameters and x):

from mpmath import mp

mp.dps = 700
a = mp.mpf(1e200)
x = mp.mpf(0.5)
log_beta = 2 * mp.loggamma(a) - mp.loggamma(2 * a)
log_density = 2 * (a - 1) * mp.log(x) - log_beta
print(mp.nstr(log_density, 17))  # 230.37929153703981

Both shape parameters are positive, so the input is inside the documented domain. The problem comes from the formula (a-1)*log x + (b-1)*log1p(-x) - logBeta a b used in logDensity. Here a = b = 1e200 and x = 0.5. So the (a-1) and (b-1) terms sum to 2*(1e200-1)*log 0.5 = -1.3862943611e200, and logBeta 1e200 1e200 (from the allStirling path in math-functions) is -1.3862943611e200 too. Both are negative, but their true difference is about 230 (the log density). In Double they cancel out, so the subtraction is 0.0 and can't recover a single digit. Would it be possible to fix this? (combining the (a-1)*log x and logBeta terms analytically before evaluating would do it, so they don't cancel out)

Thanks in advance!

Ngôn ngữ chính
Haskell
Star
309
Fork
77
Chỉ số merge pull request
Không có pull request nào được merge trong 30 ngày

Chuẩn bị môi trường

Dự án này không cung cấp dev container, Dockerfile hay hướng dẫn đóng góp, nên bạn cần tự thiết lập môi trường: hãy bắt đầu từ README và xem hướng dẫn đóng góp lần đầu của chúng tôi để biết các bước chung.

Bắt đầu từ đâu

  1. Đọc hết issue, rồi đọc hướng dẫn đóng góp của dự án.
  2. Bình luận trên issue rằng bạn sẽ nhận — tránh hai người làm cùng một việc.
  3. Fork repository và làm thay đổi trên một nhánh.
  4. Mở pull request có tham chiếu số hiệu của issue.

Issue khác của haskell/statistics

Tất cả issue của haskell/statistics

Issue tương tự

Thêm issue về Haskell

Nhận issue mới trong hộp thư của bạn

Bản tóm tắt ngắn những issue GitHub phù hợp với người mới.