ocpy package

Submodules

ocpy.custom_types module

ocpy.data module

class ocpy.data.Data(minimum_time, minimum_time_error=None, weights=None, minimum_type=None, labels=None)[source]

Bases: DataModel

Container for eclipse minimum timing data.

The Data class stores and manages observational times of minima (e.g., eclipse timings of binary stars) together with optional metadata such as timing uncertainties, observational weights, minimum type (primary/secondary), and labels.

Internally the data are stored in a pandas.DataFrame with standardized column names. The class provides utilities for safely manipulating the dataset, including:

  • Filling or computing timing uncertainties and weights

  • Loading data from files

  • Computing O–C (Observed minus Calculated) values

  • Grouping or merging datasets

The class is designed to behave like a lightweight table object while preserving domain-specific semantics required for O–C analysis of eclipsing binaries.

Parameters:
  • minimum_time (list-like) – Times of observed minima (typically in Julian Date or BJD). This field is required.

  • minimum_time_error (list-like, optional) – Uncertainties associated with each minimum time. If provided, the length must match minimum_time.

  • weights (list-like, optional) – Weights assigned to each observation.

  • minimum_type (BinarySeq, optional) – Indicator of minimum type (e.g., primary or secondary eclipse). Accepted representations may include integers, strings (e.g., "I", "II", "primary", "secondary"), or other binary encodings.

  • labels (list-like, optional) – Optional labels or identifiers for each observation.

Raises:
  • ValueError – If minimum_time is None.

  • LengthCheckError – If provided sequences do not match the length of minimum_time.

Notes

All input sequences are normalized internally so that their lengths match the length of minimum_time. Missing values are automatically expanded or filled using utilities from ocpy.utils.

The underlying storage is a pandas.DataFrame with the following standard columns:

  • minimum_time

  • minimum_time_error

  • weights

  • minimum_type

  • labels

Most modification methods return a new `Data` instance rather than mutating the existing object.

Examples

Create a dataset with minimum times:

>>> from ocpy import Data
>>> d = Data(minimum_time=[2450000.1, 2450001.3, 2450002.5])

With uncertainties and types:

>>> d = Data(
...     minimum_time=[2450000.1, 2450001.3],
...     minimum_time_error=[0.0002, 0.0003],
...     minimum_type=["I", "II"]
... )

Access the underlying table:

>>> d.data
classmethod from_file(file, columns=None)[source]

Create a Data object from a tabular file.

This method reads timing data from a CSV or Excel file and converts it into a Data instance. Column names in the file can optionally be mapped to the standardized column names used by the Data class.

Parameters:
  • file (str or pathlib.Path) –

    Path to the input file. Supported formats are:

    • .csv

    • .xls

    • .xlsx

  • columns (dict of str to str, optional) –

    Mapping between the standardized Data column names and the column names present in the file.

    Two mapping styles are accepted:

    1. Standard → file column name

      Example:

      >>> columns = {"minimum_time": "BJD", "minimum_time_error": "err"}
      
    2. File column name → standard

      Example:

      >>> columns = {"BJD": "minimum_time", "err": "minimum_time_error"}
      

    The standard column names recognized by the Data class are:

    • minimum_time

    • minimum_time_error

    • weights

    • minimum_type

    • labels

Returns:

A new Data instance containing the imported observations.

Return type:

Data

Raises:
  • ValueError – If the file format is not supported.

  • ValueError – If the file does not contain a column corresponding to minimum_time.

Notes

The minimum_time column is required to construct a valid dataset. All other columns are optional and will be set to None if not present in the file.

Internally, the file is loaded using pandas.read_csv() or pandas.read_excel().

Examples

Load a dataset from a CSV file:

>>> d = Data.from_file("minima.csv")

Load a file with custom column names:

>>> d = Data.from_file(
...     "observations.csv",
...     columns={"BJD": "minimum_time", "err": "minimum_time_error"}
... )

Load an Excel file:

>>> d = Data.from_file("minima.xlsx")
fill_errors(errors, override=False)[source]

Fill or assign timing uncertainties for the dataset.

This method returns a new Data object in which the minimum_time_error column is populated using the provided values. Existing values can optionally be preserved or replaced.

Parameters:
  • errors (list, tuple, numpy.ndarray, or float) –

    Uncertainty values to assign to minimum_time_error.

    • If a scalar is provided, it will be applied to all rows.

    • If an array-like object is provided, its length must match the number of observations in the dataset.

  • override (bool, default False) –

    If True, all existing values in the minimum_time_error column are replaced.

    If False, only entries that are currently NaN are filled, leaving existing values unchanged.

Returns:

A new Data instance with updated minimum_time_error values.

Return type:

Data

Raises:

LengthCheckError – If errors is array-like and its length does not match the number of rows in the dataset.

Notes

This method does not modify the original object. Instead, it returns a new Data instance with the updated values.

Internally, column assignment is handled by the private _assign_or_fill() method to ensure consistent behavior across similar operations.

Examples

Assign a constant uncertainty to all observations:

>>> d2 = d.fill_errors(0.0002)

Fill only missing uncertainties:

>>> d2 = d.fill_errors([0.0002, 0.0003, 0.00025])

Replace all existing uncertainties:

>>> d2 = d.fill_errors(0.0002, override=True)
fill_weights(weights, override=False)[source]

Fill or assign observational weights for the dataset.

This method returns a new Data object in which the weights column is populated using the provided values. Existing weights can optionally be preserved or replaced.

Parameters:
  • weights (list, tuple, numpy.ndarray, or float) –

    Weight values to assign to the weights column.

    • If a scalar is provided, it will be applied to all rows.

    • If an array-like object is provided, its length must match the number of observations in the dataset.

  • override (bool, default False) –

    If True, all existing values in the weights column are replaced.

    If False, only entries that are currently NaN are filled, leaving existing values unchanged.

Returns:

A new Data instance with updated weights values.

Return type:

Data

Raises:

LengthCheckError – If weights is array-like and its length does not match the number of rows in the dataset.

Notes

This method does not modify the original object. Instead, it returns a new Data instance with the updated values.

Internally, column assignment is handled by the private _assign_or_fill() method to ensure consistent behavior across similar operations.

Examples

Assign a constant weight to all observations:

>>> d2 = d.fill_weights(1.0)

Fill only missing weights:

>>> d2 = d.fill_weights([1.0, 0.5, 2.0])

Replace all existing weights:

>>> d2 = d.fill_weights(1.0, override=True)
calculate_weights(method=None, override=True)[source]

Calculate observational weights based on timing uncertainties.

This method computes weights for each observation, typically using the inverse-variance method, and returns a new Data instance with the updated weights column.

Parameters:
  • method (callable, optional) –

    A custom function to compute weights from the minimum_time_error series. It must accept a pandas.Series of errors and return a pandas.Series of weights.

    If None (default), the inverse-variance method is used:

    \[w_i = \frac{1}{\sigma_i^2}\]

    where \(\sigma_i\) is the timing uncertainty for the i-th observation.

  • override (bool, default True) – If True, existing weights values are replaced. If False, only missing entries (NaN) are filled.

Returns:

A new Data instance with updated weights.

Return type:

Data

Raises:
  • ValueError – If minimum_time_error contains NaN values.

  • ValueError – If minimum_time_error contains zero values, which would cause division by zero in the default method.

  • TypeError – If method is provided but is not callable.

Notes

  • The default inverse-variance weighting gives higher weight to observations with smaller uncertainties.

  • This method does not modify the original Data instance; it returns a new instance with updated weights.

  • Internally, column assignment uses _assign_or_fill() to respect the override flag.

Examples

Compute default inverse-variance weights:

>>> d2 = d.calculate_weights()

Compute weights with a custom method:

>>> def custom_weights(errors):
...     return 1 / errors
>>> d2 = d.calculate_weights(method=custom_weights, override=True)

Fill only missing weights without overwriting existing ones:

>>> d2 = d.calculate_weights(override=False)
calculate_oc(reference_minimum, reference_period, model_type='lmfit')[source]

Compute Observed minus Calculated (O–C) values for the dataset.

This method calculates the O–C values for each observed minimum based on a reference minimum time and period. The O–C values quantify the difference between observed and predicted timings, which is fundamental for analyzing period variations in eclipsing binaries.

reference_minimumfloat

The reference time of minimum (e.g., initial epoch) used to compute predicted minima.

reference_periodfloat

The reference orbital period of the system. This is used to compute the expected timing of each cycle.

model_typestr, default=’lmfit’

Specifies the type of O–C model to return. Supported options:

  • 'lmfit' or 'lmfit_model' – returns an OCLMFit object.

  • 'pymc' or 'pymc_model' – returns an OCPyMC object.

  • Any other string – returns a generic OC object.

OC

An instance of the appropriate O–C model class (OC, OCLMFit, or OCPyMC) containing:

  • minimum_time – observed minima

  • cycle – computed cycle numbers (integer or half-integer for secondary minima)

  • oc – O–C values

  • Additional columns from the original Data (errors, weights, labels, minimum_type)

ValueError

If the minimum_time column is missing.

  • Cycle calculation: The phase of each observation is computed as:

    \[ext{phase} =\]

rac{t - ext{reference_minimum}}{ ext{reference_period}}

The cycle number is the nearest integer to the phase.

  • Secondary minima: If minimum_type is present and indicates a secondary eclipse (e.g., “II”, “secondary”, 2), the cycle is adjusted to half-integer values:

    \[ext{cycle}_{ ext{sec}} = ext{round}( ext{phase} - 0.5) + 0.5\]
  • O–C computation: The O–C value for each observation is:

    \[ext{O–C} = t_{ ext{obs}} - ( ext{reference_minimum} + ext{cycle} imes ext{reference_period} )\]
  • The method does not modify the original `Data`. The returned object contains a copy of the original data along with the computed cycle and oc arrays.

  • The model_type determines which class is instantiated for further modeling of the O–C diagram.

Compute O–C values using the default LMFit model:

>>> oc_model = d.calculate_oc(reference_minimum=2450000.0, reference_period=1.2345)

Compute O–C values using a PyMC model:

>>> oc_model = d.calculate_oc(
...     reference_minimum=2450000.0,
...     reference_period=1.2345,
...     model_type="pymc"
... )

Access the computed O–C values:

>>> oc_model.oc
[0.0001, -0.0002, 0.0003, ...]
merge(data)[source]

Merge the current dataset with another Data object.

This method concatenates the rows of the current Data instance with those of another Data object, returning a new Data instance. Column alignment is based on column names.

Parameters:

data (Data) – Another Data instance to merge with the current dataset.

Returns:

A new Data instance containing all rows from both datasets.

Return type:

Data

Notes

  • The original datasets are not modified.

  • Missing columns in either dataset will result in NaN values in the merged dataset, following pandas’ concatenation rules.

  • Indexes are reset in the merged dataset for consistency.

Examples

>>> d1 = Data(minimum_time=[2450000.1, 2450001.2])
>>> d2 = Data(minimum_time=[2450002.3, 2450003.4])
>>> d_merged = d1.merge(d2)
>>> len(d_merged)
4
group_by(column)[source]

Split the dataset into groups based on a column.

This method groups the Data object by the values in a specified column and returns a list of new Data instances, each containing one group of rows.

Parameters:

column (str) – Name of the column to group by.

Returns:

A list of Data objects, each corresponding to one group. If the column is missing or contains only NaN values, a list with a single copy of the original dataset is returned.

Return type:

list of Data

Notes

  • Grouping is performed using pandas.DataFrame.groupby().

  • The original Data object is not modified; each group is a deep copy.

  • NaN values are treated as a separate group unless dropna=True in the internal pandas grouping.

Examples

Group a dataset by the minimum_type column:

>>> groups = d.group_by("minimum_type")
>>> len(groups)
2  # e.g., one group for primary, one for secondary minima

Access the first group:

>>> groups[0].data
minimum_time  minimum_type
2450000.1    I
2450002.3    I

If the grouping column does not exist:

>>> groups = d.group_by("nonexistent_column")
>>> len(groups)
1  # returns a single copy of the original Data

ocpy.errors module

exception ocpy.errors.LengthCheckError[source]

Bases: Exception

Raised when the length of data is nor sufficient

ocpy.model_data module

class ocpy.model_data.DataModel(minimum_time, minimum_time_error=None, weights=None, minimum_type=None, labels=None, ecorr=None, oc=None)[source]

Bases: ABC

abstractmethod classmethod from_file(file, columns=None)[source]

Read data from file

abstractmethod fill_errors(errors, override=False)[source]

Fills th errors

abstractmethod fill_weights(weights, override=False)[source]

Fills th weights

abstractmethod calculate_weights(method=None, override=True)[source]

Calculates weights using errors

abstractmethod calculate_oc(reference_minimum, reference_period, model_type='lmfit')[source]

Calculates the O-C for this Data

abstractmethod merge(data)[source]

Appends data to this DataModel

abstractmethod group_by(column)[source]

Group data by column’s data

ocpy.model_oc module

class ocpy.model_oc.ParameterModel[source]

Bases: ABC

value: float
min: float
max: float
std: float
fixed: bool = False
class ocpy.model_oc.ModelComponentModel(args=None)[source]

Bases: ABC

abstractmethod model_function()[source]

Definition of the function to fit

class ocpy.model_oc.OCModel(minimum_time, minimum_time_error=None, weights=None, minimum_type=None, labels=None, ecorr=None, oc=None)[source]

Bases: ABC

abstractmethod classmethod from_file(file, columns=None)[source]

Read data from file

abstractmethod bin(bin_count=1, bin_method=None, bin_error_method=None, bin_style=None)[source]

Bins the data and returns each a new Self

abstractmethod merge(oc)[source]

Appends oc to this oc

abstractmethod residue(coefficients)[source]

Removes the fit from current data

abstractmethod fit(model_components)[source]

Fits the given ModelComponents to the O-C

abstractmethod fit_keplerian(parameters)[source]

Makes a keplerian fit (also known as lite)

abstractmethod fit_lite(parameters)[source]

Makes a lite fit (also known as keplerian fit)

abstractmethod fit_linear(parameters)[source]

Makes a linear fit

abstractmethod fit_quadratic(parameters)[source]

Makes a quadratic fit

abstractmethod fit_sinusoidal(parameters)[source]

Makes a sinusoidal fit

ocpy.newtonian module

class ocpy.newtonian.NewtonianModel(*, integrator='ias15', dt=0.01, integrator_params=None, units=None, reference_time=0.0, t_start=None, t_end=None, stop_at_exact_time=True, escape_radius=100000.0, min_distance=0.0001, precision_integration_steps=0, integration_grid=None, central_mass=1.0, bodies=None, orbit_type='heliocentric', orbit_output_type='heliocentric', T0_ref=0.0, P_ref=1.0, compute_xyz=True, compute_orbital=True, name=None)[source]

Bases: ModelComponent

name = 'newtonian'
integrate(times, params_dict=None, *, compute_orbital=None)[source]
model_func(x, **kwargs)[source]
class ocpy.newtonian.NewtonianOp(model)[source]

Bases: Op

make_node(x, *args)[source]

Construct an Apply node that represent the application of this operation to the given inputs.

This must be implemented by sub-classes.

Returns:

node – The constructed Apply node.

Return type:

Apply

perform(node, inputs, outputs)[source]

Calculate the function on the inputs and put the variables in the output storage.

Parameters:
  • node – The symbolic Apply node that represents this computation.

  • inputs – Immutable sequence of non-symbolic/numeric inputs. These are the values of each Variable in node.inputs.

  • output_storage – List of mutable single-element lists (do not change the length of these lists). Each sub-list corresponds to value of each Variable in node.outputs. The primary purpose of this method is to set the values of these sub-lists.

Notes

The output_storage list might contain data. If an element of output_storage is not None, it has to be of the right type, for instance, for a TensorVariable, it has to be a NumPy ndarray with the right number of dimensions and the correct dtype. Its shape and stride pattern can be arbitrary. It is not guaranteed that such pre-set values were produced by a previous call to this Op.perform(); they could’ve been allocated by another Op’s perform method. An Op is free to reuse output_storage as it sees fit, or to discard it and allocate new memory.

grad(inputs, output_gradients)[source]

Construct a graph for the gradient with respect to each input variable.

Each returned Variable represents the gradient with respect to that input computed based on the symbolic gradients with respect to each output. If the output is not differentiable with respect to an input, then this method should return an instance of type NullType for that input.

Using the reverse-mode AD characterization given in [1], for a \(C = f(A, B)\) representing the function implemented by the Op and its two arguments \(A\) and \(B\), given by the Variables in inputs, the values returned by Op.grad represent the quantities \(\bar{A} \equiv \frac{\partial S_O}{A}\) and \(\bar{B}\), for some scalar output term \(S_O\) of \(C\) in

\[\operatorname{Tr}\left(\bar{C}^\top dC\right) = \operatorname{Tr}\left(\bar{A}^\top dA\right) + \operatorname{Tr}\left(\bar{B}^\top dB\right)\]
Parameters:
  • inputs – The input variables.

  • output_grads – The gradients of the output variables.

Returns:

The gradients with respect to each Variable in inputs.

Return type:

grads

References

ocpy.oc module

class ocpy.oc.Parameter(value=None, min=None, max=None, std=None, fixed=False, distribution='truncatednormal')[source]

Bases: ParameterModel

Represents a model parameter with optional bounds, uncertainty, and distribution.

Parameters:
  • value (float, optional) – The current or initial value of the parameter.

  • min (float, optional) – Lower bound for the parameter, if applicable.

  • max (float, optional) – Upper bound for the parameter, if applicable.

  • std (float, optional) – Standard deviation of the parameter, used for uncertainty.

  • fixed (bool, default False) – Whether the parameter is fixed (not free to vary) during fitting.

  • distribution (str, default 'truncatednormal') – The assumed probability distribution for Bayesian sampling (e.g., ‘truncatednormal’, ‘uniform’, ‘normal’).

Notes

  • This class inherits from ParameterModel.

  • Can be used to define parameters for any ModelComponent.

  • Supports integration with fitting routines and Bayesian inference.

Examples

Define a free parameter with initial value 1.0:

>>> p = Parameter(value=1.0)

Define a bounded parameter with uncertainty:

>>> p = Parameter(value=0.5, min=0.0, max=1.0, std=0.1)
value: float | None = None
min: float | None = None
max: float | None = None
std: float | None = None
fixed: bool | None = False
distribution: str = 'truncatednormal'
class ocpy.oc.ModelComponent(args=None)[source]

Bases: ModelComponentModel

Base class for a mathematical model component with named parameters.

This class provides the foundation for any model component used in O–C analysis or other curve-fitting tasks. It manages parameters, allows updating from inference data, and supports switching between numerical backends such as NumPy or PyMC.

params

A dictionary of parameter names to Parameter objects. These parameters define the model and can be updated during fitting.

Type:

dict of str -> Parameter

math_class

Numerical backend used for mathematical operations (default: numpy).

Type:

module

_atan2

Function used for two-argument arctangent, adjusted to the current math backend.

Type:

callable

set_math(mathmod)[source]

Switch the math backend for calculations (e.g., NumPy or PyMC).

model_function()[source]

Return the underlying model function.

update_parameters(params_dict)[source]

Update parameter values from a dictionary of {name: value}.

update_from_idata(inference_data, group='posterior', stat='median')[source]

Update parameters from a PyMC InferenceData object.

_param(v)[source]

Convert a number or existing Parameter into a Parameter object.

Notes

  • Intended to be subclassed by concrete model components such as Linear, Sinusoidal, or Keplerian.

  • Supports both classical and Bayesian fitting workflows.

  • Provides compatibility with different computational backends, including NumPy and PyMC (if available).

  • Parameter values can be updated manually or via inference results.

Examples

Create a simple linear model component:

>>> linear = Linear(a=1.0, b=0.0)
>>> linear.params
{'a': Parameter(value=1.0), 'b': Parameter(value=0.0)}

Switch to a PyMC math backend (if PyMC installed):

>>> linear.set_math(pm.math)

Update parameters manually:

>>> linear.update_parameters({'a': 2.0})
>>> linear.params['a'].value
2.0
params: Dict[str, Parameter]
math_class = <module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/oc-py/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>
set_math(mathmod)[source]

Set the mathematical backend for model computations.

This method allows switching the numerical backend used in all mathematical operations of the component, including functions like sin, cos, sqrt, and arctan2. Supports NumPy by default and PyMC if available.

Parameters:

mathmod (module) – The module providing mathematical functions (e.g., numpy, pymc.math, or pytensor.tensor).

Returns:

self – Returns the instance itself to allow method chaining.

Return type:

ModelComponent

Notes

  • The _atan2 function is updated to match the chosen backend.

  • If PyMC is available and mathmod corresponds to PyMC’s math module, _atan2 uses pm.math.arctan2 or pytensor.tensor.arctan2.

  • If the backend does not provide arctan2, it falls back to numpy.arctan2.

  • This is useful for switching between classical numerical computation and probabilistic programming frameworks.

Examples

Use NumPy as the backend (default):

>>> component.set_math(np)

Use PyMC backend if installed:

>>> import pymc as pm
>>> component.set_math(pm.math)
model_function()[source]

Return the underlying model function of the component.

The model function defines the mathematical relationship between the independent variable(s) and the output, based on the component’s parameters. It is intended to be used for evaluation, fitting, or plotting.

Returns:

The function implementing the model. The callable typically has the signature:

f(x, param1, param2, ...)

where x is the independent variable and the remaining arguments are the component’s parameters.

Return type:

callable

Notes

  • Subclasses must define the model_func method, which this property returns.

  • Useful for passing to optimization routines or plotting functions that require a callable.

Examples

Retrieve and call the model function of a linear component:

>>> linear = Linear(a=2.0, b=1.0)
>>> f = linear.model_function()
>>> f(3, linear.params['a'].value, linear.params['b'].value)
7.0
update_parameters(params_dict)[source]

Update the values of existing parameters from a dictionary.

Parameters:

params_dict (dict of str -> float) – Dictionary mapping parameter names to their new numeric values. Only parameters that already exist in the component will be updated; any unknown keys are ignored.

Returns:

self – Returns the instance itself to allow method chaining.

Return type:

ModelComponent

Notes

  • This method does not create new parameters; it only updates values of existing ones.

  • All values are cast to float before assignment.

  • Useful for manually adjusting parameters between fitting steps or after loading results from inference.

Examples

Update a linear model component’s parameters:

>>> linear = Linear(a=1.0, b=0.0)
>>> linear.update_parameters({'a': 2.5, 'b': -1.0})
>>> linear.params['a'].value
2.5
>>> linear.params['b'].value
-1.0

Attempting to update a non-existent parameter has no effect:

>>> linear.update_parameters({'c': 10.0})
>>> 'c' in linear.params
False
update_from_idata(inference_data, group='posterior', stat='median')[source]

Update model parameters from a PyMC InferenceData object.

This method extracts parameter values from a Bayesian posterior (or other group) and updates the component’s parameters accordingly. It automatically matches variables that are prefixed with the component’s name.

Parameters:
  • inference_data (arviz.InferenceData) – The InferenceData object containing posterior samples from a PyMC model.

  • group (str, default 'posterior') – The group within the InferenceData to use (e.g., ‘posterior’, ‘prior’).

  • stat (str, default 'median') – Statistic to use for updating parameters. Options include: - ‘median’: use the median of the posterior samples. - ‘mean’: use the mean of the posterior samples. - Any other value: uses the first sample as a fallback.

Returns:

self – The component instance with updated parameter values.

Return type:

ModelComponent

Notes

  • Only parameters present in self.params and matching the component’s prefix are updated.

  • Useful for propagating posterior estimates back into model components for further computation or plotting.

  • Requires PyMC and ArviZ to be installed.

Examples

Suppose idata contains posterior samples for a linear component named ‘linear’:

>>> linear = Linear(a=1.0, b=0.0)
>>> linear.update_from_idata(idata, group='posterior', stat='median')
>>> linear.params['a'].value  # updated from posterior median
2.3
>>> linear.params['b'].value
-0.1
class ocpy.oc.Linear(a=1.0, b=0.0, *, name=None)[source]

Bases: ModelComponent

Linear model component of the form f(x) = a * x + b.

Represents a simple linear relationship between an independent variable x and the output, parameterized by slope a and intercept b.

Parameters:
  • a (float or Parameter, default 1.0) – Slope of the linear model. Can be a numeric value or a Parameter object.

  • b (float or Parameter, default 0.0) – Intercept of the linear model. Can be a numeric value or a Parameter object.

  • name (str, optional) – Optional name for the component. Defaults to “linear”.

params

Dictionary of parameters: {‘a’: Parameter, ‘b’: Parameter}.

Type:

dict

name

Name of the component (default: “linear”).

Type:

str

model_func(x, a, b)[source]

Evaluate the linear model for input x using parameters a and b.

Examples

Create a linear component with default parameters:

>>> linear = Linear()
>>> linear.model_func(3, linear.params['a'].value, linear.params['b'].value)
3.0

Create a linear component with custom slope and intercept:

>>> linear = Linear(a=2.0, b=1.0)
>>> linear.model_func(3, linear.params['a'].value, linear.params['b'].value)
7.0
name = 'linear'
model_func(x, a, b)[source]

Evaluate the linear model f(x) = a * x + b.

Parameters:
  • x (float or array-like) – Independent variable(s) at which to evaluate the model.

  • a (float) – Slope of the linear function.

  • b (float) – Intercept of the linear function.

Returns:

The computed value(s) of the linear function at x.

Return type:

float or np.ndarray

Examples

Evaluate a linear component at a single point:

>>> linear = Linear(a=2.0, b=1.0)
>>> linear.model_func(3, linear.params['a'].value, linear.params['b'].value)
7.0

Evaluate at multiple points:

>>> import numpy as np
>>> x = np.array([0, 1, 2, 3])
>>> linear.model_func(x, linear.params['a'].value, linear.params['b'].value)
array([1., 3., 5., 7.])
class ocpy.oc.Quadratic(q=0.0, *, name=None)[source]

Bases: ModelComponent

Quadratic model component of the form f(x) = q * x^2.

Represents a simple quadratic relationship between an independent variable x and the output, parameterized by a single coefficient q.

Parameters:
  • q (float or Parameter, default 0.0) – Quadratic coefficient. Can be a numeric value or a Parameter object.

  • name (str, optional) – Optional name for the component. Defaults to “quadratic”.

params

Dictionary of parameters: {‘q’: Parameter}.

Type:

dict

name

Name of the component (default: “quadratic”).

Type:

str

model_func(x, q)[source]

Evaluate the quadratic model for input x using parameter q.

Examples

Create a quadratic component with default parameter:

>>> quad = Quadratic()
>>> quad.model_func(2, quad.params['q'].value)
0.0

Create a quadratic component with custom coefficient:

>>> quad = Quadratic(q=3.0)
>>> quad.model_func(2, quad.params['q'].value)
12.0
name = 'quadratic'
model_func(x, q)[source]

Evaluate the quadratic model f(x) = q * x^2.

Parameters:
  • x (float or array-like) – Independent variable(s) at which to evaluate the model.

  • q (float) – Coefficient of the quadratic term.

Returns:

The computed value(s) of the quadratic function at x.

Return type:

float or np.ndarray

Examples

Evaluate a quadratic component at a single point:

>>> quad = Quadratic(q=2.0)
>>> quad.model_func(3, quad.params['q'].value)
18.0

Evaluate at multiple points:

>>> import numpy as np
>>> x = np.array([0, 1, 2, 3])
>>> quad.model_func(x, quad.params['q'].value)
array([ 0.,  2.,  8., 18.])
class ocpy.oc.Sinusoidal(*, amp=None, P=None, name=None)[source]

Bases: ModelComponent

Sinusoidal model component of the form f(x) = amp * sin(2π * x / P).

Represents a simple periodic function with amplitude amp and period P.

Parameters:
  • amp (float or Parameter, optional) – Amplitude of the sinusoidal function. Can be a numeric value or a Parameter object.

  • P (float or Parameter, optional) – Period of the sinusoidal function. Can be a numeric value or a Parameter object.

  • name (str, optional) – Optional name for the component. Defaults to “sinusoidal”.

params

Dictionary of parameters: {‘amp’: Parameter, ‘P’: Parameter}.

Type:

dict

name

Name of the component (default: “sinusoidal”).

Type:

str

model_func(x, amp, P)[source]

Evaluate the sinusoidal model for input x using parameters amp and P.

Examples

Create a sinusoidal component with amplitude 1.0 and period 2.0:

>>> sinus = Sinusoidal(amp=1.0, P=2.0)
>>> sinus.model_func(0, sinus.params['amp'].value, sinus.params['P'].value)
0.0

Evaluate at multiple points:

>>> import numpy as np
>>> x = np.array([0, 0.5, 1.0, 1.5])
>>> sinus.model_func(x, sinus.params['amp'].value, sinus.params['P'].value)
array([0. , 1. , 0. , -1.])
name = 'sinusoidal'
model_func(x, amp, P)[source]

Evaluate the sinusoidal model f(x) = amp * sin(2π * x / P).

Parameters:
  • x (float or array-like) – Independent variable(s) at which to evaluate the model.

  • amp (float) – Amplitude of the sinusoidal function.

  • P (float) – Period of the sinusoidal function.

Returns:

The computed value(s) of the sinusoidal function at x.

Return type:

float or np.ndarray

Examples

Evaluate at a single point:

>>> sinus = Sinusoidal(amp=2.0, P=4.0)
>>> sinus.model_func(1, sinus.params['amp'].value, sinus.params['P'].value)
2.0

Evaluate at multiple points:

>>> import numpy as np
>>> x = np.array([0, 1, 2, 3, 4])
>>> sinus.model_func(x, sinus.params['amp'].value, sinus.params['P'].value)
array([ 0.,  2.,  0., -2., -0.])
class ocpy.oc.Keplerian(*, amp=None, e=0.0, omega=0.0, P=None, T0=None, name=None)[source]

Bases: ModelComponent

Keplerian orbital model component for radial velocity variations.

Represents a Keplerian orbit with parameters corresponding to the orbital amplitude, eccentricity, argument of periastron, orbital period, and time of periastron passage.

The model computes radial velocity or timing variations as a function of time based on classical Keplerian motion.

Parameters:
  • amp (float or Parameter, optional) – Amplitude of the Keplerian signal.

  • e (float or Parameter, default 0.0) – Orbital eccentricity. Must be between 0 and 1.

  • omega (float or Parameter, default 0.0) – Argument of periastron in degrees.

  • P (float or Parameter, optional) – Orbital period.

  • T0 (float or Parameter, optional) – Time of periastron passage.

  • name (str, optional) – Optional name for the component. Defaults to “keplerian”.

params

Dictionary of parameters: {‘amp’, ‘e’, ‘omega’, ‘P’, ‘T0’}.

Type:

dict

name

Name of the component (default: “keplerian”).

Type:

str

model_func(x, amp, e, omega, P, T0)[source]

Evaluate the Keplerian model at times x.

_kepler_solve(M, e, n_iter=5)[source]

Solve Kepler’s equation using iterative method.

Notes

  • The input angles are in degrees and internally converted to radians.

  • Uses Newton-Raphson iteration to solve Kepler’s equation.

  • Suitable for modeling radial velocity or timing variations in binary stars or exoplanets.

Examples

Create a Keplerian component with default parameters:

>>> kepler = Keplerian()
>>> kepler.params['e'].value
0.0

Custom parameters:

>>> kepler = Keplerian(amp=2.0, e=0.5, omega=45.0, P=10.0, T0=0.0)
>>> kepler.params['amp'].value
2.0
>>> kepler.params['P'].value
10.0
name = 'keplerian'
model_func(x, amp, e, omega, P, T0)[source]

Evaluate the Keplerian model at given times.

Computes the Keplerian signal (e.g., radial velocity or timing variations) based on orbital parameters: amplitude, eccentricity, argument of periastron, period, and time of periastron passage.

Parameters:
  • x (float or array-like) – Times at which to evaluate the model.

  • amp (float) – Amplitude of the Keplerian signal.

  • e (float) – Orbital eccentricity (0 ≤ e < 1).

  • omega (float) – Argument of periastron in degrees.

  • P (float) – Orbital period.

  • T0 (float) – Time of periastron passage.

Returns:

Keplerian model values at the specified times.

Return type:

float or np.ndarray

Notes

  • The argument of periastron omega is internally converted to radians.

  • Uses _kepler_solve to determine the eccentric anomaly.

  • Handles both scalar and array inputs for x.

  • Suitable for modeling timing variations in eclipsing binaries or radial velocity signals.

Examples

Evaluate at a single time:

>>> kepler = Keplerian(amp=1.0, e=0.1, omega=30.0, P=10.0, T0=0.0)
>>> kepler.model_func(5, 1.0, 0.1, 30.0, 10.0, 0.0)
0.85  # Example value (approximate)

Evaluate at multiple times:

>>> import numpy as np
>>> times = np.linspace(0, 10, 5)
>>> kepler.model_func(times, 1.0, 0.1, 30.0, 10.0, 0.0)
array([0.0, 0.52, 0.85, 0.85, 0.52])
class ocpy.oc.KeplerianOld(*, amp=None, e=0.0, omega=0.0, P=None, T0=None, name=None)[source]

Bases: ModelComponent

name = 'keplerian'
model_func(x, amp, e, omega, P, T0)[source]
class ocpy.oc.OC(oc, minimum_time=None, minimum_time_error=None, weights=None, minimum_type=None, labels=None, cycle=None)[source]

Bases: OCModel

Observed-minus-Calculated (O–C) data container and analysis class.

Represents timing variations of observed events (e.g., eclipses, pulsations) relative to a reference ephemeris. Provides methods for binning, merging, computing O–C values, and fitting models to timing variations.

data

Internal DataFrame storing columns: ‘minimum_time’, ‘minimum_time_error’, ‘weights’, ‘minimum_type’, ‘labels’, ‘cycle’, ‘oc’.

Type:

pd.DataFrame

from_file(file, columns=None)[source]

Create an OC instance from a CSV or Excel file.

