bayesian_blocks

astropy.stats.bayesian_blocks(t, x=None, sigma=None, fitness='events', **kwargs)[source]

Compute optimal segmentation of data with Scargle’s Bayesian Blocks

This is a flexible implementation of the Bayesian Blocks algorithm described in Scargle 2013 [1].

Parameters:
tnumpy:array_like

data times (one dimensional, length N)

xnumpy:array_like, optional

data values

sigmanumpy:array_like or python:float, optional

data errors

fitnesspython:str or object

the fitness function to use for the model. If a string, the following options are supported:

  • ‘events’ : binned or unbinned event data. Arguments are gamma, which gives the slope of the prior on the number of bins, or ncp_prior, which is \(-\ln({\tt gamma})\).

  • ‘regular_events’ : non-overlapping events measured at multiples of a fundamental tick rate, dt, which must be specified as an additional argument. Extra arguments are p0, which gives the false alarm probability to compute the prior, or gamma, which gives the slope of the prior on the number of bins, or ncp_prior, which is \(-\ln({\tt gamma})\).

  • ‘measures’ : fitness for a measured sequence with Gaussian errors. Extra arguments are p0, which gives the false alarm probability to compute the prior, or gamma, which gives the slope of the prior on the number of bins, or ncp_prior, which is \(-\ln({\tt gamma})\).

In all three cases, if more than one of p0, gamma, and ncp_prior is chosen, ncp_prior takes precedence over gamma which takes precedence over p0.

Alternatively, the fitness parameter can be an instance of FitnessFunc or a subclass thereof.

**kwargs

any additional keyword arguments will be passed to the specified FitnessFunc derived class.

Returns:
edgesndarray

array containing the (N+1) edges defining the N bins

See also

astropy.stats.histogram

compute a histogram using bayesian blocks

References

[2]

Bellman, R.E., Dreyfus, S.E., 1962. Applied Dynamic Programming. Princeton University Press, Princeton. https://press.princeton.edu/books/hardcover/9780691651873/applied-dynamic-programming

[3]

Bellman, R., Roth, R., 1969. Curve fitting by segmented straight lines. J. Amer. Statist. Assoc. 64, 1079–1084. https://www.tandfonline.com/doi/abs/10.1080/01621459.1969.10501038

Examples

Event data:

>>> t = np.random.normal(size=100)
>>> edges = bayesian_blocks(t, fitness='events', p0=0.01)

Event data with repeats:

>>> t = np.random.normal(size=100)
>>> t[80:] = t[:20]
>>> edges = bayesian_blocks(t, fitness='events', p0=0.01)

Regular event data:

>>> dt = 0.05
>>> t = dt * np.arange(1000)
>>> x = np.zeros(len(t))
>>> x[np.random.randint(0, len(t), len(t) // 10)] = 1
>>> edges = bayesian_blocks(t, x, fitness='regular_events', dt=dt)

Measured point data with errors:

>>> t = 100 * np.random.random(100)
>>> x = np.exp(-0.5 * (t - 50) ** 2)
>>> sigma = 0.1
>>> x_obs = np.random.normal(x, sigma)
>>> edges = bayesian_blocks(t, x_obs, sigma, fitness='measures')