Marco Casazza
email: marco.casazza@unimi.it

Alberto Ceselli
email: alberto.ceselli@unimi.it

Iterative exact algorithms for the Maximum Lifetime Tree in Wireless Sensors Networks

Keywords: sensor networks, energy conservation, data aggregation, spanning tree, exact algorithms

We tackle an optimization problem arising in the design of sensor networks: given a set of sensors, only one being connected to a backbone, to establish connection routes from each of them to the sink. Under a shortest path routing protocol, the set of connections form a spanning tree. Energy is required to transmit and receive data, and sensors have limited battery capacity: as soon as one sensor runs out of battery, a portion of the network is disconnected. We therefore search for the spanning tree maximizing the time elapsed before such a disconnection occurs, and therefore maintenance is required. We propose new mathematical formulations for the problem, proving and exploiting theoretical results on its combinatorial structure. On that basis, we design algorithms offering a-priori guarantees of global optimality. We undertake an extensive experimental campaign, showing our algorithms to outperform previous ones from the literature by orders of magnitude. We also identify which instance features have higher impact on network lifetime.

Introduction

A wireless sensor network is composed by a set of low-energy, battery-powered sensor nodes (sensors in the remainder). Sensor networks find many applications [START_REF] Tsiontsiou | Optimal probabilistic energy-aware routing for duty-cycled wireless sensor networks[END_REF]; for instance they are used to collect data from areas having no existing infrastructures, to a single edge node connected to a backbone (sink in the remainder).

Sensors connect one another by means of wireless links of limited range: if the data collection area is wide, the underlying communication network is sparse, preventing a direct connection of all sensors to the sink. Typical protocols from the literature assume each sensor to forward data with a shortest path philosophy, always to a single parent sensor [START_REF] Buragohain | Power aware routing for sensor databases[END_REF] [START_REF] Zhu | An exact algorithm for maximum lifetime data gathering tree without aggregation in wireless sensor networks[END_REF]. Therefore, for each sensor, a route (that is a sequence of sensors) connecting to the sink must be found, such that each sensor of the route has both its predecessor and its successor within the range of its wireless links, thereby allowing transmissions to be possible. The set of all routes origins a spanning tree rooted in the sink node.

A fundamental issue in wireless sensor networks is preserving energy: for a general survey on that matter we remind to [START_REF] Anastasi | Energy conservation in wireless sensor networks: A survey[END_REF]. In fact, sensors have clocks, making them to work in duty cycles: at each cycle the sensor receives data by the set of its children sensors, and transmits data to its parent. Since batteries replenishment is usually impossible or very expensive, the conservation of the battery energy is a critical issue. While sensors consume a negligible amount of energy when in standby, they spend most of the energy during both transmitting and receiving operations. It follows that duty control is critical especially for those sensors belonging to several routes. Researchers have proposed diverse techniques. A line of research aims at saving energy by selectively activating sensors over time (see e.g. [START_REF] Francesco Carrabs | Exact and heuristic approaches for the maximum lifetime problem in sensor networks with coverage and connectivity constraints[END_REF] and references therein). Other recent approaches [START_REF] Castaño | An exact approach to extend network lifetime in a general class of wireless sensor networks[END_REF] propose to adjust sensing ranges by scaling them up or down over time. Mathematical programming algorithms prove useful in many of these contexts.

While effective, specific techniques are not always applicable. For instance, sensing range scaling requires sensors with specialized capabilities, and turning off sensors is not feasible in all contexts (e.g. in smart homes for healthcare).

That is why another stream of research focuses on the management of network topology, and more specifically on a wise design of the routes between inner sensors and sink. A core combinatorial problem arises, that is finding the spanning tree maximizing the lifetime of the network, defined as the number of duty cycles that can be performed without making any sensor run out of battery.

In the literature such a problem is indeed known as Maximum Lifespan Tree Problem (MLTP). It has been studied in [START_REF] Buragohain | Power aware routing for sensor databases[END_REF], where the authors prove the NP-completeness of the problem and provide a routing algorithm with worst case constant factor performance guarantees. In [START_REF] Liang | A maximum lifetime algorithm for data gathering without aggregation in wireless sensor networks[END_REF] an approximation algorithm is proposed for a version of the problem where data is not aggregated during transmission. In [START_REF] Zhu | An exact algorithm for maximum lifetime data gathering tree without aggregation in wireless sensor networks[END_REF] and [START_REF] Ma | Exact algorithms for maximizing lifetime of wsns using integer linear programming[END_REF] two exact algorithms are presented to solve both the MLTP with and without the aggregation of the messages. However, both algorithms are strongly based on a technique that decomposes the network in independent and smaller networks when there exist links that, if removed, leave the network disconnected. It follows that their methodology is strongly dependent from the topology of a network.

Our contributions are the following. We propose new mathematical formulations for the problem, proving and exploiting theoretical results on the combinatorial structure of the MLTP (Section 2). On that basis, we design algorithms offering a-priori guarantees of global optimality (Section 3). We undertake an extensive experimental campaign on both datasets from the literature and new larger networks. We focus on both algorithms efficiency and on the evaluation of features which are more predictive of lifetime (Section 4). We finally collect some conclusions (Section 5).

Mathematical formulations

Let us be given the set of n sensors N = {1 . . . n}. We model the MLTP as a problem on a directed graph G = (N 0 , A), where N 0 = {0} ∪ N is a set of nodes representing the n sensors, together with a sink node 0, while A ⊆ N 0 × N 0 is a set of arcs, each representing a link between a pair of sensors within the range of communication.

For each node i ∈ N 0 we are also given a battery capacity b i > 0, that we assume infinite for the sink (b 0 = +∞). When a message is transmitted from a sensor to another, the transmitter consumes e t > 0 energy units of its battery, while the receiver consumes e r > 0. We assume that the payload of all the transmitted messages is of equal size at each duty cycle. No sensor can send or receive messages if its battery is empty, and we define as lifetime of a sensor node the number of duty cycles before its battery is fully drained.

A feasible solution to MLTP is a spanning tree T = {D 0 , D 1 , . . . , D n } rooted in node 0, where D i is the set of descendants of node i, with D 0 = N , D i = ∅ if i is a leaf, and D i = ∪ j ({j} ∪ D j ) for all children j of i otherwise. An example of a solution is depicted in Figure 1.

A feasible solution T is also optimal if it maximizes the lifetime of the tree l(T ), computed as the minimum lifetime over all nodes, that is

l(T ) = min i∈N b i (e t + e r ) • |D i | + e t . ( 1
) (a) (b)
Figure 1: Example of a solution to MLTP: in (a) a graph with 7 sensor nodes and a single sink (node 0) is given. Dotted lines represent the potential connections between sensor nodes. In (b) the set of the solid arcs represent a feasible solution defining the spanning tree that connects all the sensor nodes and the flow of data going from the sensor nodes to the sink. The MLTP can be formulated as an integer programming problem (IP ) as follows:

max min i∈N b i (e t + e r ) • y i + e t (2)
s.t.

(j,i)∈A

z ji = 1 if i = 0 0 otherwise ∀i ∈ N 0 (3)
x ij ≥ y j + 1 -|N | • (1 -z ij ) ∀(i, j) ∈ A (4)
y i = (i,j)∈A z ij ∀i ∈ N 0 (5)
z ij ∈ B ∀(i, j) ∈ A (6)
x ij ∈ N 0 ∀(i, j) ∈ A (7)
y i ∈ N 0 ∀i ∈ N 0 (8)
where variable z ij is set to 1 if arc (i, j) is selected in the tree, and 0 otherwise, variable x ij is the number of messages travelling arc (i, j) at each duty cycle, while variable y i is the number of descendants of node i in the tree, that is the number of messages collected by i at each duty cycle. The objective function (2) maximizes the minimum lifespan. Constraints (3) ensure that each node but the root has exactly one incoming arc. When z ij = 0, constraints (4) have no effect; when z ij = 1, they enforce x ij to take the number of messages going through arc (i, j). Constraints (5) force y i to be set to the number of descendants of node i. Subtours are avoided by the combination of ( 4) and [START_REF] Francesco Carrabs | Exact and heuristic approaches for the maximum lifetime problem in sensor networks with coverage and connectivity constraints[END_REF].

As stressed in [START_REF] Zhu | An exact algorithm for maximum lifetime data gathering tree without aggregation in wireless sensor networks[END_REF], such a formulation is nonlinear because the objective function include (a) a min operator, (b) a floor operator, and (c) variables at the denominator.

However, let us be given the lifetime of a tree l(T ) and the fractional lifetime of the same tree l(T ), computed as follows: l(T ) = min

i∈N b i (e t + e r ) • |D i | + e t . (9)
We observe that l(T ) = l(T ) , and that our objective function is

max ∀T l(T ) = max ∀T l(T ) . (10)
Therefore, our MLTP can be solved as the problem of maximizing the fractional lifetime, thus simplifying the mathematical programming model. Also, we can further observe that max ∀T l(T ) = min ∀T 1 l(T ) [START_REF] Dolan | Benchmarking optimization software with performance profiles[END_REF] and it follows that our objective function in Equation ( 1) is equivalent to

min ∀T max i∈N (e t + e r ) • |D i | + e t b i . ( 12
)
We therefore propose to obtain an equivalent linear formulation by introducing a continuous variable v that is the reciprocal of the lifetime, and the following set of constraints

v ≥ (e t + e r ) • y i + e t b i , ∀i ∈ N (13)
and by changing the objective function to min v. (14)

Flow formulation

We now propose an improved formulation based on network flows: let G = (N 0 , A ) be a directed graph that extends G with an additional terminal node t. In details, N 0 = N 0 ∪ {t} and A = {(j, i) : (i, j) ∈ A} ∪ {(i, t), ∀i ∈ N }.

Our flow formulation (F M ) reads as follows:

min v (15) s.t.
(i,j)∈A

x ij -

(j,i)∈A x ji =      |N |, if i = 0 -|N |, if i = t 0, otherwise ∀i ∈ N 0 (16) v ≥ ((e t + e r ) • ( (j,i)∈A x ji -1) + e t )/b i ∀i ∈ N (17) x ij ≤ 1, if j = t |N | • z ij , otherwise ∀(i, j) ∈ A (18) (j,i)∈A z ji ≤ 1 ∀i ∈ N 0 ( 19
)
x ij ≥ 0 ∀(i, j) ∈ A (20) z ij ∈ B ∀(i, j) ∈ A (21) v ≥ 0 (22)
where each variable x ij is the flow traversing arc (i, j) (i.e. the number of descendants of node j), and each variable z ij is set to 1 if any unit of flow traverses arc (i, j), and 0 otherwise. Variable v is the reciprocal of the lifetime. The objective function (15) still minimizes the reciprocal of the lifetime, therefore maximizing the lifetime. Constraints (16) ensure flow conservation. Constraints (17) force l to take the maximum reciprocal of the nodes lifetime. Constraints (18) impose that no flow can traverse (i, j) unless z ij = 1 (except for the terminal). Constraints (19) impose single source conditions.

