Skip to content

Fallback for time_change #89

@JoseKling

Description

@JoseKling

The implementation of time_change for inhomogeneous Poisson processes should be fallback for AbstractPointProcess. Fallback for integrated_ground_intensity can also de defined from ground_intensity.

Something along these lines

function integrated_ground_intensity(pp::AbstractPointProcess, h::History, a, b)
    func = (t, p) -> ground_intensity(pp, t, h) 
    prob = IntegralProblem(func, a, b);
    sol = solve(prob, QuadGKJL(); abstol = 1e-14, reltol = 1e-14);
    return sol.u;
end

function time_change(h::History, pp::AbstractPointProcess)
    edges = vcat(h.tmin, h.times, h.tmax)
    deltas = [integrated_ground_intensity(pp, h, edges[i], edges[i + 1]) for i in 1:(length(edges) - 1)]
    T = eltype(deltas)
    i_neg = findfirst(<(zero(T)), deltas)
    if i_neg !== nothing
        throw(
            DomainError(
                deltas[i_neg],
                """
                time_change: integrated intensity over [$(edges[i_neg]), $(edges[i_neg + 1])] is negative. This could mean:
                  1. Intensity function is negative.
                  2. Incorrect implementation of integrated_ground_intensity.
                  3. Floating-point jitter from adaptive quadrature on an interval where the true integral is near zero.
                """
            ),
        )
    end
    cs = cumsum(deltas)
    new_times = cs[1:(end - 1)]
    new_tmax = cs[end]
    return History(; times=new_times, tmin=zero(T), tmax=new_tmax, marks=h.marks)
end

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions