Defne Surujon

Yann Ponty

Peter Clote
email: clote@bc.edu

Small-world networks and RNA secondary structures

Let S n denote the network of all RNA secondary structures of length n, in which undirected edges exist between structures s, t such that t is obtained from s by the addition, removal or shift of a single base pair. Using context-free grammars, generating functions and complex analysis, we show that the asymptotic average degree is O(n) and that the asymptotic clustering coefficient is O(1/n), from which it follows that the family S n , n = 1, 2, 3, . . . of secondary structure networks is not small-world.

Introduction

Small-world networks, first introduced in [START_REF] Watts | Collective dynamics of 'small-world' networks[END_REF], satisfy two properties:

(1) the shortest path distance between any two nodes is "small" (intuitively, there are six degrees of separation [START_REF] Guare | Six degrees of separation: A play[END_REF] between any two persons), and (2) the average clustering coefficient is large (intuitively, friends of a person tend to be friends of each other). Smallworld networks appear to be ubiquitous in biology, sociology, and information technology; indeed, examples include the neural network of Caenorhabditis elegans [START_REF] Watts | Collective dynamics of 'small-world' networks[END_REF], metabolic networks of 43 organisms representing all three domains of life [START_REF] Jeong | The large-scale organization of metabolic networks[END_REF], the gene co-expression in S. cerevisiae [START_REF] Van Noort | The yeast coexpression network has a small-world, scale-free architecture and can be explained by a simple model[END_REF], protein folding networks, where nodes correspond to conformations (self-avoiding walks on a 2D lattice) and edges exist between nodes that are connected by an elementary move [START_REF] Scala | Small-world networks and the conformation space of a short lattice polymer chain[END_REF], and Markov state models of protein folding networks inferred by the software MSMBuilder [START_REF] Bowman | Using generalized ensemble simulations and Markov state models to identify conformational states[END_REF] from molecular dynamics folding trajectories for the protein villin Bowman and Pande (2010), etc. Additionally, [START_REF] Wuchty | Small worlds in RNA structures[END_REF] showed by exhaustive enumeration of low energy RNA secondary structures of E. coli phe-tRNA, that the corresponding network architecture displays small-world properties. For additional examples, see the excellent review of [START_REF] Albert | Statistical mechanics of complex networks[END_REF].

In this paper, we investigate asymptotic properties of degree and clustering coefficient for the ensemble of all RNA secondary structures, by using methods from algebraic combinatorics. In particular, we rigorously prove that the network of RNA secondary structures is asymptotically not small-world, although it displays strong differences from random networks.

Preliminaries

In this section, we define notions of RNA secondary structure, move sets M S 1 , M S 2 , and small-world networks. An RNA secondary structure of length n, subsequently called length n structure, is defined to be a set s of ordered pairs (i, j), with 1 ≤ i < j ≤ n, such that:

(1) There are no base triples; i.e. if (i, j), (k, ) ∈ s and {i, j} ∩ {k, } = ∅, then i = k and j = . (2) There are no pseudoknots; i.e. if (i, j), (k, ) ∈ s, then it is not the case that i < k < j < . (3) There are at least θ = 3 unpaired bases in a hairpin loop; i.e. if (i, j) ∈ s, then j -i > θ = 3. Note that base pairs are not required to be Watson-Crick or wobble pairs, as is the case for RNA molecules, such as that depicted in Figure 1a. This definition, sometime called homopolymer secondary structure, permits the combinatorial analysis we employ to show that RNA networks are not small-world.

Let S n denote the set of all length n structures. The move sets M S 1 and M S 2 , defined in [START_REF] Flamm | RNA folding at elementary step resolution[END_REF] for RNA secondary folding kinetics, describe elementary moves that transform a structure s into another structure t. Move set M S 1 [resp. M S 2 ] consists of either removing or adding [resp. removing, adding or shifting] a single base pair, provided the resulting set of base pairs constitutes a valid structure, where shift moves are depicted in Figure 2. We overload the notation S n to also denote the M S 1 network [resp. M S 2 network], whose nodes are the length n structures, where an undirected edge between structures s, t exists when t is obtained from s by a single move from M S 1 [resp. M S 2 ]. Figure 1b shows the M S 1 network (8 red edges) [resp. M S 2 network (8 red and 8 blue edges)] for length 7 structures, where there are 8 nodes, M S 1 degree 16 8 = 2 and M S 2 degree 32 8 = 4. See Clote (2015b) and [START_REF] Clote | Network Properties of the Ensemble of RNA Structures[END_REF] for dynamic programming algorithms that compute, respectively, the M S 1 and M S 2 degree for the network of secondary structures of a given RNA sequence.

Small-world networks satisfy two conditions: (1) on average, the minimum path length between any two nodes is small, (2) neighbors of a node tend to be connected to each other. The global clustering coefficient, defined in equation ( 77) of [START_REF] Newman | Random graphs with arbitrary degree distributions and their applications[END_REF], is given by

C g (G) = 3 × number of triangles number of connected triples (1)
where a triangle is a set {x, y, z} of nodes, each of which is connected by an edge, and a (connected) triple is a set {x, y, z} of nodes, such that there is an edge from x to y and an edge from x to z. Following [START_REF] Cont | Small-world graphs: characterization and alternative constructions[END_REF], the family {S n , n = 1, 2, 3, . . .} of RNA networks is small-world if the following conditions hold. (1) There is a constant c 1 ≥ 0, such that the minimum path length between any two nodes of S n is bounded above by c 1 ln n.

(2) There is a constant c 2 ≥ 0, such that the average network degree of S n is bounded above by c 2 ln n.

(3) The global clustering coefficient is bounded away from zero. By Theorem 2, the network size of S n is exponential in n. Since there are at most n/2 base pairs in any length n structure, condition (1) is satisfied for both the M S 1 and M S 2 networks of RNA structures. It is easy to see that the clustering coefficient of the M S 1 network of RNA structures is zero, so in the remainder of the paper, we concentrate on conditions (2) and ( 3) for the M S 2 RNA network. Specific properties of RNA networks critically depend on the chosen definition of neighborhood of a structure, leading to the investigation of various move sets in the RNA kinetics literature. In addition to the move sets M S 1 and M S 2 [START_REF] Flamm | RNA folding at elementary step resolution[END_REF], which latter also models defect diffusion [START_REF] Pörschke | Model calculations on the kinetics of oligonucleotide double-helix coil transitions: Evidence for a fast chain sliding reaction[END_REF], other groups have considered more general move sets that allow helix formation and disassociations [START_REF] Isambert | Modeling RNA folding paths with pseudoknots: application to hepatitis delta virus ribozyme[END_REF].

The overall method used is as follows:

(1) Give a context-free grammar that generates the set of all secondary structures, possibly containing a specific motif. (2) Use Table 1 to derive and then solve a functional relation for the complex generating function S(z), with the property that the nth Taylor coefficient of S(z), denoted [z n ]S(z), is equal to the number of length n structures, possibly containing a specific motif. (3) Determine the dominant singularity and apply complex analysis [START_REF] Flajolet | Singularity analysis of generating functions[END_REF] to obtain the asymptotic value of [z n ]S(z). For step (3), we use the Flajolet-Odlyzko Theorem, stated as Corollary 2, part (i) on page 224 of [START_REF] Flajolet | Singularity analysis of generating functions[END_REF]. Before stating the theorem, we define the dominant singularity of complex function f (z) to be the complex number ρ having smallest absolute value (or modulus) at which f (z) is not differentiable.

Theorem 1 (Flajolet and Odlyzko) Assume that f (z) has a dominant singularity at z = ρ > 0, is analytic for z = ρ satisfying |z| ≤ |ρ|, and that

lim z→ρ f (z) = K(1 -z/ρ) α .
(2)

Then, as n → ∞, if α / ∈ 0, 1, 2, ..., f n = [z n ]f (z) ∼ K Γ(-α) • n -α-1 • ρ -n
where ∼ denotes asymptotic equality and Γ denotes the Gamma function.

The plan of the paper is now as follows. In Section 3, we show that the average M S 2 degree of S n is O(n).

Expected network degree

Due to space constraints, details for the computation of the asymptotic number of secondary structures as well as for M S 1 expected degree for homopolymers cannot be given in this paper. Nevertheless, these computations can be found in Clote (2015a), from which we take the following results. Recalling the notation ∼ for asymptotic equality, we have

Theorem 2 If S(z) is the generating function for the number of secondary structures for a homopolymer, then

[z n ]S(z) ∼ 0.713121 • n -3/2 • 2.28879 n
If M S 1 degree(n) denotes the M S 1 expected network degree for a homopolymer, then