__getitem__(item)[source]

Access a column, row, or filtered OC subset.

__setitem__(key, value)[source]

Assign values to a column.

__len__()[source]

Return the number of entries.

bin(bin_count=1, bin_method=None, bin_error_method=None, bin_style=None)[source]

Bin O–C data using weighted averages.

merge(oc)[source]

Merge with another OC instance.

calculate_oc(reference_minimum, reference_period, model_type='lmfit')[source]

Compute O–C values relative to a reference ephemeris.

plot(model=None, ax=None, res_ax=None, res=True, ...)[source]

Visualize O–C data and optional model fits.

Notes

  • Supports both scalar and array-like input for all columns.

  • cycle is automatically computed or can be provided explicitly.

  • Binning can be performed with custom reducers and bin styles.

  • Integration with ModelComponent subclasses enables model fitting.

Examples

Create an OC instance from raw data:

>>> oc = OC(oc=[0.1, -0.05, 0.0],
...         minimum_time=[2450000.0, 2450001.0, 2450002.0],
...         weights=[1.0, 1.0, 1.0])

Access O–C column:

>>> oc['oc']
0    0.10
1   -0.05
2    0.00
Name: oc, dtype: float64

Bin O–C data:

>>> binned = oc.bin(bin_count=2)

Merge with another OC:

>>> merged = oc.merge(other_oc)
classmethod from_file(file, columns=None)[source]

Load O–C data from a CSV or Excel file into an OC instance.

Parameters:
  • file (str or Path) – Path to the input file. Supported formats: .csv, .xls, .xlsx.

  • columns (dict, optional) – Mapping of file column names to OC attribute names. Keys are original column names in the file; values are the corresponding attribute names (minimum_time, oc, weights, cycle, etc.).

Returns:

An OC instance containing data loaded from the file.

Return type:

OC

Raises:

ValueError – If the file type is unsupported.

Notes

  • Any columns not present in the file are filled with None.

  • Column renaming allows flexibility for files with nonstandard column names.

  • The method automatically handles array-like columns and converts them for internal storage in a pandas.DataFrame.

Examples

Load a CSV file with standard columns:

>>> oc = OC.from_file("observations.csv")

Load an Excel file with custom column names:

>>> oc = OC.from_file("observations.xlsx",
...                   columns={"Time": "minimum_time",
...                            "O-C": "oc",
...                            "Weight": "weights"})
bin(bin_count=1, bin_method=None, bin_error_method=None, bin_style=None)[source]

Bin the O–C data along the cycle axis and compute weighted averages and errors.

Parameters:
  • bin_count (int, optional) – Number of bins to create. Default is 1.

  • bin_method (callable, optional) – Function to compute binned values from array and weights: func(array: np.ndarray, weights: np.ndarray) -> float. Default is weighted mean.

  • bin_error_method (callable, optional) – Function to compute binned errors from weights: func(weights: np.ndarray) -> float. Default is 1 / sqrt(sum(weights)).

  • bin_style (callable, optional) – Function to define custom bin edges: func(df: pd.DataFrame, bin_count: int) -> np.ndarray of shape (N, 2).

Returns:

New OC instance containing binned data.

Return type:

OC

Raises:

ValueError – If cycle, oc, or weights columns are missing or contain NaNs. If bin_count is not positive.

Notes

  • Uses cycle as the x-axis and oc as the y-axis.

  • By default, bins are equal-width and weighted averages are used.

  • Supports custom binning logic via bin_style.

  • Binned rows have:
    • labels set to “Binned”

    • minimum_time and weights set to NaN

    • minimum_type set to None

Examples

Bin the data into 5 equal-width bins:

>>> binned_oc = oc.bin(bin_count=5)

Bin using a custom error function and custom bin edges:

>>> def custom_bin_style(df, n_bins):
...     return np.array([[0, 2], [2, 5], [5, 10]])
>>> binned_oc = oc.bin(bin_count=3, bin_error_method=lambda w: np.sqrt(np.sum(w)), bin_style=custom_bin_style)
merge(oc)[source]

Merge another OC instance into the current one.

Parameters:

oc (OC) – Another OC instance whose data will be appended.

Returns:

A new OC instance containing the combined data of both instances.

Return type:

OC

Notes

  • The merge is performed row-wise, preserving all columns.

  • Indexes are reset in the resulting DataFrame.

  • The original OC instances are not modified; the merge returns a new instance.

Examples

Merge two OC datasets:

>>> oc1 = OC.from_file("observations1.csv")
>>> oc2 = OC.from_file("observations2.csv")
>>> merged = oc1.merge(oc2)
>>> len(merged) == len(oc1) + len(oc2)
True
calculate_oc(reference_minimum, reference_period, model_type='lmfit')[source]

Compute O–C (Observed minus Calculated) values based on a reference ephemeris.

Parameters:
  • reference_minimum (float) – Reference time of a primary minimum (e.g., first observed eclipse).

  • reference_period (float) – Reference period for the cycle calculation.

  • model_type (str, optional) – Type of model to return: - “lmfit” : returns an OCLMFit instance if available. - other : returns an OC instance. Default is “lmfit”.

Returns:

New OC (or OCLMFit) instance with calculated cycle and oc columns.

Return type:

OC

Raises:

ValueError – If the minimum_time column is missing. If minimum_type contains unrecognized values for secondary/primary distinction.

Notes

  • The method computes the cycle number as (t - reference_minimum) / reference_period, rounded to the nearest integer.

  • If minimum_type is present, secondary minima (e.g., “II”, “sec”) are offset by 0.5 cycles.

  • Calculated O–C values are added as the oc column; cycle numbers as cycle.

  • Other columns (minimum_time_error, weights, labels, etc.) are preserved.

Examples

Compute O–C for an OC dataset:

>>> oc = OC.from_file("observations.csv")
>>> oc_calc = oc.calculate_oc(reference_minimum=2450000.0, reference_period=1.2345)
>>> oc_calc.data[['cycle', 'oc']].head()
   cycle     oc
0    0.0  0.012
1    1.0 -0.005
2    2.0  0.003
residue(coefficients)[source]

Compute the residuals of the O–C data relative to a fitted model.

Parameters:

coefficients (ModelResult) – Fitted model object (from lmfit) containing the optimized parameters used to compute model-predicted O–C values.

Returns:

New OC instance containing the residuals (observed - model) in the oc column. Other columns are preserved from the original instance.

Return type:

OC

Notes

  • This method does not modify the original OC instance; a new instance is returned.

  • Residuals are computed using the model defined by coefficients applied to the cycle column.

  • Typically used after OC.fit or similar model fitting methods.

Examples

Compute residuals after fitting a linear model:

>>> result = oc.fit_linear(a=1.0, b=0.0)  # returns a ModelResult
>>> oc_resid = oc.residue(result)
>>> oc_resid.data[['cycle', 'oc']].head()
   cycle       oc
0    0.0   0.002
1    1.0  -0.001
2    2.0   0.003
fit(functions)[source]

Fit one or more model components to the O–C data.

Parameters:

functions (ModelComponentModel or list of ModelComponentModel) – Single or multiple model components (e.g., Linear, Quadratic, Sinusoidal, Keplerian) to fit to the oc data. Each component defines its own functional form and parameters.

Returns:

Fitted model result (from lmfit) containing optimized parameters, fit statistics, and the model-predicted O–C values.

Return type:

ModelResult

Raises:
  • ValueError – If the required columns (cycle, oc, or weights) are missing or contain NaNs.

  • TypeError – If functions is not a ModelComponentModel or list thereof.

Notes

  • Supports combining multiple model components into a single composite fit.

  • Uses weighted fitting based on the weights column.

  • The resulting ModelResult can be used to evaluate residuals or plot the fitted model.

  • Internally, converts functions to a composite lmfit model if multiple components are provided.

Examples

Fit a linear and sinusoidal model to O–C data:

>>> linear = Linear(a=0.0, b=0.0)
>>> sinus = Sinusoidal(amp=0.01, P=10.0)
>>> result = oc.fit([linear, sinus])
>>> result.best_values
{'a': 0.0012, 'b': -0.0023, 'amp': 0.0098, 'P': 9.97}
fit_keplerian(*, amp=None, e=None, omega=None, P=None, T=None)[source]

Fit a Keplerian (elliptical orbit) model to the O–C data.

Parameters:
  • amp (ParameterModel, optional) – Amplitude of the Keplerian signal.

  • e (ParameterModel, optional) – Orbital eccentricity.

  • omega (ParameterModel, optional) – Argument of periastron in degrees.

  • P (ParameterModel, optional) – Orbital period.

  • T (ParameterModel, optional) – Time of periastron passage.

Returns:

Fitted Keplerian model component with optimized parameters.

Return type:

ModelComponentModel

Notes

  • Uses the Keplerian formula to model periodic variations in O–C data.

  • Supports optional ParameterModel objects for each orbital parameter; unspecified parameters can be estimated from the data.

  • The returned model component can be used to evaluate predicted O–C values or for further combination with other model components.

Examples

Fit a Keplerian model with initial guesses:

>>> kepler_model = oc.fit_keplerian(
...     amp=Parameter(value=0.01),
...     e=Parameter(value=0.1),
...     omega=Parameter(value=90),
...     P=Parameter(value=1000),
...     T=Parameter(value=2450000)
... )
>>> kepler_model.params['amp'].value
0.0098
fit_lite(*, amp=None, e=None, omega=None, P=None, T=None)[source]

Fit a simplified Keplerian model to the O–C data.

Parameters:
  • amp (ParameterModel, optional) – Amplitude of the Keplerian signal.

  • e (ParameterModel, optional) – Orbital eccentricity.

  • omega (ParameterModel, optional) – Argument of periastron in degrees.

  • P (ParameterModel, optional) – Orbital period.

  • T (ParameterModel, optional) – Time of periastron passage.

Returns:

Fitted simplified Keplerian model component with optimized parameters.

Return type:

ModelComponentModel

Notes

  • This method provides a lighter or approximate version of a full Keplerian fit.

  • Useful for initial estimates or for datasets where full Keplerian fitting is unnecessary.

  • Supports optional ParameterModel objects for each orbital parameter; unspecified parameters can be estimated from the data.

  • The returned model component can be evaluated for predicted O–C values or combined with other model components.

Examples

Fit a simplified Keplerian model:

>>> kepler_lite = oc.fit_lite(
...     amp=Parameter(value=0.01),
...     e=Parameter(value=0.05),
...     P=Parameter(value=1000),
... )
>>> kepler_lite.params['P'].value
1000.0
fit_linear(*, a=None, b=None)[source]

Fit a linear model to the O–C data.

Parameters:
  • a (ParameterModel, optional) – Slope of the linear model. If not provided, it will be estimated from the data.

  • b (ParameterModel, optional) – Intercept of the linear model. If not provided, it will be estimated from the data.

Returns:

Fitted linear model component with optimized parameters.

Return type:

ModelComponentModel

Notes

  • The linear model has the form oc = a * cycle + b.

  • Optional ParameterModel objects can be used to fix or initialize parameters.

  • The returned model component can be used to compute predicted O–C values or combined with other model components for composite fitting.

Examples

Fit a linear trend to O–C data:

>>> linear_model = oc.fit_linear(a=Parameter(value=0.001), b=Parameter(value=-0.002))
>>> linear_model.params['a'].value
0.001
>>> linear_model.params['b'].value
-0.002
fit_quadratic(*, q=None)[source]

Fit a quadratic model to the O–C data.

Parameters:

