Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Extrapolation options #356

Open
wants to merge 14 commits into
base: master
Choose a base branch
from

Conversation

SouthEndMusic
Copy link
Member

@SouthEndMusic SouthEndMusic commented Nov 13, 2024

Fixes #355.

This is what I have in mind:

using DataInterpolations
using Plots
using Random

Random.seed!(2)

t = cumsum(rand(5))
u = rand(5)

t_eval = range(-1, 3.5, length = 100)
p = scatter(t,u; label = "data")
title!("Extrapolation types")

for extrapolation_type in [:constant, :linear, :extension]
    A = QuadraticSpline(u,t; extrapolation_up = extrapolation_type, extrapolation_down = extrapolation_type)
    plot!(t_eval, A.(t_eval); label = string(extrapolation_type))
end

p

figure

Let me know whether this is a good approach before I apply it to all interpolation types, integrals, derivatives.

@SouthEndMusic SouthEndMusic marked this pull request as draft November 13, 2024 17:36
@SouthEndMusic SouthEndMusic changed the title QuadraticSpline POC Extrapolation options Nov 13, 2024
@ChrisRackauckas
Copy link
Member

This looks like a great idea. It should be an EnumX enum though instead of just symbols in order to autocomplete and spellcheck properly.

@SouthEndMusic
Copy link
Member Author

SouthEndMusic commented Nov 14, 2024

Derivative extrapolation now also looking good:

using DataInterpolations
using CairoMakie
using Random

Random.seed!(2)

t = cumsum(rand(5))
u = rand(5)

t_eval = range(-1, 3.5, length = 1000)
fig = Figure()
ax_itp = Axis(fig[1,1]; title = "(inter/extra)polation")
ax_deriv = Axis(fig[2,1]; title = "(inter/extra)polation first derivative")
ax_deriv2 = Axis(fig[3,1]; title = "(inter/extra)polation second derivative")
scatter!(ax_itp, t,u; label = "data")

for extrapolation_type in [ExtrapolationType.constant, ExtrapolationType.linear, ExtrapolationType.extension]
    A = QuadraticSpline(u,t; extrapolation_up = extrapolation_type, extrapolation_down = extrapolation_type)
    lines!(ax_itp, t_eval, A.(t_eval); label = string(extrapolation_type))
    lines!(ax_deriv, t_eval, DataInterpolations.derivative.(Ref(A), t_eval); label = string(extrapolation_type))
    lines!(ax_deriv2, t_eval, DataInterpolations.derivative.(Ref(A), t_eval, 2); label = string(extrapolation_type))
end

axislegend(ax_itp)
axislegend(ax_deriv)
axislegend(ax_deriv2)
fig

figure

@SouthEndMusic
Copy link
Member Author

And the integral:

using DataInterpolations
using CairoMakie
using Random

Random.seed!(2)

t = cumsum(rand(5))
u = rand(5)

t_eval = range(-2, 5.5, length = 1000)
fig = Figure()
ax_itp = Axis(fig[1,1]; title = "(inter/extra)polation")
ax_antideriv = Axis(fig[2,1]; title = "(inter/extra)polation antiderivative")
scatter!(ax_itp, t,u; label = "data")

for extrapolation_type in [ExtrapolationType.constant, ExtrapolationType.linear, ExtrapolationType.extension]
    A = QuadraticSpline(u,t; extrapolation_up = extrapolation_type, extrapolation_down = extrapolation_type)
    lines!(ax_itp, t_eval, A.(t_eval); label = string(extrapolation_type))
    lines!(ax_antideriv, t_eval, DataInterpolations.integral.(Ref(A), t_eval); label = string(extrapolation_type))
end

axislegend(ax_itp)
axislegend(ax_antideriv)
fig

figure

@SouthEndMusic
Copy link
Member Author

SouthEndMusic commented Nov 16, 2024

I've found myself in a refactoring rabbit hole again 😅

  • I'm refactoring integration because I found the previous implementation, where _integral gives the antiderivative unintuitive. With my refactor _integral integrates over an interval that lies between 2 consecutive data points, which also allows for some simplification / reduction of operations
  • I also want to refactor QuadraticInterpolation because the way we compute it now is unnecessarily computationally intensive and cumbersome. I wonder how many people actually use that one anyway

@sathvikbhagavan
Copy link
Member

Can you do the refactors in different PRs?

@SouthEndMusic
Copy link
Member Author

Not really, or I have to do them before this one. Is it fine by you if I wrap up what I'm doing here and make some follow-up issues?

@sathvikbhagavan
Copy link
Member

Thats alright, refactors that are not strictly needed in this one can have follow up issues/PRs.

@SouthEndMusic SouthEndMusic marked this pull request as ready for review November 17, 2024 16:01
@SouthEndMusic
Copy link
Member Author

@sathvikbhagavan sorry for the spam but just to let you know that this PR is ready for review

Copy link
Member

@sathvikbhagavan sathvikbhagavan left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A few comments

docs/src/extrapolation_methods.md Outdated Show resolved Hide resolved
docs/src/extrapolation_methods.md Show resolved Hide resolved
src/DataInterpolations.jl Outdated Show resolved Hide resolved
src/DataInterpolations.jl Outdated Show resolved Hide resolved
src/interpolation_caches.jl Show resolved Hide resolved
src/interpolation_caches.jl Show resolved Hide resolved
Copy link
Member

@sathvikbhagavan sathvikbhagavan left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Almost there. A few minor comments.

return if order == 1
_derivative(A, t, iguess)
elseif order == 2
function _extrapolate_derivative_down(A, t, order)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
function _extrapolate_derivative_down(A, t, order)
function _extrapolate_derivative_left(A, t, order)

end
end

function _extrapolate_derivative_up(A, t, order)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
function _extrapolate_derivative_up(A, t, order)
function _extrapolate_derivative_right(A, t, order)

function _extrapolate_derivative_down(A, t, order)
(; extrapolation_left) = A
typed_zero = zero(first(A.u) / one(A.t[1]))
if extrapolation_left == ExtrapolationType.none
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Slightly knit picky, but I think checking this one line above would prevent a few wasteful computation of typed_zero

return total
end

function _extrapolate_integral_down(A, t)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
function _extrapolate_integral_down(A, t)
function _extrapolate_integral_left(A, t)

end
end

function _extrapolate_integral_up(A, t)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
function _extrapolate_integral_up(A, t)
function _extrapolate_integral_right(A, t)

end
end

function _extrapolate_down(A, t)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
function _extrapolate_down(A, t)
function _extrapolate_left(A, t)

end
end

function _extrapolate_up(A, t)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
function _extrapolate_up(A, t)
function _extrapolate_right(A, t)


function _extrapolate_down(A, t)
(; extrapolation_left) = A
if extrapolation_left == ExtrapolationType.none
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same as derivative/integral comment

test_constant_extrapolation(LinearInterpolation, u, t)

for extrapolation_type in [ExtrapolationType.linear, ExtrapolationType.extension]
# Down extrapolation
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
# Down extrapolation
# Left extrapolation

@test DataInterpolations.integral(A, t_eval) == -0.5
t_eval = 3.0

# Up extrapolation
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
# Up extrapolation
# Right extrapolation

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

Configurable extrapolation behavior
3 participants