Thibaud Rahier

Sylvain Marié

A Pre-Screening Approach for Faster Bayesian Network Structure Learning

Keywords: Bayesian networks, Structure learning, Information theory, Conditional entropy, Determinism, Functional relations, Screening

Learning the structure of Bayesian networks from data is a NP-Hard problem that involves optimization over a super-exponential sized space. Still, in many real-life datasets a number of the arcs contained in the final structure correspond to strongly related pairs of variables and can be identified efficiently with information-theoretic metrics. In this work, we propose a meta-algorithm to accelerate any existing Bayesian network structure learning method. It contains an additional arc pre-screening step allowing to narrow the structure learning task down to a subset of the original variables, thus reducing the overall problem size. We conduct extensive experiments on both public benchmarks and private industrial datasets, showing that this approach enables a significant decrease in computational time and graph complexity for little to no decrease in performance score.

Introduction

Bayesian networks are probabilistic graphical models that present interest both in terms of knowledge discovery and density estimation. Learning Bayesian networks from data has been however proven to be NP-Hard by [START_REF] Chickering | Learning Bayesian networks is NP-complete[END_REF]. There has been extensive work on tackling the ambitious problem of Bayesian network structure learning (BNSL) from observational data. Algorithms fall under two main categories: constraint-based and score-based.

Constraint-based structure learning algorithms rely on testing for conditional independence relations that hold in the data in order to reconstruct a Bayesian network encoding these independence relations. The PC algorithm by [START_REF] Spirtes | Causation, prediction, and search[END_REF] was the first practical application of this idea, followed by increasingly optimized approaches such as the fast incremental association (Fast-IAMB) algorithm.

Score-based structure learning relies on the definition of a network score, then on the search for the best-scoring structure among all possible directed acyclic graphs (DAGs). The number of possible DAG structures with n nodes is of order 2 O(n 2 ) , which prevents exhaustive search when n is typically larger than 30. Most of the score-based algorithms used in practice therefore rely on heuristics, as the original approach from [START_REF] Cooper | A Bayesian method for the induction of probabilistic networks from data[END_REF] which assumes a prior ordering of the variables is known, or [START_REF] Bouckaert | Bayesian belief networks: from inference to construction[END_REF] who proposed to search through the structure space using greedy hill climbing with random restarts. Since these first algorithms, different approaches have been proposed, increasingly pushing the limits of state-of-the art in BNSL: some based on the search for an optimal ordering [START_REF] Teyssier | Ordering-based search: a simple and effective algorithm for learning Bayesian networks[END_REF][START_REF] Chen | Improving Bayesian network structure learning with mutual information-based node ordering in the K2 algorithm[END_REF], others on optimizing the search task in accordance to a given score [START_REF] Scanagatta | Learning Bayesian networks with thousands of variables[END_REF][START_REF] Nie | Learning Bayesian networks with bounded tree-width via guided search[END_REF], using integer programming [START_REF] Cussens | Bayesian network learning with cutting planes[END_REF] or bio-inspired optimization heuristics (Kareem andOkur, 2019, 2021). See [START_REF] Scutari | Who learns better bayesian network structures: Accuracy and speed of structure learning algorithms[END_REF] for a recent review.

Meanwhile, data itself may contain determinism, for example in the fields of cancer risk identification [START_REF] De Morais | Handling almostdeterministic relationships in constraint-based Bayesian network discovery: Application to cancer risk factor identification[END_REF] or nuclear safety [START_REF] Mabrouk | An efficient Bayesian network structure learning algorithm in the presence of deterministic relations[END_REF]. Data is also increasingly collected and generated by software systems whether in social networks, smart buildings, smart grids, smart cities or the internet of things (IoT) in general [START_REF] Koo | An internet-of-things (iot) system development and implementation for bathroom safety enhancement[END_REF]. These systems in their vast majority rely on relational data models or semantic data models [START_REF] El Kaed | Building management insights driven by a multi-system semantic representation approach[END_REF] where the same entity may be described with several attributes. It is therefore now common to find deterministically related variables in datasets. Determinism has been shown to interfere with Bayesian network structure learning, notably constraint-based methods as mentioned by [START_REF] Luo | Learning Bayesian networks in semi-deterministic systems[END_REF].

In this paper, we first remind the background of Bayesian network structure learning (Section 2) and bring forward the following contributions:

we state some theoretical results bridging the gap between the notion of determinism and Bayesian network scoring (Section 3), we propose and study the complexity of the quasi-determinism screening BNSL (qds-BNSL) meta-algorithm, whose aim is to accelerate any existing BNSL algorithm by reducing the learning problem to a subset of the original variables via the detection of strong arcs (Section 4), we conduct experiments on both public benchmarks and private industrial datasets, demonstrating empirically that our meta-algorithm indeed accelerates the overall BNSL procedure with very low performance loss and also leads to sparser and therefore more interpretable graphs (Section 5).

2 Bayesian network structure learning

Bayesian networks

Let X = (X 1 , . . . , X n ) be a n-tuple of categorical random variables with respective value sets V al(X 1 ), . . . , V al(X n ). The distribution of X is denoted by,

∀ x = (x 1 , . . . , x n ) ∈ V al(X), p(x) = P (X 1 = x 1 , . . . , X n = x n ).
For I ⊂ J1, nK, we define X I = {X i } i∈I , and the notation p(•) and p(•|•) is extended to the marginals and conditionals of any subset of variables:

∀(x I , x J ) ∈ V al(X I∪J ), p(x I |x J ) = P (X I = x I |X J = x J ).
Moreover, we suppose that D is a dataset containing M i.i.d. instances of (X 1 , . . . , X n ). All quantities empirically computed from D will be written with a . D exponent (e.g. p D refers to the empirical distribution with respect to D). Finally, D I refers to the restriction of D to the observations of X I .

A Bayesian network is an object B = (G, θ) where (1) G = (V, A) is a directed acyclic graph (DAG) structure with V the set of nodes and A ⊂ V × V the set of arcs. We suppose V = J1, nK where each node i ∈ V is associated with the random variable X i , and π

G (i) = {j ∈ V s.t. (j, i) ∈ A} is the set of i's parents in G 4 and (2) θ = {θ i } i∈V is a set of parameters. Each θ i defines the local conditional distribution of X i given its parents in the graph, P (X i |X π(i) ). More precisely, θ i = {θ xi|x π(i) } where for i ∈ V, x i ∈ V al(X i ) and x π(i) ∈ V al(X π(i) ), θ xi|x π(i) = p(x i |x π(i) ).
A Bayesian network B = (G, θ) encodes the following factorization of the distribution of X:

for x = (x 1 , . . . , x n ) ∈ V al(X), p(x) = n i=1 p(x i |x π G (i) ) = n i=1 θ xi|x π G (i) .
Such a factorization notably implies that each variable is independent of its non-descendants given its parents.

Score-based approach to Bayesian network structure learning

For a given scoring function s : DAG V → R, where DAG V is the set of all possible DAG structures with node set V , score-based BNSL aims at solving the following combinatorial optimization problem:

G * ∈ argmax G∈DAG V s(G).
(1)

It can be shown that 2 1) where |V | = n. There are therefore 2 O(n 2 ) possible DAG structures containing n nodes [START_REF] Koller | Probabilistic graphical models: principles and techniques[END_REF]: the size of DAG V is said to be super-exponential in |V |. Most scoring functions used in practice are based on the likelihood function. The most straightforward being the max log-likelihood sore, that we now present.

