Alexis Breust

Christophe Chabot

Jean-Guillaume Dumas

Laurent Fousse

Pascal Giorgi
email: pascal.giorgi@lirmm.fr.

Recursive double-size fixed precision arithmetic

We propose a new fixed precision arithmetic package called RecInt. It uses a recursive double-size data-structure. Contrary to arbitrary precision packages like GMP, that create vectors of words on the heap, RecInt large integers are created on the stack. The space allocated for these integers is a power of two and arithmetic is performed modulo that power. Operations are thus easily implemented recursively by a divide and conquer strategy. Among those, we show that this packages is particularly well adapted to Newton-Raphson like iterations or Montgomery reduction. Recursivity is implemented via doubling algorithms on templated data-types. The idea is to extend machine word functionality to any power of two and to use template partial specialization to adapt the implemented routines to some specific sizes and thresholds. The main target precision is for cryptographic sizes, that is up to several tens of machine words. Preliminary experiments show that good performance can be attained when comparing to the state of art GMP library: it can be several order of magnitude faster when used with very few machine words. This package is now integrated within the Givaro C++ library and has been used for efficient exact linear algebra computations.

Introduction

Mathematical computations that needs integers above machine word precision are compelled to rely on a third party library. Among arbitrary precision libraries for integers, GMP [START_REF] Granlund | The GNU multiple precision arithmetic library[END_REF] (or its fork MPIRwww.mpir.org) is the most renown and efficient. The underlying structure of GMP/MPIR integers is based on an array of machine word integers that are accessed through a pointer. For instance, a 256-bits integers a can been represented by a dimension four array [a 0 , a 1 , a 2 , a 3 ] such that a = a 0 + a 1 2 64 + a 2 2 128 + a 1 2 192 . Additionally to the pointer, GMP/MPIR stores two integers that represent respectively the number of allocated words and the number of used words, ensuring the dynamic of the precision. Indeed, such a structure is well designed to efficiently handle arbitrary precision but can be too costly when one knows in advance the targeted precision. Furthermore, when dealing with multiple integers at the same time, i.e. in matrix computation, access through pointers breaks cache mechanisms and penalizes performances. Hopefully, one can access to a low level API (mpn level) to both handle fixed precision and multiple integer without suffering from an heavy interface. However, this approach let the memory management to the user and cannot be incorporated to code that have been designed at a higher level. The purpose of our work is to provide an alternative to the GMP library for fixed precision integers that allow the flexibility of a high level library while still being efficient. In particular, our approach is to design very simple codes that extend naturally machine word to other powers of two.

A common and efficient way to compute over prime finite fields is to first perform the operation over the integers and map the result back to the field via modular reduction. In this specific case, the use of arbitrary precision is not relevant since precision is fixed by the cardinality of the finite field. Indeed, doubling the precision is often sufficient. To further optimize this approach, one can use precompiled code to tackle specific range of modulus. This approach has been proposed in the MPFQ library for cryptographic size moduli [START_REF] Gaudry | The mpFq library and implementing curve-based key exchanges[END_REF]. While being very efficient to perform scalar computation over finite field, i.e. modular exponentiation, the design of this library does not allow to use lazy modular reduction that proved very efficient for linear algebra [START_REF] Dumas | Dense linear algebra over prime fields[END_REF].

Besides fixed precision integers, our package also provides prime fields support that can be easily embedded in high level code, while offering very good performance, in particular for matrix computations.

Functionality

Fixed precision arithmetic : extending the word size

With fixed precision integers, the maximum number of bit is given in advance and all arithmetic operations are done within this precision. In particular, if K is the maximum bitsize (i.e. the precision) then all computations are done modulo 2 K . This corresponds to forgetting the carry for addition/subtraction or to getting the lowest K-bits part of the integer products.

Fixing K in advance facilitates the storage of integers through static array of size K/64 on 64-bits architecture. In order to mimic word-size integers and to fully use the memory, we focus only on supporting integers of size a power of two: K = 2 k , with k >= 6. This assumption will allow us to provide a simple recursive data structure that eases manipulation and implementation of most arithmetic operations as explained in Section 3.1.

Fixed precision types

Our C++ package is devoted to provide functionality for fixed integer types. As for native integers, we provide signed and unsigned types. Furthermore, we extend this by providing a modular integer type that allows, e.g., computation over prime finite fields as a standard type. We denote by ruint<k> the type of our unsigned integer at precision 2 k (for instance we also define: using ruint128 = ruint<7>; using ruint256 = ruint<8>; etc.). Then the other types (signed and modular) are implemented with this structure. Our fixed precision integers are integrated in the Givaro library and available at https://github.com/linbox-team/givaro. Signed type. The recursive signed integer is defined as a rint<k> type. The data structure is just a ruint<k> and operations are performed via unsigned operations thanks to a two's complement representation.

Modular integer type.

A special type rmint<k> is provided for modular operations. The data structure is composed of an unsigned value, ruint<k>, as well as a global ruint<k> for the modulus that is shared by all modular elements. Modular elements modulo p are represented in the range [0, p -1] and all operations guarantee that values are always reduced to this range.

Integer operations

Our package provides basic integer operations defined both as mathematical operator (×, +, -, /) or with their functional variants (mul,add,sub,div). A simple example of using our package for doing vector dot product with 256-bit integers is given below: ruint<8> *A,*B,res=0; ... for(size_t i=0;i<n;i++) res+=A[i]*B[i]; Elements of type ruint<8> can be easily converted to rint<8> for computing with signed integer, or to rmint<8> for computing modulo a 256-bit prime. Several other functions are available (gcd, axpy, exponentiation, comparison, shift) as well as some internal functions that eases algorithm implementation, e.g. addition with carry, internal access to the structure, etc.

Two types of modular reduction. One can easily switch between two type of modular reduction. In particular, one can choose between classic modular reduction using Euclidean division or the one using Montgomery's method [START_REF] Montgomery | Modular multiplication without trial division[END_REF], described in Section 3.2. The choice of representation is done at compile time, by a template constant: rmint<k>, or rmint<K, MG INACTIVE> gives you the classical one, while rmint<K, MG ACTIVE> uses Montgomery's reduction.

Underlying theory and technical contribution

One drawback of arbitrary precision integers is to use a dynamic structure to handle the variation of integers value. There, defining contiguous structure of such integers leads to break cache prefetching when accessing data. Generally, one alternative is then to use a low level API where dynamics of the integers is delegated to the programmer. Unfortunately, this often becomes incompatible with the design of code at a high level, and then may not allow to re-use existing code. In the following section, we present our prototype designing a simple data structure preserving good performances.

Template recursive data structure

Using specificities of fixed precision integers, we can represent an integer of 2 k bits as two integers of 2 k-1 bits and therefore use a recursive representation. ruint<k> is our recursive unsigned integer of 2 k bits. The base case ruint<6> is called a limb and it corresponds to a native 64-bits integer. With these conventions, a ruint<k> can store any value between 0 and 2 2 k -1. In order to ensure a static contiguous storage, we chose to use a template recursive data structure with partial specialization. Note that the compiler will always completely unroll this structure at compiled time.

