Title: LEARNING PROCEDURAL ABSTRACTIONS AND EVAL-UATING DISCRETE LATENT TEMPORAL STRUCTURE

Abstract: Clustering methods and latent variable models are often used as tools for pattern mining and discovery of latent structure in time-series data. In this work, we consider the problem of learning procedural abstractions from possibly highdimensional observational sequences, such as video demonstrations. Given a dataset of time-series, the goal is to identify the latent sequence of steps common to them and label each time-series with the temporal extent of these procedural steps. We introduce a hierarchical Bayesian model called PRISM that models the realization of a common procedure across multiple time-series, and can recover procedural abstractions without supervision. We also bring to light two characteristics ignored by traditional evaluation criteria when evaluating latent temporal labelings (temporal clusterings) -segment structure, and repeated structure -and develop new metrics tailored to their evaluation. We demonstrate that our metrics improve interpretability and ease of analysis for evaluation on benchmark time-series datasets. Results on benchmark and video datasets indicate that PRISM outperforms standard sequence models as well as state-of-the-art techniques in identifying procedural abstractions.

1. INTRODUCTION:

A fundamental problem in machine learning is the discovery of structure in unsupervised data. In this work we are particularly interested in uncovering the latent structure in procedural datapotentially high-dimensional observational time-series resulting from some latent sequence of events. For example, a video showing how to change a tire might involve jacking up the car, removing one nut from the wheel, removing a second nut, etc. There exists an enormous wealth of videos of procedures available and inferring the latent structure of such videos could be useful in a huge range of applications, from supporting search queries ("find all segments where someone jacks up the car") to adding to nascent work on learning from observation (Borsa et al., 2017; Stadie et al., 2017; Torabi et al., 2018) in which robots may be taught merely by observing a procedure performed in a video. Procedural learning is related to activity recognition and the broader field of latent temporal dynamical learning, but focuses on the simpler but still ubiquitous setting where the latent activity sequence is a fixed procedure.

As we started to develop methods for performing latent procedure inference from temporal data, we considered how to evaluate our resulting methods. One important evaluation protocol for unsupervised learning methods is external evaluation, where predicted labelings of the discrete variables are compared to ground-truth labels (Rosenberg & Hirschberg, 2007) . But the precise criteria for this evaluation is critical. Though there are an enormous number of methods for modeling and inferring latent structure in dynamical systems, including Hidden Markov models, dynamic Bayesian networks, and temporal clustering for time-series data e.g. Fox et al. (2008b; 2014) ; Krishnan et al. (2015) ; Linderman et al. (2017) ; Zhou et al. (2008; 2013) ; Vidal & Favaro (2014) , the external evaluation for such approaches is typically done using the same clustering metrics that are used for non-temporal data: metrics like the normalized mutual information (NMI) (Strehl & Ghosh, 2002) or by computing cluster correspondences using the Munkres algorithm (Munkres, 1957) . Unfortunately such metrics disregard temporal information and can therefore make it hard to assess important considerations that arise in latent temporal structure extraction, including in procedural inference. Figure 1 : Procedural learning. Given an input set of temporal sequences all demonstrating the same procedure we wish to infer the latent structure of that procedure. Here we show the results of our method, PRISM, on inferring a cpr procedure from a set of videos in the INRIA Instructional videos dataset. The ground truth temporal clustering on the upper left represents a human labeling of the activities involved in a particular cpr video where color denotes the same activity. For example, the activity denoted by the cream color occurred 4 times across the sequence.

To address the limitations of prior evaluation criteria for the external evaluation of discrete latent temporal structure, our first contribution is the creation of new criteria. We identify at least two key aspects of temporal latent structure discovery that are not captured well with existing non-temporal metrics: segment structure, the temporal durations between transitions from one discrete latent state to another, and repeated structure, where a particular latent variable may appear in multiple segments. Figure1 shows an example latent sequence associated with a person performing cardiopulmonary resuscitation (cpr) -the colors denote the activity type (such as blowing air into the lungs) and the temporal duration of a fixed activity type is shown as a continuous block of color. Our new metrics can evaluate both segment structure and repeated structure, and we later demonstrate how these metrics may be of broader interest in latent temporal dynamical learning, as they can illustrate and illuminate different strengths and weaknesses of existing work on common benchmark datasets.

Equipped with new measures for evaluating success, we return to our key goal, inferring the latent procedure observed in a set of high dimensional demonstrations. Past expressive models like Hidden Markov models can struggle in this setting, since the resulting procedure may not be Markov: while the setting can be forced to appear Markov by employing the history as part of the state, this can reduce sample efficiency (by artificially increasing the latent space size). Additionally assuming the system is Markov when it is not can mean inferring stochastic latent dynamical structure when it is in fact deterministic. Recent work (Sener & Yao, 2018) tackles a closely related setting, but does not allow repeated primitives, which is a key aspect of many procedures. Instead we introduce a new generative Bayesian model (PRISM) for identifying procedural abstractions without supervision. Results on datasets for procedural tasks -including surgical demonstrations and how-to videosshow that our method can learn procedural abstractions effectively compared with state-of-the-art methods and popular baselines, and outperform prior work on both our new metrics and classic metrics when repeated structure is present.

2. DEFINITIONS AND SETTING:

Temporal Clustering. A time-series X i = (x i 1 , x i 2 , . . . , x i mi ), i ∈ [n] is a sequence of feature vectors (or items) of length m i with x i t ∈ R d . A temporal clustering C = (c 1 , . . . , c m

) is a sequence of cluster labels with c t ∈ Γ C where Γ C is a set of cluster labels. Temporal clusterings map each item in the feature trajectory X to a cluster label in Γ C . We use C i as shorthand to denote the temporal clustering generated by some method for X i . The ground truth temporal clustering for X i will be denoted by G i , with corresponding ground truth label set Γ G . For simplicity, we define

C a:b = (c a , c a+1 , . . . , c b ) where 1 ≤ a ≤ b ≤ m are time indices. For α ∈ Γ C , let C[α] = {t|c t = α, 1 ≤ t ≤

m} be the set of time indices in C whose cluster labels equal α. C[α] is the cluster or partition (in the standard sense) corresponding to the label α. Temporal clusterings are distinct from time-series segmentation (Chung et al., 2004) and changepoint detection (Killick et al., 2012) , which consider only identifying boundaries of different dynamic regimes in the time-series.