M S 1 degree(n) ∼ 0.473475 • n
Define the grammar G to consist of the terminal symbols ( , •, ) , , , , nonterminal symbols S, T , S, R, θ, with start symbol S. Shift moves are represented in the grammar by one of the three expressions: , , , as depicted in Figure 2. In particular, represents the right shift depicted in Figure 2a (ignoring possible intervening structure), where base pair (x, y) is transformed to (x, y ) for x < y < y; alternatively, the can represent the shift (x, y) to (x, y ) for x < y < y , as depicted in Figure 2b. The expression can represent the left shift depicted in Figure 2c, where base pair (x, y) is transformed to (x , y) for x < x < y; alternatively, can represent the shift (x, y) to (x , y) for x < x < y, as depicted in Figure 2d. The expression can represent the right-to-left shift depicted in Figure 2e, where base pair (x, y) is transformed to (y , x) for y < x < y; alternatively, can represent the shift (x, y) to (y, x ) for x < y < x , as depicted in Figure 2f. The grammar G allows us to count the number of secondary structures, that additionally contain a unique occurrence of exactly one of the three expressions: , , . Since two shift moves correspond to each of the previous three expressions, it follows that the total number of M S 2 -M S 1 (shift-only) moves, summed over all structures for a homopolymer of length n with θ = 1, is equal to 2[z n ]S † (z).

The production rules of grammar G are as follows:

S → S • | ( S ) | S ( S ) | S ( R ) | T T → R | S R | R S | S R S | R | S R | S R | S S R | R R | S R R S → • | S • | ( R ) | S ( R ) R → θ | R • | ( R ) | S ( R ) θ → • • • (3)
The nonterminal S is responsible for generating all secondary structures of length greater than or equal to 1. In contrast, the nonterminal S is responsible for generating all wellbalanced expressions of length greater than or equal to 1, that involve exactly one of the three expressions: , , . To that end, the nonterminal T is responsible for generating all such expressions, in which the rightmost symbol is either or , but not • or ) . By induction on length of sequence generated, one can show that G is an nonambiguous context-free grammar that generates all secondary structures having a unique occurrence of one of , , . As mentioned before, 2 times the number of such expressions of length n is equal to the number of M S 2 -M S 1 edges in the network of secondary structures.

As explained in [START_REF] Lorenz | Asymptotics of RNA shapes[END_REF] and [START_REF] Flajolet | Analytic Combinatorics[END_REF], it is possible to automatically transform the previous production rules into equations that relate the corresponding generating functions, where we denote generating functions of S(z), T , S(z), R(z) by the same symbols used for the corresponding nonterminals S, T , S, R. This technique is known in the literature as DSV methodology [START_REF] Lorenz | Asymptotics of RNA shapes[END_REF], or as the symbolic method Flajolet and Sedgewick (2009) -see Table 1. In this fashion, we obtain the following:

S = z S + z 2 S + z 2 S S + z 2 R S + T T = 2z 3 R + 4z 3 RS + 2z 3 RS 2 + z 3 R 2 + z 3 SR 2 S = z + zS + z 2 R + z 2 RS R = θ + zR + z 2 R + z 2 RS θ = z 3
and by eliminating all variables except S and z, we use Mathematica to obtain the quadratic equation in S having two solutions, for which the only solution analytic at 0 is the following:

S(z) = S = A + B √ P C (4)
where

P = 1 -2z -z 2 + z 4 + 3z 6 + 2z 7 + z 8 A = 3 -15z + 23z 2 -9z 3 -z 4 -9z 5 + 23z 6 -25z 7 + 7z 8 -z 9 + 6z 10 - 8z 11 + 2z 12 + 2z 13 + 2z 14 B = -3 + 12z -14z 2 + 4z 3 + 5z 5 -10z 6 + 8z 7 -2z 10 C = 2(-z 3 + 3z 4 -z 5 -z 6 -z 7 + z 8 - 3z 9 + z 10 + z 11 + z 12 )
The dominant singularity ρ of S(z) in equation ( 4) is the complex number having smallest absolute value (or modulus) at which S(z) is not differentiable. For the functions in this paper, the dominant singularity will always be the (complex) root of polynomial P under the radical, having smallest modulus -since the square root function is not differentiable over the complex numbers at zero. Letting F (z) = B √ P C and noting that the dominant singularity ρ = 0.436911, a calculation shows that

lim z→ρ F (z) = lim z→ρ B • √ P • (1 -z/ρ) 1/2 C • (1 -z/ρ) P = P 1 -z/ρ = 1 + 0.288795z -0.339007z 2 - 0.775919z 3 -0.775919z 4 -1.775919z 5 - 1.064714z 6 -0.436911z 7 C = C 1 -z/ρ = -2z 3 + 1.422410z 4 + 1.255605z 5 + 0.873822z 6 + 2z 8 -1.422410z 9 - 1.255605z 1 0 -0.873822z 1 1 and so lim z→ρ F (z) = 0.684877 • lim z→ρ (1 -z/ρ) -1/2 = 0.684877 • lim z→ρ (1 -z/0.436911) -1/2
Taking α = -1/2 in the Flajolet-Odlyzko Theorem [START_REF] Flajolet | Singularity analysis of generating functions[END_REF], we obtain:

[z n ] F (z) ∼ 0.684877 Γ(1/2) • n -1/2 • 1 ρ n = 0.3864 • n -1/2 • 2.28879 n
By Theorem 2 the asymptotic number of secondary structures for a homopolymer when θ = 3 is 0.713121 • n -3/2 • 2.28879 n , and so we have the following result.

Theorem 3 The asymptotic M S 2 -M S 1 degree of S n is 2[z n ] F (z) [z n ]S(z) ∼ 0.772801 • n -1/2 • 2.28879 n 0.713121 • n -3/2 • 2.28879 n = 1.083688 • n
Adding the asymptotic values from Theorem 2 and Theorem 3, we determine the M S 2 degree.

Corollary 4 The asymptotic M S 2 degree for the network S n of RNA structures is 1.557164• n.

Using a Taylor series expansion at zero for the functions used to determine both the M S 1 and M S 2 -M S 1 degree, we have verified that the numerical results for S n are identical with those independently computed by the dynamic programming C-implementations described in [START_REF] Clote | Expected degree for RNA secondary structure networks[END_REF] and [START_REF] Clote | Network Properties of the Ensemble of RNA Structures[END_REF]. We also note that the current approach is much simpler than the program in [START_REF] Clote | Network Properties of the Ensemble of RNA Structures[END_REF], although the latter is more general, since it computes the M S 2 degree for any user-specified RNA sequence. Using well-known methods, these asymptotic results can be extended from homopolymers to RNA sequences with Watson-Crick and wobble base pairs by using a "stickiness model", which stipulates the probability p that any two positions can form a base pair, defined by

p = 2 (p A p U + p G p U + p G p C ) (5)
where p A , p C , p G , p U are user-specified nucleotide relative frequencies. Since we consider shifts, we need an additional stickiness parameter q, which specifies the probability that a shift can occur between three randomly selected positions in which one position is fixed, defined by

q = p A p 2 U + p C p 2 G + p G p 2 C + p 2 U + 2p C p U + p U p 2 A + p 2 G + 2p A p G (6)
By including stickiness parameters into our computations, we obtained values presented in Table 2, which shows asymptotic M S 1 and M S 2 degrees for a number of classes of RNA.

4 Asymptotic M S 2 clustering coefficient Section 4.1 describes a grammar to count the number of triangles for S n with respect to M S 2 moves, while Section 4.2 describes a grammar to count two particular triples.

Counting triangles

Let G be the grammar with terminal symbols ( , ) , •, , nonterminal symbols S , S 1 , . . . , S 8 , S, R, X, θ, start symbol S and the following production rules:

S → S 1 | S 2 | S 3 | S 4 | S 5 | S 6 | S 7 | S 8 S → • | S • | ( R ) | S ( R ) R → θ | R • | ( R ) | S ( R ) X → λ | R θ → • • •
where λ denotes the empty word, and S 1 , . . . , S 8 are specified in the following 8 exhaustive and mutually exclusive cases. Note that S 1 , . . . , S 3 generate structures containing type A triangles, while S 4 , . . . , S 8 generate structures containing type B triangles.

Rule 1

The following productions generate all secondary structures s, such that for x < y < z, it is the case that s∪{(x, y)} and s∪{(y, z)} are also secondary structures, hence form a triangle:

S 1 → S 1 • | ( S 1 ) | S ( S 1 ) | S 1 ( R ) | X R R
with corresponding DSV equations

S 1 = zS 1 + z 2 S 1 + z 2 SS 1 + z 2 RS 1 + Xz 3 R 2 Rule 2
The following productions generate all secondary structures s, such that for x < y < z, it is the case that s∪{(x, y)} and s∪{(x, z)} are also secondary structures, hence form a triangle:

S 2 → S 2 • | ( S 2 ) | S ( S 2 ) | S 2 ( R ) | X R X
with corresponding DSV equations

S 2 = zS 2 + z 2 S 2 + z 2 SS 2 + z 2 RS 2 + X 2 z 3 R Rule 3
The following productions generate all secondary structures s, such that for x < y < z, it is the case that s∪{(x, z)} and s∪{(y, z)} are also secondary structures, hence form a triangle:

S 3 → S 3 • | ( S 3 ) | S ( S 3 ) | S 3 ( R ) | X X R
with corresponding DSV equations

S 3 = zS 3 + z 2 S 3 + z 2 SS 3 + z 2 RS 3 + X 2 z 3 R Rule 4
The following productions generate all secondary structures s, such that for x < y < z, it is the case that s ∪ {(x, y)}, s ∪ {(x, z)} and s ∪ {(x, w)} are also secondary structures, hence the latter form a triangle:

S 4 → S 4 • | ( S 4 ) | S ( S 4 ) | S 4 ( R ) | X R X X
with corresponding DSV equations

S 4 = zS 4 + z 2 S 4 + z 2 SS 4 + z 2 RS 4 + X 3 Rz 4
Rule 5

For x < y < z < w, let s 1 = (x, w), s 2 = (y, w), s 3 = (z, w). The following productions generate all secondary structures s, such that for x < y < z, it is the case that s ∪ {(x, w)}, s ∪ {(y, w)} and s ∪ {(z, w)} are also secondary structures, hence the latter form a triangle:

S 5 → S 5 • | ( S 5 ) | S ( S 5 ) | S 5 ( R ) | X X X R
with corresponding DSV equations

S 5 = zS 5 + z 2 S 5 + z 2 SS 5 + z 2 RS 5 + X 3 z 4 R Rule 6
For x < y < z < w, the following productions generate all secondary structures s, such that for x < y < z, it is the case that s ∪ {(x, y)}, s ∪ {(y, z)} and s ∪ {(y, w)} are also secondary structures, hence the latter form a triangle:

S 6 → S 6 • | ( S 6 ) | S ( S 6 ) | S 6 ( R ) | X X R R
with corresponding DSV equations

S 6 = zS 6 + z 2 S 6 + z 2 SS 6 + z 2 RS 6 + X 2 z 4 R 2 Rule 7
For x < y < z < w, the following productions generate all secondary structures s, such that for x < y < z, it is the case that s ∪ {(x, z)}, s ∪ {(y, z)} and s ∪ {(z, w)} are also secondary structures, hence the latter form a triangle:

S 7 → S 7 • | ( S 7 ) | S ( S 7 ) | S 7 ( R ) | X R R X
with corresponding DSV equations

S 7 = zS 7 + z 2 S 7 + z 2 SS 7 + z 2 RS 7 + X 2 z 4 R 2

Rule 8 bis

The following productions generate all secondary structures s, such that for x < y < z, it is the case that s ∪ {(x, z)}, s ∪ {(x, y)} and s ∪ {(y, z)} are also secondary structures, hence the latter form a triangle. This grammar is identical to that in rule 1 above, with the exception that S 1 is replaced by S 8 . Let S (z) denote the generating function for the number of structures containing a unique triangle motif, where triA(z) [resp. triB(z)] is the generating function for the collection of structures containing a unique occurrence of type A [type B] triangle, as treated in rules 1-3 [resp. rules 4-8]. We obtain the following compact form for the DSV equations for the grammar G that generates all structures containing a triangle:

S = triA + triB triA = triA • z + X • z • triA • z+ triA • z • R • z + X • z • R • z • R • z+ X • z • R • z • X • z + X • z • X • z • R • z triB = triB • z + X • z • triB • z+ triB • z • R • z + X 3 z 4 R + X 3 z 4 R+ X 2 z 4 R 2 + X 2 z 4 R 2 + Xz 3 R 2
Using Mathematica, we determine the following.

[z n ]S (z) = 0.870311 • 2.28879 n • n -1/2
By Theorem 2, the asymptotic number of secondary structures is 0.713121 • n -3/2 • 2.28879 n , and so we have the following result.

Theorem 5 The asymptotic average number of triangles per structure is

[z n ]S (z) [z n ]S(z) ∼ 0.870331 • n -1/2 • 2.28879 n 0.713121 • n -3/2 • 2.28879 n ∼ 1.220453 • n

Counting triples

In this section, we describe a grammar for two particular triples. Let G be the grammar having terminal symbols •, ( , ) , [ , ] , nonterminal symbols S ‡ , S † , S, R, X, θ, start symbol S ‡ , and productions given in equation ( 7) below together with the following:

S → • | S • | X ( R ) R → θ | R • | X ( R ) X → λ | R θ → • • • Triple with motif [ ] [ ] or [ [ ] ]
The following grammar generates all secondary structures s that have two special base pairs (i, j) and (x, y), designated by [ ] , which are either sequential or nested. For each structure s, which contains a unique occurrence of the sequential motif [ ] [ ] or of the nested motif [ [ ] ] , we must count four possible triples: (1) {s 1 , s 2 , s 3 }, where s 1 = s -{(i, j), (x, y)}, s 2 = s -{(i, j)}, s 3 = s -{(x, y)}. (2) {s 1 , s 2 , s 3 }, where s 1 = s, s 2 = s -{(i, j)}, s 3 = s -{(x, y)}. (3) {s 1 , s 2 , s 3 }, where s 1 = s -{(i, j)}, s 2 = s -{(i, j), (x, y)}, s 3 = s.

(4) {s 1 , s 2 , s 3 }, where s 1 = s -{(x, y)}, s 2 = s -{(i, j), (x, y)}, s 3 = s. For this reason, we multiply by 4 the asymptotic number of structures generated by the following grammar G. The grammar G has terminal symbols •, ( , ) , [ , ] , nonterminal symbols S ‡ , S † , S, R, X, θ, start symbol S ‡ , and the following production rules.

S ‡ → S ‡ • | ( S ‡ ) | S ( S ‡ ) | S ‡ ( R ) | [ S † ] | S [ S † ] | S † [ R ] | S † ( S † ) S † → S † • | ( S † ) | S ( S † ) | S † ( R ) | [ R ] | S [ R ] (7)
When applying the Flajolet-Odlyzko Theorem in the current case, we have ρ = 0.436911 and α = -3/2. A computation shows that

lim z→ρ S ‡ (z) = 0.0177098 (1 -z/ρ)) -3/2 [z n ]S ‡ (z) ∼ 0.0199834 • n 1/2 • 2.28879 n [z n ]S ‡ (z) [z n ]S(z) ∼ 0.0199834 • n 1/2 • 2.28879 n 0.713121 • n -3/2 • 2.28879 n ∼ 0.0280225 • n 2
As mentioned, the number of triples contributed in the current case is 4 times the last value. Thus the expected number of triples involving a structure

containing [ ] [ ] or [ [ ] ] is 4 • 0.0280224 • n 2 = 0.1120896 • n 2 .
Theorem 6 The asymptotic average number of triples per structure, for the triples described in this section, is

4[z n ]S ‡ (z) [z n ]S(z) ∼ 0.11209 • n 2
From Theorems 5 and 6, we obtain an upper bound for the global clustering coefficient, defined in equation ( 1).

Theorem 7 (Bound on global clustering coefficient)

C g (G) = 3 × number of triangles number of connected triples = O 1 n
and hence the family S n , n = 1, 2, 3, . . . of RNA secondary structures is not small-world.

Discussion

In this paper, we have used methods from algebraic combinatorics [START_REF] Flajolet | Analytic Combinatorics[END_REF] to determine the asymptotic average degree and asymptotic clustering coefficient of the M S 2 network S n of RNA secondary structures. Since the clustering coefficient is not bounded away from zero, it follows that the family S n , n = 1, 2, 3, . . ., of networks is not small-world. Our rigorous result differs from conclusions drawn from computer simulations of [START_REF] Bowman | Protein folded states are kinetic hubs[END_REF]; [START_REF] Wuchty | Small worlds in RNA structures[END_REF], which suggest that molecular folding networks are small-world. However, this paradoxical result is possible, since the notion of finite smallworld network is not precisely defined due to absence of an exact bound for both average path length between any two nodes and for the clustering coefficient. of secondary structures for move sets MS 1 and MS 2 for different values of threshold θ and base pair and triple stickiness parameters p and q, defined in Equations ( 5) and ( 6) respectively. Five models are considered: hp -homopolymer model with p = q = 1; wc -Watson-Crick pairing model with uniform compositional frequency p A = p C = p G = p U = 1 4 , hence p = 1 4 , and q = 0.0625); wcw -Watson-Crick and wobble pairing model with uniform compositional frequency, hence p = 3 8 , and q = 0.15625; wc tRNA -Watson-Crick base pairing model based on the compositional frequencies (p A = 1288 4534 , p U = 1029 4534 , p G = 1223 4534 , p U = 994 4534 ) observed in family RF00005 of 4,534 tRNAs in the Rfam 12.0 database [START_REF] Nawrocki | Rfam 12.0: updates to the RNA families database[END_REF], hence p = 0.259427, q = 0.066394; wcw tRNA -Watson-Crick and Wobble base pairing model using compositional frequency of RF00005, so that p = 0.377699, q = 0.158476. x y′ y (a) (x, y) → (x, y )

Type of nonterminal Generating function

A → B | C A(z) = B(z) + C(z) A → B C A(z) = B(z)C(z) A → t A(z) = z A → ε A(z) = 1
x y′ y of motif 1, . . . , 12. Note that S2 is not the same as R 2 -the former nonterminal S2 simply generates all secondary structures of length 2 or greater.

(b) (x, y) → (x, y ) x x′ y (c) (x, y) → (x , y) x x′ y (d) (x, y) → (x , y) x y′ y (e) (x, y) → (y , x) x x′ y (f) (x, y) → (y, x )
Λ B → Λ B • | X ( Λ B ) | Λ B ( R ) | R 1 | • • • | R 12 (type B triples) R 1 → X R R X R 2 → X X R X R 3 → X R S2 | X R S2 | X R S S R 4 → X S2 | X S2 S2 | X S2 S | X S S2 | X S S S R 5 → X X R X R 6 → X X R X R 7 → X R X R R 8 → X R X R R 9 → X R R R R 10 → X X R R R 11 → X R X R R 12 → X S2 R | X S S R | X S2 R

A.3 Type C and D triples

The nonterminal symbol triples generates all secondary structures containing a unique occurrence of one of the 25 type C or D triple motifs, as depicted in Figure 9 together with the following, where we write Λ C in place of Λ CD for reasons of brevity , where s 2 , s 3 are obtained from s 1 by adding a base pair with the property that there is a shift move from s 2 to s 3 . This type of triangle is described in rules 1,2,3 in Section 4.1. For instance, the motif for triangle rule 1 indicates that there is a base pair (x, y) ∈ s 2 that can be shifted to the base pair (y, z) ∈ s 3 ; similarly, the motif for rule 2 [resp. 3] indicates that there is a base pair (x, y) ∈ s 2 [resp. (x, z) ∈ s 2 ] that can be shifted to the base pair (x, z) ∈ s 3 [resp. (y, z) ∈ s 3 ]. Figure 6: Type B triangles are constituted by structures s 1 , s 2 , s 3 , where s 2 , s 3 are obtained from s 1 by a shift with the property that there is a shift move from s 2 to s 3 . This type of triangle is described in rules 4-8 in Section 4.1. For instance, the motif in triangle rule 4 indicates that there is a base pair (x, y) ∈ s 1 can be shifted to base pair (x, z) ∈ s 2 and (x, w) ∈ s 3 . The other panels have analogous meanings. Figure 7: Type A triples are constituted by structures s 1 , s 2 , s 3 , where s 2 , s 3 are obtained from s 1 by a base pair addition with the property that there is a no M S 2 move between s 2 and s 3 . Type A triples are given in rules 1-3 in Section 4.2, are are represented in panels 1,2,3,4 of this figure. Panel 1 [resp. 2] indicates that structures s 2 , s 3 can be obtained from structure s 1 by the addition of disjoint base pairs (u, v) (x, y) which are not nested, i.e. ( ) ( ) , [resp. which are nested, i.e. ( ( ) ) ]. Panel 3 indicates that structures s 2 , s 3 can be obtained from structure s 1 by the addition of disjoint base pairs (u, v) (x, y) which would form a pseudoknot if added simultaneously to s 1 , i.e. ( [ ) ] . Panel 4, which is identical to panel 1 of Figure 5, represents a non-triangular connected triple only in the non-polymer case. This panel indicates that structures s 2 , s 3 can be obtained from structure s 1 by the addition of non-disjoint base pairs (u, v) (v, w) which share a base. To each triple motif that has 2 extended connected components (see text) there corresponds a quadrilateral, as shown in the panel with label 2-component motif, where edges are labeled by a/r for base pair addition/removal. To each corner of the quadrilateral, there corresponds a unique triple -thus, panels 1 and 2 actually each represent 4 triples. It follows that the average number of triples per structure for type A(1) and type A(2) triples must be multiplied by 4. The same remark holds for type C,D triples in in Figure 9, which have 2 extended connected components. Figure 8: Type B triples are of the form s 1 , s 2 , s 3 , where s 1 , s 2 are connected by a shift, as are s 1 , s 3 , but s 2 , s 3 are not connected by any M S 2 move. Note that all motifs are connected.

Λ B = R1 + R2 + R3 + R4 + R5 + R6 + R7 + R8 + R9 + R10 + R11 + R12+ Λ B z + Xz 2 Λ B + Λ B z 2 R Λ C 1 = XzS † zXz Λ C 2 = Xz 5 (S2 + S2 2 + 2SS2 + S 3 ) Λ C 367 = Xz∆ A z Λ C 4 = X 4 z 5 R Λ C 5 = XzRzS † z Λ C 8 = XzRzS † z Λ C 9 = X 2 z 5 (2S2 + S 2 ) Λ C 10 = X 4 z 5 R Λ C 11 = X 2 z 3 S † Λ C 12 = X 2 z 5 R(2S2 + S 2 ) Λ C 13 = X 2 z 5 R(2S2 + S 2 ) Λ C 14 = Xz 5 R(S2 + S2 2 + 2SS2 + S 3 ) Λ C 15 = X 2 z 5 (S2 + S2 2 + 2SS2 + S 3 ) Λ C disc = ∆ A zRz + ∆ A zS † z + (S † zRzRz + S † zXzRz + S † zRzXz)(1 + X 2 z 2 ) Λ C 20 = Xz 3 RS † Λ C 21 = Xz 3 RS † Λ C 22 = Xz 5 R(S2 + S2 2 + 2SS2 + S 3 ) Λ C 23 = X 2 z 5 R(2S2 + S 2 ) Λ C = Λ C • z + X • z • Λ C • z + Λ C • z • R • z + 4 • Λ C 1 + Λ C 2 + 4 • Λ C 367 + Λ C 4 + 4 • Λ C 5 + 4 • Λ C 8 + Λ C 9 + Λ C 10 + 4 • Λ C 11 + Λ C 12 + Λ C 13 + Λ C 14 + Λ C 15 + Λ C disc + 4 • Λ C 20 + 4 • Λ C 21 + Λ C 22 + Λ C 23 ∆ A = ∆ A • z + X • z • ∆ A • z + ∆ A • z • R • z + X • z • R • z • R • z+ X • z • R • z • X • z + X • z • X • z • R • z ∆ B = ∆ B • z + X • z • ∆ B • z + ∆ B • z • R • z + X 3 z 4 R + X 3 z 4 R+ X 2 z 4 R 2 + X 2 z 4 R 2 + Xz 3 R 2 ∆ = ∆ A + ∆ B Λ = Λ A + 4 • Λ B + Λ C + 3 • ∆
To each type B triple motif, there corresponds a quadrilateral, shown in panel with label 2-component motif, where edges are labeled by a/r for base pair addition/removal or by s for base pair shift. For each of these 12 motifs, there is a unique base pair that can shift to the remaining two base pairs -for instance, in motif 1, the base pairs are (x, y), (y, z) and (x, u), for x < y < z < u, where (x, y) can be shifted to each of (y, z) and (x, u). This uniquely defined base pair should be located in the quadrilateral diagonally opposite the corner labeled by emptyset. Since each corner of the quadrilateral corresponds to one of 4 triples associated with the motif, it follows that the average number of type B triples per structure must be multiplied by 4. 9,10,12,13,14,15,22,23. The 15 motifs having 2 extended connected component are: 1, 3, 5, 6, 7, 8, 11, 16, 17, 18, 19, 20, 21, 24, 25. Of the latter, the six motifs 16,17,18,19,24,25 are disconnected successive, and the nine motifs 1, 3,5,6,7,8,11,20,21 are disconnected nested. To each motif, which has 2 extended connected components, there actually correspond 4 triples, as explained in the caption to

Figure 1 :

 1 Figure 1: (a) Consensus secondary structure of the type III hammerhead ribozyme from Peach Latent Mosaic Viroid (PLMVd) AJ005312.1/282-335 (isolate LS35, variant ls16b), taken from Rfam Gardner et al. (2011) family RF00008. (b) Network for size 7 homopolymer with θ = 3, having 8 nodes and 8 red M S 1 edges (base pair addition or removal), 8 blue M S 2 -M S 1 edges (base pair shift), hence a total of 16 M S 2 edges. It follows that M S 1 degree is 16 8 = 2, while MS2 is 32 8 = 4.

Figure 2 :

 2 Figure 2: Illustration of possible shift moves, where each subcaption indicates the terminal symbols involved in the corresponding production rule.

Figure 3 :

 3 Figure 3: (a) M S 2 -degree distribution for the 106,633 secondary structures for a 20-nt homopolymer with θ = 3 (green shaded curve), with Poisson distribution of the same mean. (b) M S 2 -degree distribution for the 32 nt selenocystein insertion (SECIS) element fruA with sequence CCUCGAGGGG AACCCGAAAG GGACCCGAGA GG.

Figure 4 :Figure 5 :

 45 Figure4: Complete listing of all possible moves for triangles and non-triangular connected triples with a designated first structure, up to equivalence. Consider, for instance, the triangle T (not shown) in which s 1 → s 2 by a shift, s 1 → s 3 by base pair removal, and s 3 → s 2 by a base pair addition. Then T is equivalent to a triangle of type A, where node 1 is occupied by s 3 , node 2 resp 3 is occupied by s 1 resp. s 2 . Other instances of triangles the reader may consider are analogously equivalent to a triangle of type A or B.

Figure 9 :

 9 Figure9: Type C and D triples are of the form s 1 , s 2 , s 3 , where s 1 , s 2 are connected by a shift, and s 3 is obtained from s 1 by a base pair addition, but s 2 , s 3 are not connected by any M S 2 move. Triples may have either one or two extended connected components (see text); in the former case, the asymptotic expected number of triples is O(n), while in the latter case the expected expected number of triples is O(n 2 ), as is the case for the type A triple motifs in panels 1 and 2 of Figure7. The 10 motifs having 1 extended connected component are: 2, 4, 9, 10, 12, 13, 14, 15, 22, 23. The 15 motifs having 2 extended connected component are: 1,3, 5, 6, 7, 8, 11, 16, 17, 18, 19, 20, 21, 24, 25. Of the latter, the six motifs16,17,18,19,24,25 are disconnected successive, and the nine motifs 1,3,5,6,7,8,11,20,21 are disconnected nested. To each motif, which has 2 extended connected components, there actually correspond 4 triples, as explained in the caption to Figure7. For instance, to the motif 16, given by undirected edges (1, 2), (3, 4), (3, 5), there correspond four structures s 1 , s 2 , s 3 , s 4 , where (3, 4) ∈ s 1 , (3, 5) ∈ s 2 , (1, 2), (3, 4) ∈ s 3 , (1, 2), (3, 5) ∈ s 4 , with the following four triples: (1)s 1 → s 2 , s 1 → s 3 (type C triple); (2) s 2 → s 3 , s 2 → s 1 (type C triple); (3) s 3 → s 4 , s 3 → s 1 (type D triple); (4) s 4 → s 3 , s 4 → s 2 (type D triple). This is analogous to the situation summarized in the panel in Figure7with label 2-component motif. Note that motifs 3,6,7 have a type A triangle contained within an outer designated base pair [ ] .

Figure 7 .

 7 Figure9: Type C and D triples are of the form s 1 , s 2 , s 3 , where s 1 , s 2 are connected by a shift, and s 3 is obtained from s 1 by a base pair addition, but s 2 , s 3 are not connected by any M S 2 move. Triples may have either one or two extended connected components (see text); in the former case, the asymptotic expected number of triples is O(n), while in the latter case the expected expected number of triples is O(n 2 ), as is the case for the type A triple motifs in panels 1 and 2 of Figure7. The 10 motifs having 1 extended connected component are: 2, 4, 9, 10, 12, 13, 14, 15, 22, 23. The 15 motifs having 2 extended connected component are: 1,3, 5, 6, 7, 8, 11, 16, 17, 18, 19, 20, 21, 24, 25. Of the latter, the six motifs16,17,18,19,24,25 are disconnected successive, and the nine motifs 1,3,5,6,7,8,11,20,21 are disconnected nested. To each motif, which has 2 extended connected components, there actually correspond 4 triples, as explained in the caption to Figure7. For instance, to the motif 16, given by undirected edges (1, 2), (3, 4), (3, 5), there correspond four structures s 1 , s 2 , s 3 , s 4 , where (3, 4) ∈ s 1 , (3, 5) ∈ s 2 , (1, 2), (3, 4) ∈ s 3 , (1, 2), (3, 5) ∈ s 4 , with the following four triples: (1)s 1 → s 2 , s 1 → s 3 (type C triple); (2) s 2 → s 3 , s 2 → s 1 (type C triple); (3) s 3 → s 4 , s 3 → s 1 (type D triple); (4) s 4 → s 3 , s 4 → s 2 (type D triple). This is analogous to the situation summarized in the panel in Figure7with label 2-component motif. Note that motifs 3,6,7 have a type A triangle contained within an outer designated base pair [ ] .

Table 1 :

 1 Translation between context-free grammars and generating functions. Here, G = (V, Σ, S, R) is a given context-free grammar, A, B, C are any nonterminal symbols in V , and t is a terminal symbol in Σ. The generating functions for the languages L(A), L(B), L(C) are respectively denoted by A(z), B(z), C(z).

	Move set θ	hp	wc	wcw	wc tRNA wcw tRNA
	MS 1	1 0.55279 0.42265 0.46158 0.42157	0.46021
	MS 1	3 0.47348 0.35130 0.38531 0.35038	0.38408
	MS 2 -MS 1 1 1.44721 0.57735 0.84796 0.57985	0.85119
	MS 2 -MS 1 3 1.08369 0.42908 0.31409 0.21550	0.63055
	MS 2	1 2.00000 1.00000 1.30954 1.00142	1.31140
	MS 2	3 1.55717 0.43527 0.32336 0.22162	0.63971

Table 2 :

 2 The asymptotic expected degree of the network S θ,p n

  . Nonterminals Λ C 1 , Λ C 2 , Λ C 367 , Λ C 4 , Λ C 5 , Λ C 8 , Λ C 9 , Λ C 10 , Λ C 11 , Λ C 12 , Λ C 13 , Λ C 14 , Λ C 15 , Λ C disc , Λ C 20 , Λ C 21 , Λ C 22 , Λ C 23generate respectively the collection of secondary structures containing a unique occurrence of motif 1, . . . , 25, whereby nonterminal Λ C 367 is the same rule for motif 3,6 and 7, and nonterminal Λ C disc is the same rule for all disconnected successive motifsi.e. the 6 type C and D triple motifs16, 17, 18, 19, 24, 25 which have exactly 2 extended connected components of successive form [ ] [ ] .

Acknowledgments

This research was supported in part by National Science Foundation grant DBI-1262439 to PC and the French/Austrian RNALands project (ANR-14-CE34-0011 and FWF-I-1804-N28) to YP. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the authors and do not necessarily reflect the views of the National Science Foundation.

A Supplement: Complete grammar

In this supplement, we list the complete grammar for all triangles and connected triples, together with the corresponding functional equations. We do not provide any detailed explanation for the complete listing of all possible triples or the complete grammar, since our intent is to provide suggestive supplementary figures and a general orientation for the reader wishing to work through the Mathematica code, available upon request, for our computation of the asymptotic clustering coefficient value 20.7728/n.

A.1 Triangles and type A triples

Nonterminal S generates all non-empty secondary structures; R [resp. S2] generates all secondary structures of length at least 3 [resp. 2]. The non-terminal S † [resp. S ‡ ] generates all secondary structures containing a unique occurrence of a motif that has exactly 1 [resp. 2] extended connected components. Thus S ‡ generates all structures having a unique occurrence of motif 1 or 2 for type A triple, as shown in Figure 7.

A.2 Type B triples

The nonterminal symbol Λ B generates all secondary structures containing a unique occurrence of one of the 12 type B triple motifs, as depicted in Figure 8. Nonterminals R 1 , . . . , R 12 generate respectively the collection of secondary structures containing a unique occurrence This gives rise to the following functional equations for all connected triples, both nontriangular triples, as well as 3 triples associated with each triangle. Since rules R 1 , . . . , R 12 are written as R1, . . . , R12, we write S2 in place of R2, which latter had been defined by


Alasdair Newson

Andrés Almansa

Yann Gousseau

Patrick Pérez

TEMPORAL FILTERING OF LINE SCRATCH DETECTIONS IN DEGRADED FILMS

Keywords: Film defects, line scratches, affine dominant motion, temporal filtering

The film defect known as the line scratch is difficult to restore automatically due to the large number of false alarms present in scratch detection algorithms. In this paper, an algorithm for dealing with these false alarms is proposed. Validating true scratches, which is the approach generally proposed in the literature, is a difficult task since scratch characteristics are hard to determine, making tracking these defects problematic. Instead, we eliminate false alarms by analysing their compatibility with a global motion estimation. We compare our algorithm with two other scratch detection methods from the literature. Experiments show that our algorithm outperforms these two, and that the proposed temporal filtering greatly improves precision while maintaining high recall.

INTRODUCTION

Recent high-profile projects, such as the restoration of George Melies's "Voyage dans la Lune" (1902) 1 reflect interest in the subject of digital film restoration. The digitalisation and transfer of old films into high definition formats implies that high quality restoration is required.

Common defects in films include dust/dirt, blotches, flicker and line scratches. In this paper, we concentrate on the last of these defects, the line scratch. Line scratches present certain singular characteristics, such as temporal persistence, which require special treatment. Scratches are usually detected by identifying straight, vertical lines which are darker or brighter than their surroundings [START_REF] Kokaram | Detection and removal of line scratches in degraded motion picture sequences[END_REF][START_REF] Bruni | Fast removal of line scratches in old movies[END_REF]. In this paper we present a method to deal with the difficult problem of false alarms due to the line scratch detection process. These false alarms are mostly due to thin vertical scene structures which closely resemble line scratches. Contrary to previous work [START_REF] Joyeux | Tracking and map reconstruction of line scratches in degraded motion pictures[END_REF][START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF], rather than try to validate true scratches, whose characteristics are difficult to determine, we explicitly reject false alarms. Using dominant affine motion estimation, we transform the problem of determining the false alarms' trajectories into the simpler problem of detecting straight lines in a binary detection map. Once these trajectories are detected, the false alarms belonging to the trajectories may be eliminated.

