Marco Corneli
email: marco.corneli@malix.univ-paris1.fr

Pierre Latouche
email: pierre.latouche@univ-paris1.fr

Fabrice Rossi
email: fabrice.rossi@apiacoa.org

Multiple change

Keywords: Dynamic networks, non homogeneous Poisson point processes, stochastic block model, variational EM, PELT

published or not. The documents may come

Introduction

The development of communication infrastructures has led to an unprecedented collection of data stored in the form of interactions between units of interest. For instance, interactions can correspond to email exchanges between employees of a company. They can also be characterized by posts between individuals via social media (Twitter, Facebook, Linkedin) or text messaging. In practice, interactions occur at specific time points, which are often recorded. Thus, if actor i interacts with actor j, at time η, the corresponding interaction can be represented by a triple (i, j, η). The set of all triples can then be used to build a dynamic graph where each actor is associated with a node, and an edge between two nodes (i, j) is present at time η if the corresponding interaction (i, j, ν) is recorded in the data set.

Statistical models for dynamic graphs are usually discrete in time, i.e. predefined time intervals are considered and interactions during those time intervals are aggregated to obtain snapshots. In the binary case, two nodes of a graph snapshot are connected if an interaction between them occurred during the corresponding time frame. The number of interactions between the nodes during the time frame can also be recorded. Most of these models consider the stochastic block model (SBM) [START_REF] Wang | Stochastic blockmodels for directed graphs[END_REF][START_REF] Nowicki | Estimation and prediction for stochastic blockstructures[END_REF] and extend it to the dynamic framework. SBM assumes that vertices are spread in latent clusters and that the probability for two nodes to connect depends on their respective clusters.

No assumption is made on the connection probabilities such that different types of clusters of nodes can be taken into account. In its original form, SBM cannot deal with dynamic interactions. [START_REF] Yang | Detecting communities and their evolutions in dynamic social networks a bayesian approach[END_REF] proposed a dynamic version of SBM by allowing the cluster of each node to switch at time t + 1 depending on its current state at time t. The switching probabilities are all characterized by a transition matrix. The alternative approach of [START_REF] Xu | Dynamic stochastic blockmodels: Statistical models for timeevolving networks[END_REF] characterizes the temporal changes through a state space model and relies on the Kalman filter along with the Rauch-Tung-Striebel smoother, for inference. Contrary to [START_REF] Yang | Detecting communities and their evolutions in dynamic social networks a bayesian approach[END_REF], the edge probabilities are seen as time varying parameters in [START_REF] Xu | Dynamic stochastic blockmodels: Statistical models for timeevolving networks[END_REF]. The work of [START_REF] Yang | Detecting communities and their evolutions in dynamic social networks a bayesian approach[END_REF] was generalized by [START_REF] Matias | Statistical clustering of temporal networks through a dynamic stochastic block model[END_REF] to deal with other types of edges. In their paper, they also showed that it was not possible to let both the connectivity parameters and cluster memberships vary over time without incurring into identifiability issues. Other static variants of the SBM model have also been adapted to the dynamic context. For instance, [START_REF] Xing | A state-space mixed membership blockmodel for dynamic network tomography[END_REF], [START_REF] Ho | Evolving cluster mixed-membership blockmodel for timeevolving networks[END_REF] and [START_REF] Kim | Nonparametric multi-group membership model for dynamic networks[END_REF] extended the mixed membership SBM of [START_REF] Airoldi | Mixed membership stochastic blockmodels[END_REF] in order to look for overlapping clusters of nodes, through time. Moreover, the dynamic random subgraph model (DRSM) [START_REF] Zreik | The dynamic random subgraph model for the clustering of evolving networks[END_REF] was built upon the RSM model [START_REF] Jernite | The random subgraph model for the analysis of an acclesiastical network in merovingian gaul[END_REF] to uncover clusters within subgraphs provided a priori. Note that the popular static latent position model of [START_REF] Hoff | Latent space approaches to social network analysis[END_REF] was also extended by [START_REF] Sarkar | Dynamic social network analysis using latent space models[END_REF] and [START_REF] Friel | Interlocking directorates in irish companies using a latent space model for bipartite networks[END_REF] to deal with dynamic interactions.

The scope of application of the models for dynamic graphs we have discussed so far is limited. Indeed, aggregating the data leads to a loss of information and the choice of the time intervals used to build the snapshots has a strong impact on the results [START_REF] Matias | Statistical clustering of temporal networks through a dynamic stochastic block model[END_REF]. In order to deal with dynamic interactions on a continuous time frame, a natural choice is to consider point processes. Thus, [START_REF] Matias | Statistical clustering of temporal networks through a dynamic stochastic block model[END_REF] relied on the so called (doubly stochastic) non homogenous Poisson point processes (NHPPP). Following a SBM like approach, nodes are assumed to belong to hidden clusters. Each pair of nodes is then associated to a NHPPP whose intensity function depends on the respective clusters. A variational expectation maximization (VEM) algorithm is finally employed to non parametrically estimate these functions and to uncover the clusters. This work is partially related to [START_REF] Dubois | Stochastic blockmodelling of relational event dynamics[END_REF] who relied on a parametric form for the intensity functions, which depend on the past network history and other predefined statistics. An alternative approach was proposed by Corneli et al (2016b) where the intensity functions of the Poisson point processes are assumed to be piecewise constant on predefined time intervals, each time interval belonging to an hidden time cluster. In their model, the value of each intensity function at time t not only depends on the clusters of nodes, but also on the corresponding time cluster. A greedy search algorithm was proposed to uncover the clusters.

In this paper, we extend the work of [START_REF] Matias | Statistical clustering of temporal networks through a dynamic stochastic block model[END_REF] to simultaneously uncover clusters of nodes sharing connection profiles, and to look for adjacent time intervals on which the connectivity patterns between pairs of clusters are stationary. In practice, considering dynamic interactions over a continuous time interval, we assume the intensity functions of the NHPPP to depend on the hidden node clusters and to be piecewise constant. Moreover, they are assumed to share D -1 common discontinuities whose location and number are unknown. These discontinuities induce a segmentation of the entire time interval over which the interactions are observed. In order to perform inference, a VEM algorithm is derived. We show that the V-M step can be tackled relying on a multiple changepoint detection tecnique for univariate time series, the pruned exact linear time (PELT) [START_REF] Killick | Optimal detection of changepoints with a linear computational cost[END_REF] method, which we adapted to our framework. Finally, the number of clusters of vertices is estimated using a Bayesian information criterion (BIC) involving variational approximations.