Segment. A segment is a contiguous sequence of identical cluster labels in C. Each cluster in C may be split across multiple segments. We represent each segment as a pair of time indices (a, b), where a represents the start time-index of the segment while b the end time-index (both included). Formally, we let The procedure captures the primitives (as tokens) present in the temporal clustering and the sequence in which they occur. We say that C exhibits repeated structure if atleast one primitive is reused in C. This corresponds to having multiple tokens of the same cluster label in P(C). The procedure (A, B, C, A, B) contains repeated structure since both A and B are reused. We also define a weight function W that counts the number of running duplicates in a sequence e.g. W(A, A, A, B, B, B, C, C, C, A, A, B) = (3, 3, 3, 2, 1). Notice that W computes the length of each segment in the sequence. Weights in W(C) thus have a one-to-one correspondence with tokens in P(C), equaling the length of each token in the procedure. Taken together W(C), P(C) can be used to recover C. Lastly, we let H(P 1 , P 2 , W 1 , W 2 ) be a similarity measure between 2 temporal clusterings represented by P 1 , W 1 and P 2 , W 2 respectively. 

S(C, α) = {(a, b)|C a:b = (c a = α, c a+1 = α, . . . , c b = α), 1 ≤ a ≤ b ≤ m} be a function

3. EVALUATING LATENT TEMPORAL STRUCTURE:

External clustering evaluation is the process of comparing a clustering to a ground-truth labeling in order to evaluate it. Unfortunately, as we discuss below, using traditional clustering evaluation criteria is problematic for evaluating temporal clusterings.

Existing evaluation criteria. The use of standard clustering metrics for evaluation is prevalent in prior work that attempts to evaluate latent temporal structure (Zhou et al., 2008; 2013; Fox et al., 2008b; 2009; 2014; Krishnan et al., 2015; Sener & Yao, 2018; Hoai & De la Torre, 2012) . For a temporal clustering C the widely used purity (Rosenberg & Hirschberg, 2007) metric penalizes the presence of items from different ground truth labels in the same cluster. This is a desirable property that a clustering evaluation criterion should have. Purity is defined as purity (Rosenberg & Hirschberg, 2007) , which captures a similar idea in a different way:

= 1 m α∈Γ C max β∈Γ G |G[β] ∩ C[α]|, A related metric is homogeneity

homogeneity = 1 -H(Γ G |Γ C ) H(Γ G )

where H(Γ G ) is the entropy of the clustering G when partitioned by the labels Γ G , H(Γ G |Γ C ) is the conditional entropy of G given C. Intuitively the conditional entropy term looks inside every cluster in C, and checks the entropy of the items inside in terms of the ground truth labels. The conditional entropy will be low (and homogeneity high) if every cluster in C contains items of only a single ground truth label. Another important metric is completeness (Rosenberg & Hirschberg, 2007) , which prefers clusterings where all items from a ground-truth label lie in the same cluster i.e. the ground truth label is not split across clusters. Completeness is defined identically to homogeneity, but with the clusterings swapped G ←→ C.

Both completeness and homogeneity are considered important criteria for clustering. Maximizing either at the expense of the other is bad -high homogeneity and low completeness imply a clustering that is too fine-grained (e.g. every item in its own cluster), while low homogeneity and high completeness imply a clustering that is too coarse (e.g. a single cluster containing all items). Therefore, metrics that balance these 2 criterion often perform best (Rosenberg & Hirschberg, 2007) normalized mutual information (NMI) (Strehl & Ghosh, 2002) and V-measure (Rosenberg & Hirschberg, 2007) are examples of such metrics.

Another popular suite of metrics uses the Munkres algorithm (Zhou et al., 2008) -a one-to-one correspondence between clusters and ground-truth labels is computed so that any classification criterion (commonly accuracy) under this correspondence is maximized.

Characterizing latent temporal structure. Two major characteristics distinguish the evaluation of temporal clusterings from the standard clustering case.

• Repeated structure. The presence of repeated structure in ground-truth -ground-truth labels that are split across several segments in the trajectory. Ideally, we want that this repeated structure is identified correctly e.g. in Fig 3 identifying that the segments A 0:2 , A 4:6 are repetitions of the same primitive. • Segment structure. The locations and sizes of the segments in ground-truth e.g. in Fig 3 we would like temporal clusterings that switch between primitives at exactly the time points t = 2, 4 with the same 3 segment structure.

Neither repeated or segment structure are seen in standard clusterings, since they treat data as a bag of items.

Limitations of existing criteria. Existing metrics cannot score repeated or segment structure because they convert temporal clusterings to standard clusterings -a contingency matrix (Rosenberg & Hirschberg, 2007 ) -before evaluation. This erases the temporal information important to score these ideas. In fact, qualitatively distinct temporal clusterings can be mapped to the same contingency matrix, and receive the same score, which is clearly undesirable.

We 

[A] -C 2 into C 2 [X] and C 2 [Z] and C 3 into C 3 [X] and C 3 [Z].

In addition, C 2 also misses the repeated structure in A while C 3 does not. C 3 captures this repeated structure as XZ, with X 0:1 Z 1:2 corresponding to A 0:2 and X 4:5 Z 5:6 corresponding to A 4:6 . However, it does so at the cost of creating more segments than C 2 -ignoring labels, C 2 has a segmentation exactly like G. Intuitively, we prefer C 1 > C 3 > C 2 if we care more about repeated structure, and We now discuss new evaluation criteria that can systematically evaluate and score both repeated and segment structure.

C 1 > C 2 > C 3 , if

3.1. EVALUATING REPEATED STRUCTURE:

We describe a new evaluation criteria, called the repeated structure score (RSS) below. Algorithmically, the calculation of the RSS requires the following steps, 1. Pick a ground-truth label β ∈ Γ G and find the set of segments S(G, β) marked with β in G. These segments all exhibit repeated structure in G.

2. Find C a:b for every (a, b) ∈ S(G, β). This is the predicted temporal clustering restricted to a particular segment (a, b) in G.

3. For every such pair C a1:b1 , C a2:b2 , use a scoring function H to compute H(C a1:b1 , C a2:b2 ).

4. Aggregate these scores across all such pairs for every ground-truth label

β ∈ Γ G , RSS = β∈Γ G (a1,b1),(a2,b2)∈S(G,β) H(C a1:b1 , C a2:b2 ) 2 β∈Γ G |S(G, β)| (a,b)∈S(G,β) (b -a + 1)