Symmetry reduction

In order to improve our formulation, we study additional valid inequalities that help to reduce symmetries of feasible solutions.

First, it is easy to prove that:

Observation 1. For any feasible instance of MLTP it always exists an optimal solution where z 0i = 1 for all nodes i ∈ N such that (0, i) ∈ A .

In fact, connecting a sensor node i directly to the sink 0 avoids additional battery consumption on the parent of i, since 0 has an infinite battery. An example is provided in Figure 2.

Second, when sensor nodes can share a parent, we have that:

Observation 2. For any feasible instance of MLTP, given a clique of three nodes i, j, k ∈ N 0 such that (k, i), (k, j), (j, i) ∈ A , there always exists an optimal solution where

z kj + z ji ≤ 1. ( 23
)
Proof. Let us be given an optimal solution where z kj = z ji = 1, meaning that k is parent of j and j is parent of i in the spanning tree that is solution to MLTP. Such a solution implies that all the messages of i will pass though both j and k, consuming the batteries of both nodes. However, since arc (k, i) exists and k can be parent of i, we can always obtain a solution with z kj = z ki = 1 and z ji = 0, in such a way that the battery consumption of both nodes i and k is unaltered, while the one of node j is reduced. As none of the batteries consumption increase, the solution obtained is still optimal. An explanatory example is depicted in Figure 3. Observation 2 can be generalized into:

Observation 3. For any feasible instance of MLTP, given three nodes i, j, k ∈ N 0 , if there exists a path p = {(k, m 1 ), (m 1 , m 2 ), . . . , (m |p|-1 , j)} from k to j in G and (k, i), (j, i) ∈ A , there always exists an optimal solution where

z ji + (m,m )∈p z m,m ≤ |p|. ( 24
)
Proof. The proof follows directly from the prof of Observation 2: if all the arcs travelled in path p are selected in a solution, then connecting i to j would increase the battery consumption of all the nodes in path p, including node k. Instead, connecting i directly to k would decrease the battery consumption of all sensor nodes in p (excluding k).

An example is represented in Figure 4. Finally, let δ -(i) be the set of ingoing arcs in node i of graph G , we first prove that:

Observation 4. For any feasible instance of MLTP, given two sensor nodes i, j ∈ N sharing the same neighbourhood, that is

δ -(i) \ {j} = δ -(j) \ {i}, (25)
then there always exists an optimal solution where

z ij = z ji = 0. ( 26
)
Proof. In fact, if (25) is verified, then for any node k belonging to the shared neighbourhood of i and j, we can apply Observation 2. This means that given a solution where j is parent of i (resp. i is parent of j) we can transform such a solution into one where i is connected directly to the parent of j (resp. j is connected directly to the parent of i) reducing the battery consumption of j (resp. i). and6). In general, node 5 could have either the same parent as 4 or a different one in their neighbourhood, but we can always forbid to select arcs (4, 5) and [START_REF] Francesco Carrabs | Exact and heuristic approaches for the maximum lifetime problem in sensor networks with coverage and connectivity constraints[END_REF][START_REF] Anastasi | Energy conservation in wireless sensor networks: A survey[END_REF] in a solution as in (b). This restriction is stronger than the inequality of Observation 2.

An example is depicted in Figure 5. Furthermore, we can prove that:

Observation 5. For any feasible instance of MLTP, given an arc (i, j) ∈ A and a node k ∈ δ -(i) \ (δ -(j) ∪ {j}), there always exists an optimal solution where either

z ij = 0 or z ij = z ki = 1.
In other words, if the a node i is parent of node j, then the parent of i must be a node outside the neighbourhood of j, otherwise we know from Observation 2 that we can reduce the battery consumption on i by connecting j to k directly. Therefore, it follows that the following inequalities

z ij ≤ k∈δ -(i)\δ -(j)\{j} z ki , ∀(i, j) ∈ A (27)
are valid for MLTP. An example is represented in Figure 6.

Exact algorithm

In this section we propose three approaches that make use of model (FM) to solve MLTP with optimality guarantees.

Our approaches are inspired by the one in [START_REF] Ma | Exact algorithms for maximizing lifetime of wsns using integer linear programming[END_REF], where the authors perform a binary search over all the feasible lifetimes and check whether their model if each lifetime allows a feasible solution.

However, our approaches differ in four aspects: first, to solve feasibility subproblems, we use the mathematical programming formulation (FM) which is stronger than those from the literature; second, we integrate in the formulation our inequalities to further reduce the search space; third, we use the value of the lifetime set at each iteration to tighten the value of the coefficients in constraints (18), thereby incrementally strengthening the formulation; fourth, we design an enhanced search strategy.

First, let us rewrite constraints (18) as follows:

x ij ≤ 1, if j = t (d j + 1) • z ij , otherwise , ∀(i, j) ∈ A (28)
where d j is an upper bound to the number of descendants of node j that is initially set to n in constraints (18). Let d = {d 1 , d 2 , . . . , d n } be the vector of the descendant upper bounds used as big-M coefficients in Equation (28), and let F M d be the formulation F M using those values in d as coefficients in Equation (28). We remark that:

Remark 1. Given two flow models F M d and F M d such that d ≤ d , then F M d is a relaxation of F M d .
In fact, any feasible solution to F M d is feasible to F M d , but the opposite may be false, and an optimal solution to F M d is at least as good as an optimal solution to F M d . Now, let µ(i, l(T )) be the function computing the maximum number of descendants that a node i can have to allow a lifetime l(T ), that is:

µ(i, l(T )) = b i /l(T ) -e t e t + e r . ( 29
)
We prove that:

Theorem 1. If we are given a solution T * d , which is optimal for model F M d when setting

d i = µ(i, l(T * d )), (30)
then T * d is optimal for any other model F M d with d < d .

Proof. By way of contradiction, let us assume that it exists an optimal solution

T * d to F M d that is better than T * d , that is l(T * d ) > l(T * d ).
Then, for each node i, T * d would allow a number of descendants that is at most the one of