An alternative model to Poisson processes based ones has been proposed in [START_REF] Guigourès | A triclustering approach for time evolving graphs[END_REF][START_REF] Guigourès | Discovering patterns in time-varying graphs: a triclustering approach[END_REF]. As the model proposed in the present paper, this approach does not aggregate the temporal network prior analyzing it. It also produces both clusters of nodes and a segmentation of time. It is therefore a natural reference model in our context. By being based on the non parametric MODL approach [START_REF] Boullé | Data grid models for preparation and modeling in supervised learning, Microtome Casteigts A, Flocchini P, Quattrociocchi W, Santoro N[END_REF], the triclustering technique of [START_REF] Guigourès | A triclustering approach for time evolving graphs[END_REF][START_REF] Guigourès | Discovering patterns in time-varying graphs: a triclustering approach[END_REF] can handle non Poisson distributed interaction counts but it is somewhat blind to intensity changes as will be illustrated in Section 4.1.2. In particular, time segments cannot be detected when the type of connectivity structure is persistent through time but subject to parallel shifts in the interaction intensity levels.

The paper is organized as follows. In Section 2, we introduce the model and notations. In Section 3, we derive the VEM algorithm along with the model selection procedure. Finally, in Section 4, some experiments on simulated and real data are carried out to assess the proposed methodology.

A generative model for continuous time dynamic graphs

We consider interaction data with time-stamps and assume that a set of N actors {1, . . . , N } is provided. Those actors are assumed to interact (possibly repeatedly) at arbitrary times, that are recorded. A flat (i.e. non graphical) representation of those interactions is a finite set D = {(i m , j m , ν m )} 1≤m≤M , subset of {1, . . . , N } 2 × R + , in which a triple (i, j, ν) represents an interaction between actors i and j, at time ν. The temporal period under study is the interval [0, T ]. Without loss of generality, we assume that D is sorted on the time variable. In addition, each interaction time is supposed to be unique (see Section 2.3). Then 0 < ν 1 < ν 2 < . . . < ν M < T 1 .

In order to simplify the presentation, self interactions and directed interactions are not considered in this paper. Extensions to those cases are straightforward. The proposed approach makes the implicit assumption that interaction time spans do not play a significant roles in the general behaviour of the actors. Lifting this assumption and modelling explicitly the duration of the interactions is out of the scope of this paper.

Continuous time dynamic graph

Those interaction data can be seen as a dynamic (or time evolving) graph. The actor set corresponds to the set of nodes of the graph. Then each interaction (i, j, ν) can be seen as an edge in the graph, connecting i and j, and decorated by the time ν. This graph is in general a multiple graph in the sense that two nodes i and j can be connected by multiple edges (with different times). Another way of seeing the interaction data as a graph would be to follow [START_REF] Guigourès | Discovering patterns in time-varying graphs: a triclustering approach[END_REF] and to use the general framework of Casteigts et al (2012). In both cases, the most important aspect is that the proposed model is not based on a time series of graphs but rather on a fixed graph whose edges have a temporal dimension. This latter approach is more general as it allows one to model time as a continuous variable rather than a discrete one.

Node/actor clusters

Following the Stochastic Block Model (SBM) generative scheme [START_REF] Wang | Stochastic blockmodels for directed graphs[END_REF][START_REF] Nowicki | Estimation and prediction for stochastic blockstructures[END_REF], we assume that the connectivity in a dynamic graph is explained by hidden roles: each node/actor belongs to a cluster of nodes/actors and the way actors interact depends only on their respective clusters. In addition, cluster memberships are assumed to be fixed through time. Conversely, as we shall see in the following section, the connectivity parameters are allowed to change. These two modeling assumptions respect the identifiability principle for dynamic SBMs provided in [START_REF] Matias | Statistical clustering of temporal networks through a dynamic stochastic block model[END_REF] who showed that the connectivity parameters and group memberships should not both be allowed to change over time.

Thus, each node i is associated to a random variable Z i sampled from a multinomial distribution

P(Z i = k) = π k , ∀k ∈ {1, . . . , K},
where K is the number of clusters and

K k=1 π k = 1.
1 Interaction times are assumed to be strictly positive and ν M +1 is assumed to fall outside the time interval [0, T ]: ν M +1 ≥ T a.s.

We denote Z = {Z 1 , . . . , Z N } the set of all latent variables. Note that the 0-1 notation Z i = (Z i1 , . . . , Z iK ) with Z ik = 1 if node i belongs to the k-th cluster, 0 otherwise, will be used interchangeably with the notation where Z i ∈ {1, . . . , K}, when no confusion occurs. Finally, the vector of cluster proportions is denoted π = (π 1 , . . . , π K ).

Time-stamped interactions as point processes

Let us consider two fixed nodes, i and j. We denote M (i,j) the number of edges between i and j, in other words the number of distinct times ν such that (i, j, ν) ∈ D. Without loss of generality, those interaction times can be sorted into the following list A (i,j) := {ν

(i,j) 1 , . . . , ν (i,j) M (i,j) }, (1)
with ν

(i,j) 1 < ν (i,j) 2 < • • • < ν (i,j)
M (i,j) . In probabilistic terms, A (i,j) can be seen as a point process. As such a point process takes values in [0, T ], it is naturally associated to a counting process {M (i,j) (t)} t∈[0,T ] . The random variable M (i,j) (t) counts the number of interactions, between i and j, that happened before (or exactly at) t, i.e.

M (i,j) (t) = A (i,j) ∩]0, t] ,
where |S| denotes the cardinal of the set S. A simple yet flexible generative model for the interaction times in A (i,j) is the nonhomogeneous Poisson point process (NHPPP). Such a process is characterized by an intensity function κ (i,j) (•), positive and integrable on [0, T ]. Denoting

κ(i,j) (t) = t 0 κ (i,j) (s)ds, t ≤ T,
assuming that A (i,j) is generated by a NHPPP with intensity function κ (i,j) means that for all t ∈ [0, T ], M (i,j) (t) follows a Poisson distribution with parameter κij (t).

As pointed out in the previous section, the proposed model follows the SBM rationale in which the way actors interact depends only on their hidden role (a.k.a. cluster). Combined with the NHPPP assumption, this leads to the following assumptions:

1. given Z, the A (i,j) for all i > j are generated by N (N -1)/2 independent nonhomogeneous Poisson point processes with intensity functions {κ (i,j) (•)} i>j ; 2. there are

K(K + 1)/2 positive integrable functions λ = {λ kg (•)} k,g defined on [0, T ] such that κ (i,j) (t) = λ Z i Z j (t) for all t ∈ [0, T ].
With those assumptions, the conditional likelihood of a set of interactions A (i,j) between two nodes i and j is given by [START_REF] Daley | An introduction to the theory of point processes: volume I: Elementary Theory and Methods[END_REF] p(A (i,j)

|Z i = k, Z j = g, λ kg ) = M (i,j) m=1 λ kg (ν (i,j) m ) exp (-Λ kg (T )) , (2)
with