n(n-1) 2 ≤ |DAG V | ≤ 2 n(n-
The max log-likelihood (MLL) score Let l D (θ) = log(p θ (D)) be the log-likelihood of the set of parameters θ given the dataset D. For a given DAG structure G ∈ DAG V , we define the MLL score of G with respect to D as:

s M LL D (G) = max θ∈Θ G l D (θ).
where Θ G is the set of all θ's such that B = (G, θ) is a well defined Bayesian network. The MLL score is very straightforward and intuitive, but it favors denser structures: if G 1 = (V, A 1 ) and G 2 = (V, A 2 ) are two graph structures such that A 1 ⊂ A 2 , we can show that:

s M LL D (G 1 ) ≤ s M LL D (G 2
). This problem is generally solved by using a score that induces a goodness-of-fit versus complexity tradeoff, such as BIC [START_REF] Schwarz | Estimating the dimension of a model[END_REF], which is a penalized version of the MLL score, or BDe [START_REF] Heckerman | Learning Bayesian networks: The combination of knowledge and statistical data[END_REF], which is derived from the marginalization of the likelihood, which implicitly penalizes the model's complexity through a Dirichlet prior on the parameters. In this paper, we will use the BDe score to evaluate a BN structure's quality, as it is done in several papers as [START_REF] Teyssier | Ordering-based search: a simple and effective algorithm for learning Bayesian networks[END_REF], or [START_REF] Nie | Learning Bayesian networks with bounded tree-width via guided search[END_REF]. This score is known to be a good indicator of generalization performance, which is what we are aiming to optimize. In the following sections, we propose to look for such a solution by constructing sparse graphs with minimal MLL.

3 Determinism and Bayesian networks

Definitions

We propose the following definitions of determinism and deterministic DAGs using the notion of conditional entropy. In this paper, determinism will always be meant empirically, with respect to a dataset D.

Definition 1 Determinism wrt D Given a dataset D containing observations of variables X i and X j , the relationship

X i → X j is deterministic with respect to D iff H D (X i |X j ) = 0, where H D (X i |X j ) = - xi,xj p D (x i , x j ) log(p D (x i |x j ))
is the empirical conditional Shannon entropy.

It is straightforward to prove that Definition 1 relates to a common and intuitive perception of determinism, as the one presented by [START_REF] Luo | Learning Bayesian networks in semi-deterministic systems[END_REF]. Indeed,

H D (X i |X j ) = 0 ⇔ ∀x j ∈ V al(X j ), there exists a unique x i ∈ V al(X i ) s.t. p D (x i |x j ) = 1.
This definition is naturally extended to X I and

X J for I, J ⊂ V , i.e. X I → X J is deterministic with respect to D iff H D (X J |X I ) = 0. Definition 2 Deterministic DAG wrt D G ∈ DAG V is said to be deterministic with respect to D iff ∀i ∈ V s.t. π G (i) = ∅, X π G (i) → X i is deterministic wrt D.

Deterministic trees and MLL score

We first recall a lemma that relates the MLL score presented in Section 2 to the notion of empirical conditional entropy. This result is well known and notably stated by [START_REF] Koller | Probabilistic graphical models: principles and techniques[END_REF].

Lemma 1 For G ∈ DAG V associated with variables X 1 , . . . , X n observed in a dataset D, s M LL D (G) = -M n i=1 H D (X i |X π(i) )
where by convention

H D (X i |X ∅ ) = H D (X i ).
The next proposition follows then straightforwardly. We remind that a tree is a DAG in which each node has exactly one parent, except for the root node which has none.

Proposition 1 If T is a deterministic tree with respect to D, then T is a solution of (1):

s M LL D (T ) = max G∈DAG V s M LL D (G).
It is worth noticing that complete DAGs also maximize the MLL score. The main interest of Proposition 1 resides in the fact that, under the (strong) assumption that a deterministic tree exists, we are able to explicit a sparse solution of ( 1), with n -1 arcs instead of n(n-1) 2 for a complete DAG.

Deterministic forests and the MLL score

The deterministic tree assumption of Proposition 1 is very restrictive. In this section, it is extended to deterministic forests, defined as follows:

Definition 3 Deterministic forest wrt D F ∈ DAG V is said to be a deterministic forest with respect to D iff F = p k=1 T k ,
where T 1 , . . . , T p are p disjoint deterministic trees wrt D V T 1 , . . . , D V Tp respectively and s.t.

p k=1 V T k = V .
In the expression

p k=1 T k , ∪ is the canonical union for graphs: G ∪ G = (V G ∪ V G , A G ∪ A G ). For a given deterministic forest F with respect to D, we define R(F ) = {i ∈ V | π F (i) = ∅} the set of F 's

roots (the union of the roots of each of its trees).

Proposition 2 Suppose F is a deterministic forest wrt D. Let G * R(F ) be a solution of the BNSL optimization problem (1) for X R(F ) and the MLL score i.e.

s M LL D R(F ) (G * R(F ) ) = max G∈DAG R(F ) s M LL D R(F ) (G).
Then,

G * = F ∪ G * R(F ) is a solution of (1) for X, i.e. s M LL D (G * ) = max G∈DAG V s M LL D (G).
As opposed to Proposition 1, the assumptions of Proposition 2 are always formally verified: if there is no determinism in the dataset D, then R(F ) = V , and every tree T k is formed of a single root node. In that case, solving problem (1) for G * R(F ) is the same as solving it for G * . Of course, we are interested in the case where R(F ) < n, as this enables us to focus on a smaller structure learning problem while still having the guarantee to learn the optimal Bayesian network with regards to the MLL score.

As seen in Section 2, the main issue with the MLL score is that it favors complete graphs. However, a deterministic forest F containing p trees is very sparse (n -p arcs), and even if the graph

G * R(F ) is dense, the graph G * = F ∪ G * R(F )
may still satisfy sparsity conditions.

4 Structure learning with quasi-determinism screening

Quasi-determinism

When it comes to BNSL algorithms, even heuristics are computationally intensive. We would like to use the theoretical results presented in Section 3 to simplify the structure learning problem.

Our idea is to narrow the structure learning problem down to a subset of the original variables: the roots of a deterministic forest, in order to significantly decrease the overall computation time. This is what we call determinism screening.

However, one does not always observe real empirical determinism, although there are very strong relationships between some of the variables. We therefore propose to relax the aforementioned determinism screening to quasi-determinism screening, where quasi is meant with respect to a parameter : we talk about -quasi-determinism.

There are several ways to measure how close a relationship is from deterministic. [START_REF] Huhtala | Tane: An efficient algorithm for discovering functional and approximate dependencies[END_REF] consider the minimum number of observations that must be dropped from the data for the relationship to be deterministic. Since we are in a score-maximization context, we will rather use as a threshold on the empirical conditional entropy. The following definition is the natural generalization of Definition 1.

Definition 4 -quasi-determinism ( -qd) Given a dataset D containing observations of variables X i and X j , the relation- ship X i → X j is -qd wrt D iff H D (X j |X i ) ≤ .
It has been seen in Proposition 2 that a deterministic forest is the subgraph of an optimal DAG with respect to the MLL score, while still satisfying sparsity conditions. Such a forest is therefore very promising with regards to the fitcomplexity tradeoff (typically evaluated by scores such as BDe or BIC).

Combining this intuition with the -qd criteria presented in Definition 4, we propose the quasi-determinism screening approach to BNSL, defined in the next subsections. An alternate definition with a relative is proposed in Section 6.

Quasi-determinism screening algorithm

Algorithm 1 details how to find the simplest -qd forest F from a dataset D and a threshold . Here simplest refers to the complexity in terms of number of parameters, which in the case of categorical variables corresponds to the number of states: for each variable X i that has at least one -qd parent, we select the one that has the lowest number of states from the set of all potential -qd parents π (i) (line 11).

This algorithm takes for input D (a dataset containing M i.i.d realizations of X) and (a threshold for quasi-determinism). The routine on lines 4-9 makes sure no Algorithm 1 Quasi-determinism screening (qds)

Input: D , 1: Compute empirical conditional entropy matrix

H D = H D (X i |X j ) 1≤i,j≤n 2: for i = 1 to n do 3: compute π (i) = {j ∈ J1, nK \ {i} | H D ij ≤ } 4: for i = 1 to n do 5: if ∃j ∈ π (i) s.t. i ∈ π (j) then 6: if H D ij ≤ H D ji then 7: π (j) ← π (j) \ {i} 8: else 9: π (i) ← π (i) \ {j} 10: for i = 1 to n do 11: π * (i) ← argmin j∈π (i)
|V al(X j )| (select one index in case of tie)

12: Compute forest F = (V F , A F ) where V F = J1, nK and A F = {(π * (i), i) | i ∈ J1, nK s.t. π * (i) = ∅} Output: F
cycle is introduced by the screening phase, and guarantees the next proposition holds (ensuring that Algorithm 1 is indeed well defined):

Proposition 3 For any rightful input D and , the output of Algorithm 1 is a forest (i.e. a directed acyclic graph with at most one parent per node).

Learning Bayesian networks using quasi-determinism screening

We now present Algorithm 2, which uses quasi-determinism screening to accelerate Bayesian network structure learning. This algorithm takes the following input: D (a dataset containing M realizations of X), (a threshold for quasideterminism

s M LL D (G * ) ≥ max G∈DAG V s M LL D (G) -M n .
In practice, this bound is not very tight and this result therefore has few applicative potential. However, it shows that:

s M LL D (G * ) -→ →0 max G∈DAG V s M LL D (G).
In other words, → s M LL D (G * ) is continuous in 0, and Proposition 4 generalizes Proposition 2.

Algorithm 2 is promising, notably if for small we have |R(F )| significantly smaller than n. In that case, ref-BNSL, that only has to be run on D R(F ) , can be expected to be much faster and more accurate than if it is run on the entire dataset D.

Complexity analysis

Complexity of baseline BNSL learning algorithms

The number of possible DAG structures being super exponential in the number of nodes, BNSL algorithms do not entirely explore the structure space but use smart caching and pruning methods to have a good performance & computation time trade-off.

Let ref-BNSL be a reference Bayesian network structure learning algorithm and C ref (M, n) be its complexity. C ref (M, n) should typically be thought of as linear in M and exponential, or at least high degree polynomial, in n for the best algorithms.

Complexity of Algorithm 1

We have the following decomposition of the complexity of Algorithm 1:

1. Lines 1-3: O(M n 2 ). Computation of H D : we need counts for every couple (X i , X j ) for i < j (each time going through D), which implies M n(n-1) 2 operations. 2. lines 4-9: O(n 2 ). Going through H D once. 3. lines 10-12: O(n 2 ). Going through H D once.

Overall one has that

C Alg1 (M, n) = O(M n 2 ).
Complexity of Algorithm 2 For a given dataset D, we define:

∀ ≥ 0, n r ( ) = |R(F )|.
The function n r (•), associates to ≥ 0 the number of roots of the forest F returned by Algorithm 1. The complexity of Algorithm 2 then decomposes as:

1. Line 1: O(M n 2 ). Run of Algorithm 1. 2. Lines 2-4: C ref (M, n r ( )). Run of ref-BNSL on reduced dataset D R(F ) with n r ( ) columns. This yields C Alg2 (M, n) = O(M n 2 ) + C ref (M, n r ( )).
We are interested in how much it differs from C ref (M, n), which depends mainly on:

how n r ( ) compares to n, how C ref (M, n) varies with respect to n.

C ref (M, n) is known to be typically exponential in n for the best exact structure learning algorithms, as those presented by [START_REF] Silander | A simple approach for finding the globally optimal Bayesian network structure[END_REF] or [START_REF] Cussens | Bayesian network learning with cutting planes[END_REF], and it is expected to be significantly larger than O(M n 2 ) for high-performing heuristics. We therefore expect an important decrease in computational time when running Algorithm 2 compared to its baseline version, as long as n r ( ) is sufficiently smaller than n. In the next section, we run a reference structure learning algorithm and Algorithm 2 on benchmark datasets in order to confirm this intuition.

Experiments

Experimental setup

Data Table 1 summarizes the data used in our experiments. We considered the largest open-source categorical datasets among those presented5 by [START_REF] Davis | Bottom-up learning of Markov network structure[END_REF] and available on the UCI repository [START_REF] Dheeru | UCI machine learning repository[END_REF]: 20 newsgroup, adult, book, covertype, kddcup 2000, msnbc, msweb, plants, reuters-52 and uscensus. Moreover, as it was done by [START_REF] Scanagatta | Learning treewidth-bounded Bayesian networks with thousands of variables[END_REF], we chose the largest Bayesian networks available in the literature6 , for each of which we simulated 10000 observations: andes, hailfinder, hepar 2, link, munin 1-4, pathfinder and win95pts.

We also include two industrial datasets containing descriptive metadata on which we have privileged access, priv-metadata 1 and priv-metadata 2. Programming details and choice of ref-BNSL Most of the code associated with this project was done in R, enabling an optimal exploitation of the bnlearn package from [START_REF] Scutari | Learning Bayesian networks with the bnlearn R package[END_REF], which is a very good reference among open-source packages dealing with Bayesian networks structure learning.

We need a BSNL algorithm to obtain a baseline performance. After carefully evaluating several algorithms implemented in the bnlearn package, we chose to use Greedy Hill Climbing with random restarts and a tabu list, as it consistently outperformed other built-in algorithms both in time and score, in addition to being also used as a benchmark algorithm in the literature, notably by [START_REF] Teyssier | Ordering-based search: a simple and effective algorithm for learning Bayesian networks[END_REF]. In this section, we refer to this algorithm as ref-BNSL.

Choice of for qds-BNSL An approach to choosing in the case of the qds-BNSL algorithm is to pick values for n r ( ), and manually find the corresponding values for . For a given dataset and ρ ∈ [0, 1], we define ρ = n -1 r ( ρn ). In other words, ρ is the value of for which the number of roots of the qd forest F represents a proportion ρ of the total number of variables (more details in Rahier ( 2018)). The computation of ρ is not problematic: once H D is computed and stored, evaluating n r ( ) is done in constant time, and finding one of n r (•)'s quantiles is doable in O(log(n)) operations (dichotomy), which is negligible compared to the overall complexity of the screening. In the case of the priv-metadata datasets, choosing = 0 leads to a dramatic decrease of the number of variables that are considered by the baseline algorithm, since these datasets contain several truly deterministic relationships by design.

Algorithm evaluation

The algorithms are evaluated using 3 axes of performance:

-BDe score of Section 2 with a uniform prior and equivalent sample size (ESS) equal to 5, inspired from [START_REF] Teyssier | Ordering-based search: a simple and effective algorithm for learning Bayesian networks[END_REF] and referred to as BDeu. -Number of arcs of the learned Bayesian network.

The BDeu score naturally penalizes overly complex models (in terms of number of parameters), it is however interesting to look at the number of arcs, as it is a straightforward way to evaluate how complex a Bayesian network appears to a human expert (and thus how interpretable this structure is). -Computing time t run (all algorithms were run on the same machine).

It is essential to remark that ref-BNSL is used both to obtain a baseline performance and inside qds-BNSL. In both cases, it is run with the same settings until convergence. The comparison of computing times is thus fair.

We present the obtained results for our selected baseline algorithm ref-BNSL, and 3 versions of qds-BNSL. For each dataset, we selected ∈ { 0.9 , 0.75 , 0.5 }), corresponding to a restriction of ref-BNSL to 90%, 75% and 50% of the original variables respectively (for the priv-metadata datasets, these three choices of are merged into the single choice = 0, which results in a decrease of more than 50% of the original variables).

The results are shown in Table 2, one group of columns per evaluation criterion, and each value is the median of 10 runs with different seeds. In each table, the median value of the criterion is displayed for ref-BNSL (ref ), and the relative difference is displayed for the three versions of qds-BNSL we consider (qds 0.9 , qds 0.75 and qds 0.5 ).

Results

Score It appears in Table 2 that the decrease in BDeu score is smaller than 5% for all the considered datasets when 90% of the variables remain after the -186 -0.0 -0.0 -3.9 8550 -7.9 -13 -39 903 +0.0 +0.0 -8.5 pathfinder -27 -0.7 -0.7 -4.9 231 -14 -35 -69 161 -4.3 -8.7 -24 win95pts -9.2 +0.1 -1.1 -9.2 132 -6.0 -31 -69 115 +0.0 -0.9 -12 priv-meta1 -8.72 +1.1 +1.1 +1.1 13794 -99 -99 -99 70 -41 -41 -41 priv-meta2 -8.72 +12.5 +12.5 +12.5 4346 -99 -99 -99 102 -59 -59 -59 pre-screening (qds 0.9 ), and for most of them when 75% of the variables remain (qds 0.75 ). This is also observed with 0.5 for datasets that contain a lot of very strong pairwise relationships as kddcup, msweb, or munin 2-4. For priv-metadata datasets, our approach increases the score (slightly for priv-metadata1 and of more than 12% for priv-metadata2).

Computing time Table 2 shows a significant decrease in computational time for qds-BNSL, which is all the more important as is large. In the best cases, we have both a very small decrease in BDeu score, and an important decrease in computational time. We suspect that this is also due to the presence of many strong pairwise relationships. For example, the algorithm qds-BNSL with = 0.5 is 55% faster for msweb, and 54% for munin 3, while implying only around 1% decrease in score compared to ref-BNSL. If we allow a 5% score decrease, qds-BNSL can be up to 70% faster (20 newgroups, book, msnbc, kddcup, hepar2, pathfinder). On the industrial datasets this computational time decrease is astonishing: it is of more than 2 orders of magnitude in the case of priv-metadata 1 and more than 3 orders of magnitude for priv-metadata 27 . These results confirm the complexity analysis of the previous section, in which we supposed that the screening phase had a very small computational cost compared to the standard structure learning phase.

Complexity As showed by Table 2, Bayesian networks learned with qds-BNSL are consistently less complex than those learned with ref-BNSL. Several graphs learned with qds 0.5 are more than 30% sparser while still scoring less than 5% below the baseline algorithm: 20 newsgroups, book, kddcup 2000, msnbc, msweb and hepar 2. They provide an interesting example of the sparsity induced by qds-BNSL. After the qd 0.5 -screening phase, half of the variables (corresponding to the nodes in white) are considered to be sufficiently explained by V 1. They are therefore not taken into account by ref-BNSL, which is run only on the variables corresponding to the nodes in gray (more details in [START_REF] Rahier | Screening strong pairwise relationships for fast Bayesian network structure learning 2nd Italian-French Statistics Seminar-IFSS[END_REF]). In the msnbc case, this learning problem restriction implies only a small decrease in the final graph's generalization performance (as seen in BDeu scores), while being 7 times faster to compute and enabling a significantly better readability.

In this processed version of the msnbc dataset [START_REF] Davis | Bottom-up learning of Markov network structure[END_REF], each variable contains a binary information regarding the visit of a given page from the msnbc.com website8 . The Bayesian network displayed in Figure 2 shows in a compact way the influence between the different variables. For instance, we see that visits of the website's pages corresponding to nodes in white (e.g. 'weather' (V 8), 'health' (V 9) or 'business' (V 11)) are importantly influenced by whether the user has also visited the frontpage (V 1). For example, learned parameters show that a user who did not visit the website's frontpage (V 1) is about 10 times more likely to have visited the website's 'summary' page (V 13) than a user who did visit the frontpage. Such information is much harder to read from the graph learned with ref-BNSL displayed in Figure 1 right.

Concluding remarks

We have seen that, both in theory and in practice, the quasi-determinism screening approach enables a decrease in computational time and complexity for a small decrease in graph score. This tradeoff is all the more advantageous as there actually are strong pairwise relationships in the data, that can be detected during the screening phase, thus enabling a decrease in the number of variables to be considered by the baseline structure learning algorithm during the second phase of Algorithm 2. Optimal cases for this meta-algorithm take place when n r ( ) is significantly smaller than n for reasonably small compared to the variable's entropies. Among benchmark datasets this is reasonably frequent (e.g 20 newsgroup, msnbc, munin2-4, webkb), and we argue it is extremely frequent among industrial datasets, as we have shown in our priv-metadata datasets which are only a small sample of the kind of datasets in which we can find very strong (even completely deterministic) relations.

Besides, we still have potential to improve the qds-BNSL meta-algorithm, by paralellizing the computation of H D , and implementing it in C instead of R.

Our main research perspectives are (1) to understand how one can anticipate how good the score/computation time/complexity trade-off can be before running any algorithm all the way through, saving us from running qds-BNSL on datasets in which there are no strong pairwise relationships to be detected, (2) find a principled way to choose and (3) tighten the bound of Proposition 4 and generalize it to the BDeu score.

In another directions, we have some insights on ways to generalize our quasideterminism screening idea. The proof of Proposition 2 suggests that the result still holds when F is any kind of deterministic DAG (and not only a forest). We could therefore use techniques that detect determinism in a broader sense than only pairwise, to make the screening more efficient. For this purpose we could take inspiration from papers of the knowledge discovery in databases (KDD) community, as [START_REF] Huhtala | Tane: An efficient algorithm for discovering functional and approximate dependencies[END_REF], or more recently [START_REF] Papenbrock | Functional dependency discovery: An experimental evaluation of seven algorithms[END_REF] who evaluate functional dependencies discovery methods. We also could broaden our definition of quasi-determinism: instead of considering the information-theoretic quantity H D (X|Y ) to describe the strength of the relationship Y → X, one could choose H D (X|Y )

H D (X) , which represents the proportion of X's entropy that is explained by Y . Moreover, H D (X|Y ) H D (X) ≤ can be rewritten as M I D (X,Y )

H(X)

≥ 1 -, which gives another insight to quasi-determinism screening: for a given variable X, this comes down to finding a variable Y such that M I D (X, Y ) is high. This is connected to the original idea of [START_REF] Chow | Approximating discrete probability distributions with dependence trees[END_REF], and later [START_REF] Cheng | Learning belief networks from data: An information theory based approach[END_REF], for whom pairwise empirical mutual information is central. This alternate definition of -quasi-determinism does not change the algorithms and complexity considerations described in Section 4.

  ) and ref-BNSL (a BNSL baseline algorithm, taking for input a dataset, and returning a Bayesian network structure). The extension of the def-Algorithm 2 Bayesian network structure learning with quasi deterministic screening (qds-BNSL) Input: D, , ref-BNSL 1: Compute F by running Algorithm 1 with input Dand 2: Identify R(F ) = {i ∈ J1, nK | π F (i) = ∅}, the set of F 's roots. 3: Compute G * R(F ) by running ref-BNSL on D R(F ) 4: G * ← F ∪ G * R(F )Output: G * inition of determinism to quasi-determinism (Definition 3) prevents us to have 'hard' guarantees as those presented in Proposition 2. However, we are able to explicit bounds for the MLL score of a graph G * returned by Algorithm 2, as stated in the following Proposition. Proposition 4 Let , D and ref-BNSL be rightful input to Algorithm 2, and G * the associated output. Then, if ref-BNSL is exact (i.e. always returns an optimal solution) with respect to the M LL score, we have the following lower bound for s M LL D (G * ):

  8 7093 -20 -22 -44 879 +0.0 +0.0 -13 munin3 -165 +0.0 +0.0 -1.1 11558 -37 -29 -54 898 +0.0 +0.0 -7.8 munin4