q (ParameterModel, optional) – Quadratic coefficient. If not provided, it will be estimated from the data.

Returns:

Fitted quadratic model component with optimized parameter(s).

Return type:

ModelComponentModel

Notes

  • The quadratic model has the form oc = q * cycle^2.

  • Optional ParameterModel can be used to fix or initialize the quadratic coefficient.

  • The returned model component can be used to compute predicted O–C values or combined with other model components for composite fitting.

Examples

Fit a quadratic trend to O–C data:

>>> quad_model = oc.fit_quadratic(q=Parameter(value=0.0001))
>>> quad_model.params['q'].value
0.0001
fit_sinusoidal(*, amp=None, P=None)[source]

Fit a sinusoidal model to the O–C data.

Parameters:
  • amp (ParameterModel, optional) – Amplitude of the sinusoidal signal. If not provided, it will be estimated from the data.

  • P (ParameterModel, optional) – Period of the sinusoidal signal. If not provided, it will be estimated from the data.

Returns:

Fitted sinusoidal model component with optimized parameters.

Return type:

ModelComponentModel

Notes

  • The sinusoidal model has the form oc = amp * sin(2 * pi * cycle / P).

  • Optional ParameterModel objects can be used to fix or initialize parameters.

  • The returned model component can be evaluated for predicted O–C values or combined with other model components in a composite fit.

Examples

Fit a sinusoidal variation to O–C data:

>>> sinus_model = oc.fit_sinusoidal(amp=Parameter(value=0.01), P=Parameter(value=1000))
>>> sinus_model.params['amp'].value
0.01
>>> sinus_model.params['P'].value
1000
fit_parabola(*, q=None, a=None, b=None)[source]

Fit a parabolic model to the O–C data.

Parameters:
  • q (ParameterModel, optional) – Quadratic coefficient. If not provided, it will be estimated from the data.

  • a (ParameterModel, optional) – Linear coefficient. If not provided, it will be estimated from the data.

  • b (ParameterModel, optional) – Constant term (intercept). If not provided, it will be estimated from the data.

Returns:

Fitted parabolic model component with optimized parameters.

Return type:

ModelComponentModel

Notes

  • The parabolic model has the form oc = q * cycle^2 + a * cycle + b.

  • Optional ParameterModel objects can be used to fix or initialize coefficients.

  • The returned model component can be used to compute predicted O–C values or combined with other model components for composite fitting.

Examples

Fit a parabolic trend to O–C data:

>>> parab_model = oc.fit_parabola(q=Parameter(value=0.0001),
...                               a=Parameter(value=0.001),
...                               b=Parameter(value=-0.002))
>>> parab_model.params['q'].value
0.0001
>>> parab_model.params['a'].value
0.001
>>> parab_model.params['b'].value
-0.002
plot(model=None, *, ax=None, res_ax=None, res=True, title=None, x_col='cycle', y_col='oc', fig_size=(10, 7), plot_kwargs=None, extension_factor=0.05, model_components=None)[source]

Plot the O–C data along with optional model predictions and residuals.

Parameters:
  • model (InferenceData, ModelResult, or list of ModelComponent, optional) – Model or model components to overlay on the O–C plot.

  • ax (matplotlib.axes.Axes, optional) – Axes object for the main plot. If None, a new figure and axes are created.

  • res_ax (matplotlib.axes.Axes, optional) – Axes object for residuals plot. Only used if res=True.

  • res (bool, default True) – If True, show residuals between the data and model on a separate subplot.

  • title (str, optional) – Title for the plot.

  • x_col (str, default "cycle") – Column name in the data to use for the x-axis.

  • y_col (str, default "oc") – Column name in the data to use for the y-axis.

  • fig_size (tuple, default (10, 7)) – Figure size if a new figure is created.

  • plot_kwargs (dict, optional) – Additional keyword arguments passed to the plotting function (e.g., marker, linestyle, color).

  • extension_factor (float, default 0.05) – Fractional extension for axis limits to improve visual spacing.

  • model_components (list, optional) – List of individual model components to display separately on the plot.

Returns:

If res=False, returns the main Axes object. If res=True, returns a dictionary with keys ‘main_ax’ and ‘res_ax’ containing the corresponding Axes objects.

Return type:

matplotlib.axes.Axes or dict

Notes

  • This method uses the Plot.plot function from the visualization module.

  • Supports overlaying multiple model components or posterior predictive samples.

  • Automatically handles binning or weighting if present in the data.

  • Residuals plot shows data minus model values for quick assessment of fit quality.

Examples

Plot O–C data with a fitted sinusoidal model:

>>> sinus_model = oc.fit_sinusoidal(amp=Parameter(value=0.01), P=Parameter(value=1000))
>>> oc.plot(model=sinus_model, title="O–C with Sinusoidal Fit")

Plot data with residuals on a separate subplot:

>>> oc.plot(model=sinus_model, res=True, fig_size=(12, 8))

ocpy.oc_lmfit module

class ocpy.oc_lmfit.OCLMFit(oc, minimum_time=None, minimum_time_error=None, weights=None, minimum_type=None, labels=None, cycle=None)[source]

Bases: OC

O–C data handler with fitting capabilities using lmfit.

This class extends OC by providing methods to fit O–C data to linear, quadratic, sinusoidal, Keplerian, and combined models using lmfit.Model. It supports parameter constraints, fixed parameters, and custom initial values. Fitted models can be used to compute residuals or generate new O–C predictions.

math

Mathematical module used in component evaluations. Defaults to numpy.

Type:

module

Notes

  • Each fit method (fit_linear, fit_quadratic, fit_sinusoidal, etc.) returns an lmfit.ModelResult with the best-fit parameters.

  • Parameters can be passed as Parameter instances, floats, or None. None defaults are used.

  • Residuals can be computed using residue() to generate a new OCLMFit object.

  • Designed for O–C analysis of periodic phenomena in astronomy, such as eclipsing binaries.

Examples

Fit a linear model to the O–C data:

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_linear(a=0.0, b=0.0)

Fit a Keplerian model with custom initial parameters:

>>> result = oc_fit.fit_keplerian(amp=0.01, e=0.1, P=1000.0, T0=2450000.0)

Compute residuals after fitting:

>>> oc_resid = oc_fit.residue(result)
math = <module 'numpy' from '/home/docs/checkouts/readthedocs.org/user_builds/oc-py/envs/latest/lib/python3.12/site-packages/numpy/__init__.py'>
fit(model_components, *, nan_policy='raise', method='leastsq', **kwargs)[source]

Fit the O–C data using one or more model components via lmfit.

This method combines one or more ModelComponentModel instances into a single composite lmfit model, applies parameter constraints, and fits the O–C data (cycle vs oc) according to the provided weights.

Parameters:
  • model_components (list of ModelComponentModel) – A list of model components (e.g., Linear, Quadratic, Keplerian, Sinusoidal) to fit simultaneously.

  • nan_policy ({'raise', 'omit', 'propagate'}, optional) – How to handle NaN values in the data. Default is ‘raise’.

  • method (str, optional) – The fitting algorithm to use. Passed directly to lmfit.Model.fit. Default is ‘leastsq’.

  • **kwargs – Additional keyword arguments forwarded to lmfit.Model.fit.

Returns:

The lmfit.ModelResult object containing best-fit parameters, uncertainties, and fit statistics.

Return type:

ModelResult

Raises:
  • ValueError – If the weights column contains NaN values.

  • TypeError – If any model component is not compatible with lmfit.

Notes

  • Each component can have initial parameter values, bounds (min, max), and fixed status.

  • The method automatically generates unique prefixes for each component to avoid parameter name collisions.

  • Fitted parameters are applied to the respective ModelComponentModel instances internally.

  • Supports multiple components of the same type with unique prefixes.

Examples

Fit a linear and sinusoidal model simultaneously:

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> linear_comp = Linear(a=0.0, b=0.0)
>>> sinusoidal_comp = Sinusoidal(amp=0.01, P=1000)
>>> result = oc_fit.fit([linear_comp, sinusoidal_comp], method='leastsq')

Access fitted parameters:

>>> result.params['linear_a'].value
>>> result.params['sinusoidal_amp'].value
residue(coefficients, *, x_col='cycle', y_col='oc')[source]

Compute the residual O–C values after subtracting a fitted model.

This method evaluates the fitted model at the observed cycles and subtracts the model values from the observed O–C data, returning a new OCLMFit instance containing only the residuals.

Parameters:
  • coefficients (ModelResult) – The fitted lmfit.ModelResult object containing the best-fit parameters for the model.

  • x_col (str, optional) – Name of the column in self.data to use as the independent variable. Default is ‘cycle’.

  • y_col (str, optional) – Name of the column in self.data to use as the dependent variable (O–C values). Default is ‘oc’.

Returns:

A new OCLMFit instance with the same metadata as the original object, but with the oc column replaced by residuals (observed minus fitted values).

Return type:

OCLMFit

Raises:

ValueError – If x_col or y_col is not present in the data.

Notes

  • The returned object preserves all other columns (minimum_time, weights, labels, etc.) from the original dataset.

  • Useful for iterative fitting, diagnostics, or residual analysis.

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_linear(a=0.0, b=0.0)
>>> residuals = oc_fit.residue(result)
>>> residuals.data.head()
fit_linear(*, a=None, b=None, **kwargs)[source]

Fit a linear model to the O–C data using least-squares optimization.

The model has the form:

y = a * x + b

where x is the cycle number and y is the O–C value. Both slope (a) and intercept (b) can be provided as fixed values, free parameters, or left to be estimated from the data.

Parameters:
  • a (Parameter, float, or None, optional) – Initial guess or fixed value for the slope. If None, defaults to 0.0.

  • b (Parameter, float, or None, optional) – Initial guess or fixed value for the intercept. If None, defaults to 0.0.

  • **kwargs – Additional keyword arguments are passed to the underlying fit method, e.g., method, nan_policy, etc.

Returns:

An lmfit.ModelResult object containing the optimized parameters, fit statistics, and methods to evaluate the model or extract residuals.

Return type:

ModelResult

Raises:

ValueError – If the O–C data contains NaN values in the weights column, which are required for weighted fitting.

Notes

  • This method internally wraps the a and b values as Parameter objects if they are not already, allowing fixed or free fitting.

  • Useful for detecting long-term trends in O–C diagrams.

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_linear(a=0.0, b=0.0)
>>> print(result.best_values)
{'a': 0.0012, 'b': -0.03}
>>> residuals = oc_fit.residue(result)
fit_quadratic(*, q=None, **kwargs)[source]

Fit a quadratic model to the O–C data using least-squares optimization.

The model has the form:

y = q * x^2

where x is the cycle number and y is the O–C value. The quadratic coefficient q can be provided as a fixed value, a free parameter, or left to be estimated from the data.

Parameters:
  • q (Parameter, float, or None, optional) – Initial guess or fixed value for the quadratic coefficient. If None, defaults to 0.0.

  • **kwargs – Additional keyword arguments are passed to the underlying fit method, e.g., method, nan_policy, etc.

Returns:

An lmfit.ModelResult object containing the optimized parameter q, fit statistics, and methods to evaluate the model or extract residuals.

Return type:

ModelResult

Raises:

ValueError – If the O–C data contains NaN values in the weights column, which are required for weighted fitting.

Notes

  • This method internally wraps q as a Parameter object if it is not already, allowing fixed or free fitting.

  • Useful for detecting parabolic trends in O–C diagrams, e.g., due to period changes.

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_quadratic(q=1e-6)
>>> print(result.best_values)
{'q': 2.3e-6}
>>> residuals = oc_fit.residue(result)
fit_parabola(*, q=None, a=None, b=None, **kwargs)[source]