Λ kg (t) := t 0 λ kg (s)ds.
Notice that as the A (i,j) are assumed to be generated by Poisson point processes, all the interaction times are distinct (almost surely). This justifies the assumption used at the beginning of the present section.

The data set D can be seen as the union of all the A (i,j) , with the added information of the interacting pairs at each interaction time: the pair (i m , j m ) corresponds to the nodes actually having an interaction at time ν m . Combining equation ( 2) applied to all those processes with the conditional independence assumption between them, we obtain the following complete data likelihood p(D, Z|λ, π) = p(D|Z, λ)p(Z|π)

= exp   - N j>i Λ Z i Z j (T )   M m=1 λ Z im Z jm (ν m ) N i=1 π Z i .
(3)

Common change points

The generative model proposed above assumes some structure on the node set but none on the time interval. While this approach can be interesting in some settings (Corneli et al, 2016a), it reduces to some extent interpretation possibilities.

In particular, such a flexible model does not emphasize brutal changes in the interaction structure between clusters of nodes.

In order to emphasize those aspects, we assume that the intensity functions of the NHPPPs are piecewise constant and that the D intervals on which they are constant are shared between all functions. In other words, we assume that there are D -1 change points

0 = η 0 < η 1 < • • • < η D-1 < η D = T, (4)
such that for all 0