Figure 1 Fig. 1 .

 11 Figure 1 displays two Bayesian networks learned on the msnbc dataset.

Table 1 .

 1 Datasets presentation

	name	short name	n	M
	20 newsgroups 20ng	930 11293
	adult	adult	125 36631
	book	book	500	8700
	covertype	covertype	84 30000
	kddcup 2000	kddcup	64 180092
	msnbc	msnbc	17 291326
	msweb	msweb	294 29441
	plants	plants	69 17412
	reuters 52	r52	941	6532
	uscensus	uscensus	68 2458285
	andes	andes	223 10000
	hailfinder	hailfinder	56 10000
	hepar 2	hepar2	70 10000
	link	link	724 10000
	munin 1	munin1	186 10000
	munin 2	munin2	1003 10000
	munin 3	munin3	1041 10000
	munin 4	munin4	1038 10000
	pathfinder	pathfinder	109 10000
	windows 95 pts win95pts	76 10000
	priv-metadata 1 priv-meta1	43	1000
	priv-metadata 2 priv-meta2	41	1000

Table 2 .

 2 For algorithms ref , q0.5, q0.75 and q0.5, and benchmark datasets, we display the Bayesian Network's (1) BDeu score averaged by observation, (2) learning time (including prescreening) and (3) number of arcs. Every result that corresponds to a BDeu score less than 5% smaller than ref-BNSL's score is boldfaced.

		BDeu score		Computation time	Number of arcs
	dataset	ref q0.9 q0.75 q0.5 ref q0.9 q0.75 q0.5 ref q0.9 q0.75 q0.5
		(%)	(%)	(%)	(s) (%) (%) (%) (nb) (%) (%) (%)
	20ng	-143 -0.7 -2.1 -4.8 21495 -1.6 -43 -73 3136 -4.5 -15 -32
	adult	-13 -0.2 -0.1 -4.0 102 -6.6 -22 -61 371 +3.2 +7.0 -14
	book	-35 -0.8 -1.7 -4.6 7600 -24 -40 -71 2196 -11 -19 -40
	covertype	-14 -0.2 -1.2 -12 565 -6.8 -33 -71 337 -0.9 -11 -38
	kddcup	-2.4 -0.3 -1.0 -3.8 2167 -11 -33 -74 285 -5.3 -19 -39
	msnbc	-6.2 -0.1 -2.6 -4.6 252 -21 -61 -86 102 -7.8 -33 -64
	msweb	-9.8 +0.0 -0.1 -1.0 4701 -6.3 -9.9 -55 1264 -2.5 -3.6 -35
	plants	-13 -2.6 -7.6 -21 455 -47 -62 -84 320 -6.2 -18 -42
	r52	-95 -0.8 -2.0 -6.1 18630 -14 -38 -77 2713 -3.6 -9.1 -25
	uscensus	-23 -0.3 -1.8 -10 21782 -0.4 -32 -78 220 -10 -20 -38
	andes	-93 -0.5 -6.2 -17 898 -2.2 -27 -70 336 -0.9 -7.1 -23
	hailfinder	-50 -0.1 -2.7 -10	46 -5.3 -17 -55 64 -1.6 +6.2 -16
	hepar2