Our C++ template structure ruint<> is given hereafter together with an example for 256-bits integers, where High and Low correspond respectively to the leading and trailing 2 k-1 -bits of the 2 k -bit integer. Thanks to the recursive structure of ruint<k>, the implementation of arithmetic operations can be done by a recursive approach. For instance, addition can be done by two recursive calls with a carry propagation in between. Of course, the base case is mapped to the corresponding word-size integer operation. One major interest of such approach is that the compiler will unroll all the recursive calls leading to reduced control flow overhead and better instruction scheduling. As shown in Figure 1 (left), such approach leads to better performance against the GMP library for the addition of two integers (when the result is stored in one of the operand) up to 1024 bits. Multiplication is handled via a naive approach (or via Karatsuba's method), still using a recursive implementation, as shown in the following code: Still, without using any specific assembly code, contrary to GMP, our approach allows to have good speed-up for small integers up to 256-bits, as shown in Figure 1 (right). Remark: We use Karatsuba's method only for more than 1024 bits. Montgomery's modular reduction. Classical modular reduction uses integer division to map the result back into the desired integer range. Besides good theoretical complexity [START_REF] Brent | Modern Computer Arithmetic[END_REF], this approach is not efficient in practice since hardware divisions are costly. If there exists a radix β such that divisions with β are inexpensive and gcd(β, M ) = 1, Montgomery gives in [START_REF] Montgomery | Modular multiplication without trial division[END_REF] a method for modular reduction without trial division. The idea of Montgomery is to compute

Cβ -1 mod M instead of computing C mod M . As shown in [2, Algorithm 2.7],
the fast version of this reduction requires to compute integer products modulo β and integer division by β. Assuming β to be a power of 2, this provides an algorithm without any division.

When the integer C to be reduced has a precision 2 k+1 and the modulo has a precision of 2 k , i.e. β = 2 2 k , this boils down to two multiplications at precision 2 k . The design of our recursive integers is naturally compliant to such method and implementation is almost straightforward. From our first experiment, this straightforward implementation gives 1.5 speedup against GMP for 128-bits integer but becomes not competitive for larger modulus. Inversion modulo 2 2 k . On of the main operation in the setup of Montgomery reduction is the computation of the inverse of M modulo β = 2 2 k . The classical algorithm for this is to use the extended Euclidean algorithm. But with the special structure of B, it is better to use a Newton-Raphson iteration [START_REF] Arazi | On calculating multiplicative inverses modulo 2 m[END_REF][START_REF] Dumas | On Newton-Raphson iteration for multiplicative inverses modulo prime powers[END_REF], doubling the precision at each iteration. Thus here also our recursive structure is perfectly suited to this kind of algorithm, as shown in Figure 2 (left). Note that GMP uses an extended gcd. Recursive division. For the integer division, our structure is also well suited to a recursive algorithm. Thus, we use that of [START_REF] Burnikel | Fast recursive division[END_REF] which uses two sub-algorithms dividing respectively 2 digits by 1 digit and 3 halves by 2. They allow then a recursive division of an s-digits integer by a r-digits integer with complexity O(rs log(3)-1 + r log(s)). Recursive shift. One operation is slightly more complicated on a recursive structure: shifting. At each recursive level, five cases have to be explored: no shift; shifting the lower part on the high part, in parts or completely; shifting more than half the word; shifting more than the whole word. For instance, if the shift d is such that 0

< d < 2 2 k-1 , then b = a << d satisfies: b.High = (a.High d) ⊕ a.Low (2 k-1 -d) and b.Low = (a.Low d).
4 Application

Dense linear algebra: Freivalds certificate

For dense linear algebra, a basic building block is the matrix-vector multiplication. We show in Figure 2 (right) a multi-precision version of a modular matrixvector multi-precision. There, we use the FFLAS-FFPACK package [START_REF] Dumas | Dense linear algebra over prime fields[END_REF], version 2.2.1 (http://linalg.org/projects/fflas-ffpack) and just change the underlying representation from GMP to RecInt. This possibility of easily changing the base representation can also be extremely important when certifying the results. In verified computing, a client (the Verifier) will check the results provided by a server (typically a cloud, the Prover). In dense linear algebra the basic tools for this is Freivalds' certificate for matrix-multiplication [START_REF] Freivalds | Fast probabilistic algorithms[END_REF]. It uses matrix-vector multiplications with random vectors to check that C = AB, via Cv -A(Bv) ? == 0 on a random projection. With this certificate it is then possible to check any fast dense linear algebra computations with any precision [12, § 5]. These certificates have a double goal. First this is a way to improve the confidence in a result, and this is even more the case if the underlying data structure is different for the computation and for the certificate. Second it provides a way to check outsourced computations. This is economically viable only if the Verifier's time is faster than the Prover's time. We show in Figure 2 (right) that RecInt makes it possible to gain on all dense linear algebra verification, when operations require a few machine words.

Sparse Linear algebra: prime power rank

The Smith normal form is useful, in topology for instance, for computing the homology of a simplicial complex over the integers. There, the involved matrices are quite often sparse and computations reduce to computing the rank of these matrices modulo some prime powers [8, § 5.1]. We report in Table 1 on timings using the LinBox library [START_REF] Dumas | LinBox: A generic library for exact linear algebra[END_REF], version 1.4.1 (https://github.com/linbox-team/ linbox). The last invariant factor of the matrix S16.231x231 is 2 63 , while that of S22.1002x1002 is 2 85 , so both need some extra precision to be computed. But as this factor is not known in advance, a strategy is to double the precision, until the rank modulo the prime power reaches the integer rank of the matrix. For the considered matrices we would stop at 256 bits, but we see in Table 1 that the strategy can be faster until 512 bits. As a comparison we also give timings modulo 3 40 to show that this is not specific to the characteristic 2.

Towards an FPGA implementation

We present here the first attempts towards an implementation on a FPGA. To build a programmable hardware (http://shiva.minalogic.net/), we had to provide basic arithmetic libraries to be used in an Elliptic curve based encryption scheme. We chose to use a dedicated software transforming C++ source into VHDL called GAUT [START_REF] Coussy | High-Level Synthesis: From Algorithm to Digital Circuit, chapter GAUT: A High-Level Synthesis Tool for DSP Applications[END_REF]. The creation of a VHDL program can be split in the following steps: compilation of the C++ source and creation of the corresponding graph; compilation of the library containing the required operations; synthesis of the VHDL program and estimation of performance. Table 2 shows some simulations of a modular exponentiation on a Virtex 5. We made the output flow vary in order to check the effect on the required size on the FPGA. We notice that required size can be significantly reduced if we accept a lower output flow. One nice point is that, due to the simple recursive structure of RecInt, these results have been obtained without significant modifications on the C++ source and automatically transformed. They are not optimal but rather promising, with no significant work load.

Conclusion and Perspective

Our RecInt package is a first attempt to provide high level API for fixed precision integers that are usable out of the box. Thus it is easy to switch from native integer types and to still provide good performance when compared to standard libraries. To further improve the performance, we need to introduce SIMD vectorized instructions in the design of our specializations and remove as much as possible conditional jumps. Modular multiplication should also benefit from either Barret's method or Montgomery-Svoboda's algorithm [2, §2.4.2].

6 RecIntFig. 1 .

 61 Fig. 1. Comparing integer addition (left) and multiplication (right) with GMP library

6 RecIntFig. 2 .

 62 Fig. 2. Newton-Raphson iteration for the inverse modulo a prime power (left) and Modular matrix-vector multiplication (right)

Table 1 .

 1 Rank modulo a prime power http://hpac.imag.fr/Matrices/Tsuchioka

	Matrix	mod int64 t	6	7	RecInt 8	9	10	GMP
	S16.231x231 3 40	-0.01s 0.05s 0.12s 0.37s 1.28s	0.19s
	S16.231x231 2 63 0.03s 0.03s 0.09s 0.29s 1.15s 4.64s	1.85s
	S16.231x231 2 64	-	-0.09s 0.30s 1.15s 4.64s	1.85s
	S22.1002x1002 3 40	-0.34s 1.04s 2.60s 7.67s 25.84s	6.32s
	S22.1002x1002 2 63	-2.60s 7.00s 23.10s 88.24s 356.79s 154.22s
	S22.1002x1002 2 86	-	-7.37s 24.17s 90.30s 357.44s 190.36s

Table 2 .

 2 FPGA area of modular exponentiation for different output flows, withRecInt integers from 128 to 512 bits

	128 bits	op/s 7812 15625 31250 62500 125000 250000 Flipflops 3040 6100 7008 9538 16617 32199
	256 bits	op/s Flipflops 3618 4391 6923 7654 8776 14542 976 1953 3906 7812 15625 31250
	512 bits	op/s Flipflops 3553 3553 3553 4704 5729 7458 61 122 244 488 976 1953

This material is based on work supported in part by the Agence Nationale pour la Recherche under Grant ANR-11-BS02-013 HPAC.


Nawal Benabbou

Patrice Perny

Solving Multi-agent Knapsack Problems Using Incremental Approval Voting

In this paper, we study approval voting for multi-agent knapsack problems under incomplete preference information. The agents consider the same set of feasible knapsacks, implicitly defined by a budget constraint, but they possibly diverge in the utilities they attach to items. Individual utilities being difficult to assess precisely and to compare, we collect approval statements on knapsacks from the agents with the aim of determining the optimal solutions by approval voting. We first propose a search procedure based on mixedinteger programming to explore the space of utilities compatible with the known part of preferences in order to determine or approximate the set of possible approval winners. Then, we propose an incremental procedure combining preference elicitation and search in order to determine the set of approval winners without requiring the full elicitation of the agents' preferences. Finally, the practical efficiency of these procedures is illustrated by various numerical tests.

INTRODUCTION

Collective decision making on a combinatorial domain appears in various contexts such as investment planning, resource allocation or group configuration. Due to strategic aspects often surrounding group decision-making and the possible divergences in individual values, developing formal methods and tools for modeling preferences and solving multi-agent combinatorial optimization problems is a critical issue. This has motivated a lot of work in the recent years, in the field of computational social choice [9]. We focus here on the multi-agent knapsack problem which consists of determining, given a finite set of items, a subset of maximal utility under a budget constraint. This is a standard example of combinatorial problem with many potential applications such as project selection, portfolio management or committee election, see e.g. [START_REF] Lu | Budgeted social choice: From consensus to personalized decision making[END_REF][START_REF] Elkind | Properties of multiwinner voting rules[END_REF][START_REF] Oren | Online (budgeted) social choice[END_REF][START_REF] Skowron | Finding a collective set of items: From proportional multirepresentation to group recommendation[END_REF] for examples of recent contributions in AI.

In combinatorial optimization problems, the agents cannot be expected to provide extensive preference models. Compact representations are needed to handle individual and collective preferences. Usually, in knapsack problems, preference over subsets of items are represented by additive utility functions. More precisely, the utility of a subset of items for an agent is defined as the sum of the utilities of its elements. Individual utilities are difficult to assess, especially on a combinatorial domain. Although the elicitation task is simplified when utilities are decomposable, elicitation methods based on systematic pairwise comparisons are practically unfeasible due to the large amount of feasible subsets and their implicit definition. Hence we are interested in designing incremental elicitation procedures, in which preference queries are selected iteratively, to be as informative as possible at every step, so as to progressively reduce the set of admissible utility profiles until the set of optimal knapsacks can be determined. This approach has been successfully used in AI for additive utility elicitation on explicit sets [START_REF] Chajewska | Making rational decisions using adaptive utility elicitation[END_REF][START_REF] Wang | Incremental Utility Elicitation with the Minimax Regret Decision Criterion[END_REF][START_REF] Boutilier | Constraint-based Optimization and Utility Elicitation using the Minimax Decision Criterion[END_REF], but also on combinatorial solution spaces [START_REF] Gelain | Elicitation strategies for soft constraint problems with missing preferences: Properties, algorithms and experimental studies[END_REF][START_REF] Benabbou | Incremental Weight Elicitation for Multiobjective State Space Search[END_REF][START_REF] Benabbou | On possibly optimal tradeoffs in multicriteria spanning tree problems[END_REF].

Remark that, even if numerical representations of individual preferences are accessible under the form of utility functions, they are generally constructed independently for each agent. Hence, it is unlikely that such representations allow the welfare of individuals to be compared. In this context, the definition of a social utility as the sum of individual utilities (utilitarism), for example, would be meaningless. Assuming that utilities of items are expressed on the same scale or that utilities are normalized would not be sufficient to overcome the problem, as shown by the following example:

Example 1. Consider a multi-agent knapsack problem involving 3 items and 2 agents with utilities:

u 1 = u 1 1 x1 + u 1 2 x2 + u 1 3 x3 and u 2 = u 2 1 x1 + u 2 2 x2 + u 2 3
x3 to be maximized under the constraint x1 + x2 + x3 ≤ 2, where xi ∈ {0, 1}, i = 1, 2, 3, are the decision variables and u i j represents the utility of item j for agent i. This problem could appear to elect a committee of size 2, given 3 candidates and 2 voters. Assume that individual preference orders over committees have been elicited, and are equal to {2, 3} 1 {1, 3} 1 {1, 2} and {1, 2} 2 {1, 3} 2 {2, 3} for agents 1 and 2 respectively. One possible numerical representation of these preferences (using the same utility scale for the two agents) is given by: (u 1 1 , u 1 2 , u 1 3 ) = (1, 2, 4) and (u 2 1 , u 2 2 , u 2 3 ) = (4, 2, 1) which leads to the following utilities for solutions of size 2:

{1, 2} {2, 3} {1, 3} u 1 3 6 5 u 2 6 3 5
Note that these individual values are consistent with preference orders 1 and 2. Now, if we are utilitarian, we could be tempted to deduce that {1, 3} is the optimal knapsack because it maximizes the total utility (5+5 = 10). However, such a conclusion would be meaningless, it is only due to the particular numerical representation chosen for individual utilities. Let us change the initial numerical scale by replacing numbers (1, 2, 4) by (0, [START_REF] Benabbou | Incremental Weight Elicitation for Multiobjective State Space Search[END_REF][START_REF] Benabbou | On possibly optimal tradeoffs in multicriteria spanning tree problems[END_REF]. In this case we obtain two new utility functions characterized by

(u 1 1 , u 1 2 , u 1 3 ) = (0, 3, 4) and (u 2 1 , u 2 2 , u 2
3 ) = (4, 3, 0) which leads to the following utilities for solutions of size 2:

{1, 2} {2, 3} {1, 3} u 1 3 7 4 u 2 7 3 4
Note that these new individual values are still consistent with preference orders 1 and 2. Yet, with the same utilitarian principle, we should admit now that {1, 3} is the least preferred knapsack. Therefore, choosing the solution maximizing the sum of individual utilities would not be a good procedure. It would merely be a consequence of arbitrary choices of numerical representations of preferences rather than a robust conclusion derived from the observed preference profile. Note that the problem persists if we normalize u 1 and u 2 utilities to obtain value 1 for set {1, 2, 3}. Other examples could be found for other aggregators (e.g. the minimum for an egalitarian aggregation).

In order to be able to compare the solutions of a knapsack problem when individual utility scales are not commensurate and/or not rich enough to allow the construction of a social utility, it seems natural to resort to a voting rule. The main advantage of a voting rule is indeed to perform an ordinal aggregation procedure. There is no need to know how the welfare of individuals should be compared, we only need to elicit individual preference orders. Some recent contributions consider voting rules in a very different perspective, for example, by investigating their ability to approximate the winner with respect to the utilitarian criterion [START_REF] Ariel | The distortion of cardinal preferences in voting[END_REF][START_REF] Caragiannis | Voting almost maximizes social welfare despite limited communication[END_REF][START_REF] Boutilier | Optimal social choice functions: A utilitarian view[END_REF].

Preference elicitation can be performed incrementally so as to determine the winner with a reduced amount of preference queries. Various incremental elicitation procedures have been proposed and studied in the context of single-winner elections with incomplete preferences [START_REF] Kalech | Practical voting rules with partial information[END_REF][START_REF] Lu | Robust approximation and incremental elicitation in voting protocols[END_REF][START_REF] Lihi | Reaching a joint decision with minimal elicitation of voter preferences[END_REF]. In this setting, several contributions study the determination of possible and necessary winners from a partial preference profile, e.g., [START_REF] Konczak | Voting procedures with incomplete preferences[END_REF][START_REF] Xia | Determining possible and necessary winners given partial orders[END_REF][START_REF] Lang | Winner determination in voting trees with incomplete preferences and weighted votes[END_REF][START_REF] Ding | Voting with partial information: what questions to ask?[END_REF], when the set of candidates is defined explicitly.

In knapsack problems however, solutions are numerous and defined implicitly, which is an additional challenge for the winner determination. This explains the current interest for incremental voting procedures on combinatorial domains and the purpose of this paper. Our aim here is to propose an incremental voting rule in which individual preferences are progressively revealed until a collective decision can be made, and to apply this procedure on the multi-agent knapsack problem, taking advantage of the fact that individual preferences are representable by additive utilities.

It is important to note that implementing a voting rule is not in contradiction with the representation of individual values by utilities. Individual utility functions are indeed seen as convenient representations of individual preference orders, and their use will significantly contribute to relieve the preference elicitation burden. In the standard knapsack problem, due to the linearity of preferences, the set of all preference orders compatible with a given partial order can be characterized by a convex polyhedron in the utility space. This makes it possible to resort to mathematical programming to explore all possible completions of any partially known preference profile, but also to look for possible winners, and to develop an efficient incremental elicitation procedure for the determination of all winners, as it will be seen later in the paper.

Implementing a voting rule for the knapsack problem with a partially specified preference profile could be related to multi-winner voting rules studied in the field of computational social choice. Most approaches recently proposed for multi-winner elections assume that individual preferences over items are sufficient to explain preference over subsets because they derive satisfaction from their most preferred candidate (see e.g., [START_REF] John R Chamberlin | Representative deliberations and representative decisions: Proportional representation and the Borda rule[END_REF][START_REF] Burt | Fully proportional representation[END_REF][START_REF] Richard | Proportional representation broadening the options[END_REF][START_REF] Ariel D Procaccia | On the complexity of achieving proportional representation[END_REF][START_REF] Meir | Complexity of strategic behavior in multi-winner elections[END_REF][START_REF] Lu | Budgeted social choice: From consensus to personalized decision making[END_REF][START_REF] Betzler | On the computation of fully proportional representation[END_REF][START_REF] Elkind | Properties of multiwinner voting rules[END_REF][START_REF] Skowron | The complexity of fully proportional representation for single-crossing electorates[END_REF], and see [START_REF] Lu | Multi-winner social choice with incomplete preferences[END_REF] for incremental elicitation of voter preferences). This assumption is well-suited to the election of representatives. However, for any agent, it may happen that the selection of the most preferred candi-date is not sufficient to counterbalance the presence of multiple least preferred candidates in the elected committee. This limitation also applies to the selection of items in multi-agent knapsack problems.

In this paper, we focus on approval voting because this is a simple rule that can be decisive even if only a part of the preference profile is known. Approval voting is the voting method which allows each agent to approve of (vote for) as many solutions as she wishes, and the solution with the most approval votes is the winner of the election [START_REF] Brams | Approval voting[END_REF]. In approval voting, we only need to learn, for every agent, which are the approved or disapproved subsets, so as to elect a solution receiving the maximal support. We therefore investigate incremental procedures for approval voting and their application to the knapsack problem. This work differs from multi-winner approval voting [START_REF] Marc | Approval balloting for multi-winner elections[END_REF][START_REF] Aziz | Computational aspects of multiwinner approval voting[END_REF][START_REF] Aziz | Justified representation in approval-based committee voting[END_REF] which only collects approval statements over items instead of feasible subsets of items.

The paper is organized as follows: we introduce the multi-agent knapsack problem in Section 2 and study computational issues for this problem. In Section 3, we propose a search procedure for the determination of the possible winners. Then, an incremental voting procedure to determine the set of approval winners is proposed in Section 4. Finally, numerical tests are provided in Section 5.

THE GENERAL FRAMEWORK

We consider a collective decision problem where a set of agents N = {1, . . . , n} has to jointly select a set of items (e.g., candidates, projects, objects) in a set P = {1, . . . , p}. Any subset of items can be represented by a solution vector x = (x1, . . . , xp) ∈ {0, 1} p where xj = 1 if item j is in the subset and xj = 0 otherwise. Some linear constraints on variables xj, j ∈ P , are imposed to define the admissible solution vectors. For instance, one may want to impose cardinality constraints to control the size of the subset and/or to ensure gender parity in the elected committee; there may also exist budget constraints (e.g., when the decision is subject to a maximum total cost) or capacity constraints as in knapsack problems, making some subsets of items unfeasible. For the simplicity of the presentation, we will only consider the standard knapsack constraint of the form j∈P wjxj ≤ W where wj is the (positive) weight of item j and W is a positive value representing the maximum total weight; the set of feasible solutions will be denoted by X in the sequel. Our purpose and the algorithms proposed in the paper also apply when additional (linear) feasibility constraints are considered.

We assume that the preferences of agent i, i ∈ N , can be represented by a function u i : {0, 1} p → R measuring the overall utility of any solution. Hence, given two solutions x, y representing two subsets of items, x is at least as good as y for agent i ∈ N whenever u i (x) ≥ u i (y). Here u i (x) = j∈P u i j xj, where u i j ∈ R represents the utility of item j for agent i. The profile (u 1 , . . . , u n ) of utility functions will be denoted by u. Note also that u i , as a numerical representation of a preference order, is generally not unique and any transform preserving inequalities of type u i (x) ≥ u i (y) for all solutions x, y could be considered as well.

Nevertheless, numerical representations of individual preferences by utility functions u i , i ∈ N , can be used in approval voting. Under the assumption that individual preferences are represented by utility functions u i , i ∈ N , a solution x is approved by agent i if and only if u i (x) ≥ δ i where δ i ∈ R is an approval (or utility) threshold that separates approved and non-approved solutions. The profile of thresholds (δ 1 , . . . , δ n ) will be denoted by δ in the sequel. The pair (u, δ) characterizes approved and non-approved solutions for all agents and enables the computation of approval scores for any feasible solution x. This approval score is given by f (x, u, δ) = |{i ∈ N, u i (x) ≥ δ i }|, and the winner of the election is a feasible solution maximizing this score (various tie-breaking rules can be considered). In this paper, we consider the approval multi-agent knapsack problem defined as follows:

APPROVAL MULTI-AGENT KNAPSACK PROBLEM (AMKP) Input: A finite set P of items; a positive integer W ; for each j ∈ P , a weight wj; a finite set N of agents; a positive integer K; for each i ∈ N , an approval threshold δ i , for each j ∈ P , a utility value u i j . Question: Is there a subset X ⊆ P such that j∈X wj ≤ W and |{i ∈ N, j∈X u i j ≥ δ i }| ≥ K?

Proposition 1. AMKP is NP-complete.
Proof. The proof is quite straightforward due to a simple reduction from the knapsack decision problem. Testing the existence of an admissible knapsack having a utility greater or equal to a given value K is indeed equivalent to solving an instance of the AMKP problem involving a single agent with the same utilities over items, an approval threshold equal to K and with K = 1.

There are obvious tractable cases when W is a constant and considering integer weights, for example, committee election problems (wj = 1 for all j ∈ P ) such that the committee size is a constant (W ) specified explicitly. Indeed, we have: Proposition 2. When W is constant and wj ∈ N\{0} for all j ∈ P , AMKP is in P.

Proof. Since wj is a non-zero positive integer for all j ∈ P , then we know that all feasible knapsacks necessarily include at most W items. The number of knapsacks of size at most W is equal to:

W k=0 p k = W k=0 p! k!(p -k)!
This number is obviously polynomial in p when W is a constant. Therefore, the result can be easily obtained by considering the following naive procedure: for all sets X ⊆ P of size at most W , we test whether j∈X wj ≤ W (feasibility condition), and if the test succeeds, we compute |{i ∈ N, j∈X u i j ≥ δ i }| (approval score). This procedure is polynomial in p and n when W is a constant.

Proposition 1 shows that finding the knapsack maximizing the approval score is NP-hard in the general case. Moreover, wellknown pseudo-polynomial solution methods proposed for the standard knapsack problem (based on dynamic programming) are not easily transposable to the approval winner determination problem, as shown in Example 2.

Example 2. Consider a collective decision problem where 2 agents have to choose 2 representatives from a pool of 3 candidates, i.e. N = {1, 2} and P = {1, 2, 3}. Assume that utilities and approval thresholds are the following:

u 1 1 u 1 2 u 1 3 δ 1 u 2 1 u 2 2 u 2 3 δ 2 0.5 0.3 0.2 0.5 0.1 0.5 0.3 0.7 In this case, we have f ({1}, u, δ) = 1 > 0 = f ({2}, u, δ) but f ({1, 3}, u, δ) = 1 < 2 = f ({2, 3}, u, δ).
We observe a preference reversal because {1} is preferred to {2} whereas {2, 3} is preferred to {1, 3}. Thus, preferences induced by the approval score are not additive with respect to union with disjoint items. This precludes to construct the optimal knapsack from optimal subsets of items.

Nevertheless, the winner can be obtained by solving the following mixed integer program (MIP1):

max i∈N a i s.t.                j∈P u i j xj -δ i ≥ M (a i -1), ∀i ∈ N j∈P wjxj ≤ W xj ∈ {0, 1}, a i ∈ {0, 1}, ∀i ∈ N, ∀j ∈ P
In this program, M is a constant greater than max{δ i -u i (x), i ∈ N } that allows the introduction of boolean variables a i which will be equal to 1 if and only if agent i approves solution x. Moreover, the second Equation is the knapsack constraint.

However, in practice, the full elicitation of individual utilities and approval thresholds is too expensive. Usually, we can observe simple preference statements of type "I prefer solution x to solution y", and in the case of approval voting, "I approve solution x", or "I don't approve solution y". These preference statements enable to restrict the sets of possible utility functions but generally do not allow to derive a precise utility function for each agent (see Example 1). Instead, uncertainty sets representing all possible utility functions compatible with the preference information obtained so far must be considered. The same observation applies to approval thresholds. Under utility uncertainty, we study now the determination of possible winners.

POSSIBLE WINNERS DETERMINATION

In this section, we propose an algorithm that enables to compute the set of possible approval winners given some partial knowledge of the agents' preferences. More precisely, the input of the algorithm consists of three sets of preference information for each i ∈ N : a set A i (resp. Āi ) of solutions that are known to be approved (resp. not approved) by agent i, and a set P i of pairs (y, z) such that solution y is known to be preferred to solution z by agent i. The elements of P i are not explicit approval statements, but they can be used to derive new positive or negative approval statements from those included in A i and Āi .

For each agent i ∈ N , let U i (resp. ∆ i ) denotes the set of utility functions (resp. approval thresholds) compatible with the available preference statements A i , Āi and P i . Formally, (U i , ∆ i ) is the set of all pairs (u i , δ i ) such that:

∀y ∈ A i , u i (y) ≥ δ i ; ∀y ∈ Āi , u i (y) < δ i ; ∀(y, z) ∈ P i , u i (y) ≥ u i (z)
where u i is of the form u i (x) = j∈P u i j xj and δ i ∈ R. Let U (resp. ∆) be the cartesian product U 1 × . . . × U n (resp. ∆ 1 × . . . × ∆ n ).
Given such uncertainty sets, the set of possible approval winners is defined as follows: Definition 1. The set P W (X , U, ∆) of possible approval winners is the set of all solutions x ∈ X that maximize the approval score for some utility profile u ∈ U and some approval threshold vector δ ∈ ∆. More formally:

P W (X , U, ∆) = u∈U,δ∈∆ arg max x∈X f (x, u, δ).
Recall that Example 2 shows that standard dynamic programming procedures cannot be used to determine the approval winners when utilities and approval thresholds are known. This difficulty remains when utilities and/or approval thresholds are partially known.

The branch and bound approach is the most commonly used tool for solving NP-hard optimization problems. We propose here a branch and bound procedure to compute the set PW(X , U, ∆), where nodes of the search tree represent partial instances of the decision variable vector x = (x1, . . . , xp). More precisely, each node η of the tree is characterized by a pair (P 0 η , P 1 η ) where P k η = {j ∈ P, xj = k}, k = 0, 1. Let Pη = P \ (P 0 η ∪ P 1 η ) denote the set of all undecided variables at node η. Thus, each node η is associated with a region of the solution space as follows: solution x = (x1, . . . , xp) is attached to node η if and only if xj = k for all j ∈ P k η , k = 0, 1. The set of feasible solutions attached to node η is denoted by Sη hereafter. The main features of our search procedure are the following:

Initialization. Using a heuristic, a branch and bound procedure determines some feasible solutions before performing the search so as to define an initial bound on candidate solutions.

In order to obtain such a bounding set for the knapsack problem, denoted by S0 hereafter, we propose to initially ask the agents to rank all the items by preference order. Let ri(j) denote the rank of item j in the preference order provided by agent i. We can define the score α i (j) = p-r i (j) w j for each agent i ∈ N and each item j ∈ P representing the tradeoff achieved between preference and weight. Hence, for each agent i, a "good" solution to the knapsack problem can be obtained by a greedy algorithm selecting items one by one, by decreasing order with respect to scoring function α i , skipping elements whose weight is greater than the residual weight capacity. The resulting solution is inserted in S0 for initialization because it represents a good solution from the point of view of agent i. This process is repeated for all agents i ∈ N .

Moreover, a similar procedure is used with the average scoring function defined by α(j) = 1/n n i=1 α i (j), to complete S0 with a solution which is likely to be more consensual. This solution will be denoted by s in the sequel.

Evaluation and pruning. Let S be the set of solutions found so far (initially S = S0) and O be the current set of nodes to be explored. Our pruning rule is based on the notion of setwise max regret defined as follows: the setwise max regret SR(A, B, U, ∆) of a set A ⊆ X with respect to a set B ⊆ X is the maximal feasible approval score difference between the best solution in B and the best solution in A. More formally:

SR(A, B, U, ∆) = max u∈U,δ∈∆ max b∈B f (b, u, δ) -max a∈A f (a, u, δ)
If SR(A, B, U, ∆) < 0, then we know that B does not contain any possible approval winner; it indeed induces that, for all solutions b ∈ B, for all u ∈ U and for all δ ∈ ∆, there exists a ∈ A such that f (b, u, δ) < f (a, u, δ). Therefore, we propose to prune a node η ∈ O if the setwise max regret SR(S, Sη, U, ∆) of set S with respect to set Sη is strictly negative. Note that SR(S, Sη, U, ∆) = maxx∈S η max u∈U,δ∈∆ mins∈S{f (x, u, δ) -f (s, u, δ)}. This alternative formulation of setwise max regrets enables to compute SR(S, Sη, U, ∆) as the optimal value of the mixed-integer quadratic program (denoted by MIQPη) given in Figure 1. In this program, ξ > 0 is an arbitrary small value enabling to model strict inequalities. Equations (2e-2g) enable to restrict utility functions and approval thresholds to those compatible with the available preference information. Then, since the preferences of agent i, for any i ∈ N , are invariant by positive affine transformations jointly applied to function u i and threshold δ i , we can assume without loss of generality

max t s.t.                                                                                        t ≤ i∈N a i - i∈N a i s , ∀s ∈ S j∈P 1 η u i j + j∈Pη u i j xj -δ i ≥ M (a i -1), ∀i ∈ N j∈P u i j sj -δ i + ξ ≤ M a i s , ∀s ∈ S, ∀i ∈ N j∈Pη wjxj + j∈P 1 η wj ≤ W j∈P u i j = M, ∀i ∈ N j∈P u i j yj ≥ δ i , ∀i ∈ N, ∀y ∈ A i j∈P u i j yj ≤ δ i -ξ, ∀i ∈ N, ∀y ∈ Āi j∈P u i j yj ≥ j∈P u i j zj, ∀i ∈ N, ∀(y, z) ∈ P i xj ∈ {0, 1}, ∀i ∈ N, ∀j ∈ Pη a i ∈ {0, 1}, a i s ∈ {0, 1}, ∀i ∈ N, ∀s ∈ S u i j ≥ 0, δ i ≥ 0, ∀i ∈ N, ∀j ∈ P (2a) (2b) (2c) (2d) (2e) (2f) (2g)
Figure 1. MIQPη that utilities are positive and bounded above by a constant M > 0 (see Equation (2d)). Moreover, a i is a boolean variable that will be equal to 1 iff agent i approves solution x, and a i s is a boolean variable that will be equal to 1 iff agent i approves solution s, s ∈ S. Finally, Equation (2a) introduces variable t ∈ R representing the smallest approval score difference between solution x and a solution s, s ∈ S.

Note that constraints given in Equation (2b) include quadratic terms of type u i j xj, j ∈ Pη, since u i j are also variables of the optimization problem. In order to linearize theses constraints, we introduce positive variables v i j , i ∈ N, j ∈ Pη, representing the product u i j xj and Equation (2b) is replaced by the following constraints:

       j∈P 1 η u i j + j∈Pη v i j -δ i ≥ M (a i -1), ∀i ∈ N v i j ≤ u i j , ∀i ∈ N, ∀j ∈ Pη v i j ≤ M xj, ∀i ∈ N, ∀j ∈ Pη v i j -u i j ≥ M (xj -1), ∀i ∈ N, ∀j ∈ Pη
The resulting mixed-integer linear program will be denoted by MIPη.

Branching. Setwise max regrets SR(S, Sη, U, ∆) available for all nodes η ∈ O are also used to select the next node to be explored. More precisely, we select here a node η ∈ O which maximizes SR(S, Sη, U, ∆). This branching strategy aims to maximally improve the current solution set S. The optimal solution of MIPη indeed maximizes the gap f (x, u, δ) -maxs∈S f (s, u, δ) over all u ∈ U , all δ ∈ ∆ and all x ∈ η ∈O S η . Then, Sη is split in two by considering possible instantiations of a variable xj, j ∈ Pη chosen among the variables equal to 1 in the optimal solution of MIPη.

Filtering. As we will see in Proposition 3, the proposed Branch and Bound outputs, in general, a superset of the set of possible approval winners. To remove undesirable elements, we use a final filtering process, named FILTER hereafter, which iteratively deletes all solutions s ∈ S such that SR(S\{s }, {s }, U, ∆) < 0, using a simplified version of MIPη.

The algorithm implementing these principles is referred to AS (Approval-based Search) in the sequel and it is summarized by Algorithm 1. Proof. Since, the pruning rule only prunes nodes including no possible approval winner (by definition of setwise max regrets), we know that S is a superset of PW(X , U, ∆) at the end of the while loop. Let x be an element inserted in S such that x ∈ PW(X , U, ∆), if it exists. Since x is not a possible winner, we know that, for all u ∈ U and for all δ ∈ δ, there exists x ∈ PW(X , U, ∆) ⊆ S such that f (x, u, δ) < f (x , u, δ). Therefore, x is necessarily removed from S during the filtering (by definition of the filtering process).

Depending on the uncertainty sets (U and ∆) implicitly defined by the available preference information, possible approval winners might be too numerous to be enumerated efficiently. In order to save time, one may be interested in approximating the set of possible approval winners with performance guarantees. We propose below a variant of Algorithm AS for approximating possible winners with some guarantees on the quality of the output.

Approximation. Given a constant ε > 0, a set X ⊆ X is an (1 + ε)approximation of the set of possible winners if, for all x ∈ X , all u ∈ U and all δ ∈ ∆, there exists x ∈ X such that f (x, u, δ) ≤ (1 + ε)f (x , u, δ). In order to compute an (1 + ε)-approximation of the set of possible winners, we introduce an approximate version of the setwise max regret. The setwise max ε-regret SRε(A, B, U, ∆) of a set A ⊆ X with respect to a set B ⊆ X is defined as follows:

SRε(A,B,U,∆) = max u∈U,δ∈∆ max b∈B f (b, u, δ)-max a∈A (1+ε)f (a, u, δ)
If SRε(A, B, U, ∆) ≤ 0, then we know that, for all b ∈ B, all u ∈ U and all δ ∈ ∆, there exists a ∈ A such that f (b, u, δ) ≤ (1 + ε)f (a, u, δ). Therefore, we propose a variant of AS where a node η ∈ O is pruned if SRε(S, Sη, U, ∆) ≤ 0. This pruning rule is sharper than the previous one and is more likely to prune nodes in the search tree. The implementation is simple within Algorithm AS: computations of SR values must simply be replaced by computations of SRε values. Moreover, the value SRε is obtained using MIPη in which Equation (2a) is simply replaced by:

t ≤ i∈N a i -(1 + ε) i∈N a i s , ∀s ∈ S
The resulting algorithm is denoted by ASε in the sequel. Then, the following proposition holds.

Proposition 4. ASε returns an (1 + ε)-approximation of the set PW(X , U, ∆).

Proof. Let x be a solution that does not belong to S at the end of the search procedure. Let u ∈ U and δ ∈ ∆. We want to prove that f (x, u, δ) ≤ (1 + ε)f (s, u, δ) for some s ∈ S. If x was inserted in S at some step, then x has necessarily been deleted during the filtering of S. In that case, by definition of the filtering, we know that there exists s ∈ S such that f

(x, u, δ) ≤ f (s, u, δ); hence, f (x, u, δ) ≤ f (s, u, δ) ≤ (1 + ε)f (s, u, δ).
Assume now that x was pruned at some step without ever belonging to S. In that case, we know that, at this step, there exists s ∈ S such that f (x, u, δ) ≤ (1 + ε)f (s, u, δ). If solution s belongs to S at the end of the procedure, then we can directly infer the result. If solution s is deleted during the filtering, then we know that there exists s ∈ S such that f (s, u, δ) ≤ f (s , u, δ). Therefore, f (x, u, δ) ≤ (1 + ε)f (s, u, δ) ≤ (1 + ε)f (s , u, δ) which completes the proof (note that there is no chaining of errors).

ELICITATION FOR WINNER DETERMINATION

In Section 3, we have introduced a procedure for the determination of possible approval winners. It can be used in situations where a significant part of approval judgements is available. However, with incomplete preference information, the number of possible winners remains often very large, due to the combinatorial nature of the knapsack problem. Actually, the notion of possible winner (even its approximate version) is not sufficiently discriminating to support efficiently a collective decision making process. A more discriminating concept is the notion of necessary winner defined as follows:

Definition 2. The set N W (X , U, ∆) of necessary approval winners is the set of all solutions x ∈ X that maximize the approval score for all utility profiles u ∈ U and all approval threshold vectors δ ∈ ∆. More formally :

N W (X , U, ∆) = u∈U,δ∈∆ arg max x∈X f (x, u, δ).
However, given the uncertainty sets U and ∆, it might be the case that no necessary winner exists. In this situation, we need to collect more preference information so as to be able to make a collective decision. For this reason, we focus now on preference elicitation for the approval winners determination.

One may consider a two stage procedure that consists first in eliciting the entire set of approval judgements and then in applying the approval voting method. However, the exhaustive elicitation, prior to preference aggregation, is unfeasible due to the combinatorial nature of the problem. Collecting all approval statements would indeed require O(n2 p ) queries, for n agents and p items. We propose instead to combine preference elicitation and search so as to quickly focus the search on the relevant subsets of items and to concentrate the elicitation burden on the useful part of preferences. Our proposal is to collect approval statements from individuals during the search so as to progressively reduce the uncertainty attached to approval scores until determining the actual winners, i.e. the solutions maximizing the approval score.

Note that the number of possible winners reduces with the uncertainty about the agents' preferences. More precisely, for any U ⊆ U and any ∆ ⊆ ∆, we have PW(X , U , ∆ ) ⊆ PW(X , U, ∆). Moreover, if U reduces to a single utility profile and ∆ to a single approval threshold vector, then the possible winners become the actual winners of the election. More generally, if we consider an incremental elicitation procedure that iteratively collects preference information so as to progressively reduce sets U and ∆, we will necessarily reach a point where all possible winners are necessary winners. This will generally happen long before reducing U and ∆ to singletons. At that point, we know that the possible winners are the actual winners.

Therefore, we propose now an incremental elicitation procedure progressively reducing uncertainty sets to U and ∆ such that PW(X , U , ∆ ) = NW(X , U , ∆ ). Our incremental elicitation algorithm consists in inserting preference queries in Algorithm AS (see the previous section) so as to discriminate between the current best solutions (those that were stored in S). In practice, S is now restricted to the most approved solutions found so far. Within this set, we can arbitrarily select a representant, named the incumbent. Initially, the incumbent may be any feasible solution to the knapsack problem (e.g., solution s, see the initialization in Section 3). Each time a new solution (or a challenger) is found (line 6, Algorithm 1), it is now compared to the incumbent in terms of approval scores.

Note that, given the current uncertainty sets U and ∆, a solution x ∈ X is necessarily approved by an agent i, i ∈ N, if and only if u i (xj) ≥ δ i holds for all u ∈ U and all δ ∈ ∆. This can be checked by testing whether the optimum of the following linear program is positive: min j∈P u i j xj -δ i subject to Equations (2d-2g) and δ i ≥ 0, u i j ≥ 0 ∀j ∈ P . Similarly, testing whether a solution x ∈ X is necessarily disapproved by an agent or not can be performed using linear programming. Thus, each time a new solution is found, we can efficiently determine the agents -if they exist -that necessarily approve/disapprove it. Then, a natural query generation strategy, denoted by σ0 hereafter, consists in questioning all the other agents to know whether they approve this solution or not. The challenger will be inserted in S if it has the same approval score than the incumbent. If the challenger is strictly better than the incumbent, then S is reinitialized to include only the challenger. This strategy provides an incremental search algorithm, named Algorithm IAS (Incremental Approval-based Search) in the sequel, that determines the set of approval winners, as shown by the following proposition: Proposition 5. IAS returns the exact set of approval winners. Proof. To prove this result, it is sufficient to prove the following: for any initial uncertainty sets U and ∆, Algorithm IAS terminates with uncertainty sets U ⊆ U and ∆ ⊆ ∆ such that PW(X , U , ∆ ) = NW(X , U , ∆ ). Let s ∈ S be any solution returned by IAS. For all x ∈ X , we want to prove that we have f (s, u, δ) ≥ f (x, u, δ) for all u ∈ U and all δ ∈ ∆ at the end of the search procedure. First, whenever a node η is pruned at some step k of IAS using the current incumbent s k , we know that, for all x ∈ Sη, f (x, u, δ) ≤ f (s k , u, δ) for all u ∈ U k and all δ ∈ ∆ k where U k and ∆ k are the current uncertainty sets at step k. Since U ⊆ U k and ∆ ⊆ ∆ k , we have f (x, u, δ) ≤ f (s k , u, δ) for all u ∈ U and all δ ∈ ∆ . Then, according to strategy σ0, solution s k is such that f (s k , u, δ) ≤ f (s, u, δ) for all u ∈ U and all δ ∈ ∆ . The result is obtained by transitivity.

Algorithm IAS can be interrupted at any time to obtain the best solution found so far (the incumbent). Assume that IAS is interrupted at some step k, and let s k denote the incumbent at the end of this step. It is worth noting that the quality of s k can be measured by considering the following notion of maximum regret:

M R(s k , X , U k , ∆ k ) = max x∈X max u∈U k max δ∈∆ k {f (x, u, δ) -f (s k , u, δ)}
where U k , ∆ k and O k respectively denote U , ∆ and O at the end of step k. This maximum regret is indeed an upper bound on the gap to optimality: it represents the worst-case regret of choosing the incumbent instead of any other solution in terms of approval scores. This regret can be easily obtained using the available setwise max regrets SR({s k }, Sη, U k , ∆ k ) computed to evaluate nodes η ∈ O k during the search. More precisely, the following proposition holds: Proposition 6. At any iteration step k of IAS, we have:

M R(s k , X , U k , ∆ k ) = max η∈O k SR({s k }, Sη, U k , ∆ k ) Proof. We want to prove that maxx∈X R(s k , x, U k , ∆ k ) is equal to maxη∈O k SR({s k }, Sη, U k , ∆ k ), where R(s k , x, U k , ∆ k ) = maxu∈U k max δ∈∆ k {f (x, u, δ) -f (s k , u, δ)}. Let Y be the set of solutions defined by Y = {x ∈ Sη, η ∈ O k }.
Note that Y is not empty since Algorithm IAS has been interrupted. Then, by definition of setwise max regrets, we have:

SR({s k }, Sη, U k , ∆ k ) = max u∈U k max δ∈∆ k max x∈Sη {f (x, u, δ)-f (s k , u, δ)} = max x∈Sη max u∈U k max δ∈∆ k {f (x, u, δ)-f (s k , u, δ)} = max x∈Sη R(s k , x, U k , ∆ k ) for any η ∈ O k . Therefore, maxη∈O k SR({s k }, Sη, U k , ∆ k ) = maxx∈Y R(s k , x, U k , ∆ k ).
Moreover, by definition of the pruning rule, we know that SR({s k }, Sη, U k , ∆ k ) ≥ 0 for any η ∈ O k , and so maxx∈Y R(s k , x, U k , ∆ k ) ≥ 0. As a consequence, to establish the result, it is sufficient to prove that R(s k , x, U k , ∆ k ) ≤ 0 for any solution x ∈ X \Y . Three cases may occur:

• For x = s k , we can easily infer R(s k , x, U k , ∆ k ) = 0.
• For x ∈ S, since s k is the incumbent, we know that f (x, u, δ) ≤ f (s k , u, δ) for all u ∈ U k and δ ∈ ∆ k (by definition of strategy σ0). Therefore, we can infer R(s k , x, U k , ∆ k ) ≤ 0. • For x ∈ Sη, where η is a node that has been pruned during the search, we know that SR({s k }, Sη, U k , ∆ k ) < 0 (by definition of the pruning rule). Then, since we have just proved that

SR({s k }, Sη, U k , ∆ k ) = max x ∈Sη R(s k , x , U k , ∆ k ), we nec- essarily have R(s k , x, U k , ∆ k ) < 0. Hence maxx∈X R(s k , x, U k , ∆ k ) = maxx∈Y R(s k , x, U k , ∆ k ).
This establishes the result.

Since finding the knapsack that maximizes the approval score is NPhard, even when all functions ui and thresholds δ i are known (see Proposition 1), a Branch and Bound for this problem may obviously have bad running time (in the worst case, all solutions may be enumerated). Usually, the efficiency of a Branch and Bound method is evaluated in an empirical way, analysing instances with a large number of solutions. It is also evaluated by the quality of the returned solutions. We show below the practical efficiency of our algorithms.

Approval-based Search (AS)

In this subsection, we report numerical tests aiming at evaluating the computation times (given in seconds) of possible winners calculation using AS and ASε (linear optimizations are performed using the Gurobi solver). In these experiments, instances of the multi-agent knapsack problem with p = 12 are generated as follows: weights wj, j ∈ P, are uniformly drawn in {1, . . . , 100}. Then, capacity W is set to d × j∈P wj where d = 0.3, 0.4 and 0.5 so as to vary the number of solutions: the number of maximal (w.r.t. set inclusion) feasible subsets is approximatively equal to 500, 1000 and 2000 respectively. Moreover, to evaluate the impact of the uncertainty sets, we randomly generate q = 10, 20 preference statements per agent before running the algorithms. The computation times obtained by averaging over 50 runs for instances with n = 20 are reported in Table 1. We also report #S the average number of possible winners. As expected, we can see that computation times increase with the size of the problem, and decrease as the number of available preference statements increase. Moreover, we observe that AS0.1 is drastically faster than AS. More precisely, AS0.1 determines an 1.1approximation of possible winners in a few seconds while AS needs a few minutes on average.

Incremental Approval-based Search (IAS)

In this subsection, we report various numerical experiments aiming at evaluating the performance of IAS (presented in Section 4). In these experiments, only the preference ranking over single items is initially available for each agent. Then, answers to approval queries are simulated using approval thresholds and utility functions randomly generated with negative correlations so as to obtain difficult instances.

The first series of tests aims to evaluate the performance of its pruning rule in terms of average number of explored nodes (which corresponds to the average size of the search tree). For a basic comparison, we also report the number of nodes explored by exhaustive enumeration (this basic procedure is named S0 hereafter). Due to the possibility of collecting preference information during the search, algorithm IAS is significantly faster than AS and can solve much larger instances. For example, we report here the results obtained for p = 12, 15, 18 and 20 (with d = 0.4), which approximatively represents 1000, 8500, 60000 and 250000 maximal (w.r.t. set inclusion) feasible solutions respectively. Results obtained by averaging over 30 runs are reported in Table 2 shows that the average number of nodes explored by IAS remains quite low when increasing the number of solutions, whereas it drastically increases with an exhaustive enumeration. For example, our regret-based pruning rule reduces the average size of the search tree by a factor of 30 for instances with p = 12 (1000 feasible solutions), while reducing this number by a factor of 1600 for p = 20 (250000 feasible solutions).

The second series of experiments aims to evaluate the performance of IAS in terms of computation times (given in seconds) and average number of queries per agent (denoted by #Q hereafter). Results obtained by averaging over 50 runs are reported in Table 3 for instances involving 10, 20 and 30 agents. In Table 3, we can see that IAS is very efficient both in terms of number of queries per agent and computation times; for instance, for problems with 10 agents and 20 items (250000 feasible solutions), IAS determines the set of optimal knapsacks in less than 50 seconds with less than 15 queries per agent on average, whereas the full elicitation of approval statements would require 250000 queries per agent. Moreover, although the number of feasible solutions increases exponentially with p the number of items, we can see that the average number of queries increases very slowly. This shows that the compact numerical representation of individual preferences is well exploited by IAS during the resolution. Finally, when comparing Tables 1 and3, we observe that combining elicitation and search is more efficient than asking preference queries prior to the search; for example, with 20 agents, IAS determines the approval winners with no more than 15 queries per agent for problems with 250000 solutions, whereas asking 20 queries per agent prior to the search leaves us with 90 possible winners for smaller problems (2000 solutions) while doubling computation times. This shows that elicitation driven by resolution enables us to save both computation times and preference queries needed to determine the optimum.

p = 15 p = 18 p = 20 n time #Q time #Q time #Q 10
We focus now on the evaluation of the branching strategy of IAS. The branching strategy is of crucial importance for the efficiency of the query generation strategy proposed in Section 4. This elicitation strategy indeed consists in asking agents whether they approve or not some solutions found during the search; these solutions are obviously dependent on the branching strategy. For comparison, we also consider the interactive branch and bound procedure (named Random hereafter) that differs from IAS only on the branching strategy as follows: the next node to be explored and the next variable to be instantiated are both selected at random. For both branching strategies, we compute the maximum regret (MR) attached to the incumbent (using Proposition 6) each time an agent answers an approval query during the resolution. Recall that this regret gives an upper bound on the regret of choosing the incumbent instead of any other solution in terms of approval scores. Whenever this value equals zero, the incumbent is necessarily an approval winner. Regrets are here expressed on a normalized scale assigning value 1 to the initial maximum regret (computed before collecting any preference information) and value 0 when the maximum regret is 0. Figure 2 shows the results obtained by averaging over 30 runs. In Figure 2, we can see that the maximum regret reduces much more quickly with IAS than with Random. For instance, after 240 queries (i.e. 8 queries per agent) on average, the regret of choosing the incumbent instead of any true approval winner is under 10% of the initial regret, whereas it remains above 65% with the random branching strategy. Hence, the elicitation strategy seems to be much more informative using our regret-based branching strategy instead of the random branching strategy.

Finally, we estimate the performance of IAS as an anytime algorithm by computing the maximum regret attached to the incumbent at each iteration step of IAS. This regret indeed provides a guarantee on the gap to optimality at any step of the branch and bound procedure. For comparison, we consider here also the random branch and bound procedure (named Random). The results reported in Figure 3 are obtained by averaging over 30 runs. We observe that the maximum regret decreases much more quickly with IAS than with Random, as the number of iteration steps increases. For instance, after 200 iteration steps on average, the maximum regret is under 10% with our procedure while remaining above 70% for the random strategy. Moreover, we observe that this regret drops drastically from the very beginning of IAS. This shows that IAS provides (very quickly) a good solution before the end of the search.

CONCLUSION

We have presented a new approach for incremental approval voting on combinatorial domains, illustrated on the knapsack problem. The first specificity of this approach is to exploit compact numerical representations of individual preferences (additive utilities) to propose more efficient elicitation sequences. Thus, learning that agent i approves or not solution x is no longer an isolated preference information; it induces a constraint on the utility space possibly reducing the set of weak-orders consistent with the observed preferences. This contributes to derive implicitly other approval judgements on other knapsacks, which saves many preference queries.

The second specificity of our approach is to interleave preference elicitation and search. This makes it possible to elicit preferences on a combinatorial set implicitly defined with a twofold benefit in view of winner determination: on the one hand, working on partial instances of feasible solutions in the search tree facilitates the identification of relevant queries and relieves the elicitation burden. On the other hand, the search is earlier focused on the relevant part of the solution space due to the integration of new preference information at decisive steps of the search algorithm, which saves a significant part of the computational effort. This enables to solve large instances for which the systematic elicitation of all approval statements is not feasible.

We see several directions to extend this work. The first one is to relax the additivity of individual utilities so as to be able to model interactions between items. In this line, it seems natural to use generalized additive utilities functions (GAI) that could also be elicited incrementally [START_REF] Braziunas | Minimax regret based elicitation of generalized additive utilities[END_REF]. Another direction would be to consider alternative voting rules using possibly more information than approval statements. For example, positional scoring rules may be worth investigating under the assumption that individual preferences are representable by additive utilities. This is a challenging issue because, when preferences are only partially known, the ranges of possible ranks in individual rankings is generally too large to be decisive.

Algorithm 1 : 6 if Pη = ∅ then 7 S ← S ∪ Sη 8 else 9 13 O 17 O ← O \ {η} 18 end 19 SProposition 3 .

 167891317193 Approval-based SearchInput: S0: initial solutions; U, ∆: uncertainty sets Output: PW(X , U, ∆): the set of possible approval winners1 S ← S0 2 η ← [∅, ∅] 3 O ← {η} 4 while O = ∅ do 5 Select a node η in arg max η ∈O SR(S, S η , U, ∆) Select j ∈ Pη s.t. xj = 1 in the optimal solution of MIPη 10 Generate η 0 = [P 0 η ∪ {j}, P 1 η ] and η 1 = [P 0 η , P 1 η ∪ {j}] 11 forall η ∈ {η 0 , η 1 } do 12 if S η = ∅ and SR(S, S η , U, ∆) ≥ 0 then ← F ILT ER(S) 20 return SThis algorithm is justified by the following proposition: AS returns the set PW(X , U, ∆).

Figure 2 .

 2 Figure 2. Maximum regret attached to the incumbent with respect to the total number of queries (p = 12, d = 0.4, n = 30).

Figure 3 .

 3 Figure 3. Maximum regret attached to the incumbent at each step of the branch and bound procedures (p = 15, d = 0.5, n = 30).

Table 1 .

 1 Computations of possible winners.

			d = 0.3	d = 0.4	d = 0.5
	method	q	time	#S	time	#S	time	#S
	AS	10	26.6 19.8	71.8 29.8 614.4 111.4
	AS	20	25.2 17.4	55.7 23.9 442.3	90.3
	AS 0.1	10	1.1	1.9	2.7	3.8	13.9	8.3
	AS 0.1	20	1.0	1.6	2.3	2.7	9.4	5.0

Table 2 .

 2 Table 2 for problems involving 30 agents. Number of explored nodes.

	p	S0	IAS
	12	3024	109
	15	21468	238
	18	152415	289
	20	647834	440

Table 3 .

 3 Computations of necessary winners.

		9.2 10.5	29.8 12.2	48.7 14.3
	20	17.7	9.8	42.2 10.1 162.1 14.9
	30	29.7 10.9 154.0 11.8 225.7 13.1

ACKNOWLEDGEMENTS

We wish to thank Jérôme Lang for stimulating discussions during the preparation of this paper. This work is supported by the ANR project 14-CE24-0007-01-Cocorico-CoDec.


Stéphane Chrétien
email: stephane.chretien@univ-fcomte.fr

Franck Corset
email: franck.corset@univ-grenoble-alpes.fr

A lower bound on the expected optimal value of certain random linear programs and application to shortest paths in Directed Acyclic Graphs and reliability

Keywords: Random networks, Random linear programming, Reliability, Maintenance, Dyer-Frieze-McDiarmid's inequality, Shortest path 2010 MSC: 68R10, 52B05, 90C05, 90B25

The paper studies the expectation of the inspection time in complex aging systems, i.e. the length of the shortest path in a Directed Acyclic Graph, with random costs on egdes. We give a lower bound for this expectation.

Introduction and motivations

The random shortest path problem may be a good model for describing the time to failure of very complex systems with various degradation schemes as for instance nuclear plants. In this section, we describe our motivations for studying such random shortest path problems.

Problem statement

Consider a complex system whose n degradation states have been identified by experts. Let node 1 represent the state where the system is considered as new and let node n be the state of unacceptable degradation. All maximum paths from any node of the graph end at node n as in the figure below. The system is supposed to possibly evolve from a degradation state to any neighbor in the corresponding connected directed acyclic graph. The transition time between any two given states is assumed to follow a Weibull distribution whose parameters are estimated if the number of observations is sufficiently large.

Otherwise, it is possible to make Bayesian inference in order to combine the real data with some expert opinions.

Assume we start with a brand new system. Then, evolution of the system starts in state 1. Maintenance policies require that the system be inspected before reaching state n, i.e. unacceptable degradation. We represent this by a connected Directed Acyclic Graph (DAG) G = (V, E), where |V | = n and |E| = m. In the sequel, we assume that m ≥ n. Such examples of complex systems have been studied in [START_REF] Corset | Maintenance optimization from Bayesian networks and reliability with doubly censored data[END_REF]; [START_REF] Celeux | Designing a bayesian network for preventive maintenance from expert opinions in a rapid and reliable way[END_REF]; [START_REF] Adjerid | Performance evaluation and optimisation of industrial system in a dynamic maintenance[END_REF]. Moreover, [START_REF] Chen | Selective maintenance optimization for multi-state systems[END_REF] study the shortest path, in the maintenance optimization context, for some multi-state parallel-series systems. The problem posed in this paper is to provide a lower bound on acceptable inspection times.

Inspection times and shortest paths

In order to simplify the analysis, we assume that evolution inside the degradation graph, a Directed Acyclic Graph (DAG), proceeds following the rule that starting from one node i, the system goes to state j minimizing the transition time among neighbors of state i. Therefore, acceptable inspection times will be the times lower than the shortest path from state 1 to state n where each edge is weighted by its transition time. In general situations, we thus may ask for

• an estimator of the expected length of the shortest path from 1 to n,

• a confidence interval for the expected time path. This task is in general impossible to achieve because of the huge number of observations this should require in practice. The goal of this paper is to propose a lower bound on the expected length of the shortest path.

The paper is organized as follows. Since shortest path problems are well known to be representable as linear programs, we will address in the next section the more general problem of deriving a lower bound to the expectation of linear programs with random costs. In the third section, we specialize the study of this lower bound to an appropriate linear programming formulation of the shortest path problem. Moreover, we show that in the case of exponentially distributed random costs, the Dyer-Frieze-McDiarmid upper bound is as bad as possible.

The fourth section is devoted to the application to reliability theory as motivated by the introductory example above. In particular, the Weibull distribution is proved to satisfy the assumptions under which the proposed lower bound holds.

Random linear programs

Consider the linear program with random costs given by

z = min c t x Ax = b x ≥ 0 (1)
where c is a random vector with independent components taking values on R m + , A is a matrix in R n×m and b is a vector in R n . The expectation of c i is denoted by µ i and its variance by σ 2 i . Throughout the paper, we will assume that the components of c have a continuous density w.r.t. the Lebesgue measure.

In the sequel, we assume that A is full rank, i.e. rank(A) = n, and that the constraints of (1) define a polytope which is therefore a compact set. For any subset of {1, . . . , m}, we denote by A B the matrix whose column set is the set of columns of A indexed by B. We will also use the notation x B and c B for the vectors whose components are the components of x and c which are indexed by

B. A set of indices B is called a basis if its cardinality is n and the matrix A B is full rank. A basis is called feasible if x B = A -1 B b ≥ 0.
A feasible basis is said to be optimal if x * defined by

x * B = A -1 B b x * B c = 0
is an optimal solution of (1). Random linear programs have been investigated

recently and many impressive results have been optained in the case of i.i.d.

cost vectors. For instance, the assignment problem was investigated in Wästlund In this section, we propose a lower bound on the expected value of random linear programs in the spirit of the Dyer, Frieze and McDiarmid inequality, see [START_REF] Dyer | On linear programs with random costs[END_REF]. The Dyer-Frieze-McDiarmid inequality is a powerful tool for the analysis of some linear programming and combinatorial optimization problems with random costs, as detailed in the monograph of [START_REF] Michael | Probability theory and combinatorial optimization[END_REF].

More precisely, The Dyer-Frieze-McDiarmid bounds reads as follows.

Theorem 2.1. (Dyer-Frieze-McDiarmid) Assume that all the components of c in (1) are independent and nonnegative and there exists β ∈ (0, 1] such that

E [c i | c i ≥ h] ≥ E [c i ] + βh.
Let x be a feasible solution of (1). Let z * denote the random optimal value of (1). Then, assuming

E [c 1 ] x 1 ≥ • • • ≥ E [c m ] x m . E [z * ] ≤ β -1 m i=1 E [c i ] x i . (2)
The Weibull distributions W(η, γ) (η and γ are respectively the scale and shape parameters), has density function

f (x) = γ η x η γ-1 exp - x η γ . (3)
When the edges are Weibull distributed with shape parameters γ i in the interval [1, 2] for i = 1, . . . , n, we will see in Proposition 4.1

E[c i | c i ≥ h] ≤ E[c i ] + βh (4)
with β = 1. Note that Dyer-Frieze-McDiarmid requires the reverse inequality instead, in order to hold. We will however use this property to obtain a lower bound on the expectation of the optimal value of random linear programs in Theorem 2.1.

As in [START_REF] Dyer | On linear programs with random costs[END_REF], we will need the following result which is well known to users of the simplex algorithm.

Lemma 2.1. A necessary and sufficient condition for a feasible basis to be optimal is that

c t B c ≥ c t B (A t B ) -1 A B c . ( 5
)
Definition 2.1. For a basis B, let I B be the index set

I B = i ∈ {1, . . . , n} | c t B (A t B ) -1 A B c i ≥ 0 .
Using this result and following the same reasoning as in the proof of the Dyer, Frieze and McDiarmid inequality in [START_REF] Dyer | On linear programs with random costs[END_REF], we obtain the following proposition.

Proposition 2.1. Consider the random linear program (1) with random cost vector c having independent components and satisfying (4) with β ∈ [1, +∞).

Let B be a set of feasible bases. For each B ∈ B, p B will denote the probability that B is optimal. Let I B be the index set I B = ∩ B∈B I B . Let x be any vector satisfying the constraints of (1) and such that

x I c B = 0. (6)
Then, we have

E[z] ≥ 1 β B∈B p B E[c B ]x B . (7)
Proof. Fix a basis B and let E B be the event that B be optimal. Take any x satisfying the primal constraints. Then we have

E[z | E B ] = E[c t B (A B ) -1 b | E B ] = E[c t B (A B ) -1 (A B x B + A B c x B c ) | E B ] = E[c t B x B + c t B (A B ) -1 A B c x B c | E B ]. (8)
But using (4) together with (6), we have

E[c t B c | E B , c B ]x B c ≤ E[c t B c ]x B c + βc t B (A B ) -1 A B c x B c , and thus E[z | E B ] = E[c t B | E B ]x B + E E[c t B (A B ) -1 A B c | E B , c B ] | E B x B c ≥ E[c t B | E B ]x B + 1 β E[c t B c | E B ] -E[c t B c ] x B c = 1 -1 β E[c t B | E B ]x B + 1 β E[c t | E B ]x -1 β E[c t B c ]x B c . (9) Since β ∈ [1, +∞) we can rule out the term 1 -1 β E[c t B | E B ]x B and using the fact that E[c t B c ]x B c = E[c t ]x -E[c t B ]x B , we get E[z | E B ] ≥ 1 β E[c t | E B ]x - 1 β E[c t ]x -E[c t B ]x B . (10)
Finally take the expectation over all possible bases to obtain

E[z] ≥ 1 β B∈B p B E[c t | E B ]x -1 β B p B E[c t ]x + 1 β B∈B p B E[c t B ]x B = 1 β B p B E[c t B ]x B . (11)
The result of this proposition is not completely satisfactory since the probabilities p B that B be an optimal basis are not known. In certain cases, efficient approximations of these probabilities can be obtained using a more precise expression of p B . Since in the case where the components of the cost vector c are independent we easily get such an expression from the conditions for optimality given in Lemma 2.1. The lower bound we thus obtain is summarized in the following theorem.

Theorem 2.2. Consider the random linear program (1) with random cost vector c with independent coordinates. a. Let B be a basis for this program and for all j ∈ B and i ∈ B c , let

α ij = ((A t B ) -1 A B c ) ji .
Then, we have

p B = E[ i∈B c P (c i ≥ j∈B c j α ji | c B )]. (12)
b. Let x be any vector satisfying (6) and the constraints of (1). Then

E[z] ≥ 1 β B∈B E[ i∈B c P (c i ≥ j∈B c j α ji | c B )]E[c B ] t x B . (13)
Proof. a. Due to independence of the components of c, conditionally on the value of c j , j ∈ B, the events c i ≥ j∈B c j α ji are independent. Thus, the desired formula.

b. Combine a. with Proposition 2.1.

Remark 2.1. Condition (6) is easy to verify in practice by choosing a set B which is not too large and computing I B explicitely.

With these results in hand, we will now be able to turn to the more specialized case of random shortest paths in the next section.

Random shortest paths in connected DAGs

Linear programming formulation

The shortest path problem can be represented as an equivalent linear programming problem, as is well known [START_REF] Hoffman | A note on shortest path, assignment, and transportation problems[END_REF]. In (Papadimitriou and Steiglitz, 1998, pp. 75-79) for instance, the shortest path problem is shown to be equivalent to

min x∈R n c t 0 x subject to A 0 x = b and x ≥ 0, ( 14
)
where c 0 is the column vector whose components are the transition times on each edge, A 0 is the incidence matrix of the oriented degradation graph and b is the vector [-1, 0, • • • , 0, 1] t , encoding the fact that we start the path at node 1 and end it at node n. Recall that the incidence matrix is constructed as follows.

Its rows are indexed by the nodes of the graph while its columns are indexed by its edges with an extra column of all ones. In each column indexed by edge (i, j), set the i th component to -1, the j th component to 1 and set all other entries to zero. For instance, the incidence matrix for the graph of figure 1.1 is given by

A 0 =            -1 -1 -1 0 0 0 0 1 0 0 -1 0 0 0 0 1 0 1 -1 -1 0 0 0 1 0 1 0 -1 0 0 0 0 0 1 1           
.

Any solution vector x * to this linear program whose components are binary, i.e.

∈ {0, 1} encodes a path whose edges correspond to the nonzero components of

x * . The important property is that the matrix A 0 is totally unimodular (TUM) which means that every square submatrix has determinant equal to -1, 0 or 1.

This linear programming formulation of the problem has however a drawback: the incidence matrix is not full rank and the size of its kernel is the number of connected components of the graph, see [START_REF] Bollobás | Modern graph theory[END_REF]. On the other hand, for our results to apply recall that we need the matrix A in (1) to be full rank.

In order to remedy this problem, we introduce the extended incidence matrix A, given by

A = [A 0 | e],
where e is the vector whose components are all equal to one. For instance, the extended incidence matrix for the graph of figure 1.1 is given by

A =            -1 -1 -1 0 0 0 0 1 1 0 0 -1 0 0 0 1 0 1 0 1 -1 -1 0 1 0 0 1 0 1 0 -1 1 0 0 0 0 0 1 1 1            .
In addition, let c denote the extended cost vector [c t 0 , 0]. Then, we get the following proposition. (15)

Proof. Let x be an optimal solution of the given linear program. Then, x satisfies the Karush-Kuhn-Tucker equations which are of the form:

                                 c + I A t   u v   = 0, Ax = b x ≥ 0 u i x i = 0 for i = 1, . . . , n u i ≤ 0 v j (a t j x -b j ) = 0 for j = 1, . . . , m (16)
where a t j is the j th row of A and the vectors u ∈ R n and v ∈ R m are the Lagrange multipliers. More precisely, the multipliers that compose the vector u deal with the nonegativity constraints and the components of v deal with the others. The third equation is imposed in order to select the "active" constraints at optimality. In particular, it implies that if u n+1 = 0 we must have x n+1 = 0.

On the other hand, if u n+1 = 0, x n+1 may be positive. In what follows, we show that x n+1 is always null which will readily imply that this linear program also solves the shortest path problem.

Equations (16) determine a polyhedron in R n+1 × R n × R m . Now assume that [x t , u t , v t ] t is a corner point of this polyhedron and that u n+1 = 0. Since the matrix A is now full rank, the x-part of the corner vector satisfies Ax = b, and

x i = 0, i ∈ I x
for some index set I x of cardinality n + 1 -m. Now write these last nullity constraints Cx = 0 for some appropriate matrix C. Then by Cramer's rules, we obtain that the last coordinate x n+1 is proportional to

det   C A 0 | b   .
On the other hand, we know that the sum of the rows of A is equal to zero and the same holds for the sum of the components of b. Therefore, the determinant just above is null. Therefore x n+1 = 0 as announced. From this, it is easy to deduce that the vector of the first n components of an optimal solution to the present linear program also solves ( 14).

Using this proposition, we deduce that Theorem 2.2 applies to the random shortest path problem. We now consider the Dyer-Frieze-McDiarmid bound which gives an upper bound on the expected optimal value.

The DFM upper bound

The upper bound for the expected optimal cost of random linear programs in [START_REF] Dyer | On linear programs with random costs[END_REF] is a very nice result and a major contribution to the study of random optimization problems; see also [START_REF] Michael | Probability theory and combinatorial optimization[END_REF].

It is interesting to understand to what extent the Dyer-Frieze-McDiarmid

(DFM) bound is useful for the shortest path problem. Surprisingly, the answer is that the DFM bound is as bad as possible in this case, despite is remarkable efficiency on other standard combinatorial problems as shown in (Steele, 1997, Chapter 4). To understand why this happens, consider the deterministic problem where the random costs are replaced by their expected values.

ζ = min E[c] T x Ax = b x ≥ 0 (17)
Now, notice that whatever the distribution of the cost vector c may be, the following upper bound is immediate to obtain:

E[z] ≤ ζ. ( 18
)
The following proposition shows that the DFM bound is no better than this trivial upper bound. Proof. Take x equal to the vector minimizing ( 17). It is clear that x is a binary vector. It is also clear that the number of ones in this vector is less than the number of nodes in the graph. Then, the maximum value over all sets S of cardinality n in the right hand side term in ( 2) is obtained when S is taken to be the set of indices i for which

x i = 1. Thus i∈S E[c i ]x i is exactly the cost of x, i.e. ζ.
Thus, contrarily to intuition, replacing the random costs by their expected values is far from being a safe idea for the problem of providing efficient lower bounds to the mean inspection time.

Remark 3.1. Uniqueness in (17) is not restrictive since it occurs with probability one with continuous random costs.

Application to reliability

In this section, we address the problem of finding lower bounds to the inspection time of complex systems in reliability. As explained in Section 1 our main interest in studying random shortest paths problems relies in its possible application to the analysis of the time to failure for very complex systems. We will assume in this section that the transition times between two degradation states follows a Weibull distribution. In order to apply our previous results, we will need to study the Weibull distribution a little further.

Some properties of the Weibull distribution

Let X be a random variable with Weibull distribution W eib(η, γ), i.e. with probability density function given by

f X (t) = γ η t η γ e -( t η ) γ . ( 19
)
Then, the mean residual time to failure (MRTF) is given by

G X (h) = E[X | X ≥ h] = ηe ( h η ) γ Γ(1 + 1 γ , ( h η ) γ ), (20)
where Γ(a, h) is the incomplete gamma function defined by

Γ(a, h) = +∞ h t a-1 e -t dt. (21)
Lemma 4.1. Let X be a Weibull, denoted W eib(η, γ), distributed random variable. Then, a. the first two derivatives of the MRTF for a Weibull distributed variable X are given by

G X (h) = γ h η γ η h e ( h η ) γ Γ 1 + 1 γ , ( h η ) β -1 (22)
and

G X (h) = -γ 2 η γ h (γ-1) 1 + γ( h η ) γ + γ η (γ-1) h (γ-2) e ( h η ) γ Γ(1 + 1 γ , ( h η ) γ ) γ( h η ) γ + γ -1 . (23) Moreover, b. when γ ≥ 1 we have lim h→0 G X (h) = 0, lim h→+∞ G X (h) = 1 and ∀h ≥ 0, G X (h) < 1. ( 24
)
Proof. a. We omit the proof of the formula for the first and second derivative of G X .

b. Now, since G X is clearly continuous on R + and G X (0) = 0, we obtain the first assertion in b. Using the transformation u = t -(h/η) γ , the first of the two terms between parenthesis in (22) can we written ηe

( h η ) γ Γ 1 + 1 γ , ( h η ) γ = h +∞ 0 e -u 1 + ( η h ) γ u 1 γ du. (25)
For all u ≥ 0 we have the following Taylor expansion

(1 + ( η h ) γ u 1 γ = 1 + 1 γ η γ h γ u + o( 1 h γ ). ( 26
)
Multiplying by e -u and integrating, we obtain

+∞ 0 e -u 1 + ( η h ) γ u 1 γ du = +∞ 0 e -u du + 1 γ η γ h γ +∞ 0 e -u u du + o( 1 h γ ) = 1 + 1 γ η γ h γ + o( 1 h γ ). ( 27
)
Combining with ( 22), we deduce from a simple analysis of the second derivative that lim h→+∞ G X (h) = 1 and G X (h) ≤ 1 for all h ∈ R + .

The lower bound

In this subsection, we work out an easily computable lower bound derived from Theorem 2.2. We first have the following crucial result saying that the most commonly encountered Weibull distributions in reliability theory satisfy the main assumption of Proposition 2.1 and Theorem 2.2.

Proposition 4.1. Assume that X has distribution W eib(η, γ) with γ > 1.

Then, for all h ≥ 0, we have

E[X | X ≥ h] ≤ E[X] + h.
Proof. This is a direct consequence of lemma 4.1.

In the next theorem, we derive an explicit lower bound from Theorem 2.2 in the case of Weibull distributions.

Theorem 4.1. Consider the random linear program (1) with random cost vector c with independent components and assume that each component c i ,

1 ≤ i ≤ n follows a Weibull distribution W eib(η i , γ i )
a. Let B be a basis for this program and for all j ∈ B and i ∈ B c , let

α ij = ((A t B ) -1 A B c ) ji . Then, we have p B ≥ 1 - i∈B c 4 √ e j∈B α 2 ji η 2 j Γ 1 + 2 γj 1 2 + | j∈B α ji η j Γ 1 + 1 γj | η γi i . (28)
b. Let x be any vector satisfying (6) and the constraints of (1). Then

E[z] ≥ B∈B 1- i∈B c 4 √ e j∈B α 2 ji η 2 j Γ 1 + 2 γj 1 2 + | j∈B α ji η j Γ 1 + 1 γj | η γi i E[c B ] t x B . (29)
Proof. a. Theorem 2.2.a. gives the following formula for p B :

p B = E i∈B c exp -max 0, j∈B α ji c j η i γi , (30)
where the expectation is taken with respect to the variables c j , j ∈ B. Now since exp(-x) ≥ 1 -x, we obtain that

p B ≥ E 1 -i∈B c max 0, j∈B αjicj ηi γi , = 1 -i∈B c E max 0, j∈B αjicj ηi γi . (31)
Thus,

B ≥ 1 - i∈B c E max 0, j∈B α ji c j γi η γi i . (32)
In order to simplify the subsequent computations, we will use the crude majorization: max 0,

j∈B α ji c j ≤ | j∈B α ji c j | which gives p B ≥ 1 - i∈B c E | j∈B α ji c j | γi η γi i . (33)
Now, our next goal is to use the Kintchine inequalities in order to bound the last expression by a quantity expressed in terms of the l 2 norm which will be easier to control. For this purpose, one might want to center the random variables involved in (33) and use the triangle inequality to obtain

E | j∈B α ji c j | γi 1 γ i = E | j∈B α ji (c j -E[c j ] + E[c j ])| γi 1 γ i ≤ E | j∈B α ji (c j -E[c j ])| γi 1 γ i + E | j∈B α ji E[c j ]| γi 1 γ i = E | j∈B α ji (c j -E[c j ])| γi 1 γ i + | j∈B α ji E[c j ]| (34)
Using Jensen's inequality, a standard trick gives

E | j∈B α ji (c j -E[c j ])| γi 1 γ i = E | j∈B α ji (c j -E[c j ])| γi 1 γ i ≤ E | j∈B α ji (c j -c j )| γi 1 γ i (35)
where c j , j ∈ B are i.i.d. variables independent of c j , j ∈ B and such that c j has same distribution as c j , j ∈ B. Let j , j ∈ B be standard Rademacher ±1 random variables. Since j∈B α ji (c j -c j ) has the same distribution as j∈B α ji j (c j -c j ), we have

E | j∈B α ji (c j -E[c j ])| γi 1 γ i = E | j∈B α ji j (c j -c j )| γi 1 γ i ≤ 2E | j∈B α ji j c j | γi 1 γ i ( 36
)
where we used once again the triangle inequality. Notice that E | j∈B α ji j c j | γi (41)

1 γ i = E E | j∈B α ji j c j | γi | c j , j ∈ B 1 γ i = E E | j∈B α ji j c j | γi | c j , j ∈ B 1 γ i γi 1 γ i . (37)
Combining this result with (33), ( 34), ( 36), (37), and replacing E[c j ] = η j Γ 1 + 1 γj in (34), we finally obtain the desired result.

b. This follows from part a. and Proposition 2.1.

Conclusion and perspectives

In this paper, we derived a lower bound (Theorem 4.1) on the probability that a given path is optimal for the shortest path problem with independent arc weights with Weibull distributions. For this purpose, we used the linear programming formulation of the problem and extended the work of Dyer et al.

(1986) (Theorem 2.2).

The results presented here are of a theoretical nature. Further refinements and applications to real data will be proposed in a subsequent paper.

Figure 1 :

 1 Figure 1: Graph of degradation states with Weibull transitions.

(

  2005), Krokhmal et al. (2007), Krokhmal and Pardalos (2009) in the asymptotic regime.

Proposition 3. 1 .

 1 The shortest path problem is equivalent to the linear program min x∈R n c t x subject to Ax = b and x ≥ 0.

Proposition 3. 2 .

 2 Consider problem (1) where the random costs are assumed to be independent and exponentially distributed. Assume that the solution of (17) is unique. Then the DFM bound is equal to the optimum value ζ of the associated deterministic program (17).

  On the other hand, Khintchine's inequality givesE E | j∈B α ji j c j | γi | c j , j ∈ B where C γi is equal to √ 2eγ i in the present context. Thus, E E | j∈B α ji j c j | γi | c j , j ∈ B p → E[|X| p ] 1p is an increasing function and γ i is assumed to belong to [1, 2], we obtain the simpler boundE E | j∈B α ji j c j | γi | c j , j ∈ B j∈B α ji j c j | γi | c j , j ∈ B