where the denominator normalizes RSS to lie in [0, 1]. For step 2, we find C 0:2 = (X, X) and C 4:7 = (X, Z, Z). For step 3, we must choose an appropriate scoring function H. While many choices for H are possible, we let H(C a1:b1 , C a2:b2 ) be the total weight of the heaviest common sub-sequence or substring (Jacobson & Vo, 1992) in C a1:b1 and C a2:b2 Intuitively, H tries to match procedures P(C a1:b1 ) to P(C a2:b2 ) -finding the (heaviest) sequence of tokens common to both P(C a1:b1 ) and P(C a2:b2 ), while respecting the temporal ordering of both procedures. At one extreme, if H evaluates to 0 then no overlapping tokens exist in P(C a1:b1 ) and P(C a2:b2 ). This implies that C failed to identify repeated structure in the two segments and a score of 0 is appropriate. At the other extreme, if both segments follow identical procedures P(C a1:b1 ) = P(C a2:b2 ), they receive the maximum possible score according to H; reflecting that C identified repeated structure perfectly. H can be computed efficiently using a dynamic program in O(|P 1 ||P 2 |).

3.2. EVALUATING SEGMENT STRUCTURE:

We first introduce a new information-theoretic criteria -the label-agnostic segmentation score (LASS). The main purpose of LASS is to score how well C predicts the location of transition points between segments in G. We want to identify and score cases where C either oversegments (introduces extra segments) or undersegments (omits segments) the time-series compared to G. In any temporal clustering, the transition boundaries do not depend on the actual labeling, but only on where the labels switch. Therefore, LASS is label-agnostic and only requires the transition points between segments.

Prior work in changepoint detection and time-series segmentation (Killick et al., 2012) has attempted to quantify this using detection criteria -e.g. if a changepoint close to one in ground-truth is detected, a true positive is recorded. However, these criteria (i) require specifying the tolerance interval for detection, and (ii) don't degrade gracefully with the segmentation, which limits their applicability to temporal clusterings.

Recall that S G , S C are the set of all segments in G and C. To compute LASS, the main quantity we rely on is H (S C |S G ) -a conditional entropy over segments in C. This can be written as,

H (S C |S G ) = - (a,b)∈S G (b -a + 1) n (ã, b)∈S C a:b ( b -ã + 1) (b -a + 1) log ( b -ã + 1) (b -a + 1)

Intuitively, for each ground-truth segment (a, b) ∈ S G , we compute the (weighted) entropy of C restricted to (a, b). This entropy will be > 0 only if C creates new segments in (a, b). On the other hand, to detect under-segmentation, we can simply invert the roles of G and C by calculating

H (S G |S C ). Finally, noting that H (S C |S G ) ≤ H (S C ), we define LASS := 1 -H(S C |S G )+H(S G |S C ) H(S C )+H(S G )

. Oversegmentation (and under-segmentation) can also be measured separately by defining

LASS-O := 1 -H(S C |S G ) H(S C ) and LASS-U := 1 -H(S G |S C )

H(S G ) . The algebraic form chosen for the criteria relates it back to existing information-theoretic criteria such as the NMI (Vinh et al., 2010) .

Next, we extend the completeness metric defined earlier. While completeness is desirable for any temporal clustering, using completeness directly has an unintended drawback -it penalizes repeated structure. An example of this can be seen in Fig 3a , where C 3 is assigned the same completeness score as C 2 despite capturing repeated structure. Since the RSS already provides a systematic way of evaluating presence/absence of repeated structure, we modify the completeness score slightly. We define segmental completeness to be 1 -

H(Γ C |S G ) H(Γ C )

, where the conditioning is now on ground-truth segments rather than on ground-truth labels. Equivalently, we can view this as relabeling each segment in ground-truth to be distinct, and then computing completeness. In this way, segmental completeness serves a similar function to completeness without interfering with repeated structure scoring. Homogeneity can also be extended to define segmental homogeneity similarly.

Lastly, we define an aggregate segment structure score (SSS),

SSS := 1 - H(S C |S G ) + H(S G |S C ) + H(Γ C |S G ) + H(Γ G |S C ) H(S C ) + H(S G ) + H(Γ C ) + H(Γ G ) 3.3 A COMBINED EVALUATION CRITERIA

We have introduced several new metrics: (i) the RSS systematically scores repeated structure; (ii) the LASS scores how well the transition structure is captured; (iii) the SSS provides a unified metric for assessing segment structure. Note that all of these metrics obey the useful property of n-invariance (Rosenberg & Hirschberg, 2007; Meilȃ, 2007) -scaling the number of items (stretching the temporal clustering) does not change the metric values. Finally, we provide a single evaluation measure for temporal clusterings, the temporal structure score (TSS) that balances the RSS with the SSS.

TSS := (1 + β) • RSS • SSS β • RSS + SSS

While it is preferable to examine the constituent scores individually to determine the strengths and weaknesses of methods, a unified metric allows us to choose methods that balance these criteria. β can be tuned to change the influence of the constituent metrics according to the desired applicationwe set β = 1.0 in a problem-agnostic fashion in this paper. This choice is in line with previous work that defines compound criteria, such as Rosenberg & Hirschberg (2007) .

4. IDENTIFYING PROCEDURAL ABSTRACTIONS:

So far, we have argued the importance of evaluating repeated and segment structure correctly in temporal data. Our new metrics enable a systematic evaluation of these ideas, whereas existing metrics did not specifically evaluate or emphasize them. Given this, we now return to our original goal of identifying procedural abstractions in time-series data, which we term procedure identification. Concretely, we assume a dataset of time-series such that each demonstrates the same task or phenomenon. Crucially, we assume that all time-series share this single, latent procedure -different time-series are distinct realizations of this underlying procedure. The goal is to label the temporal structure present in each time-series while identifying the common, latent procedure. For example, this may involve identifying the common sequence of steps in a set of demonstrations that follow the same cooking recipe as well as their temporal extent in each demonstration.

Procedure identification is related to recent efforts in computer vision that relate to unsupervised activity recognition in videos (Wu et al., 2015; Yang et al., 2013; Krishnan et al., 2015; Sener & Yao, 2018) as well as prior work in latent variable modeling (Yang et al., 2013; Fox et al., 2008b; 2014; Johnson & Willsky, 2013; Linderman et al., 2017) . These papers typically assume a general setting where activities or latent tokens may be reused within or across multiple time-series but there is no onus on methods to discover this structure. Unlike prior work, procedure identification requires that methods also return a common procedure in addition to the temporal clustering of each time-series. This additional requirement places a burden on the method to explicitly, and not implicitly reason about repeated structure.