Fit a combined quadratic and linear model to the O–C data using least-squares optimization.

The model has the form:

y = q * x^2 + a * x + b

where x is the cycle number and y is the O–C value. Each parameter (q, a, b) can be provided as a fixed value, a free parameter, or left to be estimated from the data.

Parameters:
  • q (Parameter, float, or None, optional) – Initial guess or fixed value for the quadratic coefficient. Defaults to 0.0 if None.

  • a (Parameter, float, or None, optional) – Initial guess or fixed value for the linear coefficient. Defaults to 0.0 if None.

  • b (Parameter, float, or None, optional) – Initial guess or fixed value for the constant term. Defaults to 0.0 if None.

  • **kwargs – Additional keyword arguments passed to the underlying fit method, e.g., method, nan_policy, weights, etc.

Returns:

An lmfit.ModelResult object containing optimized parameters q, a, b, fit statistics, and methods to evaluate the model or extract residuals.

Return type:

ModelResult

Raises:

ValueError – If the O–C data contains NaN values in the weights column.

Notes

  • Internally wraps all input values as Parameter objects for consistent handling.

  • Useful for modeling O–C diagrams that exhibit both parabolic trends (period changes) and linear trends (systematic offsets).

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_parabola(q=1e-6, a=1e-4, b=0.0)
>>> print(result.best_values)
{'q': 2.3e-6, 'a': 1.1e-4, 'b': 0.0}
>>> residuals = oc_fit.residue(result)
fit_lite(*, amp=None, e=None, omega=None, P=None, T0=None, **kwargs)[source]

Fit a simplified Keplerian model to O–C data using least-squares optimization.

This “lite” Keplerian model describes a sinusoidal O–C variation approximating the orbital motion of a binary system or planet, with the form:

y = amp * ( (1-e^2)/(1+e*cos(true_anomaly)) * sin(true_anomaly + omega)/denom + e*sin(omega)/denom )

where true_anomaly is calculated from the mean anomaly using Kepler’s equation.

Parameters:
  • amp (Parameter, float, or None, optional) – Semi-amplitude of the O–C variation. Defaults to 1e-3 if None.

  • e (Parameter, float, or None, optional) – Orbital eccentricity (0 ≤ e < 0.95). Defaults to 0.0 if None.

  • omega (Parameter, float, or None, optional) – Argument of periastron in degrees. Defaults to 90.0 if None.

  • P (Parameter, float, or None, optional) – Orbital period. Defaults to 3000.0 if None.

  • T0 (Parameter, float, or None, optional) – Reference time of periastron passage. Defaults to 0.0 if None.

  • **kwargs – Additional keyword arguments passed to the underlying fit method, such as method, nan_policy, or weights.

Returns:

An lmfit.ModelResult object containing optimized parameters amp, e, omega, P, T0, fit statistics, and methods for evaluation or residual extraction.

Return type:

ModelResult

Raises:

ValueError – If the O–C data contains NaN values in the weights column.

Notes

  • Parameters are internally wrapped as Parameter objects if they are not already.

  • Useful for rapid fitting of periodic variations in O–C diagrams without modeling additional complexities of full Keplerian motion.

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_lite(amp=0.01, e=0.1, P=5000.0)
>>> print(result.best_values)
{'amp': 0.012, 'e': 0.098, 'omega': 90.0, 'P': 4990.2, 'T0': 0.0}
>>> residuals = oc_fit.residue(result)
fit_keplerian(*, amp=None, e=None, omega=None, P=None, T0=None, **kwargs)[source]

Fit a Keplerian model to O–C data by delegating to the fit_lite method.

This method is a convenient alias for fitting standard Keplerian variations in O–C diagrams, modeling periodic deviations due to orbital motion. Internally, it calls fit_lite with the given parameters.

Parameters:
  • amp (Parameter, float, or None, optional) – Semi-amplitude of the O–C variation. Defaults to 1e-3 if None.

  • e (Parameter, float, or None, optional) – Orbital eccentricity (0 ≤ e < 0.95). Defaults to 0.0 if None.

  • omega (Parameter, float, or None, optional) – Argument of periastron in degrees. Defaults to 90.0 if None.

  • P (Parameter, float, or None, optional) – Orbital period. Defaults to 3000.0 if None.

  • T0 (Parameter, float, or None, optional) – Reference time of periastron passage. Defaults to 0.0 if None.

  • **kwargs – Additional keyword arguments passed to fit_lite and ultimately to OCLMFit.fit, such as method, nan_policy, or weights.

Returns:

An lmfit.ModelResult object containing optimized parameters amp, e, omega, P, T0, fit statistics, and methods for evaluation or residual extraction.

Return type:

ModelResult

Raises:

ValueError – If the O–C data contains NaN values in the weights column (raised in fit_lite).

Notes

  • This method provides semantic clarity for users who want a Keplerian fit, but it is functionally identical to fit_lite.

  • Parameters are internally wrapped as Parameter objects if they are not already.

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_keplerian(amp=0.01, e=0.1, P=5000.0)
>>> print(result.best_values)
{'amp': 0.012, 'e': 0.098, 'omega': 90.0, 'P': 4990.2, 'T0': 0.0}
>>> residuals = oc_fit.residue(result)
fit_sinusoidal(*, amp=None, P=None, **kwargs)[source]

Fit a sinusoidal model to O–C data.

This method fits a simple sinusoidal variation to the observed O–C values, which can model periodic effects such as light-time effects in binary systems.

Parameters:
  • amp (Parameter, float, or None, optional) – Amplitude of the sinusoidal O–C variation. Defaults to 1e-3 if None.

  • P (Parameter, float, or None, optional) – Period of the sinusoidal variation. Defaults to 3000.0 if None.

  • **kwargs – Additional keyword arguments passed to OCLMFit.fit, such as method, nan_policy, or weights.

Returns:

An lmfit.ModelResult object containing optimized parameters amp and P, fit statistics, and methods for evaluation or residual extraction.

Return type:

ModelResult

Raises:

ValueError – If the O–C data contains NaN values in the weights column.

Notes

  • This method internally wraps numerical values into Parameter objects if they are not already.

  • Useful for detecting or modeling sinusoidal periodicities in O–C diagrams.

Examples

>>> oc_fit = OCLMFit(oc, minimum_time=times, cycle=cycles, weights=weights)
>>> result = oc_fit.fit_sinusoidal(amp=0.01, P=5000.0)
>>> print(result.best_values)
{'amp': 0.011, 'P': 4985.0}
>>> residuals = oc_fit.residue(result)

ocpy.oc_pymc module

class ocpy.oc_pymc.OCPyMC(oc, minimum_time=None, minimum_time_error=None, weights=None, minimum_type=None, labels=None, cycle=None)[source]

Bases: OC

A probabilistic O–C (Observed minus Calculated) modeling class using PyMC.

This class extends the OC base class to provide Bayesian inference for O–C data, allowing modeling with various components such as linear, quadratic, sinusoidal, and Keplerian functions. It leverages PyMC and ArviZ for probabilistic sampling, posterior analysis, and visualization.

Parameters:

OC (Inherited from) – Accepts the same initialization parameters as OC, such as cycle counts, O–C values, errors, weights, labels, and minimum types.

math

Set to pymc.math to allow model components to use PyMC tensor operations.

Type:

module

fit(model_components, \*, draws, tune, chains, cores, target_accept, random_seed, progressbar, return_model, \*\*kwargs)[source]

Fit the provided model components to the O–C data using Bayesian sampling. Returns an InferenceData object with posterior samples.

clean(inference_data, drop_chains=0, filter_outliers=True, iqr_multiplier=4.0)[source]

Post-process InferenceData to drop chains or filter extreme outliers.

residue(inference_data, \*, x_col="cycle", y_col="oc")[source]

Return a new OCPyMC instance containing residuals between observed and model-predicted O–C values.

fit_linear(*, a=None, b=None, cores=None, \*\*kwargs)[source]

Fit a linear component (a * x + b) to the data.

fit_quadratic(*, q=None, cores=None, \*\*kwargs)[source]

Fit a quadratic component (q * x^2) to the data.

fit_parabola(*, q=None, a=None, b=None, cores=None, \*\*kwargs)[source]

Fit a combination of quadratic and linear components to the data.

fit_sinusoidal(*, amp=None, P=None, cores=None, \*\*kwargs)[source]

Fit a sinusoidal component (amp * sin(2π x / P)) to the data.

fit_keplerian(*, amp=None, e=None, omega=None, P=None, T0=None, name=None, cores=None, \*\*kwargs)[source]

Fit a Keplerian (orbital) component to the data.

fit_lite(**kwargs)[source]

Alias for fit_keplerian; provides a simplified interface for fitting a single Keplerian component.

corner(inference_data, cornerstyle="corner", units=None, \*\*kwargs)[source]

Plot a corner plot of posterior distributions using Plot.plot_corner.

trace(inference_data, \*\*kwargs)[source]

Plot trace plots for posterior chains using Plot.plot_trace.

Notes

  • This class is designed for O–C analyses in eclipsing binaries, exoplanet timing, and other periodic astrophysical phenomena.

  • Supports both deterministic and probabilistic parameters, including truncated distributions for constrained parameters.

  • Works seamlessly with multiple model components, automatically managing parameter naming and aggregation.

math = <module 'pymc.math' from '/home/docs/checkouts/readthedocs.org/user_builds/oc-py/envs/latest/lib/python3.12/site-packages/pymc/math.py'>
fit(model_components, *, draws=2000, tune=2000, chains=4, cores=None, target_accept=None, random_seed=None, progressbar=True, return_model=False, **kwargs)[source]

Fit a model to the observed O–C (observed minus calculated) data using Bayesian inference.

This method constructs a PyMC probabilistic model using the specified model_components, applies priors based on their parameters, and samples from the posterior distribution.

Parameters:
  • model_components (list of ModelComponent) – List of model components to fit (e.g., Linear, Quadratic, Sinusoidal, Keplerian). Each component must define a model_func and optional params.

  • draws (int, optional) – Number of posterior samples to draw after tuning (default: 2000).

  • tune (int, optional) – Number of tuning (burn-in) steps (default: 2000).

  • chains (int, optional) – Number of Markov chains to run (default: 4).

  • cores (int or None, optional) – Number of CPU cores to use. Defaults to number of chains if None.

  • target_accept (float or None, optional) – Target acceptance probability for NUTS sampler (default: None).

  • random_seed (int or None, optional) – Random seed for reproducibility.

  • progressbar (bool, optional) – Whether to show the sampling progress bar (default: True).

  • return_model (bool, optional) – If True, return the PyMC Model object without sampling.

  • **kwargs – Additional keyword arguments passed to pm.sample() or sampler step methods.

Returns:

If return_model is False, returns an ArviZ InferenceData object containing posterior samples, posterior predictive samples, and deterministic variables for model evaluation. If return_model is True, returns the PyMC Model object without sampling.

Return type:

arviz.InferenceData or pm.Model

Notes

  • All model components’ parameters are converted to Parameter objects using _to_param.

  • Observational uncertainties are taken from minimum_time_error; NaNs will raise a ValueError.

  • Deterministic variables y_model and, if possible, y_model_dense are created for visualization.

  • Components marked as _expensive bypass dense evaluation to avoid heavy computations.

clean(inference_data, drop_chains=0, filter_outliers=True, iqr_multiplier=4.0)[source]

