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

A different CARMA kernel implementation #90

Merged
merged 44 commits into from
Jan 4, 2024
Merged

A different CARMA kernel implementation #90

merged 44 commits into from
Jan 4, 2024

Conversation

ywx649999311
Copy link
Contributor

Hi @dfm,

I added my implementation of the CARMA kernel (I named it lower case carma), which basically follows the suggestion in your celerite paper and breaks a CARMA model into a combination of real exponential kernels and celerite kernels. I try to follow the notation in Kelly+14 CARMA paper but absorbed its sigma term into the beta parameters.

I have also added a bunch of tests to check the CARMA and carma kernels using the equivalent exponential and celerite kernels in test_quasisep.py. In test_kalman.py, I added a few more CARMA specifications that will break your tests there.

To note that the current implementation will only work well for CARMA(2,1) or lower. To fit any higher order CARMA models, we need a separate constructor that will sample in the root space rather than the CARMA space (see discussion in section 3.3 in the Kelly+14 paper). I am working on migrating my EzTao code here and will keep developing on this branch.

src/tinygp/kernels/quasisep.py Outdated Show resolved Hide resolved
src/tinygp/kernels/quasisep.py Outdated Show resolved Hide resolved
@dfm
Copy link
Owner

dfm commented May 2, 2022

Thanks for this! I added a few inline comments to get us started.

It would also be great to add some tests that demonstrate the issues from #89, that get fixed by this PR.

To note that the current implementation will only work well for CARMA(2,1) or lower. To fit any higher order CARMA models, we need a separate constructor that will sample in the root space rather than the CARMA space (see discussion in section 3.3 in the Kelly+14 paper). I am working on migrating my EzTao code here and will keep developing on this branch.

Is there any reason to not just always operate in that parameterization? A limit of 2,1 seems pretty restrictive!

@ywx649999311
Copy link
Contributor Author

It would also be great to add some tests that demonstrate the issues from #89, that get fixed by this PR.

Good point. I will work on that.

Is there any reason to not just always operate in that parameterization? A limit of 2,1 seems pretty restrictive!

That is because not any combination of CARMA parameters returns a valid model, that is stationary and in the so-called 'minimum phase'. An -inf likelihood will (likely) be returned if the input model is not valid. This also relates to the requirement for the log input of alphas and betas. After all, we need the roots of the autoregressive and moving-average polynomials (e.g., Equation 30 of the Kelly+14 paper) to have negative real parts in order to generate a 'valid' model.

@dfm
Copy link
Owner

dfm commented May 3, 2022

That is because not any combination of CARMA parameters returns a valid model, that is stationary and in the so-called 'minimum phase'. An -inf likelihood will (likely) be returned if the input model is not valid. This also relates to the requirement for the log input of alphas and betas. After all, we need the roots of the autoregressive and moving-average polynomials (e.g., Equation 30 of the Kelly+14 paper) to have negative real parts in order to generate a 'valid' model.

The approach I would take to this is to clearly document (probably by adding a tutorial, as well as clear docstrings) different approaches for parameterizing valid models.

@dfm
Copy link
Owner

dfm commented May 18, 2022

This is on my radar, but I've been totally swamped recently. Thanks for your patience!

@ywx649999311
Copy link
Contributor Author

This is on my radar, but I've been totally swamped recently. Thanks for your patience!

No problem, me too. I am not done making changes yet, will bug you once I am ready. Thanks!

@ywx649999311
Copy link
Contributor Author

Hi, please give a look when you get time. I removed the old implementation and tests that I created for it. I think there might be a need to make a tutorial to explain the two different way to parameterize the model, and I can work on it once this is merged.

@dfm
Copy link
Owner

dfm commented Oct 30, 2022

@ywx649999311 — I'm so sorry about taking so long to get back to this!! I've just fixed the merge conflicts and fixed some typos so you should pull before making any more changes!

@dfm
Copy link
Owner

dfm commented Oct 30, 2022

You might also consider rebasing onto main, but that's definitely not required.

Scratch that - I did the rebase. Sorry if it's rude to force push directly to your branch! I'll try not to do it again...

More comments incoming!

Copy link
Owner

@dfm dfm left a comment

Choose a reason for hiding this comment

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

I've left a few more small comments, but I think this is looking great!! Thanks again for submitting this, and thanks for your patience.

src/tinygp/kernels/quasisep.py Outdated Show resolved Hide resolved
src/tinygp/kernels/quasisep.py Outdated Show resolved Hide resolved
The roots can be re-parameterized as the coefficients of a product
of quadratic equations each with the second-order term set to 1. The
input for this constructor are said coefficients. See Equation 30 in
the paper linked above for a reference.
Copy link
Owner

Choose a reason for hiding this comment

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

It's probably better to directly link to the paper here, instead of referencing to above.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

OK, fixed it.

src/tinygp/kernels/quasisep.py Outdated Show resolved Hide resolved
@dfm
Copy link
Owner

dfm commented Oct 31, 2022

@ywx649999311 — Besides my small comments above, I'm pretty much ready to merge this. I'd love it if you could add a little more detail to all the docstrings just because I'm not super familiar with all the details of this model and I want to make sure I can maintain it :D

@dhuppenkothen
Copy link
Contributor

Just wanted to say I've been using this branch to fit a CARMA(2,1) model, and it seems more stable than the current CARMA implementation in TinyGP (where I get a lot of NaNs and my attempts at sampling a posterior have had mixed results).