The exponent G may be dropped for clarity when the referred graph is obvious from context

http://alchemy.cs.washington.edu/papers/davis10a/

http://www.bnlearn.com/bnrepository/

These results were so extreme that they could not be fully captured in the 'percentage' format choice of Table2.

8 more details: http://archive.ics.uci.edu/ml/machine-learning-databases/ msnbc-mld/msnbc.data.html


Francisco Durán

Gwen Salaün

Optimization of BPMN Processes via Automated Refactoring

Business process optimization has become a strategic aspect of companies' management due to the potential of cost reduction and throughput improvement. There are several ways to achieve process optimization, depending on the level of expressiveness of the processes at hand. In this paper, we focus on processes described using BPMN, but also including an explicit description of execution time and resources associated with tasks. We propose a refactoring procedure whose final goal is to reduce the total execution time of the process given as input. Such a procedure relies on refactoring operations that reorganize the tasks in the process by taking into account the resources used by those tasks. This process refactoring technique is fully automated by a tool that we implemented and applied on several examples for validation purposes.

Introduction

Context. Business process optimization is a strategic activity in organizations because of its potential to increase profit margins and reduce operational costs. Optimization is however a difficult task to be achieved manually since several parameters should be taken into account (execution times, resources, costs, etc.). These parameters are not systematically included in existing languages used for modelling and managing business processes. Moreover, optimization requires a high level of expertise that not all users have. Automated techniques are thus required to optimize a given process for certain criteria of interest.

