domain: sdv info

text:
M Mcevoy

M O'donovan

J P Murphy

Luc Delaby

Cultivar influences milk production of grazing dairy cows

come   L'archive ouverte pluridisciplinaire

Introduction

Perennial ryegrass is considered one of the most important forage grass species used in temperate ruminant animal production. In such regions, milk production can be achieved in an economic and cost effective manner by maximizing the proportion of fresh grass in the diet of the animal. Grass based ruminant production systems tend to target calving or lambing to coincide with the onset of grass growth and thus, maximize animal production from grazed grass. Plant breeding has largely focused on increasing annual DM yield. While this is an important trait, it fails to define the ultimate value of a grass cultivar, which is its potential to improve animal performance (Reed, 1994). The objective of this study was to quantify the difference in milk performance of lactating dairy cows offered 4 different perennial ryegrass cultivars.

Material and Methods

Four cultivars of perennial ryegrass were sown as monocultures in a split plot design. The cultivars included 2 tetraploids: Bealey (heading date 24 th May) and Astonenergy (31 st May) and 2 diploids: Spelga (22 nd May) and Abermagic 28 th May). Cultivars were sown in 2009. In March 2011, 24 cows were blocked and randomly assigned to one of four groups. Each group was then assigned to one of the four cultivars in a 4×4 latin square design. Two latin squares were completed between March to end of July. The first latin square (spring period) began in mid-March and ran for 56 d, divided into 14 day periods. Animals were offered 16 kg DM herbage and 1 kg DM concentrate cow -1 d -1 . In May, a further 4 cows were blocked and assigned randomly to one of the four groups and the second 4×4 latin square was completed over 84 days, divided into 21 day periods (summer period). Concentrate was removed and animals were offered 17 kg herbage DM cow -1 d -1 . Area was adjusted daily to ensure correct allocation of herbage. The final 5 days of each period was used for measurement, with the preceding days used to adapt cows to the cultivars. Milk yield was recorded daily, milk composition was determined from a morning and evening sample for 3 consecutive days. Data were analysed using Proc Mixed in SAS. The model included latin square, period and cultivar. Cow was included as a random effect.

Results and Discussion

Pre-grazing sward height was similar across cultivars (9.0 cm). Table 1 presents the milk performance results of the study. There was no interaction between cultivar and season for any parameter, however there was a tendency towards an interaction for protein composition (P=0.07). Cultivar had a significant effect (P<0.001) on milk yield, Bealey and Astonenergy had higher milk yields (28.7 kg cow -1 d -1 ), than Spelga and Abermagic (-1.4 kg cow -1 d -1 ). Cultivar did effect milk fat (P<0.05) and protein (P<0.001) content. Milk fat was higher for Bealey and Abermagic (43.4 g/kg) than for Astonenergy and Spelga (41.7 g/kg). Daily milk solids yield was higher for Bealey and Astonenergy (2.17 kg cow -1 d -1 ), than for Spelga and Abermagic. There was a significant effect of season on milk yield and composition (P<0.001). Milk yield, fat, lactose content and milk solids was higher in spring than summer, as a result of stage of lactation of the animal, and milk protein was higher in the summer than spring due to the quality of the herbage on offer.

Conclusions

Results from this study show that choice of ryegrass cultivar can influence the milk performance of spring calving dairy cows. The differences observed would go unnoticed if these cultivars were only evaluated in cutting plots. Therefore, these results support the need for testing grass cultivars under animal production trials. Such trials would be beneficial to farmers as it would present them with information on potential improvements in production which could be achieved when selecting cultivars for reseeding.

Table 1 .

 1 Differences in milk performance of grazing dairy cows offered perennial ryegrass monocultures in the spring

	and summer period
								P-value
		Bealey	Astonenergy	Spelga	Abermagic	SE	Cultivar	Season
	Milk yield (kg/d)	28.5 a	28.8 a	27.2 b	27.3 b	0.47	0.001	0.001
	Milk fat (g/kg)	43.3 a	41.7 b	41.6 b	43.4 a	0.82	0.05	0.001
	Milk protein (g/kg)	33.4 a	33.4 a	32.8 b	33.6 a	0.33	0.001	0.001
	Milk lactose (g/kg)	46.7	46.6	46.9	46.7	0.21	NS	0.001
	Milk solids (kg/d)	2.18 a	2.15 a	2.02 b	2.08 b	0.04	0.001	0.001
	SE=standard error

Acknowledgements

The authors wish to acknowledge Moorepark farm staff.

Reed, K.F. M. 1994. New Zealand J Agric Res. 37: 227-286.


domain: math info

text:
Jimmy Vandel
email: jvandel@toulouse.inra.fr

Brigitte Mangin
email: mangin@toulouse.inra.fr

Simon De Givry
email: degivry@toulouse.inra.fr

New Local Move Operators for Learning the Structure of Bayesian Networks

We propose new local move operators incorporated into a score-based stochastic greedy search algorithm to efficiently escape from local optima in the search space of directed acyclic graphs. We extend the classical set of arc addition, arc deletion, and arc reversal operators with a new operator replacing or swapping one parent to another for a given node, i.e. combining two elementary operations (arc addition and arc deletion) in one move. The old and new operators are further extended by doing more operations in a move in order to overcome the acyclicity constraint of Bayesian networks. These extra operations are temporally performed in the space of directed cyclic graphs. At the end acyclicity is restored and newly defined operators actually lead to a gain in graph score. Our experimental results on standard Bayesian networks and challenging gene regulatory networks show large BDeu score and recall value improvements compared to state-of-the-art structure learning algorithms when the sample size is small.

Introduction

