API Documentation

Priors(key=dist, ...)

A group of zero or more priors passed by keyword arguments. The priors must be univariate distributions from the Distributions.jl package.

Derived(key=func, ...)

A group of zero or more functions that resolve to a parameter for the model. Derived parameters must be functions accepting one argument if applied to a System, or two arguments if applied to a Planet. Must always return the same value for the same input (be a pure function) and support autodifferentiation.

AstrometryLikelihood(
    (epoch = 5000, ra = -505.7637580573554, dec = -66.92982418533026, σ_ra = 10, σ_dec = 10, cor=0),
    (epoch = 5050, ra = -505.7637580573554, dec = -66.92982418533026, σ_ra = 10, σ_dec = 10, cor=0),
    (epoch = 5100, ra = -505.7637580573554, dec = -66.92982418533026, σ_ra = 10, σ_dec = 10, cor=0),
)

Represents a likelihood function of relative astometry between a host star and a secondary body. :epoch is a required column, in addition to either :ra, :dec, :σ_ra, :σ_dec or :pa, :sep, :σ_pa, :σ_sep. All units are in milliarcseconds or radians as appropriate.

In addition to the example above, any Tables.jl compatible source can be provided.

RadialVelocityLikelihood(
    (;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),
)

Represents a likelihood function of relative astometry between a host star and a secondary body. :epoch (mjd), :rv (km/s), and :σ_rv (km/s), and :inst_idx 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.

PhotometryLikelihood(
    (band = :Z, phot=15.0, σ_phot=3.),
    (band = :J, phot=13.5, σ_phot=0.5),
    (band = :L, phot=11.0, σ_phot=1.0)
)

A likelihood for comparing measured photometry points in one or more filter bands to data (provided here). Requires the :band, :phot', and:σ_phot` columns. Can be provided with any Tables.jl compatible data source.

HGCALikelihood(;gaia_id=1234)

Load proper motion anomaly data from the HIPPARCOS-GAIA Catalog of Accelerations (Brandt et al) for a star with catalog id gaia_id. The resulting velocities are in mas/yr and have the long term trend between HIPPARCOS and GAIA already subtracted out. e.g. we would expect 0 pma if there is no companion.

ObsPriorAstromONeil2019(astrometry_likelihood, period_prior)

Given a an astrometry likelihood (AstrometryLikelihood), apply the "observable based priors" of K. O'Neil 2019 "Improving Orbit Estimates for Incomplete Orbits with a New Approach to Priors: with Applications from Black Holes to Planets".

This prior correction is only correct if you supply Uniform priors on all Campbell orbital parameters and a Uniform prior on Period (not semi-major axis). This period prior has a significant impact in the fit and recommendations for its range were not published in the original paper.

Examples

astrom_like = AstrometryLikelihood(astrom_table)

# Apply observable based priors ontop of our uniform Campbell priors:
obs_prior = ObsPriorAstromONeil2019(astrom_like)

# The astrometry lieklihood object is passed as a first parameter
# since the obserable-based priors depend on the observation 
# epochs.

@planet b Visual{KepOrbit} begin
    # Instead of a prior on sma
    # a ~ Uniform(0.001, 10000)

    # Put a prior on period:
	P ~ Uniform(0.001, 2000) # yrs
    a = cbrt(system.M * b.P^2)

    # Keep sine prior on inclination
    i ~ Sine()

    # Rest are uniform
    e ~ Uniform(0.0, 1.0)
    ω ~ UniformCircular()
    Ω ~ UniformCircular()
    τ ~ UniformCircular(1.0)
end astrom_like obs_prior

Sine()

A custom univariate distribution. The pdf is a sine function defined between 0 and π. This is a common prior distribution used when fitting orbits to astrometry.

The full Distributions.jl interface is not yet defined for this distribution, but the following methods work: pdf, logpdf, minimum, maximum, insupport, mean, var, cdf, quantile

Planet([derived,] priors, [astrometry,], name=:symbol)

A planet (or substellar companion) part of a model. Must be constructed with a block of priors, and optionally additional derived parameters and/or sastrometry. name must be a symbol, e.g. :b.

System([derived,] priors, [images,] [propermotionanom,] planets..., name=:symbol)

Construct a model of a system. Must be constructed with a block of priors, and optionally additional derived parameters. You may provide ProperMotionAnomLikelihood() and/or Images() of the system. Finally, planet models are listed last. name must be a symbol e.g. :HD56441.