T * d . In fact, if l(T * d ) > l(T * d ) we have that b i /l(T * d ) -e t e t + e r < b i /l(T * d ) -e t e t + e r , ∀i ∈ N (31) b i /l(T * d ) -e t e t + e r ≤ b i /l(T * d ) -e t e t + e r , ∀i ∈ N (32) µ(i, l(T * d )) ≤ µ(i, l(T * d )), ∀i ∈ N. (33)
However, we also know that µ(i, l(T

d i = µ(i, l(T * d )), ( 34
)
then T * d is optimal for FM and therefore for MLTP. Also, we observe that: Observation 6. If we are given two lifetime values L and L , with L < L , and two vectors d and d such that

d i = µ(i, L), d i = µ(i, L ) (35)
then if F M d is feasible and F M d is infeasible, an optimal solution T * must exists such that L > l(T * ) ≥ L.

And therefore, we can prove that:

Theorem 2. If we are given a lifetime L and two vectors d and d such that

d i = µ(i, L -1), d i = µ(i, L ) (36)
then if F M d is feasible and F M d is infeasible, any feasible solution for F M d is also optimal for MLTP.

Proof. First, if F M d is infeasible, any other F M d with d < d is infeasible too. Therefore, is no solution can be found with lifetime L , no solution can be found with lifetime higher than L , since that would mean to decrease the number of descendants for at least one node. Because of that, if a solution with lifetime L -1 exists, then such solution must be also optimal. However, any feasible solution for F M d provides a lifetime that is at least L -1, because d values are computed starting from L -1. Thus, any solution that is feasible for F M d is also optimal for F M d and because of Corollary 1 is also optimal for MLTP.

Sequential search

In our first approach we exploit Theorem 2 to design the following iterative algorithm.

Let us be given a sequence L of all the nodes lifetime, that is

L = { b i /(e t + (e t + e r ) • k) : ∀i ∈ N, k = 0 . . . n} (37)
and let us assume that L is ordered from the biggest lifetime to the lowest. We have that

L = {L 1 , L 2 , . . . , L, . . . , L * 1 , . . . } (38)
where L is an upper bound to the value of the optimal lifetime L * . Our algorithm works as follows: we compute an upper bound L, for example by solving the continuous relaxation of F M , and then we compute a vector d such that d i = µ(i, L). We then solve F M d : if a feasible solution can be found, such a solution is also optimal to MLTP because of Theorem 2, as no solution with higher lifetime than L exists, and the algorithm ends. Otherwise lifetime L is too optimistic and we enter a search process that will lead to a new upper bound L, smaller than L.

Such an upper bound is found by moving to a new value in L such that L > L ≥ L * . When moving to a new lifetime, at least one value of vector d must increase. Therefore, we select L such that

L = max i∈N b i e t + (e t + e r ) • (d i + 1) , (39)
that is the best lifetime we can obtain if we allow one more descendant for at least one node. In such a way we never skip a lifetime value during our search. Then, we solve again model F M d with d i = µ(i, L) and start again the process until the model is feasible. The pseudocode of our Sequential Search Algorithm (SSA) is sketched in Algorithm 1.

Algorithm 1 Sequential Search Algorithm L ← solve continuous relaxation of F M do d i ← µ(i, L) for all i ∈ N L * ← solve F M d until a feasible solution is found if problem is feasible then return L * is optimal else L ← max i∈N bi et+(et+er)•(di+1)
end if while problem is infeasible

Binary search

In our second approach we perform a binary search in the sequence of L values as in [START_REF] Ma | Exact algorithms for maximizing lifetime of wsns using integer linear programming[END_REF] but using our FM formulation to test the feasibility and starting our search from different bounds. First, let us be given our sequence L such that

L = {L 1 , L 2 , . . . , L, . . . , L * 1 , . . . , L, . . . } (40)
where L and L are upper and lower bounds to the value of the optimal lifetime L * , respectively. At each iteration of the algorithm we compute L = (L + L)/2 and consequently compute vector d as d i = µ(i, L). We then verify if formulation F M d allows a feasible solution: we then repeat the process after either changing the value of L or L accordingly with the result of our validation.

The pseudocode of our Binary Search Algorithm (BSA) is sketched in Algorithm 2.

Algorithm 2 Binary Search Algorithm

L ← solve continuous relaxation of F M L ← find a primal bound loop L ← (L + L)/2 d i ← µ(i, L) for all i ∈ N L * ← solve F M d until a feasible solution is found if problem is feasible then if L = L then return L * is optimal end if L ← L else L ← L -1 end if end loop
We remark that both an upper and a lower bound to L * could be found, for example, by solving F M to optimality by means of generic branch-and-bound, but stopping at the root node of the branching tree.

Hybrid Search

In our last approach we alternate the Sequential Search algorithm and the Binary Search algorithm in such a way that a variant of the first is run to find a feasible solution whose value is near the upper bound, and then the latter is run to close the gap.

The Hybrid Search works as follows: first we find an upper bound L, set vector d and solve F M d as in the Sequential Search algorithm. Once the model is solved, if a feasible solution is found then it is also optimal. Otherwise we look for a lower bound L: similarly to the Sequential Search Algorithm we search L by allowing more descendants than the ones tested at the last iteration. Given an integer positive parameter ∆, we compute L as:

L = max i∈N b i e t + (e t + e r ) • (d i + ∆) . (41)
and we iterate. When a feasible solution is found, it is used to set values L and L for initializing the Binary Search Algorithm. The pseudocode of our Hybrid Search Algorithm (HSA) is sketched in Algorithm 3.

Algorithm 3 Hybrid Search Algorithm L ← solve continuous relaxation of F M L ← L search ← SEQU EN T IAL loop if search = SEQU EN T IAL then d i ← µ(i, L) for all i ∈ N L * ← solve F M d until a feasible solution is found if problem is feasible then search ← BIN ARY else L ← L -1 L ← max i∈N bi et+(et+er)•(di+∆) end if else L ← (L + L)/2 d i ← µ(i, L) for all i ∈ N L * ← solve F M d until a feasible solution is found if problem is feasible then if L == L then return L * is optimal end if L ← L else L ← L -1 end if end if end loop
We remark that if ∆ = 1 the Hybrid Algorithm actually performs a Sequential Search. After preliminary results we found to be profitable to set ∆ = 4.

Experimental analysis and conclusions

We implemented our models and algorithms in Python 3, and we used the branch-and-cut algorithms of IBM ILOG CPLEX [START_REF] Development | IBM ILOG CPLEX optimization studio: CPLEX user's manual -version 12 release 6[END_REF] version 12.6.3 to solve all the mathematical programming formulations.

As a benchmark we considered two different datasets: the first dataset, called ZTD, is kindly provided by the authors of [START_REF] Zhu | A Branch and Bound Algorithm for Building Optimal Data Gathering Tree in Wireless Sensor Networks -Computational Results[END_REF], and it consists of 300 instances having number of nodes n that goes from 30 to 100 in steps of 5. Nodes are dispatched in a 100 × 100 square, and each of them has a battery capacity b i that is randomly drawn in the range [START_REF] Tsiontsiou | Optimal probabilistic energy-aware routing for duty-cycled wireless sensor networks[END_REF][START_REF] Zhu | A Branch and Bound Algorithm for Building Optimal Data Gathering Tree in Wireless Sensor Networks -Computational Results[END_REF]. The energy consumption e t and e r are set to 6.66e -4 and 3.33e -4, respectively, and the communication radius of two nodes is set to 20. We further produced a second, more challenging, dataset, called CCD. That consists of 540 instances with number of nodes n ∈ {50, 75, 100}; similarly to dataset ZTD the nodes of each instance of CCD are distributed in a 100 × 100 square, but their battery capacities are randomly drawn in the range [γ, 10], with γ = 0, 5, 8, 9, 9.5, 10. Furthermore, we set e t = 6.66e -4 and e r = 3.33e -4 as in ZTD, but we increased the radius of communication to 25, 30, and 35. This settings led to more dense graphs, which in turn yield more difficult instances.

A summary of the features of the two datasets is reported in Table 2.

All tests have been conducted on a PC equipped with AMD Ryzen 1950X 3.4GHz CPU and 32 GB of memory. In CPLEX we kept the default settings although we set single thread execution. A time limit of 600s of computing time has been set for each run.

Comparison of models

At first we evaluate the performance of our formulations IP , F M , and F M with the additional inequalities described in Subsection 2.2 that we call F M I.

In Figure 7 we show the average computing time required by the three formulations to solve MLTP with optimality guarantees. The results show that the time required by F M and F M I is negligible for instances with up to 60 nodes, while IP struggles even for instances of limited size. We can also observe that our inequalities have a substantial impact on the performance of our models: not only the average computing times of F M I are the smallest, but F M I is the only model that solves all the instances with up to 85 nodes.

Figure 8 shows the number of instances solved to optimality by the three models. We observe that at 45 nodes the IP formulation solves only half of the instances, and that such a number decreases rapidly as the size of the instances increases. F M and F M I instead always solve more than half of the instances for each size of the network; in particular, F M I succeeds in solving all instances but 5, whose number of nodes is higher than 85.

Finally, Figure 9 shows the performance profile [START_REF] Dolan | Benchmarking optimization software with performance profiles[END_REF] of IP , F M and F M I. In details, for each hypothetical running time limit value (x axis), the charts report the fraction of instances in the dataset solved to proven optimality within that time limit (y axis). That is, the steepest increasing lines, the better. It is easy to observe that the gap between the two flow based formulations and IP is substantial.

A more detailed summary of the performance of these three models is reported in Table 3. Each row of such table corresponds at a different size of the networks, for which we report the average number of arcs, and the number of instances solved to optimality and the average computing times of the three formulations. However, in Table 3 we only report the average computing time of those instances that are solved to optimality only. A "-" symbol in the time column signal that the model of the corresponding column could not solve any instance.

Comparison against literature

In a second run of experiments we compared the results obtained by our F M I model against the state of the art methods from the literature. In particular, we considered as a benchmark the best results of improved versions of [START_REF] Ma | Exact algorithms for maximizing lifetime of wsns using integer linear programming[END_REF] and [START_REF] Zhu | An exact algorithm for maximum lifetime data gathering tree without aggregation in wireless sensor networks[END_REF]. Although the details of their new algorithms have not been published (nor disclosed) yet, the authors kindly sent us the link to their repository where their latest computational results can be found [START_REF] Zhu | A Branch and Bound Algorithm for Building Optimal Data Gathering Tree in Wireless Sensor Networks -Computational Results[END_REF]. These results are available on instances with up to 85 nodes only.

In our test for each instance of dataset ZTD we always considered the best result among the ones obtained and reported in their repository by any of their approaches. In the following we will refer to their best result as Literature Best (LiB).

In Figure 10 we report the average computing time of both our F M I and the benchmark LiB when the size of the network increases. We observe that computing times are negligible for instances with up to 50 nodes. Even for instances of bigger size, the average computing time of F M I never grows over 15s, while LiB struggles to solve instances with 80 and 85 nodes. In fact, we observe from Figure 11 that the number of instances solved to optimality by LiB is almost halved when the size of the network reaches 80 nodes.

In Table 4 we show a detailed summary of the number of instances solved to optimality by F M I and LiB, and the average computing time. Similarly to Table 3 the average time reported is computed considering only those instances solved to optimality: that is why the computing time of LiB looks smaller than the one reported in Figure 10. Restricting to the instances solved to proven optimality by both F M I and LiB, we found F M I to be always faster, except on small instances which could be solved in less than 1s by both methods.

The performance profile in Figure 12 show the substantial gap that exists between the two approaches.

Iterative exact algorithms

In a third run of tests we compared the performance of our iterative algorithms SSA, BSA, and HSA when solving instances of dataset ZTD.

From Figure 13 we observe that the three algorithms solved all of the instances of dataset ZTD within the time limit. The average computing times were smaller than 10 seconds. Also, each instance was solved by BSA and HSA in less than 10s: SSA was the only algorithm that took more than one minute to solve an instance.

From Figure 14 BSA seems to dominate the other two algorithms, although in Table 5, where we report the average number of iterations, the number of instances solved to optimality, and the average computing time, we observe that the average computing time of HSA is usually close and sometimes better than BSA. Also, the average numbers of iterations are close meaning that the main difference between the performance of the two approaches is the time required by CPLEX either to detect the infeasibility of the problem or to find a feasible solution.

SSA performances seem to be far from the the ones of BSA and HSA, as the former requires many more iterations to solve MLTP.

Increasing graph density

In a fourth run of tests we compared the performance of our algorithms SSA, BSA, and HSA on dataset CCD. Even though such a dataset was built to be more challenging, our algorithms were still able to solve all the instances, although requiring higher computing times.

In Figure 15 we observe that the average computing time required to solve instances in dataset CCD is higher than the one required to solve those of ZTD, even though it remains smaller than one minute on average. The performance profile in Figure 16 shows that on dataset CCD the three algorithms behave more similarly than they do on dataset ZTD.

In Table 6 we observe that even though the average computing time increases, the average number of iterations remains similar to the one obtained on dataset ZTD. This behaviour suggests that dataset CCD is more challenging because it is harder for CPLEX to detect whether an instance is feasible or not. Therefore, the performance of our algorithms might be further improved by the use of ad hoc procedures to detect the infeasibility or fast heuristics to rapidly find primal solutions.

Still concerning the computing time, although SSA generally performs worse than its two competitors, we observe that for instances generated with γ either very small or very large, SSA performs better. Those cases might be the ones where the dual bound provided by F M I is very close to the optimal value, and therefore the number of iterations required by SSA is very small. In order to investigate such a phenomenon we present, in Figures 17,18, and 19, scatter plots. Each point is an instance; its number of arcs (resp. sensor battery standar deviation) define the y coordinate (resp. x coordinate) value. The size of the point is proportional to the computing time required to solve the corresponding instance. In general, we observe that the distribution of the battery capacities influences the performance of our algorithms more than the connectivity of the networks. For instance, in Figure 17 we observe that when the standard deviation of the battery is 0, meaning that all sensors have the same battery capacity, SSA performs well. The same holds even when the standard deviation of the battery is large.

Lifetime

Finally, we study how sensor network features like size, connectivity, and battery capacities influence the lifetime of the network. First, in Figures 20, 21 and 22 we plot one point for each instance: their y coordinates are number of nodes, arcs and average node degree, respectively; their x coordinates are always average sensor battery; their size is the optimal lifetime value.

The size of a network does not significantly impact on its lifetime, which instead is highly dependent on the capacity of the batteries, as one could expect. However, Figure 22 shows that there is at least one other feature of the network that has an impact on the lifetime, that is the average degree of the sensor nodes. In fact we observe that when the average degree of the nodes increases, higher lifetime values are more frequent. A dense graph with an average degree of 15 can provide a lifetime that is higher than a sparse graph, even when the latter has 30% more battery capacity.

Conclusions

The MLTP is known in the literature for both its practical relevance and its theoretical interest. Our investigation on mathematical modelling and combinatorial properties helps in understanding its structure, clarifying key features of optimal solutions, and allowing for better resolution techniques.

In fact we could exploit these new results, and in particular the flow structure underlying the MLTP and the existence of symmetries between solutions, to design effective algorithms. Besides a strong mathematical programming formulation, inequalities reducing these symmetries, and the iterative tightening of problem coefficients, we found the search strategy to be the factor needed to complement structural results in a computationally effective framework. Through an exhaustive experimental campaign we showed that our models are more effective than previous ones from the literature, and that our additional inequalities have a substantial impact.

Overall, our algorithms experimentally outperform earlier ones from the literature by two orders of magnitude, allowing to optimize (and prove optimality) of several open instances.

We studied the correlation between network features and the lifetime that can be achieved, showing that the degree of the nodes is at least as important as the battery capacities to achieve a higher network lifetime.

We could also enlarge the datasets from the literature, including much larger networks, still obtaining fast computations, thereby proving our methods to be suitable as a basis for tools tackling real sized networks.

Figure 2 :Figure 3 :

 23 Figure 2: Example of Observation 1: in (a) sensor node 3 is a child of sensor node 1 but it could be connected directly to the sink 0. Connecting sensor node 3 to the sink as done in (b) reduces the number of descendants of sensor node 1.

Figure 4 :

 4 Figure4: Example of Observation 4: in (a) a path connects node 7 to node 1. However, node 7 can be connected directly to node 1 reducing the number of descendants for nodes 6, 5, 4, and 3 as depicted in (b). Since node 7 was already a descendant of node 1, the lifetime of node 1 is unaltered.

Figure 5 :

 5 Figure 5: Example of Observation 3: in (a) nodes 4 and 5 share the same neighbourhood (nodes3 and 6). In general, node 5 could have either the same parent as 4 or a different one in their neighbourhood, but we can always forbid to select arcs (4, 5) and[START_REF] Francesco Carrabs | Exact and heuristic approaches for the maximum lifetime problem in sensor networks with coverage and connectivity constraints[END_REF][START_REF] Anastasi | Energy conservation in wireless sensor networks: A survey[END_REF] in a solution as in (b). This restriction is stronger than the inequality of Observation 2.

Figure 6 :

 6 Figure 6: Observation 5 forbids solutions as in (a) imposing that node 4 has a parent outside the neighbourhood of 5 when arc (4, 5) is selected. Although Observation 5 and Observation 2 are equivalent when integrality conditions are satisfied, they are complementary while solving the continuous relaxation of MLTP.

Figure 7 :

 7 Figure 7: Average computing time of formulations IP , F M , and F M I on dataset ZTD with a time limit of 600s.

Figure 8 :Figure 9 :

 89 Figure 8: Number of instances solved to optimality by formulations IP , F M , and F M I on dataset ZTD with a time limit of 600s.

Figure 10 :Figure 11 :

 1011 Figure 10: Average computing time of LiB and F M I on dataset ZTD with a time limit of 600s.

Figure 12 :Figure 13 :

 1213 Figure 12: Performance profile of model F M I and the best result of LiB on dataset ZTD with a time limit of 600s.

Figure 14 :

 14 Figure 14: Performance profile of algorithms SSA, BSA and HSA on dataset ZTD with a time limit of 600s.

Figure 15 :

 15 Figure 15: Average computing time of algorithms SSA, BSA and HSA on dataset CCD with a time limit of 600s.

Figure 16 :Figure 17 :Figure 18 :Figure 19 :

 16171819 Figure16: Performance profile of algorithms SSA, BSA and HSA on dataset CCD with a time limit of 600s.

Figure 20 :

 20 Figure 20: Lifetime of networks of varying nodes and battery capacity.

Figure 21 :

 21 Figure 21: Lifetime of networks of varying arcs and battery capacity.

Table 1 :

 1 Mathematical notation.

	directed graph

N = {1, . . . , n} set of sensor nodes N 0 = N ∪ {0} set of nodes including the sink node

A ⊆ N 0 × N 0 set of arcs G = (N 0 , A)

  * d )) ≤ d i for all nodes i because of constraints (28) of model F M d . Therefore, the solution T * d would allow a number of descendants for each node i that is feasible also for F M d , and therefore T Corollary 1. If we are given a solution T * d , which is optimal with respect to model F M d when setting

	It follows that:

* d would be feasible for F M d , improving T * d for problem F M d , thereby contradicting our initial assumption.

Table 2 :

 2 Summary of datasets parameters.

	Dataset Nodes	Batteries		Transmission	Reception	en-	Communication
						energy	ergy	radius
	ZTD	30 to 100 in	in [1, 10]		6.66e -4	3.33e -4	20
		steps of 5
	CCD	50, 75, and 100	in		[γ, 10]	6.66e -4	3.33e -4	25, 30, 35
			with	γ	=
			0, 5, 8, 9, 9.5, 10

Table 3 :

 3 Summary of number of instances solved to optimality and average computing time of solved instances for formulations IP , F M , and F M I on dataset ZTD with a time limit of 600s.

Table 4 :

 4 Number of instances solved to optimality and average computing time of the best LiB approach and F M I model on dataset ZTD with a time limit of 600s.

	n avg. |A|	FMI		LiB
			#sol. time (s) #sol. time (s)
	30	105.8	20	0.5	20	0.0
	35	131.6	20	0.5	20	0.0
	40	170.9	20	0.5	20	0.0
	45	218.0	20	0.6	20	0.6
	50	266.0	20	0.6	20	0.2
	55	315.2	20	0.9	20	4.5
	60	375.6	20	1.0	20	16.0
	65	459.5	20	1.2	18	4.9
	70	507.5	20	1.4	18	23.2
	75	588.4	20	1.6	19	12.6
	80	659.6	20	14.7	12	32.5
	85	756.3	20	3.7	13	28.8

Table 5 :

 5 Average number of iterations, number of instances solved to optimality, and average computing time of algorithms BSA, HSA, and SSA on dataset ZTD with a time limit of 600s.

	n	γ avg. |A|		BSA			HSA		SSA
				#iter. #sol. time (s) #iter. #sol. time (s) #iter. #sol. time (s)
	51	0.0	541.3	1.7	30	0.4	1.4	30	0.6	2.3	0.4
	51	5.0	541.3	5.0	30	0.9	5.9	30	1.9	33.9	2.8
	51	8.0	541.3	5.8	30	1.0	6.0	30	2.0	36.4	3.1
	51	9.0	541.3	6.4	30	1.0	6.2	30	2.0	36.7	3.2
	51	9.5	541.3	7.3	30	1.1	7.1	30	2.1	33.0	2.9
	51 10.0	541.3	6.6	30	1.0	7.3	30	1.8	2.3	0.6
	76	0.0	1212.0	1.8	30	1.1	1.6	30	1.1	3.5	1.2
	76	5.0	1212.0	5.4	30	10.2	5.5	30	14.1	42.2	16.9
	76	8.0	1212.0	5.6	30	10.8	5.6	30	9.7	44.7	32.4
	76	9.0	1212.0	6.2	30	8.6	6.3	30	11.5	46.0	35.5
	76	9.5	1212.0	6.3	30	10.7	6.4	30	12.9	38.8	29.3
	76 10.0	1212.0	4.4	30	10.3	4.6	30	12.4	2.0	2.9
	101	0.0	2134.6	1.2	30	1.5	1.1	30	1.5	1.0	1.4
	101	5.0	2134.6	5.7	30	11.4	5.9	30	10.4	50.9	17.3
	101	8.0	2134.6	5.9	30	21.8	6.1	30	18.8	57.0	30.4
	101	9.0	2134.6	5.9	30	13.6	6.6	30	14.9	53.2	26.7
	101	9.5	2134.6	7.0	30	17.5	6.8	30	16.0	43.8	21.2
	101 10.0	2134.6	3.4	30	4.0	3.8	30	4.5	2.0	3.0

Table 6 :

 6 Number of instances solved to optimality, average number of iterations, and average computing time of algorithms BSA, HSA, and SSA on dataset CCD with a time limit of 600s.

	20

Acknowledgments

Partially funded by Regione Lombardia, grant agreement n. E97F17000000009, Project AD-COM, and "Piano Sostegno alla Ricerca 2019-20", Università degli Studi di Milano.

We are grateful to the authors of [10] for sharing their experimental data and computational results.


Yohan Reynaud
email: yohan.reynaud@improve-innov.com

Angélique Couvent
email: angelique.couvent@improve-innov.com

Aline Manach
email: aline.manach@agroparistech.fr

David Forest
email: david.forest@inrae.fr

Michel Lopez
email: michel.lopez@improve-innov.com

Daniel Picque
email: daniel.picque@inrae.fr

Isabelle Souchon
email: isabelle.souchon@inrae.fr

Didier Remond
email: didier.remond@inrae.fr

Didier Dupont
email: didier.dupont@inrae.fr

Didier Rémond

Food-dependent set-up of the DiDGI® dynamic in vitro system: Correlation with the porcine model for protein digestion of soya-based food