Learning the structure of Bayesian networks from fully observed data is known to be an NP-hard problem [START_REF] Chickering | Learning bayesian networks is NPcomplete[END_REF] which has received a lot of attention from researchers during the last two decades [START_REF] Daly | Learning Bayesian networks: approaches and issues[END_REF]. Due to its difficulty, heuristic methods have been widely used to learn Bayesian network structures. Two main approaches have been studied: constraint-based methods and score-based methods [START_REF] Koller | Probabilistic Graphical Models: Principles and Techniques[END_REF]. Constraint-based methods aim at satisfying as much independence present in the data as possible using conditional independence tests. Unfortunately, these methods can be sensitive to failures in individual independence tests. Score-based methods use a scoring function f to score a network structure with respect to the data. They score the whole network at once, being therefore less sensitive to individual failures. Different Bayesian and non-Bayesian scoring metrics can be used such as the Bayesian Information Criterion (BIC) [START_REF] Schwarz | Estimating the dimension of a model[END_REF] or the Bayesian Dirichlet criterion (BDeu) [START_REF] Heckerman | Learning Bayesian Networks: The Combination of Knowledge and Statistical Data[END_REF]. Score-based methods explore the space of network structures to find the highest-scoring network. This space being superexponential in the number of variables, local search methods are used such as greedy ascent search (also called hill-climbing), tabu search, simulated-annealing, and other complex metaheuristics like ant colony optimization [START_REF] Daly | Learning Bayesian Network Equivalence Classes with Ant Colony Optimization[END_REF]. In spite of its simplicity, the (repeated randomized or stochastic) greedy search method reveals to be a competitive method compared to more complex algorithms [START_REF] Gámez | Learning Bayesian networks by hill climbing: efficient methods based on progressive restriction of the neighborhood[END_REF]. Starting from an initial network structure, it performs a series of local moves until a local optimum is found. Each move selects and applies the best elementary operation(s) on the current network structure. The set of candidate neighboring structures is called the neighborhood in the sequel. A classical neighborhood is composed of single arc additions, deletions, and reversals. Using larger neighborhoods efficiently allows to find better local optima, and thus better network structures. [START_REF] Moore | Optimal reinsertion: A new search operator for accelerated and more accurate bayesian network structure learning[END_REF] proposed an optimal reinsertion of a target node by removing all its edges and reinserting it optimally. However this approach is limited to small problems only. [START_REF] Holland | Competing fusion for bayesian applications[END_REF] used a restricted form of look ahead called LAGD, combining several operations in a single move. In this paper, we follow a similar approach by focusing our local operations on a target node and combining several operations guided by the global acyclicity constraint of Bayesian networks. By doing so, we are able to exploit large neighborhoods efficiently. Other approaches use a compact representation of a set of network structures. They explore either the search space of variable orderings (an optimal structure compatible with the order being easier to find) [START_REF] Cooper | A Bayesian method for the induction of probabilistic networks from data[END_REF][START_REF] Teyssier | Ordering-based Search: A Simple and Effective Algorithm for Learning Bayesian Networks[END_REF][START_REF] Schmidt | Learning graphical model structure using L1-regularization paths[END_REF][START_REF] Alonso-Barba | Structural learning of bayesian networks using local algorithms based on the space of orderings[END_REF], or the search space of Markov-equivalent network classes like Greedy Equivalence Search (GES) [START_REF] Chickering | Optimal structure identification with greedy search[END_REF][START_REF] Nielsen | On Local Optima in Learning Bayesian Networks[END_REF][START_REF] Daly | Learning Bayesian Network Equivalence Classes with Ant Colony Optimization[END_REF].

In Section 2, we give Bayesian network background. Next, we define a specific stochastic greedy search algorithm and introduce the new local move operators in Section 3. We report experimental results in Section 4 and conclude.

2 Bayesian network structure learning A Bayesian network [START_REF] Koller | Probabilistic Graphical Models: Principles and Techniques[END_REF] denoted by B = (G, P G ) is composed of a directed acyclic graph (DAG) G = (X, E) with nodes representing p random discrete variables X = {X 1 , . . . , X p }, linked by a set of directed edges or arcs E1 , and a set of conditional probability distributions P G = {P 1 , . . . , P p } defined by the topology of the graph:

P i = (X i |Pa(X i )) where Pa(X i ) = {X j ∈ X | (X j → X i ) ∈ E} is the set of parent nodes of X i in G. A Bayesian network B represents a joint probability distribution on X such that: (X) = p i=1 (X i |Pa(X i ))
. Conditional probability distributions P G are determined by a set of parameters. Given the network structure G, and the fully observed data D, parameters can be estimated by simple counting, following the maximum likelihood principle.

Learning the structure of a Bayesian network consists in finding a DAG G maximizing the posterior distribution (G|D). We have (G|D) ∝ (D|G) (G) since (D) is independent of G. Under specific assumptions, the marginal loglikelihood log( (D|G)) can be expressed as a consistent decomposable scoring function f , such as the BIC and BDeu criteria [START_REF] Chickering | Optimal structure identification with greedy search[END_REF] :

f (D, G) = p i=1 f X i (D, G) = p i=1 f X i (D, Pa(X i )) (1)
As f is consistent, maximizing f when the sample size tends to infinity leads to select the true structure or one among the Markovequivalent structures set. A set of Bayesian networks are Markovequivalent if they imply exactly the same set or map of independence constraints among variables 2 . Next, we describe a novel greedy search method maximizing f in the space of DAGs.

Stochastic Greedy Search

We define the Stochastic Greedy Search (SGS) algorithm 3 for structural learning of Bayesian networks in Algorithm 1. It collects the best DAG found by r randomized hill climbing algorithms. Stochasticity comes from two random draws. The first one, common in the structure learning community, is due to InitGraph that returns a random DAG used by the inner greedy search loop. The second is more original. It is the random DAG drawn among the best neighbors in the neighborhood of the current DAG G ( SelectRandom at line 1 of Algorithm 1). The neighborhood of G, returned by Neighbors, is composed of the usual elementary operations on DAGs: arc addition (ADD), arc deletion (DELETE), and arc reversal (REVERSE). This classical neighborhood is denoted N ADR in the sequel. Only feasible operations are considered, which do not create cycles. In the next subsections, we are going to extend this set of operations. We empirically observe that the first random draw may be counter-productive and starting from an empty graph is often better than from a random DAG, as also observed in [START_REF] Melanc ¸on | Random generation of directed acyclic graphs[END_REF].

G * ← G ; s * ← s ; return G * ;
Proposition 1. Let D be a dataset of n records that are identically and independently sampled from some distribution (•). Let f be a locally consistent scoring function. The inner loop of the SGS algorithm returns a minimal independence map of (•) as the sample size n grows large.