In this work, we assume that a description of a business process is given using the standardized workflow modelling language BPMN. This language allows us to define the set of tasks involved in a process and the order in which they should be executed. This behavioural description of the model can be extended with information on the time each task takes to execute and an explicit description of the resources required for executing each task. As a consequence, the resulting model of the process does not only take into account behavioural aspects but also quantitative aspects.

Motivations. Processes are not built once and for all in a monolithic way. During their life time, processes have to be changed or updated for several possible reasons: addition or suppression of some specific task, improvement of the process with respect to a given criterion (e.g., overall execution time), adjustment of the process to consider a new regulation or internal directive, etc. When writing a process or when updating it as suggested before, the quality or correctness of the process has to be preserved. This might not be the case if process writing or reengineering/refactoring is achieved manually by human beings. Moreover, refactoring steps may be very difficult to apply when BPMN models also take into account quantitative aspects and when some criteria, such as resource usage or execution time, are used to guide the refactoring steps so as to generate an optimal process (e.g., the process must execute as quickly as possible). Therefore, there is a need for automated refactoring techniques in order to generate an optimal version of a process during its writing or update.

Proposal. Given a process model, we propose some optimization techniques that rely on the refactoring of the given process. By changing the structure, we aim at generating a different process whose overall execution time is reduced compared to the original process. The main idea is to increase the level of parallelism of the tasks involved in the process. The refactoring steps should however be applied with care. For instance, it is not helpful to put in parallel two tasks using a same resource since they will compete for such a resource and a resource cannot be involved in two different tasks at the same time. Moreover, it does not always make sense to put two tasks in parallel, and in some cases such causal dependencies must be preserved (e.g., some product must be packaged before its delivery).

More precisely, we propose in this paper optimization techniques based on process refactoring. Our approach takes as input a BPMN process extended with time and resources associated with tasks, and generates as output a new version of this process.

To do so, we first analyze the process to identify tasks that could be executed earlier because the resources that one of these tasks requires for its execution are available before its execution begins. Alternative processes are generated by moving these tasks backwards (closer to the initial node). Our approach works using an iterative approach. Each newly generated process is similarly analyzed, and new alternative processes are generated in the same way. Since many processes can be synthesized by our approach, there are different ways to generate and handle these new processes. We have implemented and carried out experiments with two strategies: (i) an exhaustive exploration of all possible processes generated by refactoring, and (ii) a guided exploration of the new processes by using some heuristic. All these techniques are fully automated in a tool that we implemented and which has been validated on many examples.

Organization. The rest of this paper is organized as follows. Section 2 introduces the considered subset of the BPMN notation and its extension with time and resources. Section 3 provides an overview of the different steps of our approach. Section 4 focuses on the refactoring process to change the structure of the process. Section 5 presents the tool support and some experimental results to assess the accuracy and performance of our approach. Section 6 compares our solution to related work and Section 7 concludes.

BPMN with Time and Resources