Keywords:

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Animal experiments are increasingly difficult to perform due to the economic resources they require and the ethical concerns that must be addressed. Despite this, it is still standard-indeed, recommended-for assessments of protein and amino acid digestibility to be conducted using human or pig models (FAO, 2013). In vitro systems could present suitable alternatives, but they first must be shown to mimic as much as possible real physiological conditions and must be validated using in vivo data. For example, Rozan et al. (1997) examined quantitative correlations between digestion end-points obtained in vivo and in vitro, specifically by comparing true digestibility assessed in rats with the final degree of hydrolysis measured using the pH-stat method. Other studies have gone further by comparing different in vitro protocols and methods; for example, the true digestibility of chickpea quantified in rats was found to correlate better (R 2 = 0.6785, p = 0.0640) with the degree of hydrolysis determined using O-phthalaldehyde after incubation with pepsin/pancreatin (Tavano, Neves and da Silva Júnior, 2016). Similarly, evaluations of in vitro static nitrogen or amino acid digestibility-in which samples were either filtered with a 1.2-µm filter (Millipore) (Saunders, Connor and Booth, 1973) or a 1 kDa-cut off dialysis bag (Vachon et al., 1987;Rozan et al., 1997)-have found favourable correlations with different nutritional quality indicators of proteins measured in rats. Instead, correlations between the apparent ileal digestibility of amino acids obtained from pigs and in vitro static amino acid digestibility (after filtration of the samples on a glass filter, pore size 40-90 pm; [START_REF] Cho | Evaluation of the Apparent Ileal Digestibility (AID) of Protein and Amino Acids in Nursery Diets by In vitro and In vivo Methods[END_REF] were mixed, with the coefficient of correlation variable from one amino acid to another. Proteolysis has also been examined in vitro using enzyme-linked immunosorbent assays (ELISA; [START_REF] Ménard | Validation of a new in vitro dynamic system to simulate infant digestion[END_REF] or by LC-MS/MS identification of released peptides [START_REF] Egger | The harmonized INFOGEST in vitro digestion method: From knowledge to action[END_REF]. In these studies, though, comparisons with in vivo results were globally affected by a high degree of inter-individual variability, which forced the authors to differentiate each animal in the presentation of the results. This complicated the analysis, because it was difficult to determine if a strong correlation was due to the quality of the simulation or an artefact of the high in vivo variability (Walther et al., 2019). In general, the number of possible set-ups to replicate digestion in vitro is large, and a major remaining challenge is to identify those that are most physiologically relevant [START_REF] Bohn | Correlation between in vitro and in vivo data on food digestion. What can we predict with static in vitro digestion models?[END_REF]. To address this, the INFOGEST network recently published standardised https://doi.org/10.1016/j.foodchem.2020.128276 Received 30 May 2020; Received in revised form 25 September 2020; Accepted 29 September 2020 approaches for static protocols, i.e. those that do not take into consideration the kinetics of digestion [START_REF] Brodkorb | INFOGEST static in vitro simulation of gastrointestinal food digestion[END_REF]. However, no such consensus has yet emerged for dynamic systems of digestion. Indeed, these models are even more complex than their static counterparts, and require extensive knowledge of how digestive physiology changes over time (e.g., gastric emptying kinetics or pH changes). To address this challenge, different teams have designed multiple systems that aim to replicate these dynamic conditions in vitro. Of these, the Vatier [START_REF] Nunez | Interests of the "artifical stomach" techniques to study antacid formulations: Comparison with in vivo evaluation[END_REF] and the Mainville [START_REF] Mainville | A dynamic model that simulates the human upper gastrointestinal tract for the study of probiotics[END_REF]) models consist of a succession of reactors, each connected to the next, which mimic the stomach and other compartments of the small intestine. More recently, the DiDGI® (Digesteur Dynamique Gastro-Intestinale; Ménard et al., 2014 or the ESIN (Engineered Stomach and small Intestine;[START_REF] Guerra | Development and Validation of a New Dynamic Computer-Controlled Model of the Human Stomach and Small Intestine[END_REF] systems were developed on the same basis. The most well-known system is the TIM-1 [START_REF] Minekus | The TNO Gastro-Intestinal Model (TIM)[END_REF], which also simulates contractile aspects of digestion. However, only limited data exist that assess the correlation of these in vitro systems with in vivo models, and most of these data focus on the digestion of milk proteins. To our knowledge, correlations between in vivo and in vitro digestion of plant protein have only been performed using static in vitro models and in vivo data from the rat model. However, the rat model has several drawbacks with respect to mimicking human digestion, as the eating behaviour of rats and their physiological parameters of digestion are not directly comparable to those of humans [START_REF] Davies | Physiological Parameters in Laboratory Animals and Humans[END_REF]. Instead, the gastrointestinal conditions in pigs are much more similar to those of humans, both anatomically and physiologically [START_REF] Henze | The pig as a preclinical model for predicting oral bioavailability and in vivo performance of pharmaceutical oral dosage forms : A PEARRL review[END_REF] and for this reason, the porcine model should be preferred for studies of protein digestion and digestibility (FAO, 2013).

The aim of this study was to examine the correlation between a porcine model of food-protein digestion and a dynamic in vitro gastrointestinal system, the DiDGI®. Specifically, we used data from three sets of experiments on cannulated mini-pigs (gastric, duodenal, and ileal) to (i) determine optimal parameters for the set-up of the in vitro system (pH, chyme transfer) and (ii) compare the evolution of digestion in both in vivo and in vitro models, specifically with respect to dry matter content, local apparent degree of proteolysis, and molecular weight distributions of peptides. This latter goal was accomplished by comparing the digestion of two soya-based foods (soya milk and tofu) in the DiDGI® and in mini-pigs. Ileal digestibilities were also obtained from both models and compared.

Materials and methods

Test meals

Two soya-based food matrices were tested: commercial tofu and commercial soya milk. Soya milk was UHT-(Ultra-High Temperature) packed. Compositions of foodstuffs have been previously published (Reynaud, Lopez, Riaublanc, Souchon, & Dupont, 2020;[START_REF] Binnerts | Soluble chromium indicator measured by atomic absorption in digestion experiments[END_REF]. Food matrices were supplemented with maltodextrin, sugar, and soya oil (Emile Noël, France) in order to standardise the amounts of protein (30.0 g), fat (23.1 g), and calories (980 kcal) in each diet (Table S1). The meals used for the DiDGI® experiments had the same nutrient proportions (7.5 g of protein, 5.8 g of fat and carbohydrates for a total of 245 kcal).

Dry matter content (DM) was determined by oven-drying at 105 °C overnight. Crude fat content was determined using the Randall method after 12 h incubation in 1 M HCl (ISO, 2008). Total nitrogen content was measured using the Kjeldahl method according to Thiex, Manson, Anderson, and Persson (2002). To calculate crude protein content, conversion factors of 5.36 for pea proteins, 5.50 for soya, and 5.50 for wheat proteins were used according to [START_REF] Mariotti | Converting nitrogen into protein -Beyond 6.25 and Jones' factors[END_REF]. Starch content was measured using a polarimeter (European Commission, 2009). Gross energy was calculated from macronutrient content using standard conversion factors.

In vivo digestion study on mini-pigs

All procedures were conducted in accordance with European guidelines (Directive 2010/63/EU) and were approved by the Auvergne Animal Experimentation Ethics Committee (CEMEAA) and the French Government. Yucatan mini-pigs were housed in subject pens in a ventilated room with controlled temperature (21 °C). Between sampling days, they received 500 g/d of a standard diet for minipigs, distributed in two equal portions at 8:00 and 16:00, and had free access to water. They were divided into three groups.

Gastric group

The gastric test included four mini-pigs (weight = 24.9 ± 1.2 kg, 8 months old). Three weeks before the trials began, mini-pigs were surgically fitted with a cannula (silicone rubber; 12-mm internal diameter, 17-mm outer diameter) at the greater curvature in the lateroventral region of the corpus. Digesta were sampled a minimum of 5 min before ingestion, then 20,40,60,90,120,180,270,360, and 450 min after ingestion of the meal. To each sample, 10 µl of pepstatin A solution (P5318, #027M4005V, Sigma-Aldrich, USA) were added per ml of digesta, at a concentration of 0.5 mg.ml -1 in methanol/acetic acid (9:1). Samples were then placed in liquid nitrogen and stored at -20 °C.

Duodenal group

The duodenal test involved six mini-pigs (weight = 24.7 ± 0.4 kg, 8 months old). Three weeks before the experiment, a cannula (silicone rubber; 12-mm internal diameter, 17-mm outer diameter) was fitted on mini-pigs 5-10 cm under the pyloric sphincter, as described previously by [START_REF] Bauchart | Peptides reproducibly released by in vivo digestion of beef meat and trout flesh in pigs[END_REF]. Each meal contained a Cr-EDTA solution (30 ml, 2.77 mg Cr .ml -1 ) as a liquid marker [START_REF] Binnerts | Soluble chromium indicator measured by atomic absorption in digestion experiments[END_REF], and 450 mg Ytterbium-acetate as a solid-phase marker (Siddons et al., 1985).

Duodenal sampling was performed a minimum of 5 min before ingestion and then 20, 40, 60, 90, 120, 180, 270, 360, and 450 min after ingestion of the meal. The pH of the collected duodenal content (15 ml maximum) was measured with an IFSET probe (1001-004 ISFET, Sentron, Holland), then the sample was divided into two parts. One part was freeze-dried (Pilot Freeze Dryer, Cryotec, France) for marker quantification; the other part was blocked with 1 mM final concentration of Pefabloc SC (AEBSF) (24305500, Roche, USA) diluted in distilled water. Then, the blocked mixture was vortexed for 5 s, placed in liquid nitrogen and stored at -20 °C until analysis.

The freeze-dried portions of the collected duodenal contents were also used to determine the dry matter content (DM) of duodenal samples (based on the weight difference before and after freeze-drying).

Freeze-dried duodenal contents were then ground and prepared for microwave plasma-atomic emission spectrometry (MP-AES) with successive steps of mineralisation (550 °C, 6 h), nitric acid dissolution (3 min boiling in 5% HNO 3 -(w/w)), and filtration (Whatman filter 4-7 µm). The amount of chromium-labelled indigestible content was then determined with a 4210 MP-AES apparatus (Agilent). The percentage of markers recovered was calculated by multiplying the measured content of the marker by the fixed volume sampled from the cannula (12 ml) and dividing this by the initial amount of the marker that was present in the meal.

Ileal group

The ileal test involved six mini-pigs (weight = 20.2 ± 1.5 kg, 8 months old). Three weeks before the experiment, mini-pigs were surgically fitted with a cannula (silicone rubber; 12-mm internal diameter, 17-mm outer diameter) at their distal ileum (10 cm before the ileocecal valve). Chromium oxide (Cr 2 O 3 ) was used as an indigestible transit marker, and was added up to the amount of 0.3% of meal dry matter content.

Contents were collected between 1 and 9 h after the distribution of meals in pre-weighed plastic bags (Sachet Whirl-Pak, 120 ml, 7.5 × 18.5 cm) attached to the ileal cannula. Full bags were replaced with new ones as necessary. Content was immediately transferred to aluminium dishes (previously weighed) and stored at -18 °C. Samples were pooled in 1-h pools and stored at -20 °C.

Ileal samples underwent the same treatments as the freeze-dried duodenal samples, i.e. freeze-drying, grinding, and mineralisation, to enable chromium quantification in samples using microwave plasmaatomic emission spectrometry (4210 MP-AES, Agilent). Then, a representative sample of the entire postprandial period (per animal-meal combination) was constituted using a fixed percentage (DM basis) from each batch of collected effluent (hourly), starting from the first detection of chromium.

Determination of operating parameters for the DiDGI® dynamic in vitro system

Gastric emptying equation

Gastric emptying was modelled according to Elashoff's equation [START_REF] Elashoff | Analysis of gastric emptying data[END_REF] and parameters were determined for each meal from mineral data collected through the duodenal cannula.

For each timepoint, the total amount of marker collected up to this time was summed. This function is thought to follow the equation:

= f t A ( ) (1 2 t t ( ) 1/2 (1)
where A is the maximum percentage of marker collected compared to initial marker intake; is the shape of the curve; and t 1/2 is the time at which 50% of the meal has been emptied from the compartment (here the stomach). The parameters of interest are t 1/2 and , from Elashoff's equation. Our model assumed that the flux of chyme was continuous and homogeneous. To avoid aberrant regression, samples that did not demonstrate a plateau in their kinetics were removed from the analysis; after this filtering step, at least five animals per meal and per marker were retained for analysis. Regression parameters are presented in Fig. 1A.

Intestinal transit flow

Transit of the bolus from the mouth to the ileal cannula was estimated with formula (1). In initial tests, the plateau of recovered chromium was not reached (Fig. 1B), so we setA (the maximum percentage of marker collected compared to initial marker intake) to 0.5. This can be interpreted as the hypothesis that the maximum content that could be recovered through the ileal cannula was 50%.

pH in the stomach and duodenum

We used food-dependent gastric acidification equations from a previous study (Reynaud, Buffière, et al., 2020), using the mean data of pH catheter probes. Duodenal kinetics of pH were deduced from a nonlinear regression of the pH measured in duodenal samples (Fig. S1).

Endogenous secretions

We simulated the saliva and basal gastric fluid present in the stomach when the chyme arrived by adding 12 ml of simulated saliva fluid (SSF) and 40 ml of simulated gastric fluid (SGF) to the initial meal. These solutions were prepared according to [START_REF] Minekus | A standardised static in vitro digestion method suitable for food -an international consensus[END_REF] and quantities were determined based on euthanasia data (data not shown).

These amounts of endogenous secretions came from the recommendations of INFOGEST [START_REF] Minekus | A standardised static in vitro digestion method suitable for food -an international consensus[END_REF]. Therefore, in order to maintain the same global ratio of protein:protease from one meal to another, enzymes were injected as follows:

In the stomach, the total amount of pepsin units (P6887, #SLBV3035, Sigma-Aldrich, USA) was 80 000 per gram of food protein and the speed of injection was 0.5 ml.min -1 (concentration of the stock solution: 3333 U pepsin .ml -1 ). In the duodenum, the total amount of trypsin units from pancreatin (P7545, #SLBV6830, Sigma-Aldrich, USA) was 4 000 per gram of food protein and the speed of injection was 0.25 ml.min -1 (concentration of the stock solution: 333 U trypsin .ml -1 ). Biliary acids solution (B8631, Sigma-Aldrich, USA) was secreted at 0.5 ml.min -1 . For the first 30 min, a 40 mmol.l -1 solution of biliary acids was used, then the concentration was reduced to 20 mmol.l -1 .

In vitro digestion study using DiDGI®

The dynamic three-compartment model of the gastrointestinal tract

This study used a three-compartment system instead of the original two-compartment system described in [START_REF] Ménard | Validation of a new in vitro dynamic system to simulate infant digestion[END_REF]. In the original system, the two compartments represent the stomach and the small intestine, respectively; to this set-up, Adouard and colleagues (2016) added a third compartment to split the small intestine into 2 independent sections i.e. the duodenum and the jejunum + ileum [START_REF] Adouard | Survival of cheese-ripening microorganisms in a dynamic simulator of the gastrointestinal tract[END_REF]. In our experiment, the third compartment emptied into a receptacle, the contents of which were analogous to those collected though the ileal cannula of mini-pigs.

Digestion programme

The digestion programme was designed following consultation with available in vivo data and existing literature. It was food-dependent and is summarised in Table 1. pH was regulated with 0.2 N HCl and 1 N NaOH solutions.

Post-hoc treatment of ileal content for the simulation of intestinal absorption

All contents collected from the third DiDGI® compartment corresponded to the in vivo effluent arriving at the ileal cannula; in both cases, ileal contents were collected in aluminium plates on the same sampling schedule. At least once per hour of digestion, once the chyme reached the in vitro "ileum", contents were manually stirred and a 2-ml sample was collected and blocked with 1 mM final concentration of Pefabloc SC (76307, #117M4033V, Sigma-Aldrich, USA) diluted in distilled water. Samples were vortexed for 5 s and put on ice before being stored at -20 °C. The remainder of the ileal contents was frozen at -20 °C before being freeze-dried.

All of the frozen samples from a given digestion were pooled, homogenised, and centrifuged (5 000 G, 10 min, 20 °C). Pellets and supernatants were weighed, and total nitrogen content was determined (Kjeldahl method). Supernatants were then subjected to different methods for simulating intestinal absorption, described below. From those processes, we generated values for simulated-absorbed and simulated-non-absorbed fractions of ileal content.

Values to be correlated between in vivo and in vitro experiments

Dry matter content (DM) in duodenal and ileal sections

Dry matter content was determined from the difference between sample weights before and after freeze-drying, with the exception of the duodenal contents sampled from the DiDGI®, for which the oven method (105 °C, overnight) was used.

Proteolysis in gastric and duodenal sections

Determination of the concentration of α-amino groups, an indicator of proteolysis in the stomach and the duodenum, was carried out using the ninhydrin reaction described by [START_REF] Moore | A modified ninhydrin reagent for the photometric determination of amino acids and related compounds[END_REF]. This technique was applied to gastric and duodenal samples that had been blocked immediately after they had been taken. The total nitrogen content of each sample was measured (Kjeldahl method) to calculate the "local apparent degree of proteolysis" (LADP) with the following equation, where the density of contents was considered to be 1:

= × LADP AG in supernatant mg N l Ntot in sample gN g (%) ( ( ). ) ( . 100 ) 10 1 1 2
(2)

Molecular weight determination by size-exclusion HPLC (HPSEC) of gastric samples

Gastric samples taken from mini-pigs and from the DiDGI® prior to ingestion and at 20, 60, 90, and 180 min post-ingestion were analysed with respect to the molecular weight distribution of peptides. Samples were centrifuged (5 000 G, 10 min, 20 °C), diluted 1:10 in an trifluoroacetic acid (TFA):acetonitrile:water buffer (0.01:2:10 v:v:v), and then filtered through a 0.40-μm pore-size filter. Peptides were separated according to their molecular weight with a column (Shodex Asahipak GF-310HQ, Interchyme, France) that was connected to a Waters e2695 separation module equipped with a Waters 2998 UV/Visible detector with diode array (Waters Inc., USA). Elution was performed under an isocratic 0.2 ml.min -1 flow of TFA:acetonitrile:water buffer at 37 °C. Detection was performed at 214 nm. Sample injection volume was 30 μl. The column was calibrated by injecting 30 μl of five markers of known molecular weight, each at a concentration of 2 mg.ml -1 : blue dextran (2000 kDa) and carbonic anhydrase (30 kDa), cytochrome C (12.384 kDa), vitamin B12 (1.855 kDa), and tripeptide GlyGlyGly (189 Da). The resulting calibration curve (log molecular weight vs retention time) enabled us to distinguish among 13 ranges of molecular weight: > 35 kDa, 35-18 kDa, 18-10 kDa, 10-6 kDa, 6-3 kDa, 3-2.5 kDa, 2.5-2 kDa, 2-1.5 kDa, 1.5-1 kDa, 1-0.8 kDa, 800-500 Da, 500-204 Da, < 204 Da.

Protein digestibility

Using the data obtained from the mini-pigs fitted with ileal cannulas, the apparent and true ileal digestibility of crude protein and amino acids were determined (Reynaud et al., 2021). These reported values were used here to study the correlation between in vivo-and DiDGI®-generated data.

Total nitrogen content (elemental analyser, vario ISOTOPE cube, Elementar) of both the simulated-absorbed and simulated-non-absorbed fractions of the DiDGI® ileal content were determined.

Estimates of digestibility obtained from the DiDGI® system are hence expressed as the ratio X simulated absorbed fraction /(X simulated absorbed fraction + X simulated non-absorbed fraction )

Data & statistical analysis

Data are disclosed as mean ± SD (standard deviation) for in vitro data and as mean ± SEM (standard error of the mean = SD / √number of replicates) for in vivo data or To present the data, we chose to use the SEM for in vivo data because n = 4 -6 and the variability is high. Bars on the graphics are reduce and are more readable. For in vitro data, we chose the SD because n = 3 and variability between replicates is moderate. Repeated measures ANOVA was performed on dry matter content (DM), local apparent degree of proteolysis (LADP), and molecular weight distributions to detect significant differences; the model included the fixed effects of model, food, time, model*time, model*food, food*time, and food*model*time. A post-hoc LSMEANS test Pepsin secretion 600 000 U pepsin over 360 min (0.5 ml.min -1 ) (P6887) Pancreatin secretion 30 000 U trypsin over 360 min (0.25 ml.min -1 ) (P7545) Bile secretion 40 mmol.l -1 until 30 min, then 20 mmol.l -1 (0.5 ml.min -1 ) (B8631)

was performed when interactions were significant. Differences were considered statistically significant if the p-values were < 0.05. All statistical analyses were performed using SAS (SAS University Edition, Release 3.71; SAS Institute Inc., Cary, NC, USA). Modelling of marker data was performed with XLSTAT (XLSTAT Quality 18.07, Addinsoft, France). Gastric data were analysed with a principal component analysis (PCA) carried out in XLSTAT.

Results

Effect of the model on dry matter content

Dry matter content (DM) in the duodenum compartment of the DiDGI® was compared to the DM of samples taken from mini-pigs fitted with a duodenal cannula. Data are presented in Fig. S2. In both the in vivo and in vitro models, DM values globally decreased after the intake of soya milk and tofu (from 13% to 5% and from 20% to 5%, respectively). For both soya milk and tofu, there was a significant model*time interaction: for the digestion of soya milk, DM values were significantly different between in vivo and in vitro models at 270 min, and the same was true for the digestion of tofu at 120 and 180 min. Beside this, no significant differences were detected for DM values throughout the duodenal digestion. Generally speaking, observations made from the different models at a given time point were consistent with each other, with the exception of the instances noted above.

Next, we examined the kinetic evolution of DM between samples taken from the ileal group of mini-pigs and those taken from the output of the DiDGI® (Fig. 1B). Because of the way the DiDGI® was set up, the tofu chyme arrived in the last compartment more slowly than the soya milk chyme, and the tofu-derived ileal content was collectable only starting from 300 min. A significant model*time interaction was detected for both soya milk (p = 0.0113) and tofu (p = 0.0329). The DM values of in vivo ileal samples were always lower than their corresponding in vitro ileal samples from the DiDGI®, but for the soya milk meal, these differences were statistically significant only for the samples taken at 120 min (p = 0.0006) and 180 min (p = 0.0282). For the tofu meal, instead, in vitro DM was significantly higher than in vivo DM for both sampling times available (300 min and 360 min, p < 0.0028). This result was expected, as in vivo digestion also includes absorption by the small intestine, which was not simulated during the in vitro dynamic digestion.

Effect of the model on local apparent degree of proteolysis

In the gastric phase (Fig. 2A &A'), the local apparent degree of proteolysis (LADP) was similar between in vivo and in vitro models for both meals; the only exception to this was found with the tofu meal sampled at 20 min, for which LADP differed in vitro and in vivo (p = 0.0031). Moreover, the effect of time was not significant, meaning that the kinetics were constant and LADP did not vary significantly as a function of time.

At the duodenal step of the soya milk digestion (Fig. S2B), significantly different LADP values were obtained from the DiDGI® and from mini-pigs only at the end of the assay, at 360 min (p < 0.0001). For the tofu meal (Fig. S2B), instead, the kinetics of the evolution of LADP were unchanged between in vivo and in vitro models (interaction effect: p = 0.5954). The LADP in the duodenum was determined for each in vitro replicate and animal/food combination between 20 and 180 min of digestion, and these values were plotted as a function of the dry matter content (DM) of the samples (Fig. 3). There was significantly more variation in duodenal DM (interaction food*time, p < 0.0001) with the tofu meal (5-23%) than with the soya milk meal (6-18%). For this reason, plotted values associated with the tofu meal were highly dispersed on the graphic, with a particularly visible distinction between DiDGI® replicates and the corresponding animal trials. Interestingly, at 20 and 60 min, tofu DM was similar between in vivo and in vitro assays, whereas LADP values obtained from the DiDGI® were lower than those from mini-pigs. At 120 and 180 min, the DiDGI® generated a tofu chyme that was less hydrolysed, and with a lower DM, than that found in mini-pig samples. For digestions of the soya milk meal, the largest differences between in vivo and in vitro tests were found at 20 and 180 min, at which time DM values were lower in the DiDGI® than in the mini-pig duodenum. At the end of digestion, the soya milk meal seemed to be more hydrolysed in the DiDGI® (8-11%) than in vivo (5-8%). Overall, the two meals followed two different patterns of digestion in the duodenum (interaction food*time of LADP: p < 0.0001, and of DM: p < 0.0001): for soya milk, LADP increased as digestion progressed and the DM decreased, while for tofu, LADP decreased with DM as digestion progressed.

Effect of the model on the molecular weight distribution of gastric peptides

Table S2 presents the results of statistical analyses that examined the relative proportions of gastric peptides, as determined by their molecular weight, in the outputs of the in vivo and in vitro systems at 20, 60, 120 and 180 min after tofu or soya milk meals ingestion. Fig. 4 is a graphical representation of a principal component analysis (PCA) of the proportions of gastric peptides grouped by molecular weight plus the gastric LADP. Only fractions for which at least one significant interaction was detected were included; thus, data for the > 35 kDa, 500-800 Da, and 204-500 Da fractions were removed from the PCA. The x-axis, the first principal component, labelled F1, contributed 46.8% of the observed variability. It was positively correlated with the proportions of low-molecular-weight fractions (< 2.5 kDa) and negatively correlated with the proportions of higher-molecular-weight fractions (> 3 kDa) (Fig. 4A). The y-axis, the second principal component, labelled F2, contributed 21.5% of the observed variability and was positively correlated with the proportion of the 2.5-3 kDa fraction (Fig. 4A). In Fig. 4B to 4E, each sample is plotted on this plane in order to observe the evolution of the in vivo or in vitro digestions between 20 and 180 min. The first observation to note is that the variability among animals was higher than the variability between in vitro assays. Indeed, for each food and at each time, the DiDGI® replicates (of tofu and soya milk digestions in Fig. 4B and Fig. 4C respectively) group together, whereas the replicates of digestion in each pig were quite distinct (Fig. 4D and4E). At 20 min of digestion, in vivo and in vivo data were on different side of the factorial plan: gastric peptides observed in vivo were larger (> 3 kDa) than those observed in vitro (< 2.5 kDa). Then, the distribution of gastric peptides tended to become similar along the digestion. As gastric digestion progressed, the proportion of 2-3 kDa peptides increased globally, so that most assays reached the same area on the factorial plane at 180 min. Exceptions were observed in two animals who ingested the soya milk meal. A high degree of divergence was observed at 20 min for in vivo data, but starting from 120 min, proteolysis patterns from most animals were similar to each other. No clear difference between meals was observed in vivo whereas in the in vitro digestion of the tofu and the soya milk meals, the relative proportion of the 2-3 kDa fraction appeared to increase more slowly for the soya milk meal than for the tofu meal.

In particular, proportions of three fractions-1.5-2 kDa (p = 0.0070), 1-1.5 kDa (p = 0.0003), and 0.8-1 kDa (p < 0.0001)-evolved differently over time in the two models (Fig. S4 C, D, andE). In each case, the DiDGI® yielded decreasing or morestable values than mini-pigs. In analysing the fractions, we also noted the following interesting results:

• In the 10-18 kDa fraction (Fig. S4 A), there were significant effects of the interactions model*food (p = 0.0347) and model*food*time (p = 0.0499), but not of the interactions food*time (p = 0.4915) and model*time (p = 0.1323). Indeed, the proportion of this fraction increased over time in the DiDGI® but remained constant in the mini-pig model.

• In the 2.5-3 kDa fraction (Fig. S4 B), there were significant effects of the interactions model*time (p = 0.0232) and model*food*time (p = 0.0439), but not the interactions food*time (p = 0.1674) and model*food (p = 0.9487). Proportions of the 2.5-3 kDa fraction evolved differently between models, and the effect of the model differed between the tofu meal and the soya milk meal.

• In the 1-1.5 kDa fraction (Fig. S4 D), we detected the significant interactions model*time (p = 0.0003) and food*time (p = 0.0015), but the third-order interaction, model*food*time, was not significant (p = 0.5582). A compensation phenomenon occurred, in which different kinetics were observed between the tofu and the soya milk meals in the DiDGI®, but not in mini-pigs. Furthermore, until 120 min, higher proportions of this fraction were detected in the soya milk meal in mini-pigs than in the DiDGI® or from the tofu meal in either model.

Effect of the model on ileal digestibility

Table 2 shows the in vivo and in vitro digestibility of the tofu and soya milk meals. For the in vitro calculations, the difference between endogenous and dietary nitrogen was not considered. Moreover, because no significant difference was observed between the two meals with respect to the true ileal digestibility in vivo (soya milk meal: 99.4 ± 2.4%; tofu meal:101.8 ± 8.9%), only the in vivo apparent ileal digestibility was compared with the in vitro digestibility of nitrogen. Although a significant effect was detected for food (p = 0.0365), there was no effect of model (p = 0.4207) or of the model*food interaction (p = 0.5805).

Indeed, when we performed a simple centrifugation of ileum content from the DiDGI®, and considered the supernatant as absorbed and the pellet as undigested, we obtained the same range of values for nitrogen digestibility as those obtained in vivo. The difference in the apparent ileal digestibility of nitrogen measured in vivo between the soya milk meal (71.3 ± 2.5%) and the tofu meal (56.5 ± 7.8%) was also found using the DiDGI® system, with the in vitro digestibility of nitrogen from soya milk and tofu equal to 72.7 ± 1.4% and 63.7 ± 3.5%, respectively.

As the soya milk and the tofu meals were totally digestible when we examined the in vivo true ileal digestibility, we expect that the nitrogen measured in the pellet of in vitro ileal content likely came from endogenous secretions.

Discussion

Other studies have compared in vivo data from pig experiments with the nitrogen balance resulting from in vitro digestions that included a filtration/dialysis step to simulate intestinal absorption [START_REF] Babinszky | An In-Vitro Method for Prediction of the Digestible Crude Protein Content in Pig Feeds[END_REF][START_REF] Boisen | Prediction of the apparent ileal digestibility of protein and amino acids in feedstuffs and feed mixtures for pigs by in vitro analyses[END_REF][START_REF] Minekus | Development and Validation of a Dynamic Model of the Gastrointestinal Tract[END_REF][START_REF] Cho | Evaluation of the Apparent Ileal Digestibility (AID) of Protein and Amino Acids in Nursery Diets by In vitro and In vivo Methods[END_REF]. There have also been investigations focused on individual amino acids, which have yielded contrasting results depending on the amino acid under consideration [START_REF] Boisen | Prediction of the apparent ileal digestibility of protein and amino acids in feedstuffs and feed mixtures for pigs by in vitro analyses[END_REF]Pujol and Torrallardona, 2007;[START_REF] Jezierny | In vitro prediction of standardised ileal crude protein and amino acid digestibilities in grain legumes for growing pigs[END_REF][START_REF] Cho | Evaluation of the Apparent Ileal Digestibility (AID) of Protein and Amino Acids in Nursery Diets by In vitro and In vivo Methods[END_REF]. In the literature, most of the recent examinations of the correlation between in vivo and in vitro protein digestion have used milk proteins, which are almost completely digested at the ileum stage [START_REF] Minekus | Development and Validation of a Dynamic Model of the Gastrointestinal Tract[END_REF]. For this reason, efforts have been made to compare the evolution of proteolysis throughout digestion and not only the end results [START_REF] Egger | The harmonized INFOGEST in vitro digestion method: From knowledge to action[END_REF]. Here, we aimed to build on these previous works by investigating the condition of proteins throughout the digestive process.

In the ileum, proportions of in vivo digestible nitrogen and in vitro simulated-absorbed nitrogen were similar, and the difference between the soya milk and tofu meals was conserved from one model to another. Early studies correlating measures of pH drop in vitro with the apparent or true fecal digestibility of protein in vivo were performed on a large spectrum of foodstuffs, including plant-based foods [START_REF] Hsu | A Multienzyme Technique for Estimating Protein Digestibility[END_REF]Perdersen and Eggum, 1981), but these results were shown to be highly dependent on the individual buffering capacity of a food created by its mineral and protein content [START_REF] Moughan | In-vitro Determination of Nitrogen Digestibility and Lysine Availability in Meat and Bone Meals and Comparison with in-vivo Ileal Digestibility Estimates[END_REF]. Later, Ihekoronye (1986) compared the final degree of hydrolysis in vitro, measured with TNBS, with measures of apparent fecal digestibility obtained from rats (Satterlee, Kendrick and Miller, 1977) and in doing so found similar values for textured soya protein (in vivo: 83.3%, in vitro: 83.7%) and solubilised soya flour (in vivo: 78.9%, in vitro: 80.0%). Here, instead, the apparent digestibility of our two soya-based food products was markedly lower than those previously published values: 10 points lower for soya milk, 20 points lower for tofu. Such a difference could originate from numerous factors, such as differences between the models used or in the food items studied.

In the duodenum, differences in DM kinetics might have arisen from different sources, such as inappropriate patterns of gastric emptying and/or dilution of the chyme in endogenous secretions (basal and postprandial, gastric and duodenal). In vivo, we expected a higher amount of endogenous secretion because of the contributions of mucus and hormone secretions, as well as cell desquamation, that were not considered in in vitro assays. The differences observed in LADP may have been due in part to variations in the quantity of digestive enzymes secreted (pepsin and pancreatin) and the dynamic of bolus transfer from one compartment to another. The most important factor, though, may have been changes in the local ratio of protease to protein. We had expected that the concentration of endogenous nitrogen compounds would be higher in the mini-pig duodenum than in the second compartment of the DiDGI®; however, no significant difference was observed between models with respect to the total nitrogen content of gastric samples (0.1-0.7 gN.100 g -1 ) or duodenal samples (0.1-0.4 gN.100 g -1 ). Moreover, in mini-pig experiments, the position of the duodenal cannula was set between the pylorus and the common bile/ pancreatic ducts. This meant that, in in vivo duodenal samples, pancreatic enzymes were perhaps not even present at their maximal levels. Some studies have investigated correlations of protein digestion between an in vivo model and a dynamic in vitro system by comparing transfers and pH modification through the digestive process [START_REF] Dupont | Can dynamic in vitro digestion systems mimic the physiological reality?[END_REF]. [START_REF] Minekus | Development and Validation of a Dynamic Model of the Gastrointestinal Tract[END_REF] examined the delivery and the absorption of dietary nitrogen in the different compartments of the system and compared this with marker data obtained from cannulated animals. In that study, though, protein hydrolysis was investigated only indirectly, by analysing the absorption of nitrogen [START_REF] Minekus | Development and Validation of a Dynamic Model of the Gastrointestinal Tract[END_REF]. Instead, Chiang and colleagues (2008) investigated the gastric digestibility of protein by TCA precipitation in vivo and in a dynamic in vitro system, and found a strong correlation between the two models (r = 0.97; p = 0.0001). TCA precipitation is a good first step to evaluate the state of proteins in a mixture, but here we went further by looking at the degree of hydrolysis of proteins and at the distribution of peptides by their molecular weights.

In the stomach, we were unable to precisely track the evolution of the local apparent degree of proteolysis of the tofu meal due to the high inter-individual variability in animals, the low number of mini-pigs involved, and the nature of the sampling scheme. Thus, the main comparisons that we were able to make between the tofu and soya milk meals with respect to gastric digestion were based on the molecular weight distributions of peptides. Here, the interaction food*time had a significant effect on the proportions of the 1-1.5 kDa and free amino acid fractions (< 204 kDa). Moreover, the 0.8-1 kDa and 204-500 Da fractions; in vivo, were more abundant in the soya milk meal than in the tofu meal. As final point of discussion, the tofu left the stomach more rapidly than the soya milk. It is because of the difference regarding the gastric emptying equations. However, after soya milk intake, gastric content was collectable from the DiDGI® for only 180 min, while in vivo, this content was collectable for 450 min. According to the recovered markers, at 180 min most of the soya milk chyme had emptied from the pig stomach, and content sampled after that time originated mostly from endogenous secretions.

Conclusion

In the present study, we elaborated a food-dependent digestion programme for the DiDGI® based on experimental data obtained from mini-pigs. We used this to then investigate the correlation between these two models with respect to the dry matter content in the duodenum and the ileum, the local apparent degree of proteolysis in the stomach and in the duodenum, and the molecular weight distribution of gastric peptides. Globally, no clear correlation was found at the gastric step because of the inter-individual variability of pigs and an inadequate release of pepsin in the DiDGI®. However, at a later stage of digestion, when the stomach was almost emptied, proteolysis patterns were largely consistent between models, especially for the tofu digestion. Differences were found in individual parameters in the duodenal phase, but the global kinetics were consistent from one model to another. Finally, from the ileal phase, we obtained similar values for nitrogen digestibility in both in vivo and in vitro models.

the DiDGI® dynamic in vitro system, and to Lindsay Higgins for correcting the English of the manuscript. Finally, we would like to acknowledge the contribution of Denys Durand, who performed surgery with Didier Rémond.

Fig. 1 .

 1 Fig. 1. A) Parameters of gastric emptying (MEAN ± SEM; n = 5-6). Mixed model ANOVA was performed on regression parameters (β, t1/2, A), with animal as random effect and meal as fixed effect. No significant difference was observed between soya milk and tofu meals. B) Parameters of the intestinal flow after tofu meal (broken blue curve) or soya milk meal (solid red line) ingestion (MEAN ± SEM; n = 5-6). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 2 .

 2 Fig. 2. Kinetics of proteolysis in the stomach and the duodenum, as measured by the evolution of the local apparent degree of proteolysis (LADP) in gastric (A: soya milk meal and A': tofu meal) and duodenal samples (B: soya milk and B': tofu). Data are MEAN ± SEM; n = 4 animals in gastric assays; n = 6 animals in duodenal assays. Repeated measures ANOVA was performed, with model, time, and model*time as fixed effects. Significant interactions were found for gastric LADP following tofu input (p = 0.0444) and duodenal LADP following soya milk input (p = 0.0002); a post-hoc LSMEANS test was therefore performed. Stars above curves indicate a significant difference between models at the given sampling time (p < 0.05).

Fig. 3 .

 3 Fig. 3. Values for local apparent degree of proteolysis in the duodenum were plotted against duodenal dry matter content following ingestion of A) the soya milk meal or B) the tofu meal. DiDGI® data (3 repetitions) are represented by crosses and mini-pig data (6 animals) are represented by solid squares. Colours indicate the time elapsed since ingestion (dark green: 20 min; light green: 60 min; orange: 120 min; red: 180 min). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 4 .

 4 Fig. 4. Results from animal/in vitro trials visualised using a Principal Components Analysis (PCA) of the molecular weight distributions of peptide fractions and the local apparent degree of proteolysis in the stomach. A) Graphical representation depicting the circle of correlation and contribution of each variable to the x-and yaxes. B) Graphical representation of the tofu meal digestion at different times in the first compartment of the DiDGI®. C) Graphical representation of the soya milk meal digestion at different times in the first compartment of the DiDGI®. D) Graphical representations of the tofu meal digestion at different times in the stomach of mini-pigs E) Graphical representations of the soya milk meal digestion at different times in the stomach of mini-pigs. n = 4 animals (except for the tofu in vivo at 180 min where n = 3) and n = 3 in vitro assays.