The outlay of the paper is as follows. In Section 2 we briefly recall previous work. In Section 3, we present a new motion estimation-based algorithm for temporal filtering of spatial detections in order to reject false alarms. Finally, experimental validations are presented in Section 4 on a series of degraded film sequences in comparison with previous work.

PRIOR WORK

Line scratch detection algorithms can be divided into two categories: spatial and temporal. Although we concentrate here on temporal filtering, we shall mention a few spatial algorithms, which are necessary for any temporal algorithm.

The first work concerning line scratch detection was carried out by Kokaram [START_REF] Kokaram | Detection and removal of line scratches in degraded motion picture sequences[END_REF], and introduced a scratch model which is widely used in other papers (e.g., [START_REF] Bruni | Line scratches detection and restoration via light diffraction[END_REF] and [START_REF] Bruni | A generalised model for scratch detection[END_REF]). The Hough transform is used in both [START_REF] Kokaram | Detection and removal of line scratches in degraded motion picture sequences[END_REF] and [START_REF] Chishima | A method of scratch removal from old movie film using variant window by hough transform[END_REF] to detect prominent lines. Recently Newson et al. [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF] presented an algorithm which uses a contrario methods [START_REF] Desolneux | From Gestalt Theory to Image Analysis : a Probabilistic Approach[END_REF] for the robust detection line scratches. This algorithm has the advantage of describing the beginning and end points of a scratch precisely. In the present method, we will build upon the work in [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF].

Temporal approaches may be found in [START_REF] Joyeux | Detection and removal of line scratches in motion picture films[END_REF], [START_REF] Joyeux | Film line scratch removal using kalman filtering and bayesian restoration[END_REF], [START_REF] Joyeux | Tracking and map reconstruction of line scratches in degraded motion pictures[END_REF] and [START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF]. The main goal of these temporal algorithms is to validate scratches using certain temporal hypotheses. In the first papers, Joyeux et al. rely on a sinusoidal hypothesis concerning scratch motion. A more recent method by Gullu et al. which exploits temporal information may be found in [START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF]. The authors use the local matching error between frames to determine if a detection is a true scratch or a vertical structure. The main limitation of these approaches is that robust and generic hypotheses concerning the spatio-temporal nature of true scratches are difficult to define.

TEMPORAL FILTERING ALGORITHM

All spatial line scratch detection algorithms suffer from false alarms due to thin vertical structures which closely resemble line scratches. In some situations, it is practically impossible to differentiate the two without prior knowledge concerning the scene structure. Unfortunately this sort of knowledge is difficult to obtain and use. Another possibility is to use the temporal information in the image sequence. This information is useful since the physical motion of a scratch is independent of that of the scene and therefore will differ in general.

The goal of our approach is the rejection of as many false alarms as possible rather than the validation of true scratches as in [START_REF] Joyeux | Tracking and map reconstruction of line scratches in degraded motion pictures[END_REF][START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF]. False alarms exhibit temporal characteristics which are easier to identify than those of true scratches.

False alarms due to vertical scene structures will naturally display similar motion to the underlying scene, since they are part of it. Therefore, we reject any scratch detections which follow a trajectory which conforms to the dominant scene motion. We refer to this hypothesis as the motion coherency hypothesis. Further simple hypotheses may be used for filtering false detections. One example, which we shall use, is a scene cut hypothesis, which stipulates that any trajectories start or stop close to a scene cut must correspond to false alarms.

These criteria do not deal with scratches which move with the scene, stop at a scene cut, or a completely still scene with static scratches. However, such situations are impossible to resolve without prior knowledge on the nature of scratches.

Motion coherency

In order to exploit the motion coherency hypothesis, we must determine the trajectories of the detected line scratches. In the present case, rather than trying to solve the exceedingly difficult task of correctly tracking objects such as scratches (as in [START_REF] Joyeux | Tracking and map reconstruction of line scratches in degraded motion pictures[END_REF]), we shall realign the spatial detections with respect to the motion of the scene. If this is done correctly, the false detections due to vertical scene structures should appear as straight lines in a spatio-temporal representation (on an x-t axis as may be seen in Figure 1 (b)), transforming the problem into one of detecting vertical lines in a binary image. This problem is much more constrained than that of tracking true scratches, and therefore easier to solve.

We suppose that the line scratches have previously been detected using a spatial detection algorithm. The detections of this algorithm should represent the scratches as line segments in an image. In our experiments, we use our previous work of [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF] to carry out this spatial detection.

Let I T be the x-t binary detection map of the spatially detected scratches. In this detection map, each scratch detection is represented by its average column position (on the x-axis) and the index of the frame in which it is found (on the y-axis). An example of this representation may be seen in Figure 1 (a).

Let us introduce some notation. Let S represent a scratch detection. We will refer to this as a segment. Let x(S) and ỹ(S) be the average column and row indices of the segment, and t(S) denote the frame in which segment was detected.

We wish to realign the scratches with the global scene motion. We recall that the motivating idea behind this realignment is to make the false detections appear as straight lines in an x-t representation, and thus make their trajectories easier to determine.

In our work, we use the algorithm of Odobez et al. [START_REF] Odobez | Robust multiresolution estimation of parametric motion models[END_REF] to estimate a dominant affine motion. At a pixel position q = (x q , y q ), the motion vector (u(q), v(q)) is expressed as :

u(q) = c 1 + a 1 x q + a 2 y q v(q) = c 2 + a 3 x q + a 4 y q , (1)
where c 1 and c 2 are the parameters describing the constant motion components, and a 1 ...a 4 are the parameters associated with the spatially varying components of the motion. If (x, y) is a pixel in a frame t and (x , y ) is the corresponding position in frame t + 1, we have the following relationship :

x y 1 = a1 + 1 a2 c1 a3 a4 + 1 c2 0 0 1 x y 1 := Ãt,t+1 x y 1
.

(2) The estimation is carried out between each consecutive pair of frames. For any x, y and frame t we can find the corresponding coordinates (x r ,y r ) in the reference frame r, with the following relationship : 1. Recall, precision and f 1 values comparison, in percentage. We compare results to those of Bruni [START_REF] Bruni | Fast removal of line scratches in old movies[END_REF] and Gullu [START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF].

  x r y r 1   = Ã-1 r,t   x y 1   , (3)
We also compare the results of our previous work [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF] on spatial line scratch detection with those after temporal filtering.

with Ãr,t := Ãt-1,t Ãt-2,t-1 ... Ãr,r+1 .

We realign each segment S with respect to a reference frame r, which yields a realigned binary map I T defined as :

I T (x, t) = 1 if ∃S | xr (S) = x, t(S) = t 0 otherwise (5)
where xr (S) is the average column index of the segment warped to the reference frame. Note that in I T , a point (x, t) whose value is 1 may have several corresponding original segments in the frame t. The trajectories of false alarms should appear in I T as vertical lines (see Figure 1(b)). We determine the scratch trajectories using the line segment detection algorithm presented in [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF]. This algorithm is well adapted to the current task, since it is designed for binary images. It also guarantees a precise and unique description of the trajectories, using the maximality and exclusion principles as defined in [START_REF] Desolneux | From Gestalt Theory to Image Analysis : a Probabilistic Approach[END_REF]. These properties are highly desirable here, since they avoid the task of having to bundle several partial trajectories together. Before detecting the trajectories, we perform a horizontal morphological dilation, with a distance τ x , on I T (x, t). This is necessary, since the spatio-temporal trajectories shown in Figure 1(c) are discretised, and therefore the trajectories may not be precisely detected without dilation.

In our experiments, we set the τ x parameter to 3 pixels.

Let us define a trajectory set T as a set of segments whose x-t positions are within a horizontal distance of τ x of a the x-t positions of a given trajectory. We note here that the exclusion principle of the scratch detection from [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF] implies that each point in I t may belong to at most one trajectory.

Once we have determined the trajectory sets, we decide whether to reject the segments belonging to them. This decision should be based on the behaviour of the global motion in the frames spanned by the trajectory set : if the scene has hardly moved, then we cannot distinguish between a true and a false scratch, and a decision should not be made. For example, in Figure 1 (c), two short trajectories (in red) which correspond to true scratches are detected, but not rejected, since there is no significant motion at that point in the video.

We take this decision by determining the maximum horizontal motion of the scene during the trajectory's timespan. To do this, we inspect the original positions of the segments and reject the trajectory set if, for any two segments Q and R belonging to T , we have :

|x(Q) -x(R)| ≥ τ m , (6)
where τ m is a motion threshold. Since the segments Q and R are considered to belong to the scene (apart from the case of no scene motion), this criterion reflects the maximum accumulated motion of the scene throughout the trajectory. In all of our evaluations, we set the parameter τ m to 10 pixels.