BPMN 2.0 (BPMN, as a shorthand, in the rest of this paper) was published as an ISO/IEC standard [START_REF]International Standard 19510, Information technology -Business Process Model and Notation[END_REF] in 2013 and is nowadays extensively used for modeling and developing business processes. In this paper, we focus on activity diagrams including the BPMN constructs related to control-flow modeling and behavioural aspects. Beyond those constructs, execution time and resources are also associated with tasks. Figure 1 summarizes the BPMN constructs supported in this work.

Specifically, the node types event, task, and gateway, and the edge type sequence flow are considered. Start and end events are used, respectively, to initialize and terminate processes. A task represents an atomic activity that has exactly one incoming and one outgoing flow. A sequence flow describes two nodes executed one after the other in a specific execution order. A task and a flow may have a duration or delay. The timing information associated to tasks and flows is described as a literal value (a non-negative real number, possibly 0). Resources are explicitly defined at the task level. A task that requires resources can include, as part of its specification, the required resources. Information about time and resources can be used jointly for a given task. In such a case, it means that the task needs those resources to be able to execute, and once the resources are acquired, the task is going to execute for the specified duration.

Sequence flows can be of two possible types, to explicitly specify flows that must be preserved during the refactoring process. A strong flow corresponds to a causal dependency between two nodes that cannot be changed (e.g., some product must be packaged before its delivery). A weak flow corresponds to a loose connection between two nodes that may be preserved or not (e.g., the product could be delivered before the client pays for it).

Gateways are used to control the divergence and convergence of the execution flow. We consider in this work the two main kinds of gateways used in activity diagrams, namely, exclusive and parallel gateways. Gateways with one incoming branch and multiple outgoing branches are called splits, e.g., split parallel gateway. Gateways with one outgoing branch and multiple incoming branches are called merges, e.g., merge parallel gateway. An exclusive gateway chooses one out of a set of mutually exclusive alternative incoming or outgoing branches. A parallel gateway creates concurrent flows for all its outgoing branches or synchronizes concurrent flows for all its incoming branches.

Overview of the Approach

In this section, we give an overview of the different steps of our approach. We start by introducing the simulation-based analysis with which possible improvements in the structure of the process in terms of resource usage are identified. Then, we present the refactoring procedure and how different strategies are used to explore the possible solutions.

Simulation-Based Analysis

In this work, the main idea of the refactoring approach is to change the structure of the process in order to reduce its total execution time. Since optimization mostly targets process execution time, we need to compute this time for a given process. To do so, we rely on simulation-based techniques which turn out to simplify the computation of execution times in the presence of resources. Indeed, a task needs to acquire the required resources to be able to execute, and if the resources are not available the task cannot execute. The competition for resources may thus induce delays, and these delays are not easy to identify. To analize the process, we simulate it a certain number of times (this is a parameter of the approach). During these executions, some information about the execution of tasks (pending, executing, completed) and the usage of resources is stored. At the end of each execution, we store the time taken for completing the process. The average of those times allows us to compute the average execution time of the process.

After completion of the simulation, there is an analysis step which explores the simulation log for retrieving specific information. In particular, we look for specific timestamps during the simulation at which a task is in a pending state (meaning that this task is still waiting to be able to execute), and all resources required for executing this task are available. This means that this specific task could execute earlier in a process, and we will use this information to change the structure of the process by trying to move this task backwards in the process. This analysis step returns as output a set of tasks that could be executed earlier. Each task in this set verifies the aforementioned constraints (pending task but with required resources available) during a period of time of the simulation.

Refactoring Procedure

A refactoring step takes as input a BPMN process and a task that can be executed earlier in that process. Then, depending on the type of the preceding node, a specific refactoring operation is applied. Such refactorings represent changes in fragments of the process, and multiple refactoring operations may be applied consecutively. To do so, the refactoring procedure relies on the generation and exploration of possible solutions and keeps track of the visited processes to avoid recomputations.

The refactoring operations are presented in Section 4. A refactoring operation can be seen as a model transformation L ⇒ R if C, with L and R subprocess patterns and C a condition on them. Each refactoring operation results in the transformation of a process (or model) into another. Thus, a refactoring operation can be applied if its left-hand side L matches a subprocess (or submodel) of the current process P , that is, P | p = Lσ, for some match σ and position p, and the condition C is Fig. 2. Overview of the approach satisfied for that match, that is, if Cσ is evaluated to true. The application of such a refactoring rule results in the replacement of the matched subprocess P | p by the right-hand side Rσ. By repeatedly applying the available refactoring transformations on different parts of a process, new processes are generated. Strategies are thus necessary to explore all these new processes, and eventually return the optimal one.

Figure 2 depicts the input and outputs as well as the main steps of the refactoring process from a global perspective. The approach takes as input a BPMN process that includes a description of time and resources as shown in Section 2. This process is moved to the queue of processes to be explored. One process is then extracted from this queue (the original one in the first iteration) and the simulation-based analysis is carried out to identify the tasks that could be executed earlier. For each of these tasks, the corresponding refactoring operation is applied, generating a new process. Since the refactoring operation just moves the given task one step backwards, the resulting process may require additional changes. Therefore, if the resulting process has not been explored yet, and if the resulting process respects the strong dependencies defined in the original process (causal dependencies corresponding to strong flows must be maintained by the refactoring process), this new process is moved to the queue of processes to be explored. Otherwise, it is moved to the queue of processes already explored.

The refactoring procedure consists in the iterative application of refactoring operations. Since there are different ways to refactor a process, it is not possible to know at the beginning of the approach which refactoring step will lead to the best result (process with the minimal execution time). The exhaustive exploration of the search space, achieved by repeatedly attempting each possible refactoring operation on the tasks that may be moved backwards, may be time-consuming. To reduce the exploration space, we propose two different strategies:

-The first strategy consists in a bounded breadth-first search. It carries on an exhaustive exploration of all processes to be explored up to a certain bound given as parameter. For each possible task to be moved earlier in the process, the corresponding refactoring operation is applied. If the resulting process has not been analyzed before, it is placed in the to-be-explored queue. The procedure continues while there are processes in such a queue and the bound has not been reached. -The second strategy, instead of exploring the search graph by applying refactoring operations on all tasks that can be moved earlier, only expands by applying refactoring on the task closer to the start event in the BPMN process. The intuition is that by moving this task backwards, it will be placed closer to the initial event, thus reducing the number of times it may be moved in the future.

In other words, we try to move first the tasks closer to their final positions in the process. This strategy is referred as heuristic-based in the rest of this paper.

Several experiments showing the behaviour of the following exploration algorithms, including the use of these two strategies are presented in Section 5.

Refactoring Operations

This section presents a set of refactoring operations. In each of these refactoring operations, given as input a process and a task that has to be moved earlier in the process, a new process is returned as output. The refactoring operation to be applied depends on what type of node precedes the task to be moved backwards. There are actually three main cases: this node can be another task, a merge node, or a split node. Therefore, we will organize the rest of this section tackling successively these three cases. For each case, the proposed refactoring operations are presented. Before starting, we also recall that this refactoring step focuses on the process structure and on the usage of resources by tasks, but does not take into account strong flows, which is handled at another level (as explained in Section 3.2).

Task

This first case is rather straightforward. If a task can be executed earlier, and that task is preceded by another task, we check whether these two tasks are sharing some resources. If they do not share any resources, the process is transformed to execute these two tasks within a common parallel gateway, as illustrated in Figure 3.

Merge Node

If the task T to be moved earlier is preceded by a merge node, there are several possible (sub)cases. First of all, this node can be an exclusive or a parallel gateway. Second, before that merge, there may be only tasks, only other merge nodes, or a combination of both. We have specific operations for each of these cases. We will introduce these operations in the rest of this subsection.