Clean an ArviZ InferenceData object by optionally removing chains and filtering outliers.

This method allows post-processing of MCMC samples to improve data quality for visualization or further analysis. Chains with extreme deviations can be dropped, and outliers in the posterior distributions can be filtered using the interquartile range (IQR) method.

Parameters:
  • inference_data (arviz.InferenceData) – The posterior samples to clean, typically returned by OCPyMC.fit().

  • drop_chains (int, optional) – Number of chains to drop based on deviation from the median across chains (default: 0, i.e., keep all chains).

  • filter_outliers (bool, optional) – Whether to remove outliers from the posterior samples (default: True).

  • iqr_multiplier (float, optional) – Multiplier for the interquartile range (IQR) to define outlier bounds (default: 4.0). Samples outside [Q1 - multiplier*IQR, Q3 + multiplier*IQR] are removed.

Returns:

A new InferenceData object with cleaned posterior samples, while preserving deterministic variables and model component metadata (_model_components and _model_prefixes).

Return type:

arviz.InferenceData

Notes

  • Outlier filtering is applied independently to each parameter in the posterior (excluding deterministic variables like y_model, y_model_dense, y_obs, dense_x).

  • Chains are ranked by their deviation from the overall median, and the drop_chains largest deviations are removed.

  • The cleaned dataset preserves the original structure of InferenceData, including posterior, prior, and observed_data groups.

  • Any errors encountered during filtering are issued as warnings without stopping execution.

residue(inference_data, *, x_col='cycle', y_col='oc')[source]

Compute the residuals between observed O–C values and the model predictions.

This method generates a new OCPyMC instance where the O–C values are replaced by the difference between the observed values and the median model predictions from a fitted PyMC InferenceData object.

Parameters:
  • inference_data (arviz.InferenceData) – The posterior samples returned by OCPyMC.fit() containing the fitted model.

  • x_col (str, optional) – Name of the column in self.data containing the independent variable (cycle), used to evaluate the model (default: “cycle”).

  • y_col (str, optional) – Name of the column in self.data containing the observed O–C values (default: “oc”).

Returns:

A new OCPyMC instance with the same metadata as the original, but with the oc values replaced by the residuals (observed minus fitted).

Return type:

OCPyMC

Notes

  • The model predictions are computed as the median across all chains and draws in the posterior (inference_data.posterior[“y_model”]).

  • Only the oc values are updated; all other attributes (weights, labels, minimum times, etc.) are preserved.

  • This method is useful for examining the remaining structure or noise after fitting a model.

fit_linear(*, a=None, b=None, cores=None, **kwargs)[source]

Fit a linear O–C model (O–C = a * cycle + b) using PyMC Bayesian inference.

The linear model is represented as:

O–C = a * x + b

where x is the cycle number.

Parameters:
  • a (float or Parameter, optional) – Initial value or Parameter object for the slope of the line. If None, defaults to 0.0.

  • b (float or Parameter, optional) – Initial value or Parameter object for the intercept of the line. If None, defaults to 0.0.

  • cores (int, optional) – Number of CPU cores to use for sampling. If None, the number of cores defaults to the number of chains.

  • **kwargs – Additional keyword arguments passed to OCPyMC.fit(), such as draws, tune, chains, target_accept, etc.

Returns:

Posterior samples for the linear model parameters and deterministic predictions for the O–C values.

Return type:

arviz.InferenceData

Notes

  • The slope (a) and intercept (b) are treated as PyMC random variables unless provided as fixed Parameter objects with fixed=True.

  • Use residue() on the returned model to compute residuals after fitting.

  • This method is a convenience wrapper that internally creates a Linear component and calls the general fit() method.

fit_quadratic(*, q=None, cores=None, **kwargs)[source]

Fit a quadratic O–C model (O–C = q * cycle^2) using PyMC Bayesian inference.

The quadratic model is represented as:

O–C = q * x^2

where x is the cycle number.

Parameters:
  • q (float or Parameter, optional) – Initial value or Parameter object for the quadratic coefficient. If None, defaults to 0.0.

  • cores (int, optional) – Number of CPU cores to use for sampling. If None, the number of cores defaults to the number of chains.

  • **kwargs – Additional keyword arguments passed to OCPyMC.fit(), such as draws, tune, chains, target_accept, etc.

Returns:

Posterior samples for the quadratic coefficient and deterministic predictions for the O–C values.

Return type:

arviz.InferenceData

Notes

  • The coefficient q is treated as a PyMC random variable unless provided as a fixed Parameter with fixed=True.

  • Use residue() on the returned model to compute residuals after fitting.

  • This method is a convenience wrapper that internally creates a Quadratic component and calls the general fit() method.

fit_sinusoidal(*, amp=None, P=None, cores=None, **kwargs)[source]

Fit a sinusoidal O–C model (O–C = amp * sin(2π * cycle / P)) using PyMC Bayesian inference.

The sinusoidal model is represented as:

O–C = amp * sin(2π * x / P)

where x is the cycle number.

Parameters:
  • amp (float or Parameter, optional) – Initial value or Parameter object for the sinusoidal amplitude. If None, defaults to 1e-3.

  • P (float or Parameter, optional) – Initial value or Parameter object for the period of the sinusoid. If None, defaults to 1000.0.

  • cores (int, optional) – Number of CPU cores to use for sampling. If None, the number of cores defaults to the number of chains.

  • **kwargs – Additional keyword arguments passed to OCPyMC.fit(), such as draws, tune, chains, target_accept, etc.

Returns:

Posterior samples for the sinusoidal parameters and deterministic predictions for the O–C values.

Return type:

arviz.InferenceData

Notes

  • The parameters amp and P are treated as PyMC random variables unless provided as fixed Parameter objects with fixed=True.

  • Use residue() on the returned model to compute residuals after fitting.

  • This method is a convenience wrapper that internally creates a Sinusoidal component and calls the general fit() method.

fit_keplerian(*, amp=None, e=None, omega=None, P=None, T0=None, name=None, cores=None, **kwargs)[source]

Fit a Keplerian O–C model using PyMC Bayesian inference.

The Keplerian model represents the O–C variations as caused by a companion in an orbit, parameterized by:

  • amp : amplitude of the timing variation

  • e : orbital eccentricity

  • omega: argument of periastron (degrees)

  • P : orbital period

  • T0 : time of periastron passage

Parameters:
  • amp (float or Parameter, optional) – Initial value or Parameter object for the amplitude of the O–C signal. Defaults to 0.001 if None.

  • e (float or Parameter, optional) – Initial value or Parameter object for orbital eccentricity. Defaults to 0.1 if None.

  • omega (float or Parameter, optional) – Initial value or Parameter object for the argument of periastron (degrees). Defaults to 90.0 if None.

  • P (float or Parameter, optional) – Initial value or Parameter object for the orbital period. Defaults to 1000.0 if None.

  • T0 (float or Parameter, optional) – Initial value or Parameter object for the time of periastron passage. Defaults to 0.0 if None.

  • name (str, optional) – Optional name for the Keplerian component. If None, defaults to ‘keplerian1’.

  • cores (int, optional) – Number of CPU cores to use for sampling. If None, defaults to the number of chains.

  • **kwargs – Additional keyword arguments passed to OCPyMC.fit(), e.g., draws, tune, chains, target_accept, etc.

Returns:

Posterior samples for all Keplerian parameters, along with deterministic predictions for the O–C curve.

Return type:

arviz.InferenceData

Notes

  • The parameters are treated as PyMC random variables unless provided as fixed Parameter objects with fixed=True.

  • Use residue() on the returned model to compute residuals after fitting.

  • This method internally creates a Keplerian component and calls the general fit() method for Bayesian inference.

fit_lite(**kwargs)[source]

Fit a simplified Keplerian O–C model using PyMC Bayesian inference.

This is a convenience wrapper around fit_keplerian with default Keplerian parameters. It is designed for quick fitting without specifying all orbital parameters explicitly.

Parameters:

**kwargs – Keyword arguments passed directly to fit_keplerian, including: - amp : amplitude of the O–C signal - e : orbital eccentricity - omega: argument of periastron (degrees) - P : orbital period - T0 : time of periastron passage - cores: number of CPU cores to use - draws, tune, chains, target_accept, etc.

Returns:

Posterior samples for the Keplerian parameters and deterministic predictions of the O–C curve.

Return type:

arviz.InferenceData

Notes

  • This method internally calls fit_keplerian with default parameter values:

    amp=0.001, e=0.1, omega=90.0, P=1000.0, T0=0.0

  • Use residue() on the returned model to compute residuals after fitting.

  • Suitable for quick exploratory fits or as an initial guess for more complex models.

fit_parabola(*, q=None, a=None, b=None, cores=None, **kwargs)[source]

Fit a combined quadratic and linear O–C model using PyMC Bayesian inference.

This method models the O–C variations as a sum of a quadratic component (representing secular period change) and a linear component (representing a constant period offset).

Parameters:
  • q (float or Parameter, optional) – Quadratic coefficient for the parabolic term. Default is 0.0.

  • a (float or Parameter, optional) – Linear coefficient for the linear term. Default is 0.0.

  • b (float or Parameter, optional) – Constant offset term. Default is 0.0.

  • cores (int, optional) – Number of CPU cores to use for PyMC sampling. Defaults to the number of chains.

  • **kwargs – Additional keyword arguments passed to the underlying fit method, such as draws, tune, chains, target_accept, and random_seed.

Returns:

Posterior samples for the quadratic and linear coefficients, along with deterministic predictions of the O–C curve.

Return type:

arviz.InferenceData

Notes

  • The quadratic component captures long-term period changes.

  • The linear component models any overall offset in the period.

  • Use the residue() method on the returned object to compute residuals of the fit.

  • Recommended as a first approach when O–C variations show a parabolic trend.

corner(inference_data, cornerstyle='corner', units=None, **kwargs)[source]

Generate a corner (pairwise posterior) plot from PyMC inference results.

This method visualizes the posterior distributions of model parameters, showing marginal distributions along the diagonal and pairwise correlations in the off-diagonal plots. It can use either the corner library or ArviZ built-in plotting depending on the cornerstyle argument.

Parameters:
  • inference_data (arviz.InferenceData) – The PyMC sampling results, typically returned by fit(), containing posterior samples.

  • cornerstyle ({"corner", "arviz"}, default "corner") – Choose the plotting backend: - “corner”: use the corner.py library for plotting. - “arviz”: use ArviZ’s built-in pairplot functionality.

  • units (dict, optional) – A mapping from parameter names to strings describing their units. These units will be displayed in axis labels.

  • **kwargs – Additional keyword arguments passed to the underlying plotting function.

Returns:

  • matplotlib.figure.Figure – The figure object for the ‘corner’ style.

  • arviz.plot_pair axes or Figure – The ArviZ axes object when cornerstyle=’arviz’.

Return type:

Figure | plot_pair

Notes

  • Useful for visual inspection of posterior distributions and parameter correlations.

  • Recommended to use after fit() and optionally after clean() for outlier removal.

  • Units help clarify parameter scales in the plot.

trace(inference_data, **kwargs)[source]

Generate trace plots for PyMC posterior samples.

Trace plots visualize the sampled parameter values across chains, allowing assessment of convergence, mixing, and sampling behavior. Each subplot shows the parameter values over draws for one chain.

Parameters:
  • inference_data (arviz.InferenceData) – The PyMC sampling results, typically returned by fit(), containing posterior samples.

  • **kwargs – Additional keyword arguments passed to the underlying plotting function.