Formally, finding a common procedure P corresponds to finding temporal clusterings C i for each time-series, such that

P = P(C 1 ) = P(C 2 ) = • • • = P(C n ).

Since a single procedure consistent with all temporal clusterings is desired, shared primitives must be identified across all videos. This requirement also makes evaluation in this setting particularly suited for our metrics. Methods that cannot relate primitives across demonstrations (identify repeated structure) or reason about their temporal extents (identify segment structure) will be unable to recover the underlying procedure and should therefore be penalized. There may also be additional repeated structure within a time-series, for example if someone mixes cake batter, adds ingredients, and then mixes the cake batter again.

A closely related setting was tackled by the method outlined in Sener & Yao (2018) , who assume that each time-series is a permutation of a shared bag of activities. However, as they acknowledge, a major limitation of their method is that they allow each activity to occur atmost once within a time-series, which limits repeated structure. Our problem setting allows for this possibility, which is often seen in real-world datasets e.g. the JIGSAWS surgical dataset (Gao et al., 2014) contains videos of surgeons doing multiple repetitions of primitives such as needle-passing or suturing. However, we do not allow the procedure to be permuted in our problem setting, which is a useful extension that we hope to explore in future work.

4 Modeling the procedure. PRISM assumes that each time-series X i follows the same underlying procedure p = (p 1 , . . . , p s ), where s is the number of steps in the procedure. p is a vector of discrete, global latent variables in the model. This implies that each time-series will split into s segments, each segment corresponding to a single step in the procedure.

Next, let |Γ C | = K be the number of cluster labels (primitives). Each step in the procedure will correspond to one of these K labels i.e. p r ∈ [K], 1 ≤ r ≤ s. Thus, p is an ordered sequence of high-level primitives e.g. to make tea, we would boil water → get cup → pour water and so on. Formally, we assume a generation process,

λ r ∼ Dirichlet(α 1 , . . . , α K ), 1 ≤ r ≤ s p r ∼ Categorical(λ r ), 1 ≤ r ≤ s

Each token in the procedure is independently generated from a separate Dirichlet prior λ r , which allows fine-grained control over the distribution of tokens at each step. Note that we do not disallow self-transitions, which can be utilized by the model to coalesce adjacent segments.

Modeling realizations in each time-series. Different time-series can be vastly different in how they implement the procedure e.g. they can be videos of different individuals making tea. These individuals, while following the same procedure, may spend varying amounts of time on each step (or even skip steps).

In PRISM, the realization of the procedure in the i th time-series (having m i time-steps) is given by w i = (w i,1 , . . . , w i,mi ), a sorted sequence of discrete random variables. Each w i,j ∈ [s] is an index into the procedure p. Its value dictates the step being followed in the procedure at the j th time-step in the i th time-series. We assume w i is generated as follows,

γ ∼ Dirichlet(β 1 , . . . , β s ) u i,j ∼ Categorical(γ), 1 ≤ i ≤ n, 1 ≤ j ≤ m i w i = Sort(u i,1 , . . . , u i,mi ), 1 ≤ i ≤ n

Mathematically, to generate w i , we first generate an unordered bag of m i independent categorical variables u i,j and then sort them (in ascending order). Each u i,j is generated from the same Dirichlet prior γ. Intuitively, if we use the stick-breaking metaphor for the Dirichlet distribution, γ is a prior over how long each token in the procedure is -γ r is the expected relative length of p r . Thus, γ allows us to impose a prior on the relative length of each step of the procedure, in every time-series. An equivalent view into our model is that each w i is directly generated by a Multinomial distribution with m i draws on probabilities dictated by γ. However, we find that the view that we adopt yields a far more efficient inference scheme (see Appendix).

Generating observations. For every time-series, we can combine the common procedure p and the realization of the procedure w i to recover a temporal clustering z i for that time-series. We write this temporal clustering as z i = (z i,1 , . . . , z i,mi ) where z i,j ∈ [K] is a discrete primitive.

For each primitive k ∈ [K], we assume a Gaussian observation model (µ k , Σ k ) (with a Normal-Inverse Wishart prior) shared by all time-series. PRISM then generates each x i j independently from its local assignment z i,j and the appropriate observation model. We use Gibbs sampling (see Appendix) to perform posterior inference of the latent variables

({p r } r∈[s] , {u i,j } i∈[n],j∈[mi] , {µ k } k∈[K] , {Σ k } k∈[K]

) in the model.

5. EXPERIMENTS WITH EVALUATION CRITERIA:

We first compare our metrics to existing criteria by evaluating several competing methods on real datasets. We note that the methodology for validating our evaluation criteria is similar to that followed by Rosenberg & Hirschberg (2007) , with the exception that we compare criteria on real-world data with a large number of methods. We compare to homogeneity (Hom), completeness (Com), NMI since they give the best combination of performance for clustering evaluation (Rosenberg & Hirschberg, 2007) . We also compare to accuracy computed using the Munkres method since it is widely used e.g. (Sener & Yao, 2018; Zhou et al., 2008) and the adjusted Rand index (ARI) (Vinh et al., 2010) .

Methods. We use several methods to generate temporal clusterings: Hierarchical Clustering (AGG); Hidden Markov Models (HMM); Hierarchical Dirichlet Process HMMs (HDP-HMM) (Teh et al., 2005) ; Sticky HDP-HMM (SHDP-HMM) (Fox et al., 2008a) ; HDP Hidden Semi-Markov Models (HDP-HSMM) (Johnson & Willsky, 2013) ; Switching Linear Dynamical Systems (SLDS) (Fox et al., 2008b) ; HDP-SLDS (Fox et al., 2008b) ; Sticky HDP-SLDS (SHDP-SLDS) (Fox et al., 2008b) ; Ordered Subspace Clustering (OSC) (Tierney et al., 2014) ; Temporal SC (TSC) (Li et al., 2015) ; Low-Rank SC (LRSC) (Vidal & Favaro, 2014) ; Sparse SC (SSC) (Elhamifar & Vidal, 2013) and Spatial SC (SPATSC) (Guo et al., 2013) .

Datasets. We use 2 common benchmark datasets (Fox et al., 2008b; 2009; 2014; Zhou et al., 2008; 2013) . BEES consists of 6 time-series, each recording the dancing behavior of a bee. Each time-series is built up of 3 primitive movements -waggle, left-turn and right-turn -repeated several times. For BEES, we treat each sequence as a separate dataset (BEES 1 ,. . . ,BEES 6 ) and run each method on them individually, as done in Fox et al. (2008b) . Here, we focus on results on BEES 1 but results on MOCAP6, a motion-capture dataset and all BEES datasets can be found in the Appendix.

