Optimal Segmentation with
Pruned Dynamic Programming
Jeff Scargle
NASA Ames Research Center
Jeffrey.D.Scargle@nasa.gov
Thanks: Brad Jackson + Students, Jay Norris, Jim Chiang
AI @ SLAC September 26, 2017

Event Time Tags: 2 μ-secWhat is Optimal Segmentation?
DATA:• any signal measurements• sequential in time, space, energy, …
MODEL:• partition data interval into blocks• block = set of consecutive data points• model of signal in each block• block cost (aka fitness, risk, etc. )OPTIMIZATION:• total model cost = sum of block costs• optimize over all possible partitions (~2N)

What good are
Segmented Time Series
Representations? (1)Detect and characterize any and all
statistically signficiant variability
supported by the data.

Owens Valley
Radio Observatory 40-m

Integrations (32 s) << sample intervals (~100,000 s)

Bayesian Block Generalizations
Any Data Mode:• point measurements• time-tagged events• live-time intervals• categorical• gaps, uneven sampling, exposure variation, real-time• … any data mode

Any Block Shape:• constant• linear• exponential• FRED (fast rise/exponential decay) /DERF• … any block shape

Improved speed:O(N2) → ~N

Data Points… must include whatever information is necessary
to compute the objective function of a data block(Typically: time + intensity + error)• measurements at a point in time - X( ti ), σ ( ti )• measurements averaged over a time interval• event times (TTE)• inter-event intervals• live-time intervals• auxiliary information (e.g. color)• categorical (e.g. 0-1, ACGT, … )• multivariate (multi-wavelength, mulit-messenger)• mixtures of any of the above• data on the circle (e.g. angles)• higher dimension (e.g. 3D galaxy positions)•…

Integration times ~sample intervals

The best of 106020
possible partitions

The Human Genome
ATCAGAGACTCGTACCCGTGCCTGAGCAGTGCTGAAGAGGCACTCGTTTGGAAAGGGGGCCCAT
CCCCCGGGACTTCGGACACTCCTGGCTGAAGCATAACGTGTAGGCGCTCTAGAGGGCTCGCCTA
CCCGAGCCCACATACCAGCGATGATGAGATCGACCCGCACCGACGCATGGAGCGCAGGCGTGCTGCTAGTCGAGGGCACGGTCGCCCCAACAAACG
ACCGCCCA … N ~100,000,000 bp“GC islands”: G or C —> 1, A or T —> 0; two data modes:• 0-1 on fixed grid• intervals between 1’s

Bayesian Block analysis of event time series is mathematicallythe same problems as constructing (optimal) histograms!

Starting value
for fit routineMaximum (log) Likelihood Block Costs
Event Data: Cn = Nn log( Nn / Tn )Nn = number of events in block n
Tn = length of block n

Point Measurements:
Cn = ( ∑ wn xn )2 / ( ∑ wn )xn = measured value
wn = weight = 1 / σn2Total model cost: ∑ Cnvariance!

Optimization…
…

The “best” partition
from a huge numberof possibilities

...
Optimum Partition Up To This Point

...

...

Prospective Last BlockDynamic Programming Algorithm
for ii = 1:num_points

sum_num
sum_len= sum_num= sum_len+ num_vec( ii );%count+ len_vec( ii );%length% Cost of last block:
cost_last= cost_function(sum_num
sum_len, ...
);% Compute and maximize total model cost:
cost_total= opt+ cost_last[ opt(ii+1) last ] = max( cost_total);
last_change( ii ) =endlast; % Store last change-point;Pruned Dynamic Programming Algorithm
for ii = 1:num_points
unpruned( ii ) = ii; % New point unpruned until found otherwise
sum_num( unpruned ) = sum_num( unpruned ) + num_vec( ii );%countsum_len( unpruned ) = sum_len( unpruned ) + len_vec( ii );%length% Cost of last block:
cost_last( unpruned ) = cost_function(sum_num( unpruned ), ...
sum_len( unpruned ));% Compute and maximize total model cost:
cost_total( unpruned ) = opt( unpruned ) + cost_last( unpruned );[ opt(ii+1) last ] = max( cost_total( unpruned ) );last_change( ii ) = unpruned( last ); % Store last change-point% update pruning
id_unpruned = find([ cost_total(unpruned) > opt(ii+1) - ncp_prior ]);unpruned = unpruned( id_unpruned );endWhat good areSegmented Time Series
Representations? (2)Detect and characterize,
without bins or smoothing:◆ Pulses (aka “flares”)◆ Pulse shapes (including the Arrow of Time)◆ Variability index◆ Variability time scales (min, max, dist, …)◆ Transient event triggers (real-time mode)What good are
Segmented Time Series
Representations? (3)Implement:◆ Exploratory Data Analysis◆ Time series classification◆ Noise suppression◆ Visual displays◆ Data compression◆ Data adaptive histograms

Thanks