Let us first consider the case in which the merge parallel gateway is preceded by a set of tasks T1...Tn, and none of these tasks share resources with the task T (left-hand side of Figure 4). In that case, all tasks are gathered in parallel just before the merge parallel gateway (right-hand side of Figure 4). Note that a split parallel gateway is added before these tasks in order to avoid that T executes before the tasks preceding the tasks before the merge. The preceding tasks might be using the same resources.

Assume now that the tasks preceding the merge parallel gateway share some resources with task T. If they all share resources, no refactoring is possible (since adding an additional task in parallel, competing for the same resources, would not improve the process execution time). If only one task shares resources with T, then T is moved after that task but before the merge as shown in Figure 5. If there are several tasks sharing resources with T (but not all), then T is moved before the merge and right after an additional merge parallel gateway for this subset of resources. Figure 6 shows such a case where T1...Ti share some resources with T. Therefore, after the refactoring operation, T appears after T1...Ti whereas Ti+1...Tn keep executing in parallel. If what precedes the task to be moved is an exclusive gateway, and if there are only tasks before the merge, there are two cases. For each task before the merge, if task T shares some resource with that task, T is moved before the merge but after that task. If task T does not share any resource with that task, T is put in parallel with that task. Figure 7 illustrates this operation, showing that in the case of T1...Ti (shared resources), T is moved after each of them, whereas in the case of Ti+1...Tn (no shared resources), they all appear in parallel in the resulting process. Note that T appears multiple times in the resulting process, because by including T in an exclusive pattern, it has to be executed once by each existing branch to maintain the intended behaviour.

As far as cascading merges are concerned, we support cascading merge exclusive gateways (possibly finishing with a merge parallel gateway), by applying several times the patterns introduced above. However, if there are cascading merge parallel gateways or a merge parallel gateway followed by a merge exclusive gateway, refactoring is too complicated and is not applied. Let us take the example of two merge parallel gateways. If we move a task within the first one, this is ok as presented earlier. However, if there are other merge parallel gateways inside the first one, we do not know where to move that task since we do not want to execute it multiple times.

Split Node

If the task T to be moved earlier is preceded by a split parallel gateway, then T is moved before the split, whatever precedes the split (task, split or merge). Figure 8 illustrates this pattern by moving the task before the split parallel gateway.

If the task T is preceded by a split exclusive gateway, we apply refactoring only if the split is preceded with a task (T1 for instance). In that case, we need to analyze the process to look for a merge gateway corresponding to the aforementioned split gateway. If the subprocess is balanced and the corresponding merge gateway is found, we still need to look at the resources used in that part of the process. If all tasks between the split and the merge do not share any resources with T1, then the whole block is moved in parallel with T1, as illustrated in Figure 9, because we cannot dissociate the contents (tasks for instance) appearing in the same branch of an exclusive structure. If the task T is preceded by a split exclusive gateway, and the preceding node is not a task (it is another split for instance), there is no simple refactoring and we keep it as is. Note that if the exclusive split preceding T is not balanced or some task in that block uses any of the resources of T, then the refactoring operation is not applied because optimization is not possible.

Implementation and Experiments

This section presents the tool support, a case study, and some experiments. Additional details about the tool and dataset used for the experiments are available online [START_REF]Workflow Refactoring Tool and Examples -(Blinded) GitHub Repository[END_REF].

Tool

A tool implementing the above refactoring operations and their systematic application following the described procedures has been developed in Python. It works as presented in Figure 2: the user provides as input a BPMN process and obtains as output an optimal version of that process with the gain in terms of process execution time. The implementation consists of several classes as illustrated in Figure 10. The core of the approach is implemented in several methods dedicated to process simulation and analysis (of the simulation), refactoring operations, and exploration of the refactored processes by using different strategies. Processes to be explored and already explored are stored into dictionaries. Hash values for processes are used as keys and are computed using the number of nodes at distance 1, 2, 3, etc. from the start event. These numbers are then concatenated to form a key. Since we may have collisions, for each hash value, we store a list of processes in the dictionary. As for the transformation from BPMN XML to our encoding of BPMN into Python (in both directions), we take advantage of the transformation capabilities available in VBPMN [START_REF] Krishna | VBPMN: Automated Verification of BPMN Processes[END_REF][START_REF] Krishna | Checking Business Process Evolution[END_REF].

Case Study

We illustrate our approach with a process describing the business trip organization given in Figure 11. Each task has an annotation with a pair of values where the first value is the duration of the task (in days for instance) and the second value is the set of resources required to execute that task. Let us introduce the process with more details. First, the assistant fills in the required documents. Then, the travel agency is in charge of booking flight tickets, followed by the reservation of accommodation by the user in that case. Visa is then prepared and, in parallel, the user has to be vaccinated. The final part of the process is executed when the user is back from the trip. All necessary documents for reimbursement are returned by the user. Reimbursement is then completed by the financial staff. Finally, all documents are archived by the assistant. It is worth noting that there are several strong flows in the original process, before and after g2, and between task ReturnDocuments and task Reimbursement. This means that the causal dependency between these tasks is important and must be preserved by the refactoring process. To compute the refactored version of this process, we use a bounded exploration (with 300 as bound). The resulting process is given in Figure 12, and was obtained after about 130 iterations. It takes about 20 seconds to compute the resulting process. The execution time of all the tasks of the original process is 43 days whereas the new version executes in 28 days.

Let us now comment on this new version of the process. We can see that causal dependencies defined by strong flows are preserved in this process: tasks Visa and Vaccination are executed before task ReturnDocuments, and task ReturnDocuments is executed before task Reimbursement. In the first part of the process, we can see that several tasks can be executed in parallel because they all use different resources. However, tasks HotelReservation and Vaccination cannot be executed in parallel because they use the same resource user. After task ReturnDocuments, we can see that the two final tasks are executed in parallel because they use different resources.

Experiments

In this section, we show some experimental results obtained when applying our tool to different BPMN processes. The main goal of this section is to evaluate our tool in terms of performance and accuracy of the results, particularly comparing both strategies (bounded and heuristic-based exploration).

Table 1 shows some experiments on ten processes, mostly taken from the literature, e.g., BPMN processes introduced in [START_REF] Durán | Stochastic Analysis of BPMN with Time in Rewriting Logic[END_REF][START_REF] Durán | A Rewriting Logic Approach to Resource Allocation Analysis in Business Process Models[END_REF][START_REF] Durán | Verifying timed BPMN processes using Maude[END_REF], and from the VBPMN database of examples [START_REF] Krishna | VBPMN: Automated Verification of BPMN Processes[END_REF][START_REF] Krishna | Checking Business Process Evolution[END_REF]. The table first characterizes each process in terms of number of tasks, flows, gateways, and strong flows (SF). Then, the table gives the current Average Execution Time of the process (AET c ) and the Average Execution Time corresponding to the optimal process (AET b ). The optimal process corresponds Fig. 12. Trip Organization process after automated refactoring to the process with the shortest execution time. This optimal process and the corresponding execution time were built and computed manually by the authors of this paper. Then, we show the results for the bounded exploration (bound fixed to 300) and for the heuristic-based approach. For each option, we give the execution time of the final process and the computation time to obtain the result.

Let us now comment on the results shown in this table. First of all, we can see that the bounded exploration succeeds in most cases to find the optimal solution, that is, the process with the shortest execution time. For larger examples (e.g., processes 7, 9 and 10 in the table), the best solution is not found because the bound (300) was too small to explore enough processes and find the best solution. However, if we increase the bound, the best solution is eventually found. As an example, for row 7, the best solution (with 18 as execution time) is found after about 750 iterations. As for the heuristic-based strategy, there are several cases for which the best solution is not found, but in very specific cases (e.g., process 9), it may give better results. Regarding computation time, the heuristic-based strategy is much faster (a few seconds) whereas the bounded exploration takes more time because it explores possibly many solutions (300 at most here). Regarding the computation time for the bounded exploration, we can see that this time increases with the size of the process (it takes more time to simulate, analyze and refactor a larger process than a simpler one). The number of refactorings does not really impact the computation time, since a bound is used to stop the exploration of potential solutions. Finally, the number of strong flows tends to reduce the number of possible refactorings thus the computation time, because strong flows can be seen as additional constraints on the process.