Results. For each of our metrics, we first illustrate that they qualitatively capture the characteristics that they were designed for. Fig 6e shows the result of evaluating with just the RSS -we can clearly see that SpatSC is able to capture far more repeated structure than SHDP-HMM, which is unable to label segments from the same ground-truth label consistently. Lastly, we see that SpatSC is judged to be the best method on BEES 1 by TSS in Fig 6f since it has the highest RSS and a SSS of 0.81 which is close to the highest.

Improving scores on evaluation metrics is common practice to demonstrate performance gains due to a method. We now highlight cases where TSS disagrees with NMI, ARI and/or Munkres in the assessment of competing methods, since this discrepancy has significant impact on which method is picked. We find that such disagreements happen very often across all datasets. An illustrative example for BEES 1 is shown in Fig 6h : TSS and Munkres agree that SSC is better than HDP-HSMM, but NMI and ARI strongly disagree. Before diving into why this may be happening, we note that examining the RSS and SSS allows us to understand why TSS considers SSC to be better than HDP-HSMM -it captures more repeated structure. On the other hand Munkres, while it agrees with TSS, is unable to give us this insight. We find that NMI is highly sensitive to the number of clusters -SSC only uses 2 clusters, which means it is heavily penalized by NMI, in spite of capturing some repeated structure. ARI is less sensitive than NMI, but is far less interpretable.

Another example is shown in Fig 6g , where the TSS for both SLDS and HDP-HMM is comparable. Interestingly, NMI/ARI and Munkres are all highly skewed and disagree strongly on the better method.

Since HDP-HMM has more clusters, the clusters are more homogenous and Hom is much higher, pushing up NMI. On the other hand, since Munkres finds a one-to-one correspondence between clusters and ground-truth labels, having more clusters means some cannot be matched. This causes the score to strongly downweigh suchh clusterings. Overall, our metrics provide both interpretability and the ability to analyze where the method is falling short.

Lastly, we conduct a sensitivity analysis for β as follows: assume that the most appropriate value for evaluation is β , but we instead use the default value of β = 1.0. Our main idea is to compute how well the ranking of methods under β = 1.0 captures that under β . To do this, we use the Normalized Discounted Cumulative Gain (NDCG), a common criteria from information retrieval. Varying β , we see in Fig 6i that β = 1.0 outputs similar rankings to a relatively wide range of other β values, other than at the extremes. This confirms that most settings of β will be captured well by our choice of β, unless there is a strong reason to prefer either the RSS or the SSS only.

6. EXPERIMENTS WITH PRISM:

We now evaluate PRISM on complex datasets with noisy features, including video data. For PRISM, we set α = 1.0, β = 0.1 (except BEES) and we set K (number of mixture components) to equal the number of ground-truth labels for all methods. Results are averaged over 10 random seeds. Experimental details including hyperparameters can be found in the Appendix.

Before our main results, we briefly conduct a simple check on performance for the BEES dataset introduced earlier. Though this only has a single example per procedure, we find that our approach still does well (Table 1 ) compared to two simple baselines -all approaches use the same observation model with Gaussian emissions.

Leveraging a fixed procedure. We now evaluate the benefit of our approach on the JIGSAWS dataset (Gao et al., 2014) , which consists of surgeons with varying skill demonstrating common procedures, such as suturing, needle-passing and knot-tying. We use kinematic features that correspond to the grippers controlled by the surgeon for all methods. For each surgical task, we select the 10 demonstrations that correspond to the 2 expert surgeons since they adhere most closely to the prescribed procedures. Even for a single expert, there is considerable variability across the 5 demonstrations they generate, so this is far from an idealized setting. Once again, all methods leverage the same observational model with the same features. Our baselines are two standard approaches -a Gaussian mixture model (GMM), and variants of a Bayesian Hidden Markov model.

As expected, PRISM outperforms these baselines by a considerable margin. On the suturing task, performance is improve by almost double-digits on both metrics, with performance gains also visible on the other tasks. Leveraging the common procedure allows PRISM to share statistical strength across demonstrations and make consistent predictions where possible, whereas the other methods must rely more strongly on the observation model to infer the latent sequence. Without very cleanly separated features for different primitives, they can often be led astray.

We also experiment with a smaller subset of 5 demonstrations provided by a single expert on the suturing task. This subset of data has the least amount of variability in terms of the expert's procedure.

Because of this, we see that PRISM gives an improvement of ∼ 18 points on NMI and ∼ 10 points on the TSS. To illustrate PRISM's sensitivity to hyperparameters, we vary s and K in Fig 7 on this single-expert suturing task. We find that PRISM's performance is quite stable with respect to boththis holds as long as s is larger than the number of segments in ground-truth, and K is close to the true number of ground-truth clusters.

Noisy procedures with video data. Next, we consider labeling demonstrations that are very noisy in their adherence to the procedure, with a great deal of variability in both the steps being followed, their order, as well as their relative lengths. Recent work by Sener & Yao (2018) achieved state-of-the-art performance without any form of supervision on two large video datasets -Breakfast actions (Kuehne et al., 2014) and INRIA instructional videos (Alayrac et al., 2016) . We expect our method to perform worse on the Breakfast actions dataset, as it has no repeated structure, in which case our model has no advantages over this prior model, and indeed uses a simpler observation model. Perhaps encouragingly, in the Breakfast actions dataset, we find that our method only performs slightly worse than but close to Sener & Yao (2018) , achieving a Munkres score of 33.5 compared to their score of 34.6. It is likely that this can be improved by tuning the number of segments for each procedure (we set s = 20 for all procedures).

However a key limitation of their work is that they do not model repeated structure. Therefore in domains with repeated structure we expect our approach to do well. Such structure is present in some of the activities in the INRIA dataset. We compare to Sener & Yao (2018) 's results on the INRIA dataset in Table 2 . PRISM outperforms their best method on 3 out of the 5 activities despite not explicitly modeling background information, with a significant difference on cpr and jump car.

On their comparable method which does not model background, PRISM is significantly better across all procedures. This is despite PRISM using a much weaker observation model (they simultaneously learn a feature representation for the videos) and no permutation model. An example for the cpr procedure is shown in Fig 1 and illustrates that unlike Sener & Yao (2018)'s method, PRISM can successfully recover repeated structure both within and across videos. We expect that incorporating a background model in PRISM could enable our approach to yield even stronger performance.

