RV: migrate relative RV to new API

OctofitterRadialVelocity/src/rv-absolute-margin.jl

@@ -6,7 +6,30 @@
 
 using Bumper
 
-# Gaussian processes not supported.
+"""
+    MarginalizedStarAbsoluteRVLikelihood(
+        (;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),
+
+        jitter=:jitter_1,
+        instrument_name="inst name",
+    )
+
+Represents a likelihood function of relative astometry between a host star and a secondary body.
+`:epoch` (mjd), `:rv` (m/s), and `:σ_rv` (m/s)  are all required.
+In addition to the example above, any Tables.jl compatible source can be provided.
+
+The`jitter` parameters specify which variables should be read from the model for the 
+    jitter of this instrument.
+    
+Unlike `StarAbsoluteRVLikelihood`, this version analytically marginalizes over the instrument
+RV zero point. This makes it faster in most cases.
+That said, if you have a specific prior you want to apply for the RV zero point, correlations between
+zero points, hierarchical models, etc, you should use `StarAbsoluteRVLikelihood`.
+
+Additionally, a gaussian and trend fit are not supported with the analytically marginalization.
+"""
 struct MarginalizedStarAbsoluteRVLikelihood{TTable<:Table,jitter_symbol} <: Octofitter.AbstractLikelihood
     table::TTable
     instrument_name::String

OctofitterRadialVelocity/src/rv-absolute.jl

@@ -4,16 +4,25 @@ const rv_cols = (:epoch, :rv, :σ_rv)
 
 """
     StarAbsoluteRVLikelihood(
-        (;inst_idx=1, epoch=5000.0,  rv=−6.54, σ_rv=1.30),
-        (;inst_idx=1, epoch=5050.1,  rv=−3.33, σ_rv=1.09),
-        (;inst_idx=1, epoch=5100.2,  rv=7.90,  σ_rv=.11),
+        (;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),
+
+        offset=:offset_1,
+        jitter=:jitter_1,
+        instrument_name="inst name",
+
+        # Optional:
+        trend_function=(θ_system, epoch)->0.,
+        gaussian_process=nothing
     )
 
 Represents a likelihood function of relative astometry between a host star and a secondary body.
-`:epoch` (mjd), `:rv` (m/s), and `:σ_rv` (m/s), and `:inst_idx` are all required.
+`:epoch` (mjd), `:rv` (m/s), and `:σ_rv` (m/s)  are all required.
+
+The `offset` and `jitter` parameters specify which variables should be read from the model for the 
+RV zero-point and jitter of this instrument.
 
-`:inst_idx` is used to distinguish RV time series between insturments so that they may optionally
-be fit with different zero points and jitters.
 In addition to the example above, any Tables.jl compatible source can be provided.
 """
 struct StarAbsoluteRVLikelihood{TTable<:Table,GP,TF,offset_symbol,jitter_symbol} <: Octofitter.AbstractLikelihood

OctofitterRadialVelocity/src/rv-relative.jl

@@ -1,42 +1,43 @@
 """
     PlanetRelativeRVLikelihood(
-        (;inst_idx=1, epoch=5000.0,  rv=−6.54, σ_rv=1.30),
-        (;inst_idx=1, epoch=5050.1,  rv=−3.33, σ_rv=1.09),
-        (;inst_idx=1, epoch=5100.2,  rv=7.90,  σ_rv=.11),
+        (;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),
+
+        jitter=:jitter_1,
+        instrument_name="inst name",
     )
 
 Represents a likelihood function of relative astometry between a host star and a secondary body.
-`:epoch` (mjd), `:rv` (m/s), and `:σ_rv` (m/s), and `:inst_idx` are all required.
+`:epoch` (mjd), `:rv` (m/s), and `:σ_rv` (m/s) are all required.
 
-`:inst_idx` is used to distinguish RV time series between insturments so that they may optionally
-be fit with different zero points and jitters.
 In addition to the example above, any Tables.jl compatible source can be provided.
+
+The  `jitter` parameter specify which variables should be read from the model for the 
+jitter of this instrument.
 """