Related Work

This section starts with a short overview of extensions of BPMN with time and resource features, before presenting and comparing our approach with existing solutions for BPMN refactoring. Several works propose extensions of BPMN with time constructs, see, e.g., [START_REF] Arévalo | A Metamodel to Integrate Business Processes Time Perspective in BPMN 2.0[END_REF][START_REF] Gagné | Time-BPMN[END_REF]. In [START_REF] Gagné | Time-BPMN[END_REF], the authors present Time-BPMN, an extension of BPMN that allows the specification of temporal constraints and dependencies within a BPMN diagram. In [START_REF] Arévalo | A Metamodel to Integrate Business Processes Time Perspective in BPMN 2.0[END_REF], a metamodel-based approach to integrate temporal constraints and dependencies is introduced. The time aspects are specified using rules and OCL constraints capture the semantics of these rules. Our solution shares similarities with the approach proposed in [START_REF] Gagné | Time-BPMN[END_REF].

As far as resource allocation is concerned, several solutions have been proposed by the research community in the business process domain. Schömig and Rau [START_REF] Schömig | A Petri Net Approach for the Performance Analysis of Business Processes[END_REF] use colored stochastic Petri nets to specify and analyze business processes in the presence of dynamic routing, simultaneous resource allocation, forking/joining of process-control threads, and priority-based queueing. Li et al. [START_REF] Li | Performance Modeling and Analysis of Workflow[END_REF] introduce multidimensional workflow nets to model and analyze resource availability and workload. Oliveira et al. [START_REF] Oliveira | Quantitative Analysis of Resource-Constrained Business Processes[END_REF] use generalized stochastic Petri nets for correctness verification and performance evaluation of business processes. In this work, we propose to associate resources to tasks, which is a flexible solution for modelling many different situations and scenarios. Moreover, our approach does not focus on the computation of metrics but changes the process structure to actually optimize some of these metrics.

Let us now focus on existing works on process refactoring. [START_REF] Silingas | Refactoring BPMN Models: From 'Bad Smells' to Best Practices and Patterns[END_REF] presents six common mistakes made by developers when modelling with BPMN: inconsistent naming, large process diagrams, inconsistent use of gateways, inconsistent use of events, inconsistent use of loops, poor diagram layout. For each problem, the authors present best practices for avoiding these issues. As an example, the authors propose to use explicit gateways instead of using multiple incoming/outgoing sequence flows. [START_REF] Dijkman | Identifying Refactoring Opportunities in Process Model Repositories[END_REF] presents a technique for detecting refactoring opportunities in process model repositories. The technique works by first computing activity similarity and then computing three similarity scores for fragment pairs of process models. Using these similarity scores, four different kinds of refactoring opportunities can be systematically identified. As a result, the approach proposes to rename activities or to introduce subprocesses. IBUPROFEN, a business process refactoring approach based on graphs, is presented in [START_REF] Fernández-Ropero | Graph-Based Business Process Model Refactoring[END_REF][START_REF] Pérez-Castillo | Business process model refactoring applying IBUPROFEN. An industrial evaluation[END_REF]. IBUPROFEN defines a set of 10 refactoring algorithms grouped into three categories: maximization of relevant elements, fine-grained granularity reduction, and completeness. All these works mostly focus on syntactic issues and propose synctatic improvements of the process by, for instance, removing unreachable nodes or by merging consecutive gateways of the same type. They do not aim at providing any kind of optimization regarding the process being designed as we do.

In [START_REF] Kumar | Optimizing Process Model Redesign[END_REF], the authors present an approach for optimizing the redesign of process models. It is based on capturing process improvement strategies as constraints in a structural-temporal model. Each improvement strategy is represented by a binary variable. An objective function that represents a net benefit function of cost and quality is then maximized to find the best combination of process improvements that can be made to maximize the objective. The BPMN subset used in [START_REF] Kumar | Optimizing Process Model Redesign[END_REF] is very similar to the one we use in this paper. However, the approach is rather different since they compute optimal redesigns with respect to some constraints, whereas we propose refactoring patterns with respect to process execution times.

Last but not least, it is worth mentioning recent works providing support for building (optimal) processes. [START_REF] Falcone | Semi-automated Modelling of Optimized BPMN Processes[END_REF] proposes a semi-automated approach for helping non-experts in BPMN to model business processes using this notation. Alternatively, [START_REF] Ivanchikj | From Text to Visual BPMN Process Models: Design and Evaluation[END_REF] presents an approach which combines notes taking in constrained natural language with process mining to automatically produce BPMN diagrams in real-time as interview participants describe them with stories. In this work, we tackle this issue from a different angle since we assume that an existing description of the process exists and that we want to automatically optimize it by updating its structure.

Concluding Remarks

In this paper, we have focused on a version of BPMN including task durations and an explicit description of resources. We have then proposed a simulation-based approach that allows us to identify some specific tasks which are waiting for being executed but for which the required resources are available. This means that, from the point of view of the process structure, these tasks could be executed earlier in the process. We then apply some refactoring transformations to move those tasks backwards in the process structure. This approach works by successively applying these refactorings and by thus exploring the possible solutions to find the optimal one. Several strategies have been implemented and vary in their way to apply these iterations. In any case, the refactoring process completes and returns as result a process whose average execution time is lower (or equal) than the one of the original process. Note that if the original process is already optimal, that process and its corresponding execution time are returned as output. All the steps of the refactoring approach are fully automated by a tool we implemented and applied on many examples of processes for validation purposes.

Fig. 1 .

 1 Fig. 1. Supported BPMN syntax

Fig. 3 .

 3 Fig. 3. Sequence of tasks

Fig. 4 .

 4 Fig. 4. Merge parallel gateway with preceding tasks (no shared resources)

Fig. 5 .Fig. 6 .Fig. 7 .

 567 Fig. 5. Merge parallel gateway with preceding tasks (shared resources with one task)

Fig. 8 .

 8 Fig. 8. Split parallel gateway

Fig. 9 .Fig. 10 .

 910 Fig. 9. Split exclusive gateway

Fig. 11 .

 11 Fig. 11. Example of "Trip Organization" process in BPMN

Table 1 .

 1 Experimental results

	BPMN		Characteristics		Bounded Explo. Heuristic
	Proc. Tasks Flows Gateways SF AET c AET b AET Time	AET Time
	1	5	13	4	0 70	50	50	1s	50	1s
	2	5	15	6	0 40	25	25	2s	25	2s
	3	6	13 2	2	0 19	6	6	2s	9	1s
	4	8	12	2	4 43	28	28	21s	35	1s
	5	8	12	2	2 21	14	14	23s	16	1s
	6	9	10	0	6 90	50	50	31s	90	1s
	7	10	20	6	0 24	18	19	38s	20	1s
	8	10	20 2	4	3 23	13	13	67s	17	1s
	9	12	29	8	0 200 120 140	50s	120 1s
	10	15	36	10	0 260 180 220	103s	260 3s

Acknowledgments. F. Durán has been partially supported by projects PGC2018-094905-B-100 and UMA18-FEDERJA-180, and by Universidad de Málaga, Campus de Excelencia Internacional Andalucía Tech. This work was also supported by the Région Auvergne-Rhône-Alpes within the "Pack Ambition Recherche" programme.