7. CONCLUSION:

We developed PRISM, a hierarchical Bayesian model to identify latent procedures in time-series data. Results on several datasets show that PRISM can be used to learn procedural abstractions more effectively than competing methods. We also introduced new evaluation criteria for the problem of externally evaluating temporal clusterings. Our metrics address gaps in temporal clustering evaluation, and provide both interpretability and the potential for easier analysis. In future work, we hope to extend PRISM by incorporating richer observation models, and relaxing the assumption that a single strictly-ordered procedure is followed by all time-series. Links to code can be found in the Appendix. 

B EXPERIMENTS WITH EVALUATION CRITERIA:

The results presented for the evaluation criteria developed in this paper can be reproduced using code available at https://github.com/StanfordAI4HI/ICLR2019_evaluating_ discrete_temporal_structure. The temporal clusterings predicted by each method, as well as ground-truth data is available at the link.

In addition, an implementation of only the newly proposed temporal clustering evaluation criteria can be found at https://github.com/StanfordAI4HI/tclust-eval.

For completeness, we include the result of evaluation across all methods, datasets and criteria. 

D.4 NON-MARKOV SIMULATION STUDY:

We generate 10 time-series using the following simple non-Markov sequence,

AABBAACCAADDAAEEAAF F AAGGAAHHAAAAAAAA:

There are 8 distinct latent tokens in the sequence, and we generate noisy Gaussian 2-dimensional observations for each such token, using a mean of (i, i) for the i th token, and a standard-deviation of 0.35 in either dimension.

We compare the Hidden Markov model variants to PRISM and the Gaussian mixture model. All methods use the same observational model with Gaussian emissions.

We expect that the main advantage of using an HMM vis-a-vis a GMM for recovering the latent temporal clustering, is in cases where the observations are noisy. In such cases, the learned transition dynamics can override the (noisy) observational likelihood to yield the correct sequence of tokens, as long as the procedure is Markov. However, this benefit should disappear in cases where the procedure is not Markov, like the sequence constructed by us above.

In contrast to these methods, PRISM makes no Markov assumption about the underlying procedure. This allows it to learn any procedure, Markov or otherwise. In addition, by sharing this procedure across multiple time-series, PRISM can avoid being led astray by noisy observations (unlike the GMM, which relies purely on the observation to infer the local latent token).

This argument is confirmed by the results in Table 12 . PRISM outperforms all methods on both evaluation metrics, while the HMM variants suffer due to their inability to abandon the Markov assumption that is baked into the model.

:



Figures and Tables:

(a) Prediction by PRISM for a cpr video. Note how repeated structure (cream primitive) is recovered. (b) Repeated structure captured by PRISM at different points in the video (3 left images) and across videos (right image).

(a) G is the ground truth temporal clustering and C1, C2, C3 are possible temporal clusterings. (b) C1, C2, C3 when converted to standard nontemporal clusterings for use with traditional metrics.

Figure 3: Problematic cases for traditional metrics, on a simple 6-step time-series.

Figure 2: From top -timeseries, temporal clustering, segments and procedure.

Fig 2 illustrates these ideas for a temporal clustering defined over the label set Γ = {A, B, C, D}. The procedure exhibits repeated structure since B has (atleast) 2 tokens/segments corresponding to it.

use Fig 3 as an illustrative example to highlight these shortcomings. In Fig 3, C 1 perfectly matches G, while both C 2 and C 3 split the cluster G

we care more about the segmentation quality. However, all 3 temporal clusterings C 1 , C 2 , C 3 in Fig 3 are assigned a perfect purity or homogeneity, since as seen in Fig 3b, no cluster in any of the 3 contains items from more than one ground-truth label. Moreover, no traditional clustering criteria can distinguish C 2 from C 3 since their cluster composition is identical, as seen in Fig 3b.

:7) = (1, 2) and running H yields heaviest sub-sequence (X) with score 3. C5:7 (dotted) does not contribute to the score.

Figure 4: Example for repeated structure scoring.

. Fig 4b uses this choice of H to score the pair of segments we picked -we run H ((X, X), (X, Z, Z)) and it returns a score of 3.

Figure 5: PRISM model

Fig 6 visualizes the best and worst performing methods on BEES 1 for each of our metrics. In Fig 6a, SSC is worst on LASS-O since it heavily over-segments the time-series by introducing many transitions, while TSC introduces far fewer transitions. In Fig 6b, under-segmentation is penalized and we can see that AGG should be worse since it contains longer segments and far fewer transitions than SSC. Combining these ideas in Fig 6c yields that HDP-HSMM has the best LASS -observe that it contains fewer closely spaced transition points and creates a more even segmentation at the right level of granularity, explaining its score.

(a) On LASS-O: SSC (worst) = 0.62, TSC (best) = 0.84 (b) On LASS-U: Agg (worst) = 0.85, SSC (best) = 0.98 (c) On LASS: SSC (worst) = 0.76, HDP-HSMM (best) = 0.86 (d) On SSS: SSC (worst) = 0.73, HMM (best) = 0.82 (e) On RSS: SHDP-HMM (worst) = 0.39, SpatSC (best) = 0.63 (f) On TSS: SHDP-HMM (worst) = 0.52, SpatSC (best) = 0.71 (g) Comparing SLDS to HDP-HMM on BEES1.(h) Comparing SSC to HDP-HMM on BEES1.

(i) Sensitivity to β.

Figure 6: (a)-(f) visualize the best and worst methods on our metrics; (g)-(h) show examples where we disagree with traditional criteria; (i) illustrates a sensitivity analysis with respect to the TSS tradeoff parameter β. All results are on the BEES 1 dataset.

Figure 7: Effect of varying hyperparameters on model performance. Results are on the JIGSAWS suturing task. Green marks the value of the hyperparameter in ground-truth.

.1 PRISM: PROCEDURE IDENTIFICATION WITH A SEGMENTAL MIXTURE MODEL We now describe a statistical model called PRISM for the procedure identification problem. PRISM is a generative hierarchical Bayesian model that explicitly separates reasoning about the procedure from its realization in each time-series in the dataset. It utilizes this to create a mixture distribution over primitives for every observation in the dataset. Fig 5 illustrates the model.

Table content: <table/>