-struct PlanetRelativeRVLikelihood{TTable<:Table,GP} <: Octofitter.AbstractLikelihood
+struct PlanetRelativeRVLikelihood{TTable<:Table,GP,jitter_symbol} <: Octofitter.AbstractLikelihood
     table::TTable
-    instrument_names::Vector{String}
+    instrument_name::String
     gaussian_process::GP
-    function PlanetRelativeRVLikelihood(observations...; instrument_names=nothing, gaussian_process=nothing)
+    jitter_symbol::Symbol
+    function PlanetRelativeRVLikelihood(observations...; instrument_name="1", gaussian_process=nothing, jitter)
         table = Table(observations...)
         if !issubset(rv_cols, Tables.columnnames(table))
             error("Expected columns $rv_cols")
         end
 
         # We sort the data first by instrument index then by epoch to make some later code faster
-        m = maximum(table.epoch)
         ii = sortperm(table.epoch)
         table = table[ii,:]
-        if isnothing(instrument_names)
-            instrument_names = ["1"]
-        end
 
         rows = map(eachrow(table)) do row′
-            row = (;inst_idx=1, row′[1]..., rv=float(row′[1].rv[1]))
+            row = (;row′[1]..., rv=float(row′[1].rv[1]))
             return row
         end
         table = Table(rows)
 
-        return new{typeof(table),typeof(gaussian_process)}(table, instrument_names, gaussian_process)
+        return new{typeof(table),typeof(gaussian_process),jitter}(table, instrument_name, gaussian_process, jitter)
     end
 end
 PlanetRelativeRVLikelihood(observations::NamedTuple...;kwargs...) = PlanetRelativeRVLikelihood(observations; kwargs...)
@@ -49,6 +50,12 @@ function Octofitter.likeobj_from_epoch_subset(obs::PlanetRelativeRVLikelihood, o
 end
 export PlanetRelativeRVLikelihood
 
+function _getparams(::PlanetRelativeRVLikelihood{TTable,GP,jitter_symbol}, θ_planet) where {TTable,GP,jitter_symbol}
+    jitter = getproperty(θ_planet, jitter_symbol)
+    return (;jitter)
+end
+
+
 """
 Radial velocity likelihood.
 """
@@ -72,7 +79,7 @@ function Octofitter.ln_like(
     epochs = vec(rvlike.table.epoch)
     σ_rvs = vec(rvlike.table.σ_rv)
     rvs = vec(rvlike.table.rv)
-    jitter = θ_planet.jitter
+    jitter = _getparams(rvlike, θ_planet).jitter
 
     @no_escape begin
         noise_var = @alloc(T, length(rvs))

docs/src/fit-rv-rel.md

@@ -37,11 +37,14 @@ rel_rv_like = PlanetRelativeRVLikelihood(
     (epoch=7200, rv= 27083.2046462065, σ_rv=1500.0),
     (epoch=7300, rv=  9174.176455190982, σ_rv=1500.0),
     (epoch=7400, rv=-22241.45434114139, σ_rv=1500.0),
-    )
+
+    instrument_name="simulated data",
+    jitter=:gamma # name of jitter variable to use
+)
 ```
-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.
+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. Currently 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.
+The relative RV likelihood does not incorporate an instrument-specific RV offset. A jitter parameter called can still be specified in the `@planet` block, as can parameters for a gaussian process model of stellar noise. Currently 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.
 
@@ -56,7 +59,7 @@ Next, create a planet and system model, attaching the relative rv likelihood to
     P = √(b.a^3/system.M)
     tp =  b.τ*b.P*365.25 + 6000 # reference epoch for τ. Choose an MJD date near your data.
 
-    jitter ~ LogUniform(0.1, 1000) # m/s
+    gamma ~ LogUniform(0.1, 1000) # m/s
 end rel_rv_like
 
 @system ExampleSystem begin