≤ d ≤ D, 1 ≤ k ≤ K and 1 ≤ g ≤ K, λ kg (•) is constant on [η d , η d+1 [. Therefore λ kg (t) = D d=1 λ kgd 1 [η d-1 η d [ (t), ∀k, g ∈ {1, . . . , K}, (5)
where λ kgd := λ kg (η d-1 ) and 1 G (.) is the indicator function over a set G 2 . In the following, η = {η 1 , . . . , η D-1 } denotes the set of change points and λ is the (K × K × D) tensor 3 with elements λ kgd .

A crucial consequence of the piecewise constant assumption is that on an interval [η d , η d+1 [, all the Poisson point processes are homogeneous and thus do not exhibit any temporal structure. On the contrary, the intensity is allowed to change arbitrarily from one interval to the next one, allowing brutal changes in the interaction patterns between clusters. In addition, taking into account the piecewise constant assumption, allows us to simplify the complete data likelihood as shown in the following proposition.

2 Since the whole time interval [0, T ] is considered, the intensity functions can be assumed (exceptionally) left continuous in t = T . In formulas: λ kg (T ) = λ kg (η d-1 ).

3 We use the same notation to denote the set of K(K + 1)/2 intensity functions and the tensor because under our assumptions they correspond to two different views of the same object. Notice that the frontal slices of λ are symmetric K × K matrices since we are dealing with undirected graphs.

Proposition 1 Using the constraints in equation (5), the complete data log likelihood becomes

log p(D, Z|θ) = - D d=1 K k,g   λ kgd ∆ d   N j>i Z ik Z jg   -log(λ kgd )   N j>i Z ik Z jg X (d) ij     + N i=1 K k=1 Z ik log π k , (6) where 1. θ := {η, λ, π} 2. ∆ d is the size of the interval [η d-1 , η d [, 3. X (d) ij := N (i,j) (η d ) -N (i,j) (η d-1 ) is the increment of conditional Poisson point process {M (i,j) (t)} t∈[0,T ] over the segment [η d-1 , η d [, i.e.
the number of interactions that occurred between i and j during the time interval

[η d-1 , η d [. Proof See appendix A.1.

Estimation

This section focuses on the inference of the model proposed above. This involves estimating the number K of clusters and the number D of adjacent time intervals, as well as the parameter θ. Therefore, we introduce first a penalized likelihood criterion for model selection. A variational approximation for this criterion is then introduced, leading to a variational expectation maximization (VEM) algorithm. It is finally shown how to integrate an efficient change point detection algorithm in the VEM loop to estimate the piecewise constant intensities.

Penalized likelihood

Given the set of all observed interactions D, our goal is to estimate the clustering structure in the actor/node space and the segmentation structure in the temporal space. In particular, we aim at choosing the number K of clusters and the number of change points (by choosing D) as well as their location η.

A natural quality measure in this context is the observed data (integrated) log-likelihood

log p(D|K, η, D) = log λ,π p(D, λ, π|K, η, D)dλdπ . (7)
Unfortunately this marginal log likelihood does not have an analytical form so we propose to replace it with a penalized likelihood (BIC like) term

log p(D|K, η, D) = max λ,π log p(D|K, η, D, λ, π) - 1 2 C(K, D) log α, (8)
where

C(K, D) := K -1 + K(K + 1)D 2 ,
accounts for the number of model parameters. The term α in equation ( 8) is related to the number of observations and will be discussed in Section 3.4.2. In other words, we consider the following optimization problem, for inference max

K,η,D,λ,π log p(D|K, η, D, λ, π) - 1 2 C(K, D) log α. (9)
While this allows to avoid to consider directly the log-likelihood in equation ( 7), two difficulties remain: log p(D|K, η, D, λ, π) is not directly calculable and the optimization over η and D is a complex non smooth problem. To tackle these issues, two strategies are proposed. A variational approach is used in order to derive a lower bound of log p(D|K, η, D, λ, π) (see Sections 3.2 and 3.3) while a change point detection technique is considered to address the optimization over η and D (see Section 3.4).

A variational bound

The first difficulty mentioned above is that computing the log likelihood term log p(D|K, η, D, λ, π) is not tractable. Indeed, it involves summing over all the

K N possible sets Z, i.e. log p(D|K, η, D, λ, π) = log Z p(D, Z|K, η, D, λ, π) .
Moreover, p(Z|D, K, η, D, λ, π) cannot be factorized so the standard expectation maximization (EM) algorithm (A. P. [START_REF] Dempster | Maximum likelihood from incomplete data via the em algorithm[END_REF] cannot be considered for inference. For more details on these issues, see [START_REF] Daudin | A mixture model for random graphs[END_REF] in the case of the standard SBM.

Therefore, we introduce an approximate distribution q(Z) for Z and use the standard variational decomposition

log p(D|K, η, D, λ, π) = L (q(Z); K, η, D, λ, π) + KL (q(Z)||p(Z|D, K, η, D, λ, π)) ,
where

L (q(Z); K, η, D, λ, π) = Z q(Z) log p(D, Z|K, η, D, λ, π) q(Z) = E q(Z) log p(D, Z|K, η, D, λ, π) q(Z) ,
and KL denotes the Kullback-Leibler divergence between the true and approximate posterior distribution q(Z) of Z, given the data and model parameters

KL (q(Z)||p(Z|D, K, η, D, λ, π)) = -E q(Z) log p(Z|D, K, η, D, λ, π) q(Z) .
Since log p(D|K, η, D, λ, π) does not depend on the distribution q(Z), maximizing L with respect to q(Z) is equivalent to minimizing the KL divergence (over q(Z) also).

As the KL divergence is non negative, L (q(Z); K, η, D, λ, π) is obviously a lower bound of log p(D|K, η, D, λ, π) for all q(Z), thus we use the standard variational bounding method. In our case, it consists in replacing the problem (9) by max

K,η,D,λ,π,q(Z) f (q(Z), K, η, D, λ, π) , (10)
where

f (q(Z), K, η, D, λ, π) = L (q(Z); K, η, D, λ, π) - 1 2 C(K, D) log α. ( 11
)
As in [START_REF] Daudin | A mixture model for random graphs[END_REF] for the SBM, the Z distribution is approximated by a factorized distribution

q(Z) = N i=1 q(Z i ) = N i=1 K k=1 τ Z ik ik . ( 12
)
This leads to the following expression for L.

Proposition 2 If q(Z) is of the form used in equation ( 12), then

L (q(Z); K, η, D, λ, π) = - D d=1 K k,g   λ kgd ∆ d   N j>i τ ik τ jg   -log(λ kgd )   N j>i τ ik τ jg X (d) ij     + N i=1 K k=1 τ ik log π k τ ik . (13)
Proof This expression can easily be obtained by taking the expectation with respect to q(Z) of the log likelihood in equation ( 6) and by adding the following entropy term H(q) := -E q(Z) [log(q(Z))].

Variational expectation maximization

Problem ( 10) is solved relying on a VEM algorithm which optimizes the function f (q(Z); K, η, D, λ, π) which respect to η, D, λ, π, and with respect to q(Z), alternately. The number K of clusters is handled via a grid search (see Section 3.5) and it is considered fixed in the present section. Given η, D, λ and π, the optimization with respect to q(Z) is straightforward. This corresponds to the E step of the algorithm (see Section 3.3.1). Then, given q(Z), the parameters λ as well as π are optimized away, and we use a change point detection procedure to optimize the criterion with respect to η and D. The optimization with respect to η, D, λ and π corresponds to the M step of the algorithm. The E and M steps are then iterated until convergence.

Maximization with respect to τ (E step)

The E step is based on the following proposition.

Proposition 3 A first order condition for f (q(Z), K, η, D, λ, π) to be maximal with respect to q(Z) is

τ ik = π k C exp    - D d=1 K g=1   λ kgd ∆ d   N j =i τ jg   -log(λ kgd )   N j =i τ jg X (d) ij        , where C = K k=1 τ ik , ∀i ∈ {1, . . . , N }, k ∈ {1, . . . , K}.
Proof See Appendix A.2.

In the E step of the algorithm, the entries of the set τ = {τ 1 , . . . , τ N } are updated in turn until convergence of τ . This corresponds to a fixed point procedure, as in [START_REF] Daudin | A mixture model for random graphs[END_REF] for instance. We emphasize that τ ik is the (approximate) posterior probability for node i to be in cluster k, given the data and model parameters. Thus, the clustering structure uncovered by the method is encoded through τ .

Maximization with respect to π

Notice that f (q(Z), K, η, D, λ, π) can be written

f (q(Z), K, η, D, λ, π) = N i=1 K k=1 τ ik log π k τ ik + g(q(Z), K, η, D, λ),
where the function g does not depend on π. Thus, q(Z) and K being fixed, maximizing f with respect to η , D, λ, π can be done independently on π and on the other parameters. The estimated value for π (under the constraint

K k=1 π k = 1) is then πk = N i=1 τ ik N
, ∀k ∈ {1, . . . , K}.

Maximization with respect to λ

Maximizing f with respect to λ leads to the following estimates

λkgd =            N j>i τ ik τ jg X d ij ∆ d N j>i τ ik τ jg when g > k, N j =i τ ik τ jk X d ij ∆ d N j =i τ ik τ jk when g = k. ( 14
)
Notice that contrary to π, λ does depend on η and D, which are also considered in the M optimization step.

Maximization with respect to η and D

Obviously, we have max

η,D,λ,π f (q(Z), K, η, D, λ, π) = max η,D max λ,π f (q(Z), K, η, D, λ, π) .
Thus, q(Z) and K being fixed, f can be evaluated in λ and π obtained in the previous sections. In detail

max λ,π f (q(Z), K, η, D, λ, π) = D d=1 G([η d-1 , η d [)- 1 2 K(K + 1)D 2 log α+const, ( 15
)
where all the terms which do not depend on η and/or D have been absorbed into the constant const and

G([η d-1 , η d [) := - K k,g  λ kgd ∆ d   N j>i τ ik τ jg   -log( λkgd )   N j>i τ ik τ jg X (d) ij     . ( 16

Segmentation

Up to a change of sign4 , the criterion (15), to be maximized with respect to η and D, has the general form used in change point detection problems (see equation (1) in [START_REF] Killick | Optimal detection of changepoints with a linear computational cost[END_REF], for instance) and efficient algorithms used in this context can be considered for inference.

Dynamic programming

Let us first recall that those types of additive costs can be solved exactly via a form of dynamic programming [START_REF] Jackson | An algorithm for optimal partitioning of data on an interval[END_REF]. Indeed, denoting F (s, D) the maximal value of the criterion (15), when restricted to interactions that happen in ]0, s], and for at most D segments, we have, as in [START_REF] Jackson | An algorithm for optimal partitioning of data on an interval[END_REF]; Killick et al ( 2012)

F (T, D) = max η,D ≤D   D d=1 G([η d-1 , η d [) - 1 2 K(K + 1) 2 log α   , = max ζ max η ,D ≤D-1   D d=1 G([η d-1 , η d [) - 1 2 K(K + 1) 2 log α   + G([ζ, T [) - 1 2 K(K + 1) 2 log α , = max ζ F (ζ, D -1) + G([ζ, T [) - 1 2 K(K + 1) 2 log α . (17

Restriction on the change point locations

There are two limitations in equation ( 17 Therefore, one has to reduce the search space for the change points to a finite set. In practice, this corresponds to introducing a grid of points which are a priori change point candidates:

P = {(t 0 , . . . , t U )|0 = t 0 < t 1 < • • • < t U = T, U ∈ N * }.
A natural choice for P is the set of all interaction times in D, but other choices can be used, such as intermediate times between interactions (e.g., times of the form

ν m+1 +ν m 2
) or arbitrary regular grids. Notice that choosing a grid immediately solves the problem of choosing the maximal value for D: it is exactly U .

The choice of P has several consequences. Firstly, the computational cost of the dynamic programming (with or without pruning) is directly linked to U (see below for details). Secondly, P acts as a minimal time resolution constraint. Thus, a high value of U allows to pinpoint change points very precisely but at a high computational cost, and vice versa. Therefore, both computational and expert considerations should be taken into account for choosing P. If the computational load is acceptable, using the set of all interaction times offers a maximal resolution, but that might emphasize unneeded details.

The last consequence of the choice of P concerns the value of α in the penalized log likelihood (see equation ( 8)). According to the hypotheses on the generative model, we observe N (N -1)/2 conditionally independent NHPPP. Each of those processes is observed on the intervals [t u , t u+1 [ via interaction counts. Those interaction counts are independent per the general hypothesis on NHPPP. Thus, we have α = U N (N -1)/2 independent observations.

Pruned exact linear time

In this paper, we consider the pruned exact linear time (PELT) method of [START_REF] Killick | Optimal detection of changepoints with a linear computational cost[END_REF] originally developed for change point detection in univariate time series. Using the recursive decomposition in ( 17) for maximization has a cost of O(K 2 U 2 ). For more details, we refer to [START_REF] Killick | Optimal detection of changepoints with a linear computational cost[END_REF]. The key observation is that for each point t u in P the gain of setting t u as the last change point before t u , has to be computed, for all t u < t u . Fortunately, some points t u can be pruned through the optimization routine. This is the principle of PELT and in practice it allows to speed up the exploration of the segmentation space. We use the following result whose general statement and formal proof are given in [START_REF] Killick | Optimal detection of changepoints with a linear computational cost[END_REF].

Again, consider t u and t u such that t u < t u and t u is not the last change point before t u . Then

F (t u , U -1) + G([t u , t u [) - 1 2 K(K + 1) 2 log α < F (t u , U ).
Moreover if t u fulfils the condition

F (t u , U -1) + G([t u , t u [) < F (t u , U ),
then t u can never be the optimal last change point prior to t u , for all t u > t u . This statement is true for all gain functions satisfying the following condition

G([t u , t u [) + G([t u , t u [) ≥ G([t u , t u [), t u < t u < t u . (18)
Proposition 4 The condition (18) is fulfilled by the gain function G(•) we defined.

Proof See appendix A.3.

The pseudocode Algorithm 1 illustrates how the PELT algorithm works.

Selection of K and initialization clusters

As any EM like approach, the algorithm proposed for inference depends on some initializations. Notice, however, that once initial values for τ and K are provided, the other model parameters (D, θ) = (D, η, λ, π) are estimated in the maximization (M) step and those estimates are employed in the E step to obtain a better estimate of τ and so on until the criterion f (q(Z), K, η, D, λ, π) (11) no longer increases (i.e. convergence). For a fixed value of K the initialization of τ can be obtained in several ways. For example a N × N adjacency matrix can be built by aggregating all the interactions over the time interval [0, T ], hence the entry (i, j) of this adjacency matrix corresponds to M (i,j) , with the previous notations. Then, clustering algorithms like k-means, hierarchical clustering or spectral clustering can be used to get an estimate of Z. Finally the initial matrix τ is built such that τ ik is one if Z i = k, zero otherwise. Another method to initialize τ consists in applying a k-means clustering on the rows of the N × U N matrix corresponding 5 Notice that in case the minimal partition is used, X

ij is trivially equal to one for the pair (i, j) interacting at t u-1 and zero for all the other pairs. Algorithm 1 PELT for dynamic SBM Input:

A grid 0 = t 0 < t 1 < • • • < t U = T . An (N × N × U ) tensor X whose entry 5 (i, j, u) is X (u) ij
A matrix τ of variational probabilities.

The penalty α = U N (N -1)/2.

A fixed positive number of clusters K.

The gain function G(.).

Initializations: F (0) = 1 4 K(K + 1) log α, cp(0) = N U LL, R 1 = 0. for η * in 1, . . . , U do Calculate F (η * ) = max η∈R η * [F (η) + G([t η+1 , t η * [) -F (0)]. Let η =argmax η∈R η * [F (η) + G([t η+1 , t η * [) -F (0)]. Set cp(η * ) = [cp(η), η]. Set R η * +1 = {η ∈ R η * ∪ {η * } | F (η) + G([t η+1 , t η * [) ≥ F (η * )}.

end for

Output: The change points stored in cp(U ).

to the mode-1 unfolding of the tensor X (see [START_REF] Kolda | Tensor decompositions and applications[END_REF], for more details). In each experiment in Section 4, all the mentioned initialization techniques are attempted. The one leading to the highest finale value of the lower bound is finally retained.

So far, we have assumed that the number of clusters was fixed. However, in practice K is unknown and has to be inferred from the data. Again, we rely on the criterion defined in equation ( 11) which involves a penalization term over K. Recalling that the optimal number of segments is selected by the PELT procedure (see Section 3.3), the VEM algorithm described in this section is run for different values of K in {1, . . . , K max }, for some fixed K max , and the value K maximizing the criterion is retained. The pseudocode Algorithm 2 summarizes the whole estimation routine.

Experiments

Simulated datasets

Some experiments on simulated data are carried out to test the proposed approach. Our model (called hereafter PELT-Dynamic SBM) is compared with the triclustering approach proposed in [START_REF] Guigourès | A triclustering approach for time evolving graphs[END_REF][START_REF] Guigourès | Discovering patterns in time-varying graphs: a triclustering approach[END_REF]. The approach is referred to as MODL (while MODL is a more generic technique [START_REF] Boullé | Data grid models for preparation and modeling in supervised learning, Microtome Casteigts A, Flocchini P, Quattrociocchi W, Santoro N[END_REF]). As pointed out in the introduction, this method is non parametric and search for clusters of nodes and time segments. It is based on a combinatorial generative

Algorithm 2 VEM algorithm

Input:

A maximum number of clusters K max .

A set D of M interaction times.

The criterion f (•) in equation ( 11).

Initializations: Store ← vector(K max ), Pmts ← list(K max ) for K in 1, . . . , K max do τ ← Some clustering algorithm {D, θ} ← Maximization(τ, D) while f (•) increases do τ ← Expectation(D, D, θ) {D, θ} ← Maximization(τ, D) end while Store[K] ← f (K, D, τ, θ) Pmts[K] ← {τ, θ} end for K * ← argmax(Store). Output: The estimated parameters in Pmts[K * ].
model estimated via a maximum a posteriori approach. As our model, it has no user tunable parameter and is therefore fully automated.

First scenario

The experiments considered in this section are related to the simulation setup in Section III.A of [START_REF] Guigourès | A triclustering approach for time evolving graphs[END_REF]. Each network generated is made of 40 nodes, grouped in four clusters: 5 vertices are in clusters 1 and 2, 10 vertices in cluster 3 and 20 vertices in cluster 4. The time interval [0, 100] is split into four segments (I 1 = [0, 20[, I 2 = [20, 30[, I 3 = [30, 60[ and I 4 = [60, 100[) and each segment is associated with a specific interaction pattern between clusters, as illustrated in Figure 1. For each number of edges varying from 50 to 10000, 50 1. A node and a random interaction time are sampled uniformly in {1, . . . , 40} and [0, 100], respectively. The vertex is then assigned to its cluster and the interaction time to its segment. 2. If the cluster of the selected vertex is connected to one or more clusters over the considered time segment (see Figure 1), a second vertex is sampled uniformly at random in the union of these clusters, and an edge is generated. 3. The first two steps are repeated until the desired number of edges is reached. 4. Finally, 30% of edges are rewired uniformly at random. The only difference between the current setup and the one used in [START_REF] Guigourès | A triclustering approach for time evolving graphs[END_REF], is that the networks in the present paper are undirected. Notice that the generative process for those data is neither Poisson based model nor the combinatorial model used by MODL.

For estimation purposes, a regular grid P (Section 3.4.2) with unitary length time intervals is used for PELT-Dynamic SBM. Both algorithms (PELT-Dynamic SBM and MODL) are applied to the generated interaction data and results are assessed at an aggregated level (cluster numbers) and at a more refined level relying on the adjusted Rand indexes [START_REF] Rand | Objective criteria for the evaluation of clustering methods[END_REF].

In Figure 2 As it can be seen, MODL provides more accurate estimates of both K and D for a small number of edges, while PELT-Dynamic SBM needs denser networks to recover the true number of clusters and time segments. These results are confirmed in Figures 3 and4 For each number of edges, adjusted Rand indexes (ARI) are computed to asses the quality of the estimates provided for Z and η by the two models. For what concerning the change point locations, when no change point is detected the ARI is zero, conversely when η1 = 20, η2 = 30 and η3 = 60 the ARI is one. For each number of edges, two box-and-whiskers plots are produced, one for MODL (on the left hand side) and the other for our method (right hand side).

As pointed out above, the data are not generated according to our model or to MODL combinatorial one. However, our model is still parametric which might explain the less accurate estimates provided here as compared to MODL. We show in the following sections that when the data are generated with a model closer to our model, the results are quite different.

Second scenario

The graphs generated in the second simulation scenario are made of 75 nodes, grouped into two clusters and undirected interactions are simulated over the time interval [0, 10]. This interval is split into three segments

I 1 = [0, η 1 [, I 2 = [η 1 , η 2 [, I 3 = [η 2 , 10[,
where the change points η 1 as well as η 2 are set to 2.1 and 6.9, respectively. Interactions are sampled by thinning [START_REF] Lewis | Simulation of nonhomogeneous poison processes by thinning[END_REF]

λ Z i Z j (t) = 0.11 I 1 (t) + 0.21 I 2 (t) + 0.051 I 3 (t) if Z i = Z j 0.051 I 1 (t) + 0.11 I 2 (t) + 0.0251 I 3 (t) if Z i = Z j ,
for all j > i. Thus, the IFs define a persistent community structure (see [START_REF] Fortunato | Community detection in graphs[END_REF], for more details) through time in which the intensity of the interactions within clusters is twice the intensity of the interactions between clusters. The following sampling procedure is used to sampled 50 dynamic graphs:

1. Each vertex is assigned to one of the two clusters with probability 1/2. 2. Interactions between each pair of nodes (i, j) are sampled according to the NHPPP with IF λ Z i ,Z j (t).

Again, to introduce some noise in the data, 10% of edges are rewired. Considering this new setting, we aim at evaluating the clusters uncovered by our methodology along with the estimates of the change point locations η 1 and η 2 . The grid P of all observed interactions in D was considered for inference. As mentioned in Section 3.4.2, this allows to pinpoint change points more accurately. The use of such a grid is made possible here because of the limited number of interactions generated.

Regarding the clustering task, we present in Figure 5 the results for MODL and a k-means algorithm applied on the adjacency matrix where all interactions over the time interval [0, 10] are aggregated. Note that contrary to MODL and PELT-Dynamic SBM, k-means was provided with the true number of clusters. In terms of change point detection, due to the particular generative structure, MODL cannot recover any time cluster and considers the dynamic graph as stationary. Indeed in order to avoid parametric assumptions, MODL uses rank based modeling for numerical values. In the triclustering context of [START_REF] Guigourès | A triclustering approach for time evolving graphs[END_REF][START_REF] Guigourès | Discovering patterns in time-varying graphs: a triclustering approach[END_REF] this means that interaction times are replaced interaction orders. This explains why MODL is blind to the time structure in the present scenario.

Conversely PELT-Dynamic SBM always retrieves the right number of change points in the data. The change point estimates can be observed in Figure 6 and Kernel density estimates are plotted along with the true change points as red vertical lines. This illustrates the accuracy of the proposed estimation procedure and its superiority to MODL in this situation.

Third scenario

This section aims at illustrating that aggregating interactions can lead to an important loss of information. Thus, each graph generated is made of 100 nodes clustered in two groups, with 50 nodes each. Moreover, the time interval [0, 12] is split into four segments of equal size delimited by the change points η 1 = 3, η 2 = 6 and η 3 = 9. Finally, interactions are sampled by thinning according to the

Λ 11 (T ) = Λ 12 (T ) = Λ 21 (T ) = Λ 22 (T ) = 3.8.
Thus, the average number of interactions is the same for all pairs of clusters which makes clusters indistinguishable when aggregating the interactions over the time interval. As for the previous simulation scenarios, 50 dynamic graphs are generated and 10% of edges of each graph are rewired uniformly at random. MODL as well as PELT-Dynamic SBM are then used to uncover clusters of vertices and to segment the time interval. A regular grid P with unitary length time intervals is considered.

The results are presented in Figure 7 as ARIs for both the change points and the cluster memberships. In this context, PELT-Dynamic SBM provides more reliable estimates than MODL, which fails in retrieving any cluster or temporal structures, in most cases. It's a form of extreme blindness to the whole structure of the data induced by the blindness to the temporal structure. Note that the results for PELT-Dynamic SBM are similar when relying on a grid with a higher time resolution.

Summary

The scenarios studied in this section show that PELT-Dynamic SBM is able to recover both the cluster and the temporal structure of dynamic graphs without the need for a strong prior aggregation of the interactions (a minimal aggregation is used for computational reasons in some of the experiments). Contrarily to MODL, PELT-Dynamic SBM discovers structural changes that are only based on a modification of the intensity of the interaction. Thus both approaches have different use cases. In particular, the temporal structure of the data is more easily captured by our model than by MODL.

Real data

We now focus on a cycle hire usage dataset, collected in London and publicly available at http://api-portal.tfl-gov.uk/docs. It characterizes the interactions that occurred on September 9, 2015, between the Santander stations. The dynamic network considered is made of 735 nodes and 64514 undirected edges (with no self loops), collected with a minute precision over the day. One edge connecting nodes i and j at a given time corresponds to a cycle hire from station i to station j or, conversely from station j to station i. To limit the computational burden of the segmentation step, we relied on a regular grid P corresponding to 96 time intervals of 15 minutes. The PELT-Dynamic SBM was then applied several times, for different values of K from 0 to 20. The highest value of the criterion f (equation 11) was attained for K = 11 clusters and D = 5 time segments. MODL does not find any temporal structure in those data despite obvious changes in the aggregated intensities (see Figure 9). Results are presented in Figure 8. The Santander stations are plotted on a London map6 different symbols/colors corresponds to different clusters identified by the model. Interestingly (and as expected), generally nearby stations are placed in the same cluster and the geographical distance between them plays a key role. In Figure 9, an histogram of the interaction times in the whole network is provided. Two peaks are visible around 8.30 and 18.30. The five segments detected by the methodology are delimited by the vertical red lines in the figure. We observe a strong alignment between the histogram and the segments estimated. In particular, the two observed picks are clearly associated with segment 2 and 3, respectively. We now concentrate on analyzing the results for cluster 3 (identified by the symbol + in Figure 8) which is made of stations from central London. This cluster actually corresponds to a community with higher values for the IF within the cluster than for all the IFs describing the levels of interactions with other clusters. Figure 10 gives some examples of intra and inter IFs related to cluster 3. The results are presented for clusters 1 ( ) and 7 ( ) which are geographically adjacent to cluster 3, and for clusters 4 (×) and 10 (⊗), which are not. Overall, as mentioned already, the intra IF is higher. However, this figure also highlight a temporal pattern. Indeed, it appears the inter IFs for the adjacent clusters are higher in the morning and in the evening than for the rest of the day. A somehow similar pattern is observed for clusters 4 and 10 but with much lower values in general. This is coherent with cycles being hired more often to go to a station close geographically.

In order to highlight another feature of the proposed methodology, we now give some results regarding cluster 8. More specifically, Table 1 provides the aggregated interactions between clusters 7 ( ) and 8 ( * ), over the segments uncovered. In the first segment, 16 interactions occurred between vertices of cluster 8 and 47 between vertices of cluster 8 and vertices of cluster 1. Thus, cluster 8 has a disassortative connectivity pattern with less intra edges. Conversely, in segment 2, they are more intra edges (502) within cluster 8 than between clusters 8 and 1. This corresponds to a community pattern. Looking through all the segments, we can observe that the community and disassortative pattern for cluster 8 alternates through time. Thus, clustering the vertices while detecting change points in the intensity of the interactions is mandatory here because the connectivity patterns of the data set keep changing. Therefore, any method aggregating the data would miss important information present in the data. In conclusion, the uncovered clusters as well as and assumed to be stepwise, with common discontinuity points. In order to rely on the exact pruned linear time method to infer the location and number of the discontinuity points, we considered a specific inference framework. Thus, we approximated a marginal log likelihood through a penalized BIC like term, for which we derived a variational lower bound. In this context, we used a variational expectation maximization algorithm for the inference. The methodology we proposed allows to do the segmentation of the intensity functions and the clustering of the nodes simultaneously. Experiments on toy data sets and real data were used to illustrate the relevance of the method.

A Proofs

A.1 Proof of Proposition 1

Proof Notice, first, that the central factor on the r.h.s. of the equality in equation ( 3) can be written as

M m=1 λ Z im Z jm (νm) = M m=1 N j>i λ Z i Z j (νm) 1 A (i,j) (νm)
, where 1 G (•) is the indicator function on a set G and A (ij) has been defined in Equation (1). By inverting the product on the right hand side, because of the indicator function we get

N j>i M m=1 λ Z i Z j (νm) 1 A (i,j) (νm) = N j>i M (i,j) m=1 λ Z i Z j (ν (i,j) m ) = N j>i M (i,j) m=1   K k,g
(λ kg (ν

(i,j) m )) Z ik Z jg   ,
where ν (i,j) m are the interaction times in the set A (i,j) , whose cardinality is M (i,j) . Thanks to the equation ( 5), the following holds Note that the last equality employs the definition of counting process M (i,j) (t) = M (i,j) m=1 1 ]0,t] (ν (i,j) m ).

By replacing (19) into the equation (3) and using that

Λ Z i Z j (T ) = K k,g Λ kg (T )Z ik Z jg = D d=1 K k,g λ kgd ∆ d Z ik Z jg ,
it suffices to take the logarithm of the likelihood and the proposition is proven.

A.2 Proof of Proposition 3

Proof The following objective function is taken into account L (q(Z); K, η, D, λ, π)

+ N i=1 l i K k=1 τ ik -1 .
This function has to be maximized with respect to both τ and the N Lagrange multipliers l 1 , . . . , l N , introduced to take into account the normality of the lines of τ . The most difficult step consists in taking the partial derivative of the objective function with respect to τ i 0 k 0 . We first focus on those terms of L(•) depending on d Q(τ, θ) := - where the last equality comes from the symmetry (i.e. interactions are undirected) of the frontal slices of tensors X and λ. Notice that the last term in the above equation is the function inside the exponential in Proposition 3. The remaining terms of L(•), not involving d, can be differentiated straightforward. Imposing the partial derivatives of L(•) equal to zero and using the above equation leads to the following system log(τ i 0 k 0 ) = log(π k 0 ) -∂Q(τ,θ) ∂τ i 0 k 0

+ l i 0 -1 K k 0 =1 τ i 0 k 0 = 1 ∀(i 0 , k 0 ).
The solution is obtained straightforward after some manipulations and this concludes the proof.

A.3 Proof of Proposition 4

Proof The following definitions are introduced to keep the notation uncluttered

S kg : = N j>i τ ik τ jg Y [s,t[ kg : = N j>i
τ ik τ jg (M (i,j) (t) -M (i,j) (s)) ∀s < t.

Moreover, for every tu e < tu f < tu g , the following short hand notation is used ∆ e,f := tu f -tu e and similarly for ∆ f,g and ∆ e,g . Hence, we get where the first and the last equalities come from the definition of G(•). This concludes the proof.

)

  Notice that the criterion to maximize (now with respect to η and D) is a sum of independent components: each gain function G([η d-1 , η d [) applies only to interactions that take place in the time interval [η d-1 , η d [. Notice in particular that for a given d, the λkgd are obtained from the quantities X (d) ij which correspond themselves to interaction counts during the time interval [η d-1 , η d [.

)

  Notice that in those equations, η corresponds to at most D -1 change points in [0, T [, ζ is a unique change point and η corresponds to at most D -2 change points in [0, ζ[.

  ): a maximal number of change points has to be specified and the optimization over ζ (i.e. over the position of a given change point) remains an open problem. While in theory this optimization is straightforward, the key point of the recurrence in equation (17) is the possibility of memorizing F (ζ, D -1) for all values of ζ. Indeed, the dynamic programming algorithm proceeds by computing and memorizing F (ζ, 1) for all ζ, and then computes and memorizes F (ζ, 2) using F (ζ, 1), etc.

Fig. 1 :

 1 Fig. 1: Image graph for each time segment

Fig. 2 :

 2 Fig. 2: The average number of clusters and time segments detected by MODL and PELT-Dynamic SBM versus number of sampled edges

  .

Fig. 3 :

 3 Fig. 3: ARIs for the change points (η)

  Fig. 4: ARIs for the cluster memberships (Z)

Fig. 5 :

 5 Fig. 5: Boxplots over 50 simulations of the ARIs for the clustering structures obtained by k-means, MODL and PELT-Dynamic SBM for scenario 2

Fig. 7 :

 7 Fig. 7: Adjusted Rand indexes for the change points η (left hand side) and the cluster memberships Z (right hand side)

Fig. 8 :

 8 Fig. 8: The 11 clusters found by PELT-Dynamic SBM represented here with 11 different symbols/colors on the left hand side

Fig. 9 :

 9 Fig. 9: An histogram shows how frequent interactions (cycle hires) are during the day. The vertical lines correspond to the estimated change points

Fig. 10 :

 10 Fig. 10: Estimated IFs for groups (3, •). Groups 1 and 7 are geographically adjacent to cluster 3 whereas 4 and 10 are not

Z

  ik Z jg 1 [η d-1 ,η d [ ik Z jg (N (i,j) (η d )-N (i,j) (η d-1 )) kgd . (19)

  kg ∆ e,f S kg + log(λ e,g kg )Y[tu e ,tu f [ kg ∆ f,g S kg + log(λ e,g kg )Y [tu f ,tu g [ kg ∆ e,f S kg + log(λ e,f kg )Y [tu e ,tu f [ kg ∆ f,g S kg + log(λ f,g k,g )Y [tu f ,tu g [ kg   = G([tu e , tu f [) + G([tu f , tu g [),

  Kernel density estimates over 50 simulations of the change points estimated by PELT-Dynamic SBM, for scenario 2. The true values of η 1 and η 2 are given by the red vertical lines following IFsλ Z i Z j (t) = 0.051 I 1 (t) + 0.11 I 2 (t) + 0.051 I 3 (t) + 0.11 I 4 if Z i = Z j 0.11 I 1 (t) + 0.051 I 2 (t) + 0.11 I 3 (t) + 0.051 I 4 if Z i = Z j ,for all j > i and I d denotes the d-th segment. By construction, integrating the IFs over[0, 12] leads to

	Density	5 10 15 20				Density	10 20 30 40 50 60
		0					0
		1.95	2.00	2.05	2.10	2.15	6.86	6.88	6.90	6.92	6.94
				η 1					η 2
	Fig. 6:

Table 1

 1

		: Aggregated interactions for clusters 1 ( ) and 8 ( * ) during the five seg-
	ments uncovered. On the main diagonal of each table, the numbers of interactions
	within clusters are reported: for cluster 1 on the left/top, for cluster 8 on the
	right/bottom. Interactions between clusters are outside the main diagonal. Com-
	munity structure for cluster 8 is indicated in blue. Disassortative structure for
	cluster 8 is indicated in red.
	104 47	742 441	1106 368	912 419	572 128
	47	16	441 502	368	338	419 984	128 108
	0.00-6.45	6.45-9.45	9.45-17	17-19.45	19.45-0.00
	the change point locations seems to be meaningful on a ground truth basis and the
	PELT-Dynamic SBM proved to be fit to uncover interaction patterns that could
	not easily be detected by other static or dynamic clustering algorithms.

In the framework of change point detection for time series, a cost function to be minimized is associated to each time segment. In contrast, we introduce a gain function to be consistent with the problem formulation used so far. However the two definitions are equivalent since the cost function can be thought as a gain function multiplied by -1.

Map data are available from http://www.openstreetmap.org and copyrighted Open-StreetMap contributors.


est destinée au dépôt

Context

-New direct and parallel reinitialization algorithm [5] was incorporated.

-Recoloring scheme [6] was used to reduce the number of LS functions needed to represent the polycrystal

Full

  field modeling in a level set framework of Zener pinning phenomenon -discussion of classical limit mean grain size equation. Benjamin Scholtes, Dmitrii Nikolaevich Ilin, Amico Settefrati, Nathalie Bozzolo, Andrea Agnoli, Marc Bernacki To cite this version: Benjamin Scholtes, Dmitrii Nikolaevich Ilin, Amico Settefrati, Nathalie Bozzolo, Andrea Agnoli, et al.. Full field modeling in a level set framework of Zener pinning phenomenon -discussion of classical limit mean grain size equation.. Superalloys 2016, Sep 2016, Seven Springs, Pennsylvania, United States. hal-01430741 Full-field Modeling of the Zener Pinning Phenomenon in a Level Set Framework -Discussion of Classical Limiting Mean Grain Size Equation B. Scholtes a , D. Ilin a , A. Sattefrati b , N. Bozzolo a , A. Agnoli c , M. Bernacki a a Mines ParisTech, CEMEF -Centre de mise en forme des matériaux, CNRS UMR 7635, CS 10207 rue Claude Daunesse, 06904 Sophia Antipolis cedex, France b Transvalor SA, 694 avenue Docteur Maurice Donat, 06250 Mougins, France c Snecma Gennevilliers, 171 boulevard Valmy, 92702 Colombes, France Numerical model Level Set Framework: LS function ψ is defined over a domain Ω as the signed distance to the interface Γ [3]: Boundary condition at precipitate/matrix interface:

t

  -time, M -mobility, γ -interface energy, ve -velocity due to stored energy gradient 2D simulation results for grain growth: New mean field model for the limiting mean grain size Current work References α = 0°(incoherent precipitates) Particle radii rp: 0.2, 0.4, 0.6, 0.8 µm Area fraction f: 1-8 % Microstructure at the early and final stages of the simulation (rp = 0,4 µm, f = 3% ) Dependence of < Rf > on fgb Parameter K as a function of normalized particle size Expression obtained for the limiting grain size: Parameter m as a function of normalized particle size

Simulation parameters:

Pinning of grain boundaries by second phase particles is widely used to control the grain size during forming process of superalloys.

-Adaptive metric based meshing remeshing tool [4] was used.