Comparison with baselines on BEES and the JIGSAWS surgical dataset. For a dataset, values in parentheses indicate the number of time-series. Upper entries are TSS and lower entries are NMI.

Table content: <table><tr><td>Dataset</td><td/><td/><td/><td>Bees</td><td/><td/><td/><td/><td>Surgery</td><td/><td/></tr><tr><td>Method</td><td>Bees-1 (1)</td><td>Bees-2 (1)</td><td>Bees-3 (1)</td><td>Bees-4 (1)</td><td>Bees-5 (1)</td><td>Bees-6 (1)</td><td>Average</td><td>Knot-Tying (10)</td><td>Needle-Passing (10)</td><td>Suturing (10)</td><td>Suturing (5)</td></tr><tr><td>GMM</td><td>61.21 19.97</td><td>62.38 35.52</td><td>44.74 12.86</td><td>63.86 31.20</td><td>74.79 37.30</td><td>66.44 36.09</td><td>62.24 28.82</td><td>25.12 6.98</td><td>21.93 2.94</td><td>34.75 19.96</td><td>49.03 32.58</td></tr><tr><td>HMM (α = 0.1)</td><td>66.43 32.11</td><td>70.80 37.39</td><td>52.35 17.14</td><td>78.20 39.63</td><td>87.38 65.77</td><td>73.62 43.27</td><td>71.46 39.22</td><td>22.94 5.05</td><td>25.44 2.99</td><td>38.91 17.52</td><td>49.92 27.49</td></tr><tr><td>HMM (α = 1.0)</td><td>65.00 29.20</td><td>69.16 37.82</td><td>53.28 17.75</td><td>73.75 36.63</td><td>87.45 65.93</td><td>80.54 50.00</td><td>71.53 39.56</td><td>24.48 5.50</td><td>22.72 3.27</td><td>38.97 18.98</td><td>48.72 29.93</td></tr><tr><td>HMM (α = 100.0)</td><td>58.69 22.14</td><td>78.28 50.33</td><td>53.56 18.73</td><td>81.51 43.59</td><td>85.85 61.16</td><td>86.45 61.99</td><td>74.06 42.99</td><td>25.02 6.11</td><td>22.91 3.21</td><td>37.73 18.38</td><td>47.32 28.65</td></tr><tr><td>Prism</td><td>65.90 25.04</td><td>75.96 46.01</td><td>54.65 18.26</td><td>78.27 41.57</td><td>83.78 56.11</td><td>82.07 53.72</td><td>73.44 40.12</td><td>28.88 9.75</td><td>28.00 15.47</td><td>46.68 33.26</td><td>62.03 51.54</td></tr></table>

Comparison with Sener & Yao (2018) on the INRIA dataset with the Munkres score.

Table content: <table><tr><td/><td colspan="3">SENER (no background model) PRISM (no background model) SENER (background model)</td></tr><tr><td colspan="2">changing tire 25.0</td><td>31.2</td><td>33.9</td></tr><tr><td>coffee</td><td>22.1</td><td>26.9</td><td>29.0</td></tr><tr><td>cpr</td><td>18.1</td><td>28.1</td><td>24.9</td></tr><tr><td>jump car</td><td>5.8</td><td>22.4</td><td>15.0</td></tr><tr><td>repot plant</td><td>19.1</td><td>24.0</td><td>23.9</td></tr></table>

Results on BEES 1 across all evaluation criteria and methods.