Further filtering condition

Apart from the motion coherency criterion, we also stipulate that trajectories should not start or stop near a scene cut. To enforce this, we simply reject the segments belonging to any trajectory set whose beginning and ending frame indices are within a temporal distance τ c from a scene cut. We set τ c to 4 frames. We used a simple Edge Change Ratio based scene cut detection algorithm (see [START_REF] Lienhart | Reliable transition detection in videos: A survey and practitioner's guide[END_REF]) for the detection of the scene cuts. This can be replaced by other scene cut detection algorithms.

RESULTS

We now present a quantitative evaluation of our algorithm in terms of recall, precision and f 1-score. We compare our algorithm to both the classical spatial detection algorithm of Bruni et al. [START_REF] Bruni | A generalised model for scratch detection[END_REF] and the most recent temporal filtering algorithm of Gullu et al. [START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF]. We also compare the proposed algorithm with our previous work on spatial scratch detection [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF] to ensure that our algorithm induces no degradation of results. Note that in the evaluation of our work, we impose a minimum scratch length of one tenth of the video height for all of our spatial detections. Bruni's algorithm requires a single parameter (the scratch colour) which was set manually for each sequence. Gullu's algorithm, which uses Bruni's work for initial detection, requires two extra parameters apart from the scratch colour. The first parameter is a maximum search distance for a block matching algorithm. We chose a maximum distance of 7 pixels either side, which reflects the maximum motion we expect in the sequence. Secondly, we need a maximum mean absolute difference (MAD) threshold which identifies the presence of a scratch. We chose 15 grey levels for this threshold. These two thresholds are unspecified in [START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF].

We evaluate our algorithm using seven annotated videos. These were annotated by hand by defining the beginning and end points of each scratch. The videos may be found at : http://www.enst.fr/~gousseau/scratchfiltering, as well as all the results and annotations. The first three videos are very common in the scratch detection literature [START_REF] Kokaram | Motion Picture Restoration: Digital Algorithms for Artefact Suppression in Degraded Motion Picture Film and Video[END_REF], but are of limited interest for temporal filtering, as they are short or contain no motion. The last four come from the internet archive (http://www.archive.org), or courtesy of Alvaro Pardo.

Recall is the number of correct detections divided by the total number of detections. The algorithm we use for spatial detection [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF] provides a precise spatial description of the scratches. Therefore, we shall evaluate our algorithm on a pixel-wise basis, that is to say the number of correct detected pixels divided by the total number of correct pixels.

Gullu's algorithm, on the other hand, uses Bruni's spatial detection algorithm which only provides a column index to describe a scratch. Therefore, for recall purposes, we consider that the detected pixels are all those in the detected column.

In terms of precision, we must evaluate our algorithm and that of [START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF] differently. Since the output of Gullu's algorithm is a list of detected columns in each frame, it would be unfair to judge that algorithm on a pixel-wise basis. Therefore, we consider a detection to be correct if it touches at least one annotated scratch pixel. For our algorithm, we evaluate precision on a pixel-wise basis. Note that our evaluation confers a great advantage on Bruni's and Gullu's algorithms.

For evaluation, a detected pixel is "correct" if it is within a horizontal distance of two pixels from an annotated pixel.

Figure 1 presents a visual result of our algorithm. The false alarms due to ropes on the boat are correctly eliminated. Table 3.1 shows the quantitative comparisons.

Our algorithm obtains better f 1 scores on all of the examples except "Star". On the "Star" example, we are at a great disadvantage in terms of precision, as many partial scratches are present, and the algorithms of [START_REF] Bruni | A generalised model for scratch detection[END_REF][START_REF] Gullu | Scratch detection via temporal coherency analysis and removal using edge priority based interpolation[END_REF] only have to touch such scratches once to obtain a correct detection.

Gullu's algorithm may induce some regression in terms of recall, and improves the precision only slightly. This performance may be explained by the extreme sensitivity to outliers in a scratch trajectory: it is sufficient that a single MAD value be quite high to consider a trajectory to be a false alarm.

In comparison to our previous work [START_REF] Newson | Adaptive line scratch detection in degraded films[END_REF], our temporal filtering reduces recall by a very slight amount, while precision generally improves greatly. For the first three examples, there is little temporal information and therefore very little to gain in precision. On the last four sequences, however, we may see that precision increases by up to 22 percent.

CONCLUSION

In this paper we have proposed a new approach to the temporal filtering of line scratch detections in films. The algorithm is based on the observation that false detections due to vertical scene structures move with the surrounding scene, and can therefore be identified and rejected. Robust estimation of dominant affine scene motion simplifies the problem of determining scratch trajectories to a vertical line segment detection problem. The proposed temporal filtering algorithm reduces recall only slightly, while greatly improving precision when temporal information is available. We have compared our algorithm with previous work in terms of recall, precision and f 1-score, and it shows considerably improved performance.

Fig. 1 .

 1 Fig.1. Different x-t binary maps, for the "Afgrunden 2" sequence. The x axis represents the average column of a scratch, and the y axis represents the frame number (t) which the scratch belongs to. Each white point corresponds to a detected scratch. The original video may be viewed at : http://www.enst.fr/~gousseau/scratchfiltering.

Fig. 2 .

 2 Fig. 2. Detections from a frame of "Keldjian" with false alarms due to thin vertical objects removed with temporal filtering.

Partly restored by Technicolor, presented at the Cannes Film Festival

ACKNOWLEDGMENTS

The authors would like to thank Vittoria Bruni and Domenico Vitulano. The authors also acknowledge the internet archive for the content which is available online (http://www.archive.org). The authors express their thanks to Alvaro Pardo and Julieta Keldjian, for the "Keldjian" sequence. Finally, the authors thank Anil Kokaram for his image sequences, which may be found in [14].


Alexandre Le Borgne
email: alexandre.le-borgne@mines-ales.fr

David Delahaye
email: david.delahaye@lirmm.fr

Marianne Huchard
email: marianne.huchard@lirmm.fr

Christelle Urtado
email: christelle.urtado@mines-ales.fr

Sylvain Vauttier
email: sylvain.vauttier@mines-ales.fr

Recovering Three-Level Architectures from the Code of Open-Source Java Spring Projects

Keywords: Component-Based Software Engineering, Model-driven engineering, Architecture retro-engineering from code, Architecture evolution, Component reuse, Architecture reuse

Despite the well-admitted benefits of keeping design decisions as a documentation all along the lifecycle of software, many software projects have lost this information. In order to use design information to guide software maintenance and evolution, this paper proposes to retroengineer architecture descriptions from source code. The originality of this work is to target a three-leveled architecture description language which represents software specification, configuration and deployment. Retro-engineering these three levels will provide a more precise source of guidance for the maintenance of software. Targeted projects are open-source Java projects that use Spring to describe the implemented "architecture".

I. Introduction

As software systems constantly become more complex, retrieving design decisions has become an increasingly important problematics when conceptual documentation is missing. However, despite numerous researches in the field of software architecture reconstruction, few work was dedicated to extract raw ("as implemented") componentbased description. It is important, in the first place, to understand design decisions to recover architectures as they are implemented and to not perform any improvement (re-engineering tasks) altogether. Moreover, it is important to represent the software at three abstraction levels in order to trace design decisions through the whole development process. To do so, we use the Dedal [START_REF] Zhang | Architecture-centric component-based development needs a three-level ADL[END_REF], [START_REF] Mokni | A formal approach for managing component-based architecture evolution[END_REF] architecture description language (ADL) developed by our team. This paper proposes to reconstruct component-based architectures from Java Spring [START_REF] Johnson | The spring framework -reference documentation[END_REF] projects.

The remainder of this paper is organized as follows. Section II presents the background of our approach. The core of the paper is developed in Section III where component-based architecture reconstruction from Java Spring projects is explained. Section IV details the existing work in software architecture reconstruction and Section V concludes on perspectives.

II. Background

A. Dedal, a Three-Level Architecture Description Language Dedal [START_REF] Zhang | Architecture-centric component-based development needs a three-level ADL[END_REF], [START_REF] Mokni | A formal approach for managing component-based architecture evolution[END_REF] is a three-level architecture description language (ADL) designed to give a representation of the entire life cycle of architectures and a support to manage their evolution. Design is represented by the Specification level which is composed of abstract component types. Those types are called roles which means that they define the functionalities present in the components of the future software. Implementation choices are captured by the Configuration level. This architecture level is composed of concrete component classes which are realizations of DOI reference number: 10.18293/SEKE2018-140 the roles. Deployment is described in the Assembly level. This level is composed of a set of component instances that define how to tailor software for specific execution contexts.

B. Spring

Java Spring framework [START_REF] Johnson | The spring framework -reference documentation[END_REF] is widely used in industry to manage the deployment of software architectures. It provides standardized architecture deployment capabilities thanks to a container that is able to handle explicit architecture descriptors. Architecture descriptors are defined as XML files or directly embedded in the code as annotations. They are based on the concept of beans, which define the objects that the container must instantiate and connect in order to set the initial architecture up. For instance, the deployment descriptor of Figure 1a is composed of four beans: lampDesk, lampSitting, clock1 and orchestrator1.

Connections between beans are handled by the container using dependency injection to preserve decoupling. Beans only declare reference attributes corresponding to their dependencies with other beans. These dependencies are then resolved at runtime thanks to the connections defined in the deployment descriptor. For instance, the deployement descriptor in Figure 1a defines three connections between the orchestrator, the lamps and the clock beans. Those connections are defined by the property tag which corresponds to the injection of dependencies. For instance, < property name = "clock" ref = "clock1"/ > sets the clock1 bean as the clock property of the orchestrator1 bean.

As compared to raw code, Spring projects provide some explicit architectural descriptions. However, these descriptions do not capture design decisions and thus cannot be considered as abstractions of the software architecture.

III. Extracting Component-Based Software Architectures

In order to ease the extraction of information from the Spring descriptor, using a model-driven approach, a small domain specific language named SpringDSL has been developed. It consists of an implementation of the Spring XML descriptor grammar in EMF 1 . XText 2 has been used to automatically generate the corresponding EMF metamodel since we could not find an already developed metamodel. This metamodel enables thus to parse Spring XML descriptors and get all their content as concept instances.

This section discusses how each concept of Dedal is extracted from the source code and a Spring deployment description.

A. Extracting Components

Considering component extraction as a model transformation between SpringDSL and Dedal, only a simple mapping is required to extract the assembly and also a small part of the configuration models. To do so, the eclipse QVTo 3 language is used since it defines model to <bean class="AdjustableLamp" id="lampDesk" /> <bean class="AdjustableLamp" id="lampSitting" /> <bean class="Clock" id="clock1" /> <bean class="Orchestrator" id="orchestrator1" > <property name="lamps"> <set><ref bean="lampDesk" /> <ref bean="lampSitting" /></set></property> <property name="clock" ref="clock1" /> </bean> Figures 1b and1c present an extract from the Java code of the Home Automation Software (HAS) example. Code introspection enables to extract complementary information required to build higher level architecture models and more detailed component definitions.

For generating the component roles of the architecture Specification, the type hierarchy of the beans classes is analyzed, in order to extract the most generic, thus reusable, architecture model as possible. The main idea is to retrieve the abstract superclasses that are realized by the bean class corresponding to a component class. To extract the component role, the type hierarchy is traversed and the role which is picked is the most generic component role which still holds all the required interfaces that are present in the corresponding component class and which preserves the connections which exist in the Configuration. Figure 2 is the component role HomeOrchestrator role that is realized by the primitive component class Orchestrator.

B. Extracting Interfaces

Two types of interfaces are distinguished: (i) provided and (ii) required interfaces. All the methods that are provided by the beans classes must be provided into respective component interfaces. However, in order to not provide only one large interface per component, the interfaces are cut according the type hierarchy of classes. In other words, each implemented interface is mapped as a component interface and if a class does not implement an interface, a "conceptual" interface is extracted, which is composed of the public methods of the beans class, except for getters / setters that are used whether to initialize properties or to manage connections. For extracting required interfaces, The HomeOrchestrator role extracted interfaces are described in Figure 2 with their names, direction and implemented types.

C. Extracting Connections

The XML description (Figure 1a) makes it possible to start the extraction of connections between components. Indeed, thanks to the dependency injection, clients and servers of connections are identifiable. For instance, in the current example, the orchestrator1 bean contains a property which refers to clock1, so it is possible to map a new connection between orchestrator1 and clock1. This connection is propagated to the Configuration level by creating a connection between the two instantiated component classes. Following the same principle, the connection between the realized component roles is created.

Then the interfaces that are implied in connections must be set. To do so, we search among two connected components which are their matching interfaces. For matching two interfaces, their types must be equal and their direction complementary (provided with required). Thus, for instance in our case, Clock provides the IT ime interface of type T ime and, Orchestrator requires an interface of the same type (Figure 3). Then those two interfaces match and the connection presented in Figure 4 can be set. Then it is propagated following the instantiate relation between component class interfaces and component instance interface (Figure 3). Finally, connection between roles are set in the same way as the connections between component classes.

Figure 3 is the visual representation of the three-level Dedal architecture which is composed of four component instances into the Assembly (that correspond to the beans of the Spring description), the component classes of the Configuration that are instantiated by the component instances and also the component roles into the Specification which are realized by component classes. The connections between components are also represented.

IV. Related Work

This section narrows the studied approaches to the ones which extract component-based architecture descriptions and, if possible, from object-oriented code. Moreover, retro-engineering approaches which consist in simply abstracting the software artifacts for retrieving raw design decisions are differentiated from re-engineering ones which intend to re-organize the extracted information and / or the software artifacts.

In their work, Ducasse et al. [START_REF] Ducasse | Software architecture reconstruction: A process-oriented taxonomy[END_REF] defined a taxonomy for categorizing software architecture reconstruction approaches. Following this taxonomy, the goals of the discussed approach are twofold. The first goal is to improve component reuse, by extracting component-based architecture descriptions, such as MAP [START_REF] Stoermer | Map-mining architectures for product line evaluations[END_REF], PuLSE / SAVE [START_REF] Knodel | Static evaluation of software architectures[END_REF] and ROMANTIC [START_REF] Alshara | Materializing architecture recovered from objectoriented source code in component-based languages[END_REF], [START_REF] Shatnawi | Reverse engineering reusable software components from object-oriented APIs[END_REF] approaches, but targeting the Dedal [START_REF] Zhang | Architecture-centric component-based development needs a three-level ADL[END_REF], [START_REF] Mokni | A formal approach for managing component-based architecture evolution[END_REF] ADL. The second goal is to provide the foundations for managing conformance checking (Bauhaus [START_REF] Eisenbarth | Locating features in source code[END_REF], [START_REF] Christl | Equipping the reflexion method with automated clustering[END_REF], DiscoTect [START_REF] Yan | DiscoTect: a system for discovering architectures from running systems[END_REF], PuLSE / SAVE [START_REF] Knodel | Static evaluation of software architectures[END_REF]), evolution, co-evolution (PuLSE / SAVE [START_REF] Knodel | Static evaluation of software architectures[END_REF], Huang et al. [START_REF] Huang | Runtime recovery and manipulation of software architecture of component-based systems[END_REF]) and maintenance using the formal rules that have previously [START_REF] Mokni | A formal approach for managing component-based architecture evolution[END_REF] been defined in Dedal.

However none of the studied methods intends to extract raw information of how the software is implemented. Moreover, all the discussed approaches only deal with two levels of abstraction (i.e., implementation and architecture) that may not correspond to the same paradigms (code vs component-based architecture description). Indeed three component-based architecture descriptions are essential for maintaining, evolving, tracking software life-cycle since it gives a more global and direct understanding to the architect which can get an overview of the code structure by managing components. Moreover, even the approaches which seem to fit with the discussed one, either recover architecture in a semi automatic manner from execution trace of software (i.e., DiscoTect [START_REF] Yan | DiscoTect: a system for discovering architectures from running systems[END_REF]) or do not reconstruct raw architecture such as ROMANTIC [START_REF] Alshara | Materializing architecture recovered from objectoriented source code in component-based languages[END_REF], [START_REF] Shatnawi | Reverse engineering reusable software components from object-oriented APIs[END_REF] approach which performs re-engineering of the deployed architecture by clustering classes into bigger semantic components that encapsulate classes. This is why a retro-engineering approach is proposed that builds three-level component-based architecture description from structural artifacts.

V. Conclusion and Future Work

This paper introduces an approach for software architecture reconstruction, using three levels of architecture models (Assembly, Configuration, Specification). An aspect of future work will be to improve and refine the extraction of the Specification for making it more abstract.

Real Spring projects have already been identified in open-source repositories in order to perform large scale experimentations on evolution and reuse. Getting projects from open-source repositories will also allow the implementation of versioning mechanisms.

  (a) HAS Spring description public class Clock implements Time { public void run() {...} } public class AdjustableLamp extends Light { ... public void switchState(State s) {...} public void adjustInensity(int value) {...} } public class Orchestrator extends HomeOrchestratorImpl { protected void run(){...} } (b) HAS classes instantiated in Figure 1a public interface Time { public void run(); } public abstract class Light { ... public abstract void switchState(State s); } public abstract class HomeOrchestrator { Light[ ] lights; Time clock; public abstract void run(); } public class HomeOrchestratorImpl extends HomeOrchestrator { protected void run(){...} public HomeOrchestratorImpl (Light[ ] lights, Time clock) {...} public Light[ ] getLights() {...} public void setLights(Light[ ] lights) {...} public Time getClock() {...} public void setClock(Time clock) {...} } (c) HAS most abstract classes

Figure 1 :

 1 Figure 1: HAS implementation

Figure 2 :

 2 Figure 2: Extracted component class: Orchestrator

Figure 3 :Figure 4 :

 34 Figure 3: HAS generated Dedal three-level architecture