The local consistency of f ensures that adding any arc that eliminates an independence constraint that does not hold in the generative distribution (•) increases the score. Conversely, deleting any arc that results in a new independence constraint that holds in (•) also increases the score. Hence, starting from any DAG, by selecting strictly improving moves, the algorithm terminates (because of a finite number of DAGs) and returns an independence map of (•) (I-map of (•): every independence implied by the resulting DAG is verified in (•), in the worst case it can be a complete graph) which is minimal thanks to arc deletion operations and local consistency.

The main interest of our randomization approach is to simulate a search in the space of score-equivalent networks. Each greedy search moves from a randomly-selected DAG instance of a Markov-equivalence class E(G) to another randomly-selected DAG of an adjacent4 Markov-equivalence class E(G ) thanks to our SelectRandom function. It results in a stronger property: Proposition 2. Let D be a dataset of n iid fully observed samples of some faithful distribution (•). Let f be a locally consistent scoring function. SGS returns a perfect map of (•) as both the sample size n and the number of restarts r grow large.

Recall that a faithful distribution admits a unique perfect map corresponding to the optimal structure. Compared to the GES algorithm [START_REF] Chickering | Optimal structure identification with greedy search[END_REF], which offers the same optimality guarantee within a two-phase greedy search, SGS chooses the orientation of some compelled arcs5 of the true DAG at random, whereas GES waits while no v-structures impose orientation constraints. See an example in Figure 1.

Notice that neither GES nor SGS find an optimal structure in polynomial time in the worst case, even when using a constant time consistent scoring function and a faithful distribution 6 . In general, without the faithfulness assumption, learning the optimal structure is NPhard even when the sample size is large and when each node has at most k parents, for all k ≥ 3 [START_REF] Chickering | Large-Sample Learning of Bayesian Networks is NP-Hard[END_REF].

We observe in the experiments that a small number of restarts r allows to find DAGs with better scores than GES, especially when the sample size n is limited, in this case GES found a local optimum and SGS is able to find other better local optima thanks to randomization. This was also observed in [START_REF] Nielsen | On Local Optima in Learning Bayesian Networks[END_REF].

When the sample size is small the learning problem becomes more difficult: the empirical distribution may be far from a perfect map resulting in many local optima and the scoring function is no more consistent, i.e. the likelihood does not dominate the penalty term of the scoring function which increases with the parent variable domain sizes [START_REF] Chickering | Optimal structure identification with greedy search[END_REF]. In this complex situation, we propose a new operator to escape from some local optima.

SWAP operator

Consider the 3-variable example in Figure 2 with observed data D, scoring function f , and initial DAG

G 0 = {X2 → X3}. Let assume f (D, {X1 → X3}) > f (D, {X2 → X3}) > f (D, {X1 → X3, X2 →
X3}) > f (D, {X3 → X1, X2 → X3}) > f (D, {X2 → X1, X2 → X3}) > f (D, ∅).
Then G 0 is a local minimum for the classical neighborhood N ADR . Our new operator, denoted SWAP(X|Y → Z), consists in changing one parent X to another parent Y for one target node Z. This is equivalent to a simultaneous pair of ADD and DELETE operators restricted to the same target node. In our example, applying SWAP(X2|X1 → X3) corresponds to DELETE(X2 → X3),ADD(X1 → X3), resulting in the better DAG G 1 = {X1 → X3} as shown in Figure 2. The extended neighborhood using the four operators is denoted N ADRS and SGS using N ADRS (respectively N ADR ) is denoted SGS 2 (Stochastic Greedy Search with Swap) (resp. SGS 1 ) in the sequel.

A typical suboptimality problem that we observed in our experiments happens when two nodes have the same parents. Figure 3 shows such an example with four variables. Because child nodes X3 and X4 are highly correlated due to their common parents X1 and X2, the first arc usually added is either X3 → X4 or X4 → X3 especially if the conditional probability distributions of X3 and X4 given X1 and X2 are close. Then adding a first parent to X3 and furthermore a second parent to X4 and X3 give the DAG G 4 . Here deleting X3 → X4 is going to decrease the score and the same negative effect results when adding X1 → X4 because it leads to a three-parent node. For this node, the increase in likelihood may not overcome the penalization term. But doing both operations simultaneously thanks to our SWAP operator results in a better local optimum DAG.

Let p be the number of variables in the current DAG and k be the maximum number of parents per node. Assuming a sparse graph, p k, the number of SWAP operations is bounded by O(kp 2 ) in the extended neighborhood N ADRS , whereas it is bounded by O(p 2 ) for ADD and O(kp) for DELETE and REVERSE. The complexity of N ADRS is therefore in O(kp 2 ), whereas it is in O(p 2 ) for the classical neighborhood N ADR . Notice that other approaches using larger neighborhoods such as h-look ahead in l good directions (LAGD) has a worst-case complexity in O(l h-1 p 2 ) [START_REF] Holland | Competing fusion for bayesian applications[END_REF] and optimal reinsertion (OR) neighborhood is in O(2 k p k+1 ) [START_REF] Moore | Optimal reinsertion: A new search operator for accelerated and more accurate bayesian network structure learning[END_REF]. In particular, LAGD complexity does not benefit from sparsity, l h-1 being constant, whereas our approach is faster when k decreases. Moreover computing the score difference of two DAGs before and after a SWAP operation is easy to do as it remains local to a single target node thanks to score decomposition.

Another source of suboptimality comes from the global acyclicity constraint of Bayesian networks. to restore graph acyclicity (G becomes valid) and such that the true score of the final DAG is greater than the score of the original one.

In Figure 4, it is obtained by applying a second SWAP.

; if ∆ + ∆ W * (G , DELETE(U * → W * )) > 0 then L ← L ∪ {DELETE(U * → W * )} ; ∆ ← ∆ + ∆ W * (G , DELETE(U * → W * )) ; Apply DELETE(U * → W * ) to G ; else /* Choose the best swap to get a positive change */ ; (U * |V * → W * ) ← 7 argmax (U→W)∈C,V∈X' ∆ W (G , SWAP(U|V → W)) ; ∆ ← ∆ + ∆ W * (G , SWAP(U * |V * → W * )) ; if ∆ > 0 then L ← L ∪ {SWAP(U * |V * → W * )} ; Apply SWAP(U * |V * → W * ) to G ; else L ← ∅ /* Abort all local operations */ ; 8 return L ;
For that purpose, we define an extended SWAP operator, denoted SWAP * , able to break all directed cycles by performing a succession of deletion or swap operations. It can be seen as a kind of greedy descent search in the space of directed cyclic graphs, trying to remove the less important arcs or to swap them in order to compensate for their loss, until a better valid DAG is found. We use local score changes to guide the search: Algorithm 2 terminates because there are at most O(pk) directed cycles to break, assuming p nodes and k maximum number of parents per node. Arcs newly created by swap operations cannot be swapped again more than p times thanks to our restricted list of alternative candidate parent nodes X' used at line 7, initialized to all the variables at the beginning, and updated by the list of nodes present in the current cycle at line 5.

∆ Xi (G, OP) = f Xi (D, G ) - f Xi (D, G),
We apply the same approach to extend the other operators ADD * and REVERSE * , breaking cycles with deletions or swaps. Because REVERSE * (X → Y) = ADD * (Y → X), we use only ADD * . The resulting neighborhood exploiting these extended operators is denoted N A * DS * and SGS using this neighborhood is denoted SGS 3 (Stochastic Greedy Search with Successive Swaps) in the experiments.

Experimental Results

In this section, we describe a set of experiments aimed at testing the performance of SGS i algorithms compared with state-ofthe-art Bayesian network structure learning algorithms on standard Bayesian networks and challenging gene regulatory networks.

Results on Standard Bayesian Networks

We used four gold-standard networks from the Bayesian Network Repository8 whose main properties are shown in Table 1. In this table, Nodes and Arcs specify the number of nodes and arcs respectively in the DAG. The Max in-degree is the maximum number of parents of a node. The Min-max states are the minimum and maximum number of states in the domain of a variable associated to a node. Finally, Longest path is the length of the longest directed path between any pair of nodes in the DAG. 100 samples of size n = 500 and n = 5000 were generated for each network using Causal Explorer [START_REF] Aliferis | Causal Explorer: A Probabilistic Network Learning Toolkit for Biomedical Discovery[END_REF]. We compared SGS i algorithms with LAGD [START_REF] Holland | Competing fusion for bayesian applications[END_REF], available in the WEKA software [START_REF] Hall | The WEKA Data Mining Software[END_REF] and GES [START_REF] Chickering | Optimal structure identification with greedy search[END_REF] implemented in Tetrad 4.4.0 [START_REF] Scheines | The TETRAD Project: Constraint Based Aids to Causal Model Specification[END_REF]. LAGD was shown to outperform repeated hill-climbing, order-based search, tabu search, and simulated annealing in [START_REF] Salehi | An empirical comparison of the efficiency of several local search heuristics algorithms for Bayesian network structure learning[END_REF]. GES, which is considered the reference algorithm for Bayesian network structure learning, was reported having similar performance to the most recent order-based search algorithms in [START_REF] Alonso-Barba | Structural learning of bayesian networks using local algorithms based on the space of orderings[END_REF]. Recall that SGS 1 is similar to repeated hill-climbing, SGS 2 uses the SWAP operator, and SGS 3 breaks cycles by successive DELETE and SWAP operators. Experiments were performed on a 3 GHz Intel Core2 computer with 4GB running Linux 2.6.

We fixed the maximum number of parents per node at k = 5 for SGS i and LAGD. LAGD exploits a h = 2-look ahead in l = 5 maximum improving directions. GES was restricted on the number of adjacent nodes: d = 7 for Hailfinder and d = 10 for Pigs network as done in [START_REF] Alonso-Barba | Structural learning of bayesian networks using local algorithms based on the space of orderings[END_REF]. All methods were initialized by an empty graph and optimized the BDeu score with equivalent sample size α = 1 and no prior on the network structures. For each sample, we recorded the best score obtained by GES, and by r = 10 randomized greedy searches for SGS i (Algorithm 1 at line 2) as for LAGD 9 .

In order to test the statistical significance of the difference in BDeu score between two methods, we applied a non-parametric paired test, the Wilcoxon signed-rank test [START_REF] Wilcoxon | Individual Comparisons by Ranking Methods[END_REF]. Table 2 presents the test results for all the pairs of methods by using an unilateral alternative (no difference versus better), the pairwise type I error is fixed to 6.2510 -4 which corresponds to a familywise error rate of 5% with Bonferroni correction for multiple comparisons. The above results were summarized in Table 3, in a summarized Wilcoxon score for each method, obtained by summing the positive comparisons and subtracting the negative ones. The non significant comparisons were ignored.

SGS 3 was the best method for the four networks, except for Pigs network with n = 5000 which is more accurately estimated by GES. We conjecture that in this case, GES was closed to its asymptotic optimal behavior, which may be due to the smallest network features of Pigs network among all. SGS 2 improved over SGS 1 and reached the second position especially for small sample sizes, probably for reasons commented in Figure 3. LAGD got poor results, the difference with SGS 1 can be explained by a better randomization in SGS 1 (Algorithm 1 at line 1). LAGD failed on the Pigs network due to the large number of variables p = 441 that makes the exploration of 2look ahead neighborhoods infeasible in a reasonable time. GES was the worst method of this evaluation (except for Pigs) due to limited sample sizes.

Although the algorithms are designed to maximize a (BDeu) score, we generally look for a network structure as close as possible to the true one. When the sample size grows large, improving the BDeu score leads to reduce the structural Hamming distance (SHD) between the learned and the true network. We report in Table 4 the means over 100 datasets (rounded values to the nearest integer) of the missing and spurious edges without taking into account the edge orientations. The SHD is the sum of the above values. We do not report SGS 1 and SGS 2 results that are outperformed by SGS 3 . SGS 3 (resp. GES) got the best SHD in 4 (resp. 5) configurations and outperformed LAGD (which performed as SGS 3 in 1 configuration). GES performed extremely well on the Pigs network, finding the true network with 5,000 samples, whereas SGS 3 learned too many edges but recovered all the true edges (even with n = 500). The spurious edges learned by SGS 3 are exclusively due to random orientations of compelled arcs in v-structures (see Figure 1). Assuming X1 → X3 ← X2 in the true network (v-structures are very frequent in the Pigs network) and a large sample size, if during its greedy search SGS 3 adds first X1 ← X3 and X3 → X2 then it will add next a covering edge X1 → X2 or X1 ← X2 in order to find a minimal independence map (see Proposition 1). We corrected this fault by comparing for each covered v-structure the 3 possible non-covered v-structures independently from the other variables and selecting the highest score configuration. After this post-processing step on Pigs network, we had only 1 (resp. 2) spurious edges for SGS 3 with 500 (resp. 5000) samples without increasing the missing edges number. The same postprocessing had no effect on the other smaller networks.

Detailed analysis on the Alarm network

We conducted a series of algorithm analyzes on the Alarm network with a fixed sample size n = 500. Table 5 shows a typical example of the BDeu score reached by the algorithms for a single greedy search (r = 1) on a particular sample. It also reports the number of local moves Iter. (Algorithm 1 at line 1), the number of local score f Xi (D, G) computations for different parent configurations, i.e. assuming a perfect BDeu cache, and the CPU time in seconds. As expected, the number of iterations for LAGD was less than half that of other methods, since LAGD applies two local operators (ADD, DELETE, or REVERSE) at each move. The swap operations used by SGS 2 and SGS 3 slightly increased the number of moves allowing to find better local optima. Comparing CPU times is a difficult task due to the various implementation details (language choice, caching strategy, sample compilation techniques, see for instance Chapter 18.4.3 in [18]). Nevertheless SGS i algorithms are comparable. In our implementation, SGS 3 was twice slower than SGS 1 . A fairer comparison is the number of calls to new local score computations. SGS 1 needed the fewest number, whereas GES needed the most as it tends to explore larger neighborhoods in the space of Markov-equivalent networks. Notice that SGS 2 and SGS 3 required less score computations than LAGD due to graph sparsity (see Section 3.1). The greater BDeu cache size of SGS 3 compared to SGS 2 is mainly due to the increased number of local moves.

Number of restart

BDeu score

Gold SGS 3 Empty SGS 3 Random SGS 2 Empty SGS 1 Empty GES SGS 2 Random SGS 1 Random Figure 5.
Best BDeu scores, averaged on 30 Alarm samples (n = 500), found by SGS i algorithms as the number of restarts r increases and starting either from an empty (solid line) or a random graph (dashed line). Results of GES (dotted line) and BDeu score of the true network Gold (dash-dotted line) are also given. Methods are sorted from the best, SGS 3 with an empty graph, to the worst, SGS 1 with a random initial graph, at r = 1.

We further analyzed the impact on performances of the number of restarts r and the initial graph used by SGS i algorithms (see InitGraph in Algorithm 1). Figure 5 reports averaged BDeu scores on 30 Alarm samples (the length of the 95% confidence interval of the BDeu score mean was around ±30). The 30 initial random graphs, common for SGS i algorithms, are composed of 71 arcs with at most two parents per node 10 . All SGS i methods reached a better BDeu score than GES in this small sample size situation by the use of less than r = 4 restarts. SGS i methods converged quite rapidly as r increases. Only SGS 3 found a better score than the true network. Initial random graphs were counter-productive for all the methods, except for SGS 3 . It shows that starting from a random graph is useful if the available operators allow to move efficiently in the search space, which is illustrated by the following experiment. On the contrary, using randomness within the greedy selection process (see SelectRandom in Algorithm 1) was always beneficial.

Finally, we analyzed how often our new local move operators are used during the search. Table 6 shows the mean number of local moves using a specific operator averaged on r = 50 greedy searches over 30 Alarm samples (n = 500). Labels ADD+ and SWAP+ mean that at least two local operators were successively applied at a given move in order to break cycles 11 . The new local move operators are mostly used when starting from a random graph (there are more cycles to break), but when starting from an empty graph, only a few applications of the extended moves allow to significantly improve the BDeu score as previously shown in Figure 5 and Table 3.

Results on Gene Regulatory Networks

Gene regulatory network reconstruction from gene expression data using Bayesian network structure learning was first proposed in [START_REF] Friedman | Using bayesian networks to analyse expression data[END_REF].

We used simulated expression datasets of the Systems Genetics Challenge A from the international reverse enginnering competition DREAM5 [START_REF]The DREAM5 Systems Genetics Challenges[END_REF]. Genetics data were not used as they require additional modelling to be taken into account, see e.g. [START_REF] Vignes | Gene Regulatory Network Reconstruction Using Bayesian Networks, the Dantzig Selector, the Lasso and Their Meta-Analysis[END_REF]. Expression data were generated using the SysGenSIM generator [START_REF] Pinna | Simulating system genetics data with SysGenSIM[END_REF] based on ordinary differential equation simulation. Five datasets are available for three different sample sizes (n = 100, 300, and 999). The 15 datasets were obtained from different known gene networks composed of 1, 000 variables and containing directed cycles. For each sample size, the five network structures contain a different number of edges varying from ≈ 2, 000 (Net1) to more than 5, 000 (Net5). We discretized gene expression levels into 2 to 4 states using an adaptive method based on an adapted k-means algorithm and the more general framework of Gaussian mixture models as described in [START_REF] Vignes | Gene Regulatory Network Reconstruction Using Bayesian Networks, the Dantzig Selector, the Lasso and Their Meta-Analysis[END_REF].

With such large networks, we had to adapt the learning procedure of SGS i algorithms 12 . We decided to restrict their lists of candidate parents as done in [START_REF] Goldenberg | Tractable learning of large Bayes net structures from sparse data[END_REF]: we selected for each variable X the set of parents S such that each element Y of S improves the local BDeu score when it is considered as a unique parent compared to the orphan situation ( f X (D, {Y → X}) > f X (D, ∅)). This filtering process was done prior to the search 13 . In these experiments, SGS i algorithms have a maximum number of parents per node fixed at k = 5 and use r = 10 restarts. Instead of LAGD (which was too slow), we used MMHC [START_REF] Tsamardinos | The max-min hill-climbing Bayesian network structure learning algorithm[END_REF] having two steps similar to SGS i . It first selects the skeleton of the graph using mutual information measures (MMPC [START_REF] Tsamardinos | Time and sample efficient discovery of Markov blankets and direct causal relations[END_REF] filtering step) and then orientates edges in a greedy manner. We recorded the best BDeu score of 10 runs for MMHC, by randomly permuting the input variables at each run. The skeleton being built locally for each variable, MMHC can handle large networks. All the methods started from an empty graph and optimized the BDeu score with α = 1 and no prior on the network structures.

It was difficult to perform statistics with the results of this experiment. Indeed, contrary to the standard network experiment, there were no replicated samples of the same network. We decided to pool results per sample size and we performed the Wilcoxon test on the groups. With only five results per group to compare the methods, we are aware that the power of the test is very low. So, we applied a pairwise type I error of 5% and we did not try to correct for multiple comparisons, see Table 7. However, it is worth noting that SGS 3 was always the best method and that it increased the BDeu score by about 2% in average. Surprisingly, GES appeared to be better on smaller sample sizes compared to MMHC. As explained below, MMHC was penalized by its filtering process, especially on the smallest sample size, whereas GES had no restrictions on the candidate parent sets.

In these experiments, the structural Hamming distance (SHD) was not informative due to the poor results reached by all tested algorithms for such large networks. SHD, equal to the sum of spurious edges (s) and missing edges (m), was greater than the number of edges (e) of the true graph. In this situation, even the empty structure appears better. For this reason, we computed another aggregated structural quality measure based on the precision ( e-m e-m+s ) and recall ( e-m e ) measures and took the Euclidean distance to the origin to combine them ( precision 2 + recall 2 ). Contrary to SHD, a high distance indicates a better structural quality. We observed in Table 8 contrasted performances between the tested methods depending on the sample size: for n=100, MMHC got the best results, for n = 300, it was GES, and finally SGS 3 performed the best for the largest sample size. Better BDeu scores are not always synonymous with a better structural quality, the limited sample size in addition to the non faithfulness of the data (generated by ordinary differential equations with cycles) could explain this behavior. We should notice that the good optimization behavior of SGS 3 in terms of BDeu scores resulted in a better structural quality than GES and MMHC as sample size increases. Moreover, we compared the quality of the filtering process used by MMPC and SGS 3 . Table 9 reports for the first network (Net1) of each sample size, the total number of arcs kept by the filtering process and the recall value, which represents the percentage of true edges in the filter. Our first observation is the poor recall of both filtering processes, which indicates strong structural differences between the true network and the learned network even with n = 999 sample size. Our filtering approach obtained better recall values with similar compression sizes than MMPC.

Finally, we tried the same methods as in Section 4.1 on 50 small gene regulatory networks (Web50 Artificial Gene Networks with 50 nodes and 50 arcs) without doing any filtering process. For a sample size n = 500, SGS 3 was significantly better than GES and LAGD in terms of BDeu scores and slightly better in terms of Euclidean distances to the origin of the (precision, recall) values.

Conclusion

We have presented in this paper a new greedy search algorithms called SGS i exploiting stochasticity from two random draws. We have developed a new local move operator called SWAP and extended versions for ADD and SWAP operators to overcome frequent limitations of local search methods which are local maxima and cyclic situations. We compared SGS 3 using SWAP and extended operators to state-of-the-art methods and we obtained significant BDeu and recall value improvements on classical benchmarks when the sample size is small. The complexity of SGS 3 stays moderate with sparse networks. In case of large networks with many variables we applied a filtering process in preprocessing using the same criterion as for the search. This process kept more edges from the true network than mutual information-based methods with significant reduction of the search space.

In the future, we would like to test our new operators with other local search methods like tabu search.

Figure 1 .Figure 2 .

 12 Figure 1. Four adjacent Markov-equivalence classes found by GES during its first phase of edge and v-structure insertions. (a) GES and SGS start from the empty graph . (d) The true DAG is found after three moves. The orientation of X3 → X4 and X1 → X3 edges are chosen at random by SGS, whereas GES waits until its third move to decide on edge orientations based on DAG score comparisons (enforcing the v-structure X1 → X3 ← X2 as stated by the extra ADD parameter {X1}, and forbidding X1 → X3 ← X4 in its second move).

3. 2 Figure 4 .

 24 Figure 4. Applying an extended SWAP * operation breaking a cycle by an additional SWAP operation: SWAP * (X2|X7 → X3) = {SWAP(X2|X7 → X3), SWAP(X4|X5 → X6)}.

Figure 3 .

 3 Figure 3. Problems with the classical neighborhood N ADR when two nodes (here X3 and X4) have the same parents. True DAG is G 6 . Starting from the empty graph (not shown), G 4 can be locally optimum for N ADR when the sample size is small (score of G 5 lower than G 4 ), whereas in N ADRS , SWAP(X3|X1 → X4) moves directly from G 4 to G 6 .

Algorithm 2 : 5 /*

 25 SWAP * (X|Y → Z) operator. Input : operation X|Y → Z, sample D, score f , DAG G(X, E) Output : a set of local operations L L ← ∅ /* Initialize output operations to the empty set */ ; X' ← X /* Candidate parent set for future swaps */ ; G ← G /* Copy of input DAG */ ; ∆ = ∆ Z (G , SWAP(X|Y → Z)) /* Putative score improvement */ ; 3 if ∆ > 0 then L ← L ∪ {SWAP(X|Y → Z)} ; Apply SWAP(X|Y → Z) to G ; /* Repeat deletion or swap operations until no more cycles */ while ∆ > 0 ∧ (C ←NextCycle(G )) ∅ do 4 X' ← X' \ nodes(C) ; Choose the best deletion to break cycle C */ ; (U * → W * ) ← 6 argmax (U→W)∈C\{Y→Z} ∆ W (G , DELETE(U → W)) ; /* Test if the sum of local score changes is positive */

  with G the result of applying the local move operator OP ∈ {DELET E, S WAP} to G. A negative sum of local changes aborts the search. Recall that finding a minimum number of arc dele-tions in order to restore acyclicity is NP-hard, see minimum feedback arc set problem in[START_REF] Festa | Feedback Set Problems[END_REF]. We use a greedy approach instead. The pseudo-code of SWAP * is given in Algorithm 2. The local score improvement of the initial SWAP operation is evaluated at line 3. It corresponds to a putative gain on the current score. If it is positive then this operation is applied to a copy of the input DAG G, checking next if it creates some directed cycles. Each cycle is retrieved by the NextCycle function 7 and the algorithm searches for an arc deletion in this cycle with minimum deterioration of the local score at line 6. If the combined local score change of the SWAP and DELETE operations is positive then it applies the selected arc deletion and continues to test if there are no more directed cycles at line 4. If the combined local score change is negative then it tries to swap an arc of the cycle such that the combined local score change is maximized (line 7) and positive. If it fails to find such an operation then it stops breaking cycles and returns an empty operation set. Finally if it succeeds, breaking all cycles, then it returns a feasible set of SWAP and DELETE operations resulting into a new valid DAG G with a better score than G. The true score improvement is equal to ∆.

Table 1 .

 1 Bayesian network properties

		Alarm Insurance Hailfinder Pigs
	Nodes	37	27	56	441
	Arcs	46	52	66	592
	Max in-degree	4	3	4	2
	Min-max states	2-4	2-5	2-11	3-3
	Longest path	11	10	14	6

Table 2 .

 2 Wilcoxon test comparing pairs of algorithms (familywise error rate = 5%). For Method1 versus Method2, + means that Method1 is significantly better than Method2,means that Method1 is significantly worse than Method2 and ∼ means there is no significant result

			Alarm	Insurance
	Sample size	500 5 000	500 5 000
	SGS 3 vs SGS 1	+	+	+	+
	SGS 3 vs SGS 2	+	+	+	+
	SGS 2 vs SGS 1	+	∼	∼	∼
	SGS 3 vs GES	+	+	+	+
	SGS 3 vs LAGD	+	+	+	+
	SGS 2 vs GES	+	∼	+	+
	SGS 2 vs LAGD	∼	∼	+	+
	SGS 1 vs GES	+	∼	+	+
	SGS 1 vs LAGD	∼	∼	+	+
	LAGD vs GES	+	∼	+	+
		Hailfinder		Pigs
	SGS 3 vs SGS 1	∼	+	+	+
	SGS 3 vs SGS 2	∼	+	+	+
	SGS 2 vs SGS 1	∼	∼	∼	∼
	SGS 3 vs GES	+	+	+	-
	SGS 3 vs LAGD	∼	+	n/a	n/a
	SGS 2 vs GES	+	+	∼	-
	SGS 2 vs LAGD	∼	+	n/a	n/a
	SGS 1 vs GES	+	+	∼	-
	SGS 1 vs LAGD	∼	+	n/a	n/a
	LAGD vs GES	+	+	n/a	n/a

Table 3 .

 3 Summarized Wilcoxon scores (the higher the score, the better the method)

		Alarm Insurance	Hailfinder	Pigs
	SGS 3	8	8	5	4
	SGS 2	0	2	2	-3
	SGS 1	-2	2	2	-3
	LAGD	-1	-4	-1	n/a
	GES	-5	-8	-8	2

Table 4 .

 4 Number of spurious edges (+) and missing edges (-) to sum for the structural Hamming distance. Both scores are in bold when giving the best distance per configuration

				Alarm	Insurance
	Sample size		500 5 000	500 5 000
	SGS 3	+ -	8 3	6 2	4 20	2 8
	LAGD	+ -	11 4	8 2	4 20	5 11
	GES	+ -	6 5	4 2	2 23	3 12
			Hailfinder		Pigs
	SGS 3	+ -	17 24	16 13	32 0	41 0
	LAGD	+ -	21 26	20 19	n/a n/a	n/a n/a
	GES	+ -	15 24	11 22	2 7	0 0

Table 5 .

 5 Single greedy search analysis on the Alarm network (n = 500, r = 1)

		BDeu score	Iter. BDeu cache Time(s)
	SGS 3	-5490.53	66	4543	3.2
	SGS 2	-5513.89	58	4310	2.3
	SGS 1	-5541.55	55	3305	1.5
	LAGD	-5544.61	35	4782	3.2
	GES	-5659.26	72	5531	2.4

Table 6 .

 6 Operator usage depending on the initial graph

		Empty graph init.	Random graph init.
		SGS 1 SGS 2 SGS 3 SGS 1 SGS 2 SGS 3
	ADD	53	53	52	53	9.5	11.5
	DELETE	0.5	0.1	1.4	64	20	43
	REVERSE	0.9	1.1	0.7	5.3	3.4	1.9
	SWAP	0	1.1	3.2	0	47	33
	ADD+	0	0	2	0	0	14
	SWAP+	0	0	1.6	0	0	9.5

Table 7 .

 7 Wilcoxon test (error rate = 5%) for different gene network sample sizes

		100 300 999
	SGS 3 vs MMHC	+	+	+
	SGS 3 vs GES	+	+	+
	MMHC vs GES	-	∼	+

Table 8 .

 8 Euclidean distances to the origin of the (precision, recall) values. Best distances per sample size are in bold

	n		SGS 3	MMHC	GES
		Net1	0.170	0.218	0.196
		Net2	0.213	0.295	0.232
	100	Net3	0.214	0.266	0.236
		Net4	0.201	0.265	0.214
		Net5	0.206	0.247	0.243
		Net1	0.510	0.483	0.464
		Net2	0.342	0.337	0.385
	300	Net3	0.484	0.488	0.505
		Net4	0.453	0.478	0.498
		Net5	0.419	0.397	0.428
		Net1	0.578	0.537	0.549
		Net2	0.581	0.510	0.505
	999	Net3	0.454	0.441	0.484
		Net4	0.476	0.450	0.476
		Net5	0.479	0.471	0.458

Table 9 .

 9 Total size and recall of candidate parent lists for SGS3 and MMPC on the most sparse gene network Net1

			100	300	999
	BDeu test	size recall	2670 2430 5984 9% 18% 35%
	MMPC	size recall	2568 3064 3842 5% 12% 23%

In the paper, we use G = E when the set of nodes is implicit.

BIC/BDeu give the same score for Markov-equivalent DAGs.

Don't mistake for the SGS algorithm in[START_REF] Spirtes | Causality from probability[END_REF] which is a constraint-based algorithm, but ours is score-based.

Two equivalence classes E(G), E(G ) are adjacent iff G is an I-map of G or vice-versa and the number of edges in the graphs G and G differs by one.

An arc X → Y in G is compelled if that arc exists in every DAG of E(G), otherwise it is said reversible.

The number of possible sets S in the GES operator ADD(E,S ) is exponential in the maximum degree d of the current graph in the worst case and r is unbounded for SGS.

We implemented an incremental acyclicity test[START_REF] Haeupler | Faster algorithms for incremental topological ordering[END_REF] which returns a shortest directed cycle using breadth first search.

http://www.cs.huji.ac.il/site//labs/compbio/Repository/

We randomly permute the input variables at each run.

For each node, we randomly select two parents and remove a parent if it creates a directed cycle.

We found at most 4 (resp. 5) successive operations for the 1500 runs starting from an empty (resp. random) graph.

GES managed in ∼1-hour CPU time each network thanks to its better implementation of caching and heap data structure.

[START_REF] Goldenberg | Tractable learning of large Bayes net structures from sparse data[END_REF] It could also be done during the search as in[START_REF] Gámez | Learning Bayesian networks by hill climbing: efficient methods based on progressive restriction of the neighborhood[END_REF].


domain: sdv sdv.bv sde

text:
Audrey Didier

François Balfourier

About the PGR networks in France: collaboration between users and the Genetic Resource Centre on small grains cereals

Keywords:

Main question

How to act, to better promoting genetic resources present in collections to answer both to socio-economic evolution of agriculture and to the developpment of basic scientific knowledge in genetic and genomic?

Big PGR collections

Poor evaluation data

How to share responsability between actors managing GR? How to manage the collection in order to be able to propose adapted samples to different users?

Forage Crops Grain legumes

Versailles



  Diffusion of knowledges and exchange/distribution of biological resources Definition of a policy of enrichment of the collections.

European Plant Genetic ResourcesConference 2011, Wageningen, The Netherlands, April 5-7, 2011

Aknowledgement • UFS cereal section -E. Margalé (Serasem) -F. Minard (RAGT) -J. Auzanneau (Agri-Obtentions) -P. Giraudeau (Secobra) -S. Dutriez (Caussade Semences) -S. Sunderwirth (Momont) -B. Duperrier (Limagrain-Verneuil) -V. Lein (Staaten Union) -J.M Delahaye (Lemaire Deffontaines) -S. Regnault (Unisigma) -P. Lonnet (Florimond Desprez) -S. Caiveau (Syngenta) • INRA GDEC Clermont-Fd -L. Bardy (seeds multiplication, distribution) -A. Bouguennec (Triticale sp.) -E. Boulat (network coord.) -J. Bordes (data analysis) -J. Koenig (Avena sp.) -A. Lapierre (phenotyping) -F. Exbrayat (genotyping) -N. Guilhot (bioinformatics) -C. Poncet (genotyping platform) -B. Debote (experimental farm) • INRA URGI Versailles -S. Durand (database)


Gilles Martel

F Herzog

O Huguenin-Elie

Diagram of analysis plan INFLUENCES OF LOCAL POLICIES AND OPPORTUNITIES ON FARMERS' STRATEGIES AND GRASSLAND MANAGEMENT

• Landscape complexity in mountains is enhanced by the diversity of grassland management • Natural factors (slope, soil…) do not explain all the diversity of management What is the role of local policies on farm strategy and of these strategies on grassland management ?

Material and Methods

• 3 areas in the same pedoclimatic context of the folder Jura: Doubs (F), Vaud (CH), Neuchâtel (CH)

• Multifactorial analysis followed by hierarchical clustering:

 For grassland management we used data of 634 plots on:

-Date of first and last utilization -Uses: number (0 to 4), type (grazing/mowing/depending on grass growth) and season (spring/summer/autumn) -Fertilization: type (liquid manure/manure/mineral or their combinations) and number (0 to 4)

Results

• 9 strategies identified: • 16 grassland usages and their different use between strategies:

• Interactions between strategy and area/grassland management have been studied by Chi-square test

•

  Search for type of grassland management at the farm level • Evaluate floristic diversity in plots of each management Outlook • Strategies 1,2 & 3: found in the three areas • Strategies 4 & 5 in Vaud: related to mountain pasture • Grazing milk strategy in Doubs: PDO cheeses specifications • Intensive strategy in Neuchâtel: subsidies and local dynamics • Minimizing unproductive times in Switzerland: market for heifers • Big farms overrepresentation in Doubs: agricultural dynamics  For strategy analysis we used data of 33 farms on:-number and kind of animals bred and sold; -area managed; -periods of production; -multiple variables on intensity of animal production (age at first calving, animal production expectations, quantity of concentrate…)

•

  In an area: both ubiquitous and specific strategies, due to:  Local dynamics of agriculture  Specifications of PDO cheeses • Diversity of grassland management :  Inside each strategy despite preferred/avoided uses  Enhanced by a combination of strategies

	Conclusions
	Local policies and opportunities
		Influence
	• Differences in local policies and opportunities:  Protected Designation of Origin (PDO) cheeses	Farmers' strategies
	 Specific policies on mountain pasture (fertilization, stocking rate)	Modify
	 Ecological policies	Grassland managements
	(green : more used ; red : less used than expected by chi-square statistics)
	• Each strategy has under & overrepresented grassland
	management i.e.:
	 Grazing milk farms (n°6) need plots for grazing early and avoid
	fertilization with liquid manure. They don't need to stock grass
	(few mowing)

• Test the strategies effect on floristic biodiversity at farm level • Discuss with the local stakeholders of the results