Table 1

 1 Parameters of the DiDGI® digestion.

	Parameters	Soya milk								Tofu
	Time of digestion	360 min
	Temperature	37 °C
	Gastric pH equation	1.8	+	1	+	exp	5.5 (0.40(2	))	+	exp	(0.04(	t	174))	1.5	+	1	+	exp	(0.40(3	+ 4.7 t )) exp(0.02(t 9))	+	× 1.5 ln(1 exp( t / 300)) + ln(2)
	Duodenal pH equation	2.86	+	+ (1 191.26) 0.89 3.41
	Gastric emptying equation	2	t ( 73	) 1.02									2	t ( 57	) 0.98
	Fixed volume of the duodenum	50 ml
	Intestinal flow equation	2	t ( 386	) 2.44									2	t ( 1101	) 2.52

Acknowledgements

We thank our animal housing manager, Julien Hermet, and our animal housing workers, Philippe Lhoste and Yohan Delorme. We are grateful to Caroline Buffière, Benoît Cohade, Mélissa Vauris and Kevin Liebermann who participate to in vivo experiments, to Hervé Guillemin and Thomas Cattenoz for helping to develop food-dependent set-up of

Funding

The project was funded by IMPROVE SAS and INRAE. A PhD grant was funded by IMPROVE SAS through a CIFRE procedure (CIFRE n°2016/0719) Fundamental and Clinical Pharmacology, 12, 573-583. Perdersen, B., & Eggum, O. (1981). Prediction of protein digestibility by in vitro procedures based on two mullti-enzyme systems. Zeitschrift für Tierphysiologie Tierernährung und Futtermittelkunde, 45(1-5), 190-200. Pujol, S., & Torrallardona, D. (2007). Evaluation of in vitro methods to estimate the in vivo nutrient digestibility of barley in pigs. Livestock Science, 109(1-3), 186-188. https://doi.org/10.1016/j.livsci.2007.01.143. Reynaud, Y., Buffière, C., Cohade, B., Vauris, M., Liebermann, K., Hafnaoui, N., ... Rémond, D. (2021). True ileal amino acid digestibility and digestible indispensable amino acid scores (DIAASs) of plant-based protein foods. Food Chemistry, 338, Article 128020. https://doi.org/10.1016/j.foodchem.2020.128020 ISSN 0308-8146. Reynaud, Y., Buffière, C., David, J., Cohade, B., Vauris, M., Lopez, M., ... Rémond, D. (2020). Temporal changes in postprandial intragastric pH: Comparing measurement methods, food structure effects, and kinetic modelling. Food Research International. Elsevier Ltd. https://doi.org/10.1016/j.foodres.2019.108784. Reynaud, Y., Lopez, M., Riaublanc, A., Souchon, I., & Dupont, D. (2020). Hydrolysis of plant proteins at the molecular and supra-molecular scales during in vitro digestion. Food Research International. Elsevier, 134(September 2019), Article 109204. https:// doi.org/10.1016/j.foodres.2020.109204. Rozan, P., et al. (1997). In Vivo and in Vitro Digestibility of Soybean, Lupine, and Rapeseed Meal Proteins after Various Technological Processes. Journal of Agricultural and Food Chemistry, 45(5), 1762-1769. https://doi.org/10.1021/jf960723v. Satterlee, L. D., Kendrick, J. G., & Miller, G. A. (1977). Rapid in vitro assays for estimating protein quality [of food]. Nutrition Reports International. Saunders, R. M., Connor, M. A., & Booth, A. N. (1973). Measurement of Digestibility of Alfalfa Protein Concentrates by in vivo and in vitro Methods. The Journal of nutrition, 103(4), 530-535. Siddons, R. C., et al. (1985). Ytterbium acetate as a particulate-phase digesta-flow marker. The British journal of nutrition, 54(2), 509-519. https://doi.org/10.1079/ BJN19850136. Tavano, O. L., Neves, V. A., & da Silva Júnior, S. I. (2016). 'In vitro versus in vivo protein digestibility techniques for calculating PDCAAS (protein digestibility-corrected amino acid score) applied to chickpea fractions', Food Research International. Elsevier Ltd, 89, 756-763. https://doi.org/10.1016/j.foodres.2016.10.005. Thiex, N. J., Manson, H., Anderson, S., & Persson, J.Å. (2002). Determination of crude protein in animal feed, forage, grain, and oilseeds by using block digestion with a copper catalyst and steam distillation into boric acid: Collaborative study. Journal of AOAC International, 85(2), 309-317. Vachon, C., et al. (1987). Relationship between in vitro digestion of proteins and in vivo assessment of their nutritional quality. Reproduction and Nutrition Dévelopment, 27(3), 659-672. Walther, B., et al. (2019). GutSelf: Inter-Individual Variability in the Processing of Dietary Compounds by the Human Gastrointestinal Tract. Molecular Nutrition & Food Research. https://doi.org/10.1002/mnfr.201900677.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Author Contributions

Conceptualization: YR, DD, DR, IS, DP and ML, Data curation: YR, Formal analysis: YR, Funding acquisition: DD, DR, IS, DP and ML, Investigation: YR, AC, DF, and AM, Methodology: DD, IS, DP, ML, AC and YR, Supervision: DD, DR, IS, DP and ML, Visualization: YR, Writing -original draft: YR, Writing -review & editing: All authors.

Appendix A. Supplementary data

Supplementary data to this article can be found online at https:// doi.org/10.1016/j.foodchem.2020.128276.