Skip to content

Commit

Permalink
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Improve numerical stability of Cholesky invlink
Browse files Browse the repository at this point in the history
penelopeysm committed Nov 30, 2024
1 parent f52a9c5 commit 9d11ba4
Showing 2 changed files with 4 additions and 5 deletions.
2 changes: 1 addition & 1 deletion Project.toml
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
name = "Bijectors"
uuid = "76274a88-744f-5084-9051-94815aaf08c4"
version = "0.15.2"
version = "0.15.3"

[deps]
ArgCheck = "dce04be8-c92d-5529-be00-80e4d2c0e197"
7 changes: 3 additions & 4 deletions src/bijectors/corr.jl
Original file line number Diff line number Diff line change
@@ -348,13 +348,12 @@ function _inv_link_chol_lkj(Y::AbstractMatrix)
T = float(eltype(W))
logJ = zero(T)

idx = 1
@inbounds for j in 1:K
log_remainder = zero(T) # log of proportion of unit vector remaining
for i in 1:(j - 1)
z = tanh(Y[i, j])
W[i, j] = z * exp(log_remainder)
log_remainder += log1p(-z^2) / 2
log_remainder += log(2 / (exp(Y[i, j]) + exp(-Y[i, j])))
logJ += log_remainder
end
logJ += log_remainder
@@ -380,10 +379,10 @@ function _inv_link_chol_lkj(y::AbstractVector)
log_remainder = zero(T) # log of proportion of unit vector remaining
for i in 1:(j - 1)
z = tanh(y[idx])
idx += 1
W[i, j] = z * exp(log_remainder)
log_remainder += log1p(-z^2) / 2
log_remainder += log(2 / (exp(y[idx]) + exp(-y[idx])))
logJ += log_remainder
idx += 1
end
logJ += log_remainder
W[j, j] = exp(log_remainder)

0 comments on commit 9d11ba4

Please sign in to comment.