mjd("2020-01-01")

Get the modfied julian day of a date, or in general a UTC timestamp.

mjd()

Get the current modified julian day of right now.

gaia_plx(gaia_id=12123)

Get a distribution (truncated Normal) of parallax distance in mas of a source with GAIA catalog id gaia_id.

The method signature of Octofitter.hmc is as follows:

advancedhmc(
    [rng::Random.AbstractRNG],
    model::Octofitter.LogDensityModel
    target_accept::Number=0.8,
    ensemble::AbstractMCMC.AbstractMCMCEnsemble=MCMCSerial();
    adaptation,
    iterations,
    drop_warmup=true,
    max_depth=12,
    initial_samples=50_000,
    initial_parameters=nothing,
    step_size=nothing,
    verbosity=2,
)

The only required arguments are system, adaptation, and iterations. The two positional arguments are system, the model you wish to sample; and targetaccept, the acceptance rate that should be targeted during windowed adaptation. During this time, the step size and mass matrix will be adapted (see AdvancedHMC.jl for more information). The number of steps taken during adaptation is controlled by adaptation. You can prevent these samples from being dropped by pasing includeadaptation=false. The total number of posterior samples produced are given by iterations. These include the adaptation steps that may be discarded. treedepth controls the maximum tree depth of the sampler. initialparameters is an optional way to pass a starting point for the chain. If you don't pass a default position, one will be selected by drawing initial_samples from the priors. The sample with the highest posterior value will be used as the starting point.

Visual{OrbitType}(..., plx=...)

This wraps another orbit to add the parallax distance field plx, thus allowing projected quantities to be calculated. It forwards everything else to the parent orbit.

For example, the KepOrbit type supports calculating x and y positions in AU. A Visual{KepOrbit} additionally supports calculating projected right ascension and declination offsets.

sonora_photometry_interpolator(:Keck_L′, [metalicity="+0.0"])

Given a supported photometric band and [M/H] metalicity (default=solar), return a function of temperature (K) and mass (M_jup) that gives the absolute magnitude of the planet in that bandpass.

Supported bands: :MKOY, :MKOZ, :MKOJ, :MKOH, :MKOK, :MKOL′, :MKOM′, :TwoMASSJ, :TwoMASSH, :TwoMASSKs, :KeckKs, :KeckL′, :KeckMs, :SDSSg′, :SDSSr′, :SDSSi′, :SDSSz′, :IRAC36, :IRAC45, :IRAC57, :IRAC79, :WISEW1, :WISEW2, :WISEW3, :WISE_W4

Supported metalicities: "+0.0", "-0.5", "+0.5"

itp = sonora_cooling_interpolator()

Create a function mapping (ageMyr, massMjup) -> temp_K using Sonora Bobcat cooling model grids.

plotchains(
    chain, planet_key;
    N=1500,
    ii = rand(1:size(chain,1)*size(chain,3), N),
    color=length(model.system.planets) == 0 || !haskey(chain, string(planet_key)*"_a") ? nothing : string(planet_key)*"_a",
    colorbartitle=color,
    clims=nothing,
    cmap=:plasma,
    alpha=30/length(ii),
    attime=nothing,
    kwargs...,
)

Draw samples from a posterior chain for a given planet given by name planet_key and visualize them in some way. Use kind to control what plot is made. A few options: :astrometry, :radvel, :trueanom, :meananom, :eccanom, :x, :y, :z, (:x, :y), :raoff, :decoff, :pmra, :pmdec, :accra, :accdec, :radvel, :posangle, :projectedseparation. See PlanetOrbits documentation for more details.

Inputs:

• chain The chain to draw from
• planet_key Planet name in the model (symbol)
• N=1500 Number of samples to draw for the plot
• kind=nothing Specify what kind of plot to make.
• ii=... Specific row numbers to use, if you want to e.g. plot the same 100 samples in a few different plots
• color="planet_key_a" Column name to to map colors to. Semi-major axis by default but can be any column or an arbitrary array.
• colorbartitle=color Name for colourbar
• clims=nothing Tuple of colour limits (min and max)
• cmap=:plasma Colormap
• alpha=... Transparency of the lines

projectpositions(chains.planets[1], mjd("2020-02-02"))

Given the posterior for a particular planet in the model and a modified julian date(s), return ra and dec offsets in mas for each sampling in the posterior.