Fit Relative RV Data
Octofitter includes support for fitting relative radial velocity data. Currently this is only tested with a single companion. Please open an issue if you would like to fit multiple companions simultaneously.
The convention we adopt is that positive relative radial velocity is the velocity of the companion (exoplanets) minus the velocity of the host (star).
To fit relative RV data, start by creating a likelihood object:
using Octofitter
using OctofitterRadialVelocity
using CairoMakie
using Distributions
rel_rv_like = PlanetRelativeRVLikelihood(
(;epoch=5000.0, rv=−6.54, σ_rv=1.30),
(;epoch=5050.1, rv=−3.33, σ_rv=1.09),
(;epoch=5100.2, rv=7.90, σ_rv=.11),
)
PlanetRelativeRVLikelihood Table with 4 columns and 3 rows:
inst_idx epoch rv σ_rv
┌──────────────────────────────
1 │ 1 5000.0 -6.54 1.3
2 │ 1 5050.1 -3.33 1.09
3 │ 1 5100.2 7.9 0.11
The columns in the table should be . See the standard radial velocity tutorial for examples on how this data can be loaded from a CSV file.
The relative RV likelihood does not incorporate an instrument-specific RV offset. A jitter parameter called jitter
can still be specified in the @planet
block, as can parameters for a gaussian process model of stellar noise. Unlike the StarAbsoluteRVLikelihood
, only a single instrument jitter parameter is supported. If you need to model relative radial velocities from multiple instruments with different jitters, please open an issue on GitHub.
Next, create a planet and system model, attaching the relative rv likelihood to the planet. Make sure to add a jitter
parameter (optionally jitter=0) to the planet.
@planet b Visual{KepOrbit} begin
a ~ truncated(Normal(10, 4), lower=0, upper=100)
e ~ Uniform(0.0, 0.5)
i ~ Sine()
ω ~ UniformCircular()
Ω ~ UniformCircular()
θ ~ UniformCircular()
tp = θ_at_epoch_to_tperi(system,b,50420)
jitter ~ LogUniform(1, 1000) # m/s
end rel_rv_like
@system ExampleSystem begin
M ~ truncated(Normal(1.2, 0.1), lower=0)
plx ~ truncated(Normal(50.0, 0.02), lower=0)
end b
model = Octofitter.LogDensityModel(ExampleSystem)
LogDensityModel for System ExampleSystem of dimension 12 and 6 epochs with fields .ℓπcallback and .∇ℓπcallback
using Random
rng = Random.Xoshiro(1234)
chain = octofit(rng, model)
Chains MCMC chain (1000×29×1 Array{Float64, 3}):
Iterations = 1:1:1000
Number of chains = 1
Samples per chain = 1000
Wall duration = 37.85 seconds
Compute duration = 37.85 seconds
parameters = M, plx, b_a, b_e, b_i, b_ωy, b_ωx, b_Ωy, b_Ωx, b_θy, b_θx, b_jitter, b_ω, b_Ω, b_θ, b_tp
internals = n_steps, is_accept, acceptance_rate, hamiltonian_energy, hamiltonian_energy_error, max_hamiltonian_energy_error, tree_depth, numerical_error, step_size, nom_step_size, is_adapt, loglike, logpost, tree_depth, numerical_error
Summary Statistics
parameters mean std mcse ess_bulk ess_tail rh ⋯
Symbol Float64 Float64 Float64 Float64 Float64 Float ⋯
M 1.1911 0.1069 0.0055 388.8657 321.9627 1.00 ⋯
plx 50.0008 0.0194 0.0007 868.8648 493.8274 0.99 ⋯
b_a 14.1994 2.9027 0.8250 11.1044 78.5488 1.09 ⋯
b_e 0.2793 0.1498 0.0102 167.6183 115.9231 1.00 ⋯
b_i 3.0055 0.2132 0.0141 158.3943 154.3445 1.04 ⋯
b_ωy 0.0031 0.6785 0.0615 105.5329 347.8646 1.03 ⋯
b_ωx 0.1674 0.7318 0.2137 14.4038 275.1460 1.08 ⋯
b_Ωy -0.0446 0.7236 0.0661 124.0340 228.4814 1.00 ⋯
b_Ωx 0.1575 0.6921 0.0580 167.6722 540.8796 1.00 ⋯
b_θy -0.0797 0.7160 0.0588 170.8610 390.5018 1.01 ⋯
b_θx -0.0132 0.7135 0.0697 122.2277 489.3976 1.02 ⋯
b_jitter 24.8286 45.5191 9.3271 6.6731 119.4842 1.21 ⋯
b_ω 0.3801 1.7709 0.4516 20.9093 481.8734 1.07 ⋯
b_Ω 0.4291 1.8207 0.1346 212.9350 111.0761 0.99 ⋯
b_θ -0.0431 1.8932 0.1670 175.4939 414.4864 1.00 ⋯
b_tp 41327.0735 6863.8588 727.6021 90.2470 118.9133 1.03 ⋯
2 columns omitted
Quantiles
parameters 2.5% 25.0% 50.0% 75.0% 97.5%
Symbol Float64 Float64 Float64 Float64 Float64
M 0.9830 1.1150 1.1945 1.2675 1.4016
plx 49.9625 49.9878 50.0015 50.0144 50.0379
b_a 9.6008 12.1107 13.9349 15.9886 20.4772
b_e 0.0128 0.1461 0.2970 0.4143 0.4944
b_i 2.2827 3.0051 3.0788 3.1134 3.1392
b_ωy -1.0739 -0.6641 0.0893 0.5890 1.0540
b_ωx -1.0716 -0.5305 0.2886 0.8569 1.1015
b_Ωy -1.1147 -0.7493 -0.0657 0.6544 1.0739
b_Ωx -1.0526 -0.4321 0.2683 0.8006 1.0773
b_θy -1.1022 -0.7538 -0.2049 0.6345 1.0597
b_θx -1.0597 -0.7045 -0.0586 0.7162 1.0497
b_jitter 2.4507 5.1908 9.3945 23.4636 144.0943
b_ω -2.9947 -1.1238 0.7879 1.7519 3.0365
b_Ω -3.0111 -1.0369 0.7025 2.0502 3.0505
b_θ -2.9539 -1.8676 -0.0966 1.7049 3.0215
b_tp 24505.3566 36379.0640 42170.7390 47409.1067 50158.0714
octoplot(model, chain)