Table content: <table><tr><td>method</td><td colspan="10">Agg HDP-HMM HDP-HSMM HDP-SLDS HMM LRSC OSC SHDP-HMM SHDP-SLDS SLDS SSC SpatSC TSC</td></tr><tr><td>metric</td><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr><tr><td>ARI</td><td>0.28</td><td>0.27</td><td>0.22</td><td>0.14</td><td>0.36</td><td>0.10 0.39</td><td>0.23</td><td>0.14</td><td>0.17 0.06</td><td>0.39 0.30</td></tr><tr><td>Com</td><td>0.25</td><td>0.24</td><td>0.25</td><td>0.13</td><td>0.30</td><td>0.13 0.44</td><td>0.25</td><td>0.13</td><td>0.16 0.05</td><td>0.44 0.34</td></tr><tr><td>Hom</td><td>0.24</td><td>0.46</td><td>0.51</td><td>0.13</td><td>0.30</td><td>0.13 0.28</td><td>0.50</td><td>0.13</td><td>0.16 0.03</td><td>0.28 0.22</td></tr><tr><td>LASS</td><td>0.84</td><td>0.78</td><td>0.86</td><td>0.82</td><td>0.86</td><td>0.83 0.84</td><td>0.84</td><td>0.82</td><td>0.84 0.76</td><td>0.84 0.86</td></tr><tr><td colspan="2">LASS-O 0.82</td><td>0.65</td><td>0.80</td><td>0.74</td><td>0.81</td><td>0.73 0.78</td><td>0.77</td><td>0.80</td><td>0.81 0.62</td><td>0.78 0.84</td></tr><tr><td colspan="2">LASS-U 0.85</td><td>0.97</td><td>0.92</td><td>0.91</td><td>0.92</td><td>0.95 0.92</td><td>0.94</td><td>0.86</td><td>0.88 0.98</td><td>0.92 0.87</td></tr><tr><td>Munk</td><td>0.61</td><td>0.39</td><td>0.41</td><td>0.58</td><td>0.68</td><td>0.50 0.64</td><td>0.39</td><td>0.58</td><td>0.59 0.47</td><td>0.64 0.60</td></tr><tr><td>NMI</td><td>0.24</td><td>0.33</td><td>0.36</td><td>0.13</td><td>0.30</td><td>0.13 0.35</td><td>0.35</td><td>0.13</td><td>0.16 0.04</td><td>0.35 0.27</td></tr><tr><td>RSS</td><td>0.55</td><td>0.46</td><td>0.40</td><td>0.53</td><td>0.57</td><td>0.45 0.63</td><td>0.39</td><td>0.49</td><td>0.45 0.49</td><td>0.63 0.62</td></tr><tr><td>SSS</td><td>0.78</td><td>0.75</td><td>0.81</td><td>0.76</td><td>0.82</td><td>0.79 0.81</td><td>0.81</td><td>0.76</td><td>0.78 0.73</td><td>0.81 0.81</td></tr><tr><td>TSS</td><td>0.65</td><td>0.57</td><td>0.54</td><td>0.62</td><td>0.67</td><td>0.57 0.71</td><td>0.52</td><td>0.59</td><td>0.58 0.59</td><td>0.71 0.70</td></tr><tr><td>method</td><td colspan="10">Agg HDP-HMM HDP-HSMM HDP-SLDS HMM LRSC OSC SHDP-HMM SHDP-SLDS SLDS SSC SpatSC TSC</td></tr><tr><td>metric</td><td/><td/><td/><td/><td/><td/><td/><td/><td/><td/></tr><tr><td>ARI</td><td>0.26</td><td>0.32</td><td>0.28</td><td>0.14</td><td>0.61</td><td>0.19 0.44</td><td>0.25</td><td>0.12</td><td>0.10 0.05</td><td>0.44 0.30</td></tr><tr><td>Com</td><td>0.27</td><td>0.33</td><td>0.28</td><td>0.16</td><td>0.54</td><td>0.19 0.51</td><td>0.29</td><td>0.13</td><td>0.11 0.04</td><td>0.51 0.32</td></tr><tr><td>Hom</td><td>0.28</td><td>0.66</td><td>0.53</td><td>0.16</td><td>0.54</td><td>0.20 0.33</td><td>0.58</td><td>0.20</td><td>0.10 0.02</td><td>0.33 0.41</td></tr><tr><td>LASS</td><td>0.79</td><td>0.77</td><td>0.87</td><td>0.83</td><td>0.89</td><td>0.79 0.84</td><td>0.84</td><td>0.83</td><td>0.82 0.76</td><td>0.84 0.87</td></tr><tr><td colspan="2">LASS-O 0.69</td><td>0.64</td><td>0.82</td><td>0.80</td><td>0.91</td><td>0.67 0.78</td><td>0.77</td><td>0.77</td><td>0.82 0.62</td><td>0.78 0.84</td></tr><tr><td colspan="2">LASS-U 0.94</td><td>0.96</td><td>0.92</td><td>0.87</td><td>0.88</td><td>0.96 0.91</td><td>0.93</td><td>0.90</td><td>0.82 0.98</td><td>0.91 0.90</td></tr><tr><td>Munk</td><td>0.64</td><td>0.38</td><td>0.44</td><td>0.56</td><td>0.85</td><td>0.60 0.68</td><td>0.40</td><td>0.37</td><td>0.49 0.49</td><td>0.68 0.54</td></tr><tr><td>NMI</td><td>0.27</td><td>0.47</td><td>0.38</td><td>0.16</td><td>0.54</td><td>0.20 0.41</td><td>0.41</td><td>0.16</td><td>0.10 0.03</td><td>0.41 0.37</td></tr><tr><td>RSS</td><td>0.51</td><td>0.60</td><td>0.42</td><td>0.44</td><td>0.77</td><td>0.40 0.73</td><td>0.57</td><td>0.40</td><td>0.41 0.49</td><td>0.73 0.50</td></tr><tr><td>SSS</td><td>0.75</td><td>0.75</td><td>0.82</td><td>0.76</td><td>0.86</td><td>0.76 0.80</td><td>0.79</td><td>0.76</td><td>0.74 0.74</td><td>0.80 0.82</td></tr><tr><td>TSS</td><td>0.61</td><td>0.67</td><td>0.56</td><td>0.56</td><td>0.81</td><td>0.52 0.76</td><td>0.66</td><td>0.53</td><td>0.53 0.59</td><td>0.76 0.62</td></tr></table>

Results on BEES 2 across all evaluation criteria and methods.

Table content: <table><tr><td>Method</td><td>NMI</td><td>TSS</td></tr><tr><td>GMM</td><td colspan="2">75.21 73.37</td></tr><tr><td>HMM (α = 0.1)</td><td colspan="2">47.68 72.50</td></tr><tr><td>HMM (α = 1.0)</td><td colspan="2">56.45 72.39</td></tr><tr><td colspan="3">HMM (α = 100.0) 59.32 51.31</td></tr><tr><td>Prism (β = 0.1)</td><td colspan="2">79.04 82.77</td></tr></table>

Comparison of methods in their ability to extract non-Markov procedural structure on 2-dimensional simulated data.

Table content: <table/>

Peer review (iclr2019):
 This paper describes two distinct contributions: a new compound criterion for comparing a temporal clustering to a ground truth clustering and a new bayesian temporal clustering method. Globally the paper is clear and well illustrated. 1) About the new criterion: *pros: * a) as clearly pointed out by the authors, using standard non temporal clustering comparison metrics for temporal clustering evaluation is in a way "broken by design" as standard metrics disregard the very specificity of the problem. Thus the introduction of metrics that take explicitly into account time is extremely important. b) the proposed criterion combines two parts that are very important: finding the length of the stable intervals (i.e. intervals whose instants are all classified into a single cluster) and finding the sequence of labels. *cons:* a) while the criterion seems new it is also related to criteria used in the segmentation literature (see among many other https://doi.org/10.1080/01621459.2012.737745) and it would have been a good idea to discuss the relation between temporal clustering and segmentation, even briefly. b) the reliance on a tradeoff parameter in the final criterion is a major problem: how shall one chose the parameter (more on this below)? The paper does not explore the effect of modifying the parameter. c) in the experimental section, TSS is mostly compared to NMI and to optimal matching (called Munkres here). Even considering the full list of criteria in the appendix, the normalized rand index (NRI) seems to be missing. This is a major oversight as the NRI is very adapted to comparing clusterings with different number of clusters, contrarily to NMI. In addition, the authors claim that optimal matching is completely opaque and difficult to analyse, while on the contrary it gives a proper way of comparing clusters from different clusterings, enabling fine grain analysis. 2) about the new model *pros*: a) as far as I know, this is indeed a new model b) the way the model is structured emphasizes segmentation rather than temporal dependency: the so called procedure is arbitrary and no dependency is assumed from one segment to another. In descriptive analysis this is highly desirable (as opposed to say HMM which focuses on temporal dependencies). *cons* a) the way the length of the segments in the sequence are generated (with sorting) this a bit convolved. Why not generating directly those lengths? What is the distribution of those lengths under the sampling model? Is this adapted? b) I find the experimental evaluation acceptable but a bit poor. In particular, nothing is said on how a practitioner would tune the parameters. I can accept that the model will be rather insensitive to hyper-parameters alpha and beta, but I*ve serious doubt about the number of clusters, especially as the evaluation is done here in the best possible setting. In addition, the other beta parameter (of TSS) is not studied. Minor point: - do not use beta for two different things (the balance in TSS and the prior parameter in the model)