Returns:

Array of matplotlib axes objects containing the trace plots.

Return type:

matplotlib.axes.Axes

Notes

  • Useful for evaluating MCMC convergence and diagnosing sampling issues.

  • Typically used after fit() and optionally after clean() to inspect cleaned chains.

  • Can be combined with corner plots for a complete posterior overview.

ocpy.utils module

class ocpy.utils.Checker[source]

Bases: object

static length_checker(data, reference)[source]
class ocpy.utils.Fixer[source]

Bases: object

static length_fixer(data, reference)[source]
static none_to_nan(df)[source]

ocpy.visualization module

class ocpy.visualization.Plot[source]

Bases: object

static plot_data(data, *, ax=None, x_col='cycle', y_col='oc', plot_kwargs=None)[source]

Plot the raw O−C data with optional error bars and labeling.

Parameters:
  • data (OC) – The observational O−C dataset. Must have a data attribute (pandas DataFrame) containing at least columns for x_col and y_col. Optionally, it can include ‘minimum_time_error’ for y-error bars and ‘labels’ for grouping data by category.

  • ax (matplotlib.axes.Axes, optional) – Axes object on which to plot. If None, a new figure and axes are created.

  • x_col (str, default "cycle") – Name of the column to use for the x-axis.

  • y_col (str, default "oc") – Name of the column to use for the y-axis.

  • plot_kwargs (dict, optional) – Additional keyword arguments passed to matplotlib.pyplot.errorbar.

Returns:

The axes containing the plotted data.

Return type:

matplotlib.axes.Axes

Notes

  • If ‘labels’ exist in the DataFrame, points are color-coded per unique label.

  • Unlabeled points are plotted in gray.

  • If ‘minimum_time_error’ exists, it is used as y-error bars.

  • Axes are automatically labeled and a grid is added.

classmethod plot_model_pymc(inference_data, data, *, ax=None, x_col='cycle', n_points=800, sum_kwargs=None, comp_kwargs=None, plot_kwargs=None, plot_band=True, extension_factor=0.1, model_components=None)[source]

Plot a model fit to O−C data using a PyMC inference result, with optional uncertainty bands and component decomposition.

Parameters:
  • inference_data (arviz.InferenceData) – Posterior samples from a PyMC model. Expected to contain parameter variables (2D arrays with shape [chain, draw]), and optionally ‘y_model’ or ‘y_model_dense’ for precomputed model fits.

  • data (OCPyMC) – Observational O−C dataset to plot against. Must have a data attribute (pandas DataFrame) containing at least the x_col column.

  • ax (matplotlib.axes.Axes, optional) – Axes object on which to plot. If None, uses current axes or creates a new figure.

  • x_col (str, default "cycle") – Column in data to use as the x-axis.

  • n_points (int, default 800) – Number of points to use when plotting continuous model curves.

  • sum_kwargs (dict, optional) – Additional keyword arguments for the summed model line.

  • comp_kwargs (dict, optional) – Additional keyword arguments for individual model component lines.

  • plot_kwargs (dict, optional) – Additional keyword arguments for the data points plot (markers, color, etc.).

  • plot_band (bool, default True) – Whether to display a 1σ uncertainty band from posterior samples.

  • extension_factor (float, default 0.1) – Fractional extension beyond the data range for plotting the fit curve.

  • model_components (list, optional) – List of model component objects (Linear, Quadratic, Sinusoidal, Keplerian, etc.) to reconstruct and plot individual contributions.

Returns:

The axes containing the plotted data, model fit, components, and optional uncertainty band.

Return type:

matplotlib.axes.Axes

Notes

  • If ‘y_model_dense’ and ‘dense_x’ exist in inference_data, they are used for a smooth fit.

  • Otherwise, the method reconstructs the model from posterior medians of scalar parameters.

  • If model_components are provided, each component is plotted individually.

  • Uncertainty bands are computed from a subset of posterior samples (default 200 draws).

  • Automatically handles multiple components, labeling, and plotting the sum of components.

classmethod plot_model_lmfit(result, data, *, ax=None, x_col='cycle', n_points=500, plot_kwargs=None, extension_factor=0.1)[source]

Plot a model fit to O−C data using an lmfit.ModelResult object, with optional uncertainty bands.

Parameters:
  • result (lmfit.model.ModelResult) – The fitted model result returned by lmfit.Model.fit. Expected to provide eval(x=…) and optionally eval_uncertainty(x=…, sigma=1).

  • data (OCLMFit) – Observational O−C dataset to plot against. Must have a data attribute (pandas DataFrame) containing at least the x_col column.

  • ax (matplotlib.axes.Axes, optional) – Axes object on which to plot. If None, creates a new figure and axes.

  • x_col (str, default "cycle") – Column in data to use as the x-axis.

  • n_points (int, default 500) – Number of points to evaluate for a smooth model curve.

  • plot_kwargs (dict, optional) – Additional keyword arguments for the fit line (color, linewidth, label, etc.).

  • extension_factor (float, default 0.1) – Fractional extension beyond the data range for plotting the fit curve.

Returns:

The axes containing the plotted data, model fit, and optional uncertainty band.

Return type:

matplotlib.axes.Axes

Notes

  • The method evaluates the fitted model on a dense set of points across the data range, optionally extended by extension_factor.

  • If result.eval_uncertainty is available, a 1σ uncertainty band is plotted around the fit.

  • Data points are plotted separately using the plot_data method (if called externally).

classmethod plot_model_components(model_components, xline, *, ax=None, sum_kwargs=None, comp_kwargs=None, uncertainty_band=None)[source]

Plot individual model components and their sum over a specified x-range, with optional uncertainty bands.

Parameters:
  • model_components (list) – List of model component objects. Each component must have a model_func(x, **params) method and a params attribute (dict of Parameter objects or numeric values).

  • xline (np.ndarray) – Array of x-values to evaluate the component models.

  • ax (matplotlib.axes.Axes, optional) – Axes object on which to plot. If None, a new figure and axes are created.

  • sum_kwargs (dict, optional) – Keyword arguments for the summed model curve. Default color is red, linewidth 2.6, alpha 0.95.

  • comp_kwargs (dict, optional) – Keyword arguments for individual component curves. Default linewidth 1.5, alpha 0.9, linestyle ‘–‘.

  • uncertainty_band (tuple, optional) – Tuple (x_band, y_low, y_high) representing an uncertainty envelope around the sum of components. If provided, plotted as a filled area behind the curves.

Returns:

The axes containing the component curves, sum curve, and optional uncertainty band.

Return type:

matplotlib.axes.Axes

Notes

  • Each component is evaluated using its model_func and current parameter values.

  • The sum curve is drawn on top of the components, optionally with an uncertainty band.

  • Components without parameters (or with missing parameters) will raise a KeyError.

  • Useful for visualizing contributions of multiple additive model components in O−C analysis or similar time-series modeling contexts.

classmethod plot(data, model=None, *, ax=None, res_ax=None, res=True, title=None, x_col='cycle', y_col='oc', fig_size=(10, 7), plot_kwargs=None, extension_factor=0.1, model_components=None)[source]

Plot data with optional model fit and residuals.

This is a high-level plotting function that can: - Display raw O−C data points (with optional labels and error bars), - Overlay model fits from PyMC (InferenceData), lmfit (ModelResult), or a list of model components, - Display residuals below the main plot if requested.

Parameters:
  • data (OC) – The data object containing O−C measurements. Must have data attribute (pandas DataFrame) with at least columns specified by x_col and y_col.

  • model (Union[InferenceData, ModelResult, List[ModelComponent]], optional) – Model to overlay on the data: - PyMC model (InferenceData from arviz) with posterior samples, - lmfit result (ModelResult) with .eval() method, - List of component objects with .model_func and .params. If None, only the raw data is plotted.

  • ax (matplotlib.axes.Axes, optional) – Axes for the main data/fit plot. If None, a new figure is created.

  • res_ax (matplotlib.axes.Axes, optional) – Axes for residuals plot. If None and res=True, a new subplot is created.

  • res (bool, default True) – Whether to plot residuals beneath the main plot.

  • title (str, optional) – Title for the main plot.

  • x_col (str, default "cycle") – Column in data.data used for x-axis.

  • y_col (str, default "oc") – Column in data.data used for y-axis.

  • fig_size (tuple, default (10, 7)) – Figure size in inches.

  • plot_kwargs (dict, optional) – Keyword arguments passed to the main data plot (color, markers, alpha, etc.).

  • extension_factor (float, default 0.1) – Fractional extension beyond the data range for plotting model fits.

  • model_components (list, optional) – If provided, used for plotting PyMC components in plot_model_pymc.

Returns:

If res=True, returns a tuple (ax, res_ax) for main and residual plots. Otherwise, returns only ax.

Return type:

matplotlib.axes.Axes or tuple(matplotlib.axes.Axes, matplotlib.axes.Axes)

Notes

  • Automatically handles labeled and unlabeled data points.

  • Residuals are computed as y - y_model if model is provided.

  • Supports both PyMC posterior models (median and uncertainty bands) and lmfit fits.

  • Useful for O−C analysis in astronomy or any time-series with additive models.

static plot_corner(inference_data, var_names=None, cornerstyle='corner', units=None, **kwargs)[source]

Generate a corner plot (pairwise parameter plot) from PyMC/ArviZ inference data.

Parameters:
  • inference_data (arviz.InferenceData) – The posterior samples from a Bayesian model.

  • var_names (list of str, optional) – Names of parameters to include. If None, automatically selects all posterior parameters with 2 dimensions (chain, draw) excluding model output.

  • cornerstyle ({'corner', 'arviz'}, default 'corner') – Which library/style to use for the corner plot: - ‘corner’ : uses the corner package. - ‘arviz’ : uses ArviZ’s plot_pair.

  • units (dict, optional) – Mapping from parameter name to unit string, e.g., {‘P’: ‘days’}. Units are appended to axis labels.

  • **kwargs – Additional keyword arguments passed to corner.corner or arviz.plot_pair.

Returns:

  • matplotlib.figure.Figure – The figure object for the ‘corner’ style.

  • arviz.plot_pair axes or Figure – The ArviZ axes object when cornerstyle=’arviz’.

Return type:

Figure | plot_pair

Notes

  • Automatically ignores fixed parameters with negligible variation.

  • Adds median values as “truths” in the plot if using ‘corner’ and no truths provided.

  • If units are provided, labels are formatted with LaTeX, e.g., r”$P$ [days]”.

  • Supports up to 200 subplots in ArviZ style using plot.max_subplots context.

Raises:
  • ImportError – If corner library is required but not installed.

  • ValueError – If no suitable parameters are found for plotting.

static plot_trace(inference_data, var_names=None, **kwargs)[source]

Generate trace plots for posterior samples from a PyMC InferenceData object.

Trace plots show the sampled parameter values across chains and draws, allowing evaluation of convergence, mixing, and sampling behavior.

Parameters:
  • inference_data (arviz.InferenceData) – The posterior sampling results, typically returned by a PyMC fit.

  • var_names (list of str, optional) – List of parameter names to include in the trace plot. If None, all variable parameters with variation are included.

  • **kwargs – Additional keyword arguments passed to arviz.plot_trace.

Returns:

Array of matplotlib axes objects containing the trace plots.

Return type:

matplotlib.axes.Axes

Notes

  • Automatically excludes fixed parameters (near-zero variance) unless all are fixed.

  • Figures are automatically tightened with tight_layout.

  • Designed to complement plot_corner for full posterior visualization.

Module contents