Skip to content

Commit

Permalink
Remove type dependenc in uniformize() fn. Add two new functions for t…
Browse files Browse the repository at this point in the history
…urning regular matrices into regular matrices that fulfill the requirements of transition rate matrices.
  • Loading branch information
smith-garrett committed Feb 21, 2023
1 parent a6e6a18 commit 5b0040a
Show file tree
Hide file tree
Showing 2 changed files with 27 additions and 2 deletions.
1 change: 1 addition & 0 deletions src/Uniformization.jl
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@ using Distributions: Poisson, pdf, cdf, quantile

include("algorithms.jl")
export TransitionRateMatrix
export transitionratematrix, transitionratematrix!
export uniformize, standard_uniformization, discrete_observation_times, erlangization

end
28 changes: 26 additions & 2 deletions src/algorithms.jl
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,30 @@ Base.getindex(m::AbstractRateMatrix, I...) = Base.getindex(m.matrix, I...)
#Base.IndexStyle(::Type{AbstractRateMatrix{E}}) where {E,T} = IndexStyle(T)


"""
transitionratematrix(Q)
Create a copy of the matrix Q that fulfills the requirements of a transition rate matrix.
"""
function transitionratematrix(Q)
@assert issquare(Q) "Matrix must be square."
@assert off_diag_nonnegative(Q) "Matrix contains off-diagonal elements that aren't positive."
return setdiagonal!(copy(Q))
end


"""
transitionratematrix!(Q)
Convert the matrix Q to one that fulfills the requirements of a transition rate matrix.
"""
function transitionratematrix!(Q)
@assert issquare(Q) "Matrix must be square."
@assert off_diag_nonnegative(Q) "Matrix contains off-diagonal elements that aren't positive."
return setdiagonal!(Q)
end


"""
make_dtmc(Q, λ=2^10)
Expand Down Expand Up @@ -65,7 +89,7 @@ end


"""
uniformize(Q::TransitionRateMatrix, p0, k=2^10, t=0.0,
uniformize(Q, p0, k=2^10, t=0.0,
method::Function=erlangization, args...)
Approximate 𝐩(t) = exp(t𝐐)𝐩(0) using uniformization. The parameter k controls the rate of
Expand All @@ -75,7 +99,7 @@ Returns a (normalized) distribution over the states at time 𝑡.
Uses Erlangization/external uniformization by default because it seems to be the most robust
with stiff problems.
"""
function uniformize(Q::TransitionRateMatrix, p0, k=2^10, t=0.0,
function uniformize(Q, p0, k=2^10, t=0.0,
method::Function=erlangization, args...)
@assert t zero(t) "Time t must be positive."
@assert size(p0, 1) == size(Q, 1) "Initial condition p0 must be the same size as Q."
Expand Down

0 comments on commit 5b0040a

Please sign in to comment.