I would love to see this extended to higher-order CARMA processes and would be happy to help in whatever capacity I'm capable of doing so given my limited knowledge of JAX.

pre-commit-ci bot and others added 8 commits November 6, 2023 19:49
updates:
- [github.com/pre-commit/pre-commit-hooks: v4.4.0 → v4.5.0](pre-commit/pre-commit-hooks@v4.4.0...v4.5.0)
- [github.com/psf/black: 23.9.1 → 23.10.1](psf/black@23.9.1...23.10.1)
- [github.com/astral-sh/ruff-pre-commit: v0.0.292 → v0.1.4](astral-sh/ruff-pre-commit@v0.0.292...v0.1.4)

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
@ywx649999311
Copy link
Contributor Author

ywx649999311 commented Nov 8, 2023

Hi @dfm, I added more comments/docstring. Would you please see if you have further comments? Thanks!

@ywx649999311 ywx649999311 requested a review from dfm November 8, 2023 17:16
dependabot bot and others added 8 commits December 3, 2023 13:10
Bumps [pypa/gh-action-pypi-publish](https://github.com/pypa/gh-action-pypi-publish) from 1.8.10 to 1.8.11.
- [Release notes](https://github.com/pypa/gh-action-pypi-publish/releases)
- [Commits](pypa/gh-action-pypi-publish@v1.8.10...v1.8.11)

---
updated-dependencies:
- dependency-name: pypa/gh-action-pypi-publish
  dependency-type: direct:production
  update-type: version-update:semver-patch
...

Signed-off-by: dependabot[bot] <[email protected]>
Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com>
updates:
- [github.com/psf/black: 23.10.1 → 23.11.0](psf/black@23.10.1...23.11.0)
- [github.com/astral-sh/ruff-pre-commit: v0.1.4 → v0.1.6](astral-sh/ruff-pre-commit@v0.1.4...v0.1.6)

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Bumps [actions/upload-artifact](https://github.com/actions/upload-artifact) from 3 to 4.
- [Release notes](https://github.com/actions/upload-artifact/releases)
- [Commits](actions/upload-artifact@v3...v4)

---
updated-dependencies:
- dependency-name: actions/upload-artifact
  dependency-type: direct:production
  update-type: version-update:semver-major
...

Signed-off-by: dependabot[bot] <[email protected]>
Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com>
Bumps [actions/download-artifact](https://github.com/actions/download-artifact) from 3 to 4.
- [Release notes](https://github.com/actions/download-artifact/releases)
- [Commits](actions/download-artifact@v3...v4)

---
updated-dependencies:
- dependency-name: actions/download-artifact
  dependency-type: direct:production
  update-type: version-update:semver-major
...

Signed-off-by: dependabot[bot] <[email protected]>
Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com>
Bumps [actions/setup-python](https://github.com/actions/setup-python) from 4 to 5.
- [Release notes](https://github.com/actions/setup-python/releases)
- [Commits](actions/setup-python@v4...v5)

---
updated-dependencies:
- dependency-name: actions/setup-python
  dependency-type: direct:production
  update-type: version-update:semver-major
...

Signed-off-by: dependabot[bot] <[email protected]>
Co-authored-by: dependabot[bot] <49699333+dependabot[bot]@users.noreply.github.com>
updates:
- [github.com/psf/black: 23.11.0 → 23.12.1](psf/black@23.11.0...23.12.1)
- [github.com/astral-sh/ruff-pre-commit: v0.1.6 → v0.1.9](astral-sh/ruff-pre-commit@v0.1.6...v0.1.9)

Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
@dfm
Copy link
Owner

dfm commented Jan 4, 2024

@ywx649999311 — I'm so sorry for the delay with this!! I think it looks great. I just went through and made some small changes to the structure. In particular, I managed to remove the eta parameter by using a where instead that should be more stable. I'll leave some comments inline highlighting the main changes that you should check. As soon as you're happy, I'll merge. Thanks again so much for all your work on this!!!

src/tinygp/kernels/quasisep.py Show resolved Hide resolved
src/tinygp/kernels/quasisep.py Show resolved Hide resolved
src/tinygp/kernels/quasisep.py Show resolved Hide resolved
src/tinygp/kernels/quasisep.py Show resolved Hide resolved
src/tinygp/kernels/quasisep.py Show resolved Hide resolved
@ywx649999311
Copy link
Contributor Author

Hi @dfm,
I think it all looks awesome! I added comments regarding the use of [::2] in jnp.ravel(jnp.array([a, b]))[::2]. Thanks for helping with getting rid of eta. I agree the updated implementation is clearer and more stable. Please merge it when you can. I am very happy that this kernel can be used by more people!

@dfm dfm merged commit 3f49239 into dfm:main Jan 4, 2024
8 checks passed
@dfm
Copy link
Owner

dfm commented Jan 4, 2024

Merged!! 🚀

Thanks again @ywx649999311, this is such a great addition!!

@dfm
Copy link
Owner

dfm commented Jan 4, 2024

@ywx649999311 — One other thing: could you open another PR with your name, affiliation and orcid in https://github.com/dfm/tinygp/blob/main/.zenodo.json? I meant to get that as part of this PR, but I forgot. Sorry!

@ywx649999311
Copy link
Contributor Author

@dfm I created a new PR, thanks!

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.

3 participants