RVProbability

class ssapy.rvsampler.RVProbability(arc, epoch, priors=[<ssapy.rvsampler.RPrior object>, <ssapy.rvsampler.APrior object>], propagator=KeplerianPropagator(), meanpm=False, damp=-1)[source][source]

Bases: object

A class to manage MCMC sampling of orbits (parameterized by position and velocity at a given epoch) given some angular observations.

Parameters:
  • arc (QTable with one row per observation and columns:) –

    ‘ra’, ‘dec’Angle

    Observed angular position in ICRS (topocentric).

    ’rStation_GCRF’Quantity

    Position of observer in GCRF.

    ’sigma’Quantity

    Isotropic uncertainty in observed angular position.

    ’time’float or astropy.time.Time

    Time of observation. If float, then should correspond to GPS seconds; i.e., seconds since 1980-01-06 00:00:00 UTC

    The arc is a linked set of observations assumed to be of the same object.

  • epoch (float or astropy.time.Time) – The time at which the model parameters (position and velocity) are specified. It probably makes sense to make this the same as one of the observation times, but it’s not required to be. If float, then should correspond to GPS seconds; i.e., seconds since 1980-01-06 00:00:00 UTC

  • priors (list of priors, optional) – A list of class instances similar to RPrior() or APrior() to apply jointly to the posterior probability.

  • propagator (class, optional) – The propagator class to use.

  • meanpm (bool) – fit model to mean positions and proper motions rather than to individual epochs along each streak.

  • damp (float) – damp chi residuals with pseudo-Huber loss function. Default -1: do not damp.

arc[source]
epoch[source]
priors[source]
propagator[source]
lnlike(orbit)[source][source]

Calculate log likelihood of observations given orbit.

lnprob(p)[source][source]

Calculate log posterior probability of orbit parameters p given observations (up to a constant). Parameters p are [x, y, z, vx, vy, vz] at epoch.

__call__(p)[source]

Alias for lnprob(p)

Methods Summary

__call__(p)

Calculate the log posterior probability of parameters p given observations.

chi(orbit)

Calculate chi residuals for use with lmfit.

lnlike(orbit)

Calculate the log likelihood of the observations given an orbit.

lnprior(orbit)

Calculate the log prior probability of the observations given an orbit.

lnprob(p)

Calculate the log posterior probability of parameters p given observations.

Methods Documentation

__call__(p)[source]

Calculate the log posterior probability of parameters p given observations.

Parameters:

p (array (6,)) – Sampling parameters. [x, y, z, vx, vy, vz] in meters and meters per second. Understood to be the parameters at epoch.

Returns:

  • lnprob (float) – Log posterior probability.

  • lnprior (float) – Log prior probability.

chi(orbit)[source][source]

Calculate chi residuals for use with lmfit.

This is essentially the input to a chi^2 statistic, just without the squaring and without the summing over observations. I.e., it’s

[(data_i - model_i)/err_i for i in 0..nobs]

Parameters:

orbit (Orbit) – Orbit in question. The model.

Returns:

  • chilike (array_like (nobs,)) – chi residual array for likelihoods

  • chiprior (array_like (nprior,)) – chi residual array for priors

lnlike(orbit)[source][source]

Calculate the log likelihood of the observations given an orbit.

Parameters:

orbit (Orbit) – Orbit in question.

Returns:

lnlike – Log likelihood.

Return type:

float

lnprior(orbit)[source][source]

Calculate the log prior probability of the observations given an orbit.

Parameters:

orbit (Orbit) – Orbit in question.

Returns:

lnprior – Log prior probability.

Return type:

float

lnprob(p)[source][source]

Calculate the log posterior probability of parameters p given observations.

Parameters:

p (array (6,)) – Sampling parameters. [x, y, z, vx, vy, vz] in meters and meters per second. Understood to be the parameters at epoch.

Returns:

  • lnprob (float) – Log posterior probability.

  • lnprior (float) – Log prior probability.