Jean-Guillaume Dumas
email: guillaume.dumas@imag.fr.

Vincent Zucca
email: vincent.zucca@lip6.fr.

Prover efficient public verification of dense or sparse/structured matrix-vector multiplication *

With the emergence of cloud computing services, computationally weak devices (Clients) can delegate expensive tasks to more powerful entities (Servers). This raises the question of verifying a result at a lower cost than that of recomputing it. This verification can be private, between the Client and the Server, or public, when the result can be verified by any third party. We here present protocols for the verification of matrix-vector multiplications, that are secure against malicious Servers. The obtained algorithms are essentially optimal in the amortized model: the overhead for the Server is limited to a very small constant factor, even in the sparse or structured matrix case; and the computational time for the public Verifier is linear in the dimension. Our protocols combine probabilistic checks and cryptographic operations, but minimize the latter to preserve practical efficiency. Therefore our protocols are overall more than two orders of magnitude faster than existing ones.

Introduction

With the emergence of cloud computing services, computationally weak devices (Clients, such as smart phones or tablets) can delegate expensive tasks to more powerful entities (Servers). Such heavy tasks can, e.g., be cryptographic operations, image manipulation or statistical analysis of large data-sets. This raises the question of verifying a result at a lower cost than that of recomputing it. This verification can be private, between the Client and the Server, or public, when the result can be verified by any third party.

For instance within computer graphics (image compression and geometric transformation), graph theory (studying properties of large networks), big data analysis, one deals with linear transformations of large amount of data, often arranged in large matrices with large dimensions that are in the order of thousands or millions in some applications. Since a linear transformation on a vector x can be expressed by a matrixvector multiplication (with a matrix of size m×n), a weak client can use one of the protocols in the literature [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF][START_REF] Zhang | Efficient secure and verifiable outsourcing of matrix multiplications[END_REF][START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] to outsource and verify this computation in the optimal time O(m + n), i.e., linear in the input and the output size. However as these protocols use expensive cryptographic operations, such as pairings, the constants hidden in the asymptotic complexity are usually extremely large [START_REF] Walfish | Verifying computations without reexecuting them[END_REF].

In this paper, we propose an alternative protocol, achieving the same optimal behavior, but which is also practical: the overhead for the Prover is now very close to the time required to compute the matrix-vector multiplication, thus gaining two orders of magnitude with respect to the literature. Our protocol not only does this for dense matrices, but is also sensitive to any structure or sparsity of the linear transformation. For this, we first remove any quadratic operation that is not a matrix-vector multiplication (that is we use projections and rank-1 updates) and second we separate operations in the base field from cryptographic operations so as to minimize the latter.

More precisely, we first combine rank-one updates of [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF] and the projecting idea of [START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] with Freivalds' probabilistic check [START_REF] Freivalds | Fast probabilistic algorithms[END_REF]. Second, we use a novel strategy of vectorization. For instance, with a security parameter s (e.g., an s = 128-bits equivalent security), exponentiations or pairings operations usually cost about O(s 3 ) arithmetic operations. To make the whole protocol work practical, we thus reduce its cost from O s 3 mn to O µ(A) + s 3 (m + n 4/3 ) , where µ(A) < 2mn is the cost of one, potentially structured, matrix-vector multiplication. We also similarly reduce the work of the Verifier. This allows us to gain two orders of magnitude on the Prover's work and therefore on the overall costs of outsourcing, while preserving and sometimes even improving the practical efficiency of the Verifier.

Thus, after some background in Section 2, our first improvement is given in a relaxed public verification setting in Section 3 via matrix projection and probabilistic checks. Our second improvement is given in Section 4 where the verification is bootstrapped efficiently by vectorization. We then show how to combine all improvements in Section 5 in order to obtain a complete and provably secure protocol. Finally, we show in Section 6 that our novel protocol indeed induces a global overhead factor lower than 3 with respect to non verified computations. This is gaining several orders of magnitude on the Prover side with respect to previously known protocols, while keeping the Verification step an order of magnitude faster.

Background and definitions

In this paper, we want to be able to prove fast that a vector is a solution to a linear system, or equivalently that a vector is the product of another vector by a matrix. This is useful, e.g., to perform some statistical analysis on some medical data. We distinguish the matrix, a static data, from the vectors which are potentially diverse. In the following, F p will denote a prime field and we consider:

• Data: matrix A ∈ F m×n p . • Input: one or several vectors x i ∈ F n p , for i = 1..k. • Output: one or several vectors y i = Ax i ∈ F m p , for i = 1.
.k. Then, we denote by an operation performed in the exponents (for instance, for u ∈ G n and v ∈ Z n , the operation u T v actually denotes

n j=1 u[j] v[j] ).
Publicly Verifiable Computation. A publicly verifiable computation scheme, in the formal setting of [START_REF] Parno | How to delegate and verify in public: verifiable computation from attribute-based encryption[END_REF], is in fact four algorithms (KeyGen, ProbGen, Compute, Verify), where KeyGen is some (amortized) preparation of the data, ProbGen is the preparation of the input, Compute is the work of the Prover and Verify is the work of the Verifier. Usually the Verifier also executes KeyGen and ProbGen but in a more general setting these can be performed by different entities (respectively called a Preparator and a Trustee). More formally we define these algorithms as follow:

• KeyGen(1 λ , f ) → (param, EK f , V K f ): a randomized algorithm run by a Preparator, it takes as input a security parameters 1 λ and the function f to be outsourced. It outputs public parameters param which will be used by the three remaining algorithms, an evaluation key EK f and a verification key Completeness. A publicly verifiable computation scheme for a family of function F is considered to be perfectly complete (or correct) if for every function belonging to F and for every input in the function domain, an honest Prover which runs faithfully the algorithm Compute will always (with probability 1) output an encoding σ y which will pass Verify.

V K f . • ProbGen(x) → (σ x ):
Soundness. A publicly verifiable computation scheme for a family of function F is called sound when a prover cannot convince a verifier to accept a wrong result y = y except with negligible probability. More formally we evaluate the capability of an adversary A to deceive the verifier through a soundness experiment.

In this experiment, we assume that the adversary A accesses to the output of the algorithm KeyGen by calling an oracle O KeyGen with inputs 1 λ and the function to evaluate f . This oracle O KeyGen returns public parameters for the protocol param, an evaluation key EK f and a verification key V K f . Afterwards the adversary A sends its challenge input x to an oracle O ProbGen which returns σ x . Finally A outputs an encoding σ y * = σ y and runs the Verify algorithm on inputs σ y * and V K f , whether it outputs y or ⊥ the experiment has either succeeded or failed.

Definition 1. A publicly verifiable computation scheme for a family of function F is sound if and only if for any polynomially bounded adversary A and for any f in F the probability that A succeeds in the soundness experiment is negligible in the security parameter.

Adversary model. The protocol in [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF] (recalled for the sake of completeness in Appendix A) is secure against a malicious Server only. That is the Client must trust both the Preparator and the Trustee. We will stick to this model of attacker in the remaining of this paper. Otherwise some attacks can be mounted:

• Attack with a Malicious Preparator only: send A and a correctly associated W to the Server, but pretend that A is used. Then, all verifications do pass, but for y = A x and not y = Ax. • Attack with a Malicious Trustee only: there, a malicious assistant can provide a wrong VK x making the verification fail even if the Server and Preparator are honest and correct. • Attack with Malicious Server and Trustee: the Server sends any y and any z to the Trustee, who computes VK'

x [i] = e(z [i]; g 2 )/a y [i]
, that will match the verification with y and z .

Public delegatability. One can also further impose that there is not interaction between the Client and the Trustee after the Client has sent his input to the Server. Publicly verifiable protocols with this property are said to be publicly delegatable [START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF]. The protocol in [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF] does not achieve this property, but some variants in [START_REF] Zhang | Efficient secure and verifiable outsourcing of matrix multiplications[END_REF][START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] already can.

Bilinear Pairings. The protocols we present in this paper use bilinear pairings and their security is based on the co-CDH assumption, for the sake of completeness we recall hereafter these definitions.

Definition 2 (bilinear pairing).

Let G 1 , G 2 and G T be three groups of prime order p, a bilinear pairing is a map e : G 1 × G 2 → G T with the following properties:

1. bilinearity: ∀a, b ∈ F p , ∀(g 1 , g 2 ) ∈ G 1 × G 2 , e(g a 1 , g b 2 )
= e(g 1 , g 2 ) ab ; 2. non-degeneracy: if g 1 and g 2 are generators of G 1 and G 2 respectively then e(g 1 , g 2 ) is a generator of G T ; 3. computability: ∀(g 1 , g 2 ) ∈ G 1 × G 2 , there exist an efficient algorithm to compute e(g 1 , g 2 ).

Definition 3 (co-CDH assumption).

Let G 1 , G 2 and G T be three groups of prime order p, such that there exist a bilinear map e :

G 1 × G 2 → G T . Let g 1 ∈ G 1 , g 2 ∈ G 2 be
generators and a, b $ ← F p be chosen randomly. We say that the co-computational Diffie-Hellman assumption (co-CDH) holds in G 1 , if given g 1 , g 2 , g a 1 , g b 2 the probability to compute g ab 1 is negligible.

Related work. The work of [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF] introduced the idea of performing twice the computations, once in the classical setting and once on encrypted values. This enables the Client to only have to check consistency of both results. The protocol is sound under the Decision Linear and co-CDH hypothesis. Then [START_REF] Zhang | Efficient secure and verifiable outsourcing of matrix multiplications[END_REF] extended the part on matrix-vector multiplication to matrix-matrix while adding public delegatability. Finally, [START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] introduced the idea of projecting the random additional matrix and the extra-computations of the server, which allows to reduce the cost of the Verify algorithm. It also decreases the size of the verification key by a factor m. For an m×n dense matrix, the protocol in [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF] has a constant time overhead for the Prover, but this constant is on the order of cryptographic public-key operations like pairings. Similarly, the Verifier has O (mn) cryptographic public-key pre-computations and O (n) of these for the public verification. These cryptographic operations can then induce some 10 6 slow-down [START_REF] Walfish | Verifying computations without reexecuting them[END_REF] and do no improve even if the initial matrix is sparse or structured (as the rank one updates, s • t T and σ • τ T , are always dense).

A first step towards public verifiability

Freivalds' probabilistic verification of matrix multiplications [START_REF] Freivalds | Fast probabilistic algorithms[END_REF] allows for private verifiability of matrix-vector computations. This can be naturally extended in the random oracle model via Fiat-Shamir heuristic [START_REF] Fiat | How to prove yourself: Practical solutions to identification and signature problems[END_REF]. This however forces the vectors to be multiplied to be known in advance (the full details are given in Appendix B), whereas our goal is instead to obtain public verifiability with an unbounded number of vector inputs. As an upstart, we thus first present an improvement if the public verification model is slightly relaxed: in this Section, we allow the Trustee to perform some operations after the computations of the Server. We will see in next sections how to remove the need for the Trustee's intervention. For this, we combine Freivalds projection (to check that Ax i = y i , one can first precompute w T = u T A and check that w T x i = u T y i ) with Fiore & Gennaro's protocol, in order to improve the running time of both the Trustee and the Client: we let the Prover compute its projection in the group. That way most of the pairings computations of the Trustee and Client are transformed to classical operations: the improvement is from O (n) cryptographic operations to O (n) classical operations and a single cryptographic one. Further, the projection can be performed beforehand, during the precomputation phase. That way the preparation requires only one matrix-vector for the Freivalds projection and the dense part is reduced to a single vector. The cryptographic operations can still be delayed till the last check on pairings. This is shown in Figure 1.

• Preparator: secret random u ∈ F m p , t ∈ F n p , then ω T = g u T A+t T 1 ∈ G n 1 . • Preparator to Prover: A ∈ F m×n p , ω ∈ G n 1
• Preparator to Trustee: u, t in a secure channel.

• Verifier to Prover:

x i ∈ F n p • Prover to Verifier: y i ∈ F m p , ζ i ∈ G 1 such that y i = Ax i and ζ i = ω T x i . • Verifier to Trustee: x i , y i • Trustee to Verifier: h i = (u T • y i ) ∈ F p and d i = (t T • x i ) ∈ F p , then send η i = e(g 1 ; g 2 ) hi+di ∈ G T .
• Verifier public verification: e(ζ i ; g 2 )

? == η i in G T .

Verifying the dot-products by bootstrapping and vectorization

To obtain public verifiability and public delegatability, the Client should perform both dot-products, u T • y and t T • x (from now on, for the sake of simplicity, we drop the indices on x and y). But as u and t must remain secret, they will be encrypted beforehand. To speed-up the Client computation, the idea is then to let the Server perform the encrypted dot-products and to allow the Client to verify them mostly with classical operations.

For this trade-off, we use vectorization. That is, for the vectors u and y, we form another representation as

√ m× √ m matrices: U =     u 1 . . . u √ m u 1+ √ m . . . u 2 √ m . . . . . . . . . u 1+m- √ m . . . u m     and Y =     y 1 . . . y 1+m- √ m y 2 . . . y 2+m- √ m . . . . . . . . . y √ m . . . y m     . Then u T • y = T race(U Y ).
Computing with this representation is in general slower than with the direct dotproduct, O √ m 3 instead of O (m). As shown next, this can be circumvented with well-chosen left-hand sides and at least mitigated, with unbalanced dimensions.

dot-product with rank 1 left-hand side

The first case is if u is of rank 1, that is if in matrix form, u can be represented by a rank one update,

U = µ • η T for µ, η ∈ F √ m p .
Then both representations require roughly the same number of operations to perform a dot-product since then:

T race(µ • η T • Y ) = η T • Y • µ (1)
Therefore, we let the Prover compute

z T = g η T 1 Y , where z[i] = g η[j]Y [j,i] 1
, and then the Verifier can check this value via Freivalds with a random vector v:

g η T 1 (Y • v) ? == z T v.
The point is that the Verifier needs now O (m) operations to compute (Y • v), but these are just classical operations over the field. Then its remaining operations are cryptographic but there is only O ( √ m) of these. Finally, the Verifier concludes the computation of the dotproduct, still with cryptographic operations, but once again with only O ( √ m) of them. Indeed, the dot product

d = u T y = T race(U Y ) = T race(µ • η T • Y ) = η T • Y • µ is checked by e(g d 1 ; g 2 ) = e(g 1 ; g 2 ) d = √ m i=1 e(z[i]; g µ[i]
2 ) and the latter is e(g 1 ; g 2 )

µ[i]η[j]Y [j,i] = e(g 1 ; g 2 ) η T •Y •µ . In practice,
2 = Θ ( √ m): • Preparation: O (b 1 +b 2 ) in G i ; • Prover: O (m) in G i ; • Verifier: O (m) in F p , O (b 1 +b 2 ) in G i and O (b 1 ) pairings.
Proof. Correctness is ensured by Equation [START_REF] Paulo | Pairing-friendly elliptic curves of prime order[END_REF]. Soundness is given by Freivalds check. Complexity is as given in the Lemma: indeed, for the Verifier, we have: for

Y • v: O (b 2 b 1 ) = O (m) classic operations; for g η T 1 (Y v): O (b 2 ) cryptographic (

Rectangular general dot-product

Now if u is not given by a rank 1 update, one can still verify a dot-product with only O ( √ m) pairings operations but as the price of slightly more group operations as given in Figure 3.

• KeyGen(1 λ , µ, η): given the security parameter 1 λ and vectors µ ∈ F b2 p and η ∈ F b1 p such that u = µ • η t , it selects two cyclic groups G 1 and G 2 of prime order p that admit a bilinear pairing e : G 1 × G 2 → G T and generators g 1 , g 2 and g T of the three groups. Finally it outputs params

= {b 1 , b 2 , p, G 1 , G 2 , G T , e, g 1 , g 2 , g T } and EK f = (g η T 1 ) and V K f = (g µ 2 ). • ProbGen(y): from y ∈ F m p it builds Y ∈ F b1×b2 p and outputs σ x = Y . • Compute(σ x , EK f ): compute z T = g η T 1
Y and outputs σ y = (z T ). • Verify(σ y , V K f ): it starts by sampling randomly a vector v ∈ F b2 p then it checks whether z T v is equal to g η T 1 (Y v) or not. If the test passes it returns e(z[i]; g

µ[i]
2 ) and if it fails it returns ⊥.

Figure 2: Publicly delegatable protocol for the dot-product in an external group with a rank-1 left hand side.

• KeyGen(1 λ , u): given the security parameter 1 λ and vector u ∈ F m p , it selects two cyclic groups G 1 and G 2 of prime order p that admit a bilinear pairing e : G 1 × G 2 → G T and generators g 1 , g 2 and g T of the three groups and also integers b

1 , b 2 such that m = b 1 b 2 and it outputs params = {b 1 , b 2 , p, G 1 , G 2 , G T , e, g 1 , g 2 , g T }.
Then it samples a random w ∈ F b1 p , creates U ∈ F b1×b2 p from u and finally it outputs and EK f = (g U 1 ) and

V K f = (g w T •U 1 , g ω T 2 ). • ProbGen(y): from y ∈ F m p it builds Y ∈ F b2×b1 p and outputs σ x = Y . • Compute(σ x , EK f ): compute C = g U 1
Y and outputs σ y = (C).

• Verify(σ y , V K f ): it starts by sampling randomly a vector v ∈ F b1 p then it computes z = C v and checks whether e(z[i]; g Proof. Correctness is ensured by the vectorization in Equation [START_REF] Paulo | Pairing-friendly elliptic curves of prime order[END_REF]. Soundness is given by the Freivalds check. Complexity is as given in the Lemma: indeed, for the Verifier, we have: Then Each chunk can be checked with the protocol of Figure 3 and the final dot-product obtained by adding the chunks. This gives the complexity of Corollary 7.

ω[i] 2 ) is equal to e(g w T •U 1 (Y v); g 2 )
1. Y • v: O (b 2 b 1 ) = O (m) classic operations; 2. g w T U 1 (Y v): O (b 2 ) cryptographic (group) operations; 3. z = C v: O b 2 1 cryptographic (group) operations; 4. b1 i=1 e(z[i]; g w[i]
Corollary 7. There exist a protocol for the dot-product using: O n 1-a/3 cryptographic operations for the Verifier and O n 1+a/3 cryptographic operations for the Prover, for any 0 < a < 1.

Proof. We let k = n a and use n/k times the protocol of the previous point with b

1 = 3 √ k and b 2 = k 2 3
. Overall this gives n/k(k 4/3 ) for the Prover and n/k(k 2/3 ) for the Verifier.

Public delegatability via bootstrapping

To recover the public delegatability model, we use the protocol of Figure 1 but we trade back some cryptographic operations using the protocol of Figure 2 to the Verifier. With an initial matrix A ∈ F m×n p we however trade back only on the order of O ( √ m + √ n) cryptographic operations. This gives a slower verification in practice but interaction is not needed anymore. We present our full novel protocol for matrix vector product in Figure 4 (with the flow of exchanges shown in Figure 7, Appendix E). Apart from Freivalds's

• Keygen(1 λ , A): given a security parameter 1 λ and a matrix A ∈ F m×n p , it selects two cyclic groups G 1 and G 2 of prime order p that admit a bilinear pairing e : G 1 × G 2 → G T . Then it selects generators g 1 , g 2 and g T of the three groups, and integers b

1 , b 2 , c 1 , c 2 , d 1 , d 2 such that b 1 b 2 = m and c 1 c 2 = d 1 d 2 = n with b 1 + b 2 ∈ O m 1/2 , c 1 + c 2 ∈ O n 1/2 , d 1 ∈ O n 1/3 , d 2 ∈ O n 2/3
and it outputs parameters:

param = (m, n, b 1 , b 2 , c 1 , c 2 , d 1 , d 2 , p, e, G 1 , G 2 , G T , g 1 , g 2 , g T ). Then it samples randomly µ ∈ F b1 p , η ∈ F b2 p , (ρ i ) 2 i=1 ∈ F c1 p , (τ i ) 2 i=1 ∈ F c2 p , ∈ F d1 p , (γ, δ) ∈ F 2 p and v ∈ F n p with v V ∈ F d1×d2 p . Finally it computes U = µ•η T u ∈ F n p , T = ρ 1 •τ T 1 +ρ 2 •τ T 2 t T ∈ F n p , ω T = g u T •A+t T +γδv T 1 ∈ F n p
and publishes:

EK f = {A, ω T , (g τ T 1
1 ), (g

τ T 2 1 ), (g η T 1 ), (g δ•V 1 )} V K f = {(g τ T 1 1 ), (g τ T 2 1 ), (g ρ1 2 ), (g ρ2 2 ), (g η T 1 ), (g µ 2 ), (g δ T •V 1 ), (g γ 2 ), g γ 2 }. • ProbGen(x): given a vector x ∈ F n p , output σ x = x. • Compute(σ x , EK f ): given the encoded input σ x and the evaluation key EK f , it computes y = Ax, ζ = ω T x, s T i = g τ T i 1 X for i ∈ {1, 2} with x X ∈ F c2×c1 p , z T = g η T 1 Y with y Y ∈ F b2×b1 and C = g δ•V 1 X with x X ∈ F d2×d1 p
. Finally it outputs σ y = {y, ζ, (s T i ) 2 i=1 , z, C}. • Verify(σ y , V K f ): it starts by sampling randomly four vectors (v 1 , v 2 ) ∈ F c1 p , v 3 ∈ F b1 p and v 4 ∈ F d1 p . Then for i = 1, 2 it checks whether

s T i v i is equal to g τ T i 1 (X • v 1 )
or not, then it computes

D i = e(s i [j]; g ρ1[j] 2
) for both i. Then it checks whether z T v 3 is equal to

g η T 1 (Y v 3 ) or not, then it computes H = e(z[i]; g µ[i]
2 ). From there it computes ϑ = C v 4 , and it checks whether e(ϑ[i]; g

γ [i] 2 ) is equal to e(g δ T •V 1 (Xv 4 ); g γ 2 )
. If all the previous tests have passed it performs a final check which is to know if e(ζ; g 2 ) is equal to

H • D 1 • D 2 • e(T race(C); g γ 2 )
or not, if it passes it output y, but if any of the previous check has failed it outputs ⊥.

Figure 4: Proven publicly delegatable protocol for matrix-vector product checks and vectorization, we need to use a masking of the form u T A + t T (see Figure 1) indistinguishable from a random distribution, but:

1. We have to add an extra component γδv T to u T A+t T so that it is possible, when proving the reduction to co-CDH, to make up a random vector ω T = g u T A+t where the components of u T A + t T are canceled out. This component cannot be revealed to the Prover nor the Verifier, otherwise its special structure could have been taken into account by the reduction. Also this component cannot have the rank-1 update structure as its has to be a multiple of u T A + t T . Therefore only the protocol of Figure 3 can be used to check the dotproduct with g v T . 2. To be able to apply the analysis of [6, Theorem 3] while allowing fast computations with t, we use a special form for t, namely:

t T = ρ 1 τ T 1 + ρ 2 τ T 2 .
With these modifications we are able to prove the soundness of the protocol in Figure 4. Theorem 8. Let A ∈ F m×n p whose matrix-vector products costs µ(A) arithmetic operations. Protocol 4 is sound under the co-CDH assumption, perfectly complete and its number of performed operations is bounded as follows:

Preparation Prover Verifier F p µ(A)+O (m+n) µ(A) O (m+n) G i O (m+n) O m+n 4/3 O √ m+n 2/3 Pairings 0 0 O ( √ m+ √ n)
The proof of Theorem 8 is given in Appendix E.

Remark 9. Fast matrix multiplication can be used for the computation of C in the protocol of Figure 4. This decreases the O n 4/3 factor of the Prover to O n (1+ω)/3 where ω is the exponent of matrix-matrix multiplication. The currently best known exponent, given in [START_REF] François | Powers of tensors and fast matrix multiplication[END_REF], is ω ≤ 2.3728639. This immediately yields a reduced bound for the Prover of µ(A) + O m + n 1.12428797 . This together with Corollary 7 can produce a protocol with Prover complexity bounded by µ(A) + O m + n 1+0.12428797a , for any 0 < a < 1, while the Verifier complexity is O (n) classical operations and o(n) cryptographic operations.

Conclusion and experiments

We first recall in Table 1 the leading terms of the complexity bounds for our protocols and those of [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF][START_REF] Zhang | Efficient secure and verifiable outsourcing of matrix multiplications[END_REF][START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] (that is each value x in a cell is such that the actual cost is bounded by x + o(x)).

There, we denote the base field operations by • F, the cryptographic group exponentiations or pairing operations by • G, and the cost of a product of the matrix A ∈ F m×n p by a vector is µ(A). We see that our protocols are

• F + mn • G - 2mn • F + 2mn • G Trustee (ProbGen) 2(m + n) • F + 2m • G mn • F + (2m + n) • G n • G Prover (Compute) µ(A) • F + 2mn • G µ(A) • F + 2mn • G µ(A) • F + 2mn • G Verifier 2m • G 2m • G m • G Extra storage O (mn) O (mn) O (mn) Extra communications O (m) O (m) O (1) Scheme Figure 1 Figure 4 Mode Public verif. Public deleg. Preparator (KeyGen) (µ(A) + n) • F + n • G (µ(A) + m + 5n) • F + 2n • G Trustee (ProbGen) 2(m + n + 1) • F + 1 • G 0 Prover (Compute) µ(A) • F + 2n • G µ(A) • F + (2n 4/3 + m) • G Verifier 1 • G (2m + 4n) • F + (6 √ m + 2n 2/3 ) • G Extra storage O (m + n) O (n) Extra communications O (1) O n 2/3 + √ m
suitable to sparse or structured matrix-vector multiplication as they never require O (mn) operations but rather µ(A). Moreover, we see that most of the Verifier's work is now in base field operations were it was cryptographic operations for previously known protocols. As shown in Table 2 and in Figure 5, this is very useful in practice, even for dense matrices. For these experiments we compare with our own implementations of the protocols of [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF][START_REF] Zhang | Efficient secure and verifiable outsourcing of matrix multiplications[END_REF][START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] over the PBC library 1 [START_REF] Lynn | The pairing-based cryptography (PBC) library[END_REF] for the pairings and the FFLAS-FFPACK library 2 [3]

1 https://crypto.stanford.edu/pbc, version 0.5.14 for the exact linear algebra over finite fields (C++ source files are available on request via the PC and will be publicly posted on our web site if the paper is accepted). We used randomly generated dense matrices and vectors and to optimize the costs (pairings are more expensive than exponentiations), we also chose the following parameters for the vectorizations:

b 1 = √ m/10 ; b 2 = 10 √ m ; c 1 = √ n/10 ; c 2 = 10 √ n ; d 1 = n 1/3 /3 ; d 2 = 3n 2/3
. We indeed chose a type 3 pairing over a Barreto-Naehrig curve [START_REF] Paulo | Pairing-friendly elliptic curves of prime order[END_REF] based on a 256-bits prime field, which should guarantee 128 bits of security. First, with F p the 256-bits prime field 3 , G 1 is the group of F p -rational points E(F p ) with parameters: G 1 (E) : y 2 = x 3 + 6, modulo p. Second, G 2 is a subgroup of a sextic twist of E defined over F p 2 denoted E (F p 2 ) with parameters 4 : G 2 (E ) :

y 2 = x 3 + 6e, F p 2 ∼ = F p [X]/(X 2 -2), e = a 0 + a 1 X ∈ F p 2 .
The third group G T is then a subgroup of the multiplicative group of the field F p 12 . This curve is reasonably well-suited to our needs and is supported by the PBC library. We used it for our own implementation of the three protocols [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF][START_REF] Zhang | Efficient secure and verifiable outsourcing of matrix multiplications[END_REF][START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] as well as for ours, Section 5 and Figure 7.

In the first set of timings of Table 2, we also compare the latter protocols with a compiled verifiable version obtained via the Pepper software 5 [12]. This software uses a completely different strategy, namely that of compiling a C program into a verifiable one. We added the timings for n = 1000 as a comparison, but the Pepper compilation thrashed on our 64 GB machine for n ≥ 2000.

In terms of Prover time, we see that our protocols are between two to three orders of magnitude faster than existing ones (further evidence is given in Table 4, Appendix E). Moreover, overall we see that with the new protocol, the data preparation (KeyGen) is now very close to a single non-verified computation and that the work of the Prover can be less than three times that of a non-verified computation (note first, that in both Table 2 and Figure 5, the "Compute" fields include the computation of y = Ax, and, second, that the Prover overhead being asymptotically faster than the compute time, this latter overhead is rapidly amortized). Finally, only the protocol of [START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF] did exhibit a verification step faster than the computation itself for size 2000×2000 whereas, as shown in Figure 5, our protocol achieves this only from size 3000×3000. However, we see that we are competitive for larger matrices. Moreover, as shown by the asymptotics of Theorem 8, our overall performance outperforms all previously known protocols also in practice, while keeping an order of magnitude faster Verification time. j] ).

• Compute: compute y = Ax and z = W x ∈ G m 1 (that is z[i] = n j=1 W [i, j] x[
• Verify: check that e(z[i]; g 2 ) = a y[i] VK x [i], for all i = 1, . . . , m.

In this protocol, the flow of communications is as follows:

1. Preparator: secret random α ∈ F p , s ∈ F m p , t ∈ F n p , ρ ∈ F m p , τ ∈ F n p . then W = g αA+s•t T +ρ•τ T 1 2. Preparator to Prover: A ∈ F m×n p , W ∈ G m×n 1 .
3. Preparator to Trustee: s, ρ, t, τ , in a secure channel.

4. Preparator publishes and signs a = e(g 1 ; g 2 ) α ∈ G T .

5. Verifier to both Prover and Trustee: x ∈ F n p .

Trustee publishes and signs VK

x ∈ G m T such that : VK x = e(g 1 ; g 2 ) s•(t T •x)+ρ•(τ T •x) .
7. Prover to Verifier:

y ∈ F m p , z ∈ G m 1 .
8. Verifier public verification: e(z; g 2 )

? == a y VK x (component-wise in G m T ).
This protocol is sound, complete and publicly verifiable. It however uses many costly exponentiations and pairings operations that renders it inefficient in practice: even though the Client and Trustee number of operations is linear in the vector size, it takes still way longer time that just computing the matrix-vector product in itself, as shown in the experiment Section 6.

In this paper, our aim was first to adapt this protocol to the sparse/structured case, and, second, to reduce the number of cryptographic operations in order to obtain a protocol efficient in practice.

B Probabilistic verification and the random oracle model

First we recall that private verification is very fast and does not require any cryptographic routines. Then we show that this allows to obtain a very efficient protocol in the random oracle model, but for a fixed number of inputs.

B.1 Private verification

Without any recourse to cryptography, it is well known how to privately verify a matrix-vector multiplication. The idea is to use Freivalds test [START_REF] Freivalds | Fast probabilistic algorithms[END_REF], on the left, provided that multiplication by the transpose matrix is possible:

• Verifier to Prover: A, x i ∈ F n p , for i = 1, . . . , k.

• Prover to Verifier: y i ∈ F m p , for i = 1, . . . , k.

• Verifier verification: random u ∈ F m p , then w T = u T • A, and finally check, for i = 1, . . . , k, that w

T • x i ? == u T • y i in F p .
On the one hand, this protocol uses only classical arithmetic and is adaptable to sparse matrices, that is when a matrix vector product costs µ(A) operations with µ(A) < 2mn (this is the case for instance if the matrix is not structured but is sparse with µ(A)/2 < mn non-zero elements). Indeed, in the latter case, the cost for the Prover is kµ(A), where the cost for the Verifier is µ(A) + 4kn.

On the other hand, the protocol has now Freivalds probability of revealing an error in any of the y i : 1 -1/p, if F p is of cardinality p (or 1 -1/p if u is chosen in an extension of degree of F p ).

D Small fields

The protocol of Figure 1 is quite efficient. We have made experiments with randomly generated dense matrices and vectors with the PBC library 1 for the pairings and the FFLAS-FFPACK library 2 for the exact linear algebra over finite fields. For instance, it is shown in Table 3, that for a 8000×8000 matrix over a field of size 256 bits, the protocol is highly practical: first, if the base field and the group orders are of similar sizes, the verification phase is very efficient; second, the overhead of computing ζ for the server is quite negligible and third, the key generation is dominated by the computation of one matrix-vector product.

Table 3: Verification of a 8000×8000 matrix-vector multiplication with different field sizes via the protocol in Figure 1 on a single core @3.4GHz. Differently, if the base field is small, say machine word-size, then having to use cryptographic sizes for the group orders can be penalizing for the Key Generation: multiplying a small field matrix A with a large field vector u T is much slower than y = Ax with x and A small. First of all, the computations must be compatible. For this, one possibility is to ask and verify instead for y = Ax over Z and then to let the Verifier compute y mod p for himself. There, to reduce the overhead of computing u T A, one can instead select the m values of the vector u as u = αr i s j with = i √ m + j for α a randomly chosen large value and r i , s j some randomly chosen small values. Indeed then u T A can be computed by first performing (rs T )A via O ( √ m) matrix-vector computations with s (or a √ m×n √ m matrix-vector multiplication) followed by O (n √ m) multiplications by r (or a n× √ m matrix-vector multiplication) where s j and r i are small values. Then it remains only to multiply a vector of small values by α. We have traded O (mn) operations with large values for O ( √ mn √ m + n √ m) operations with small values and O (n) with large values. Now, in order for the values to remain correct over Z, the value of (u T A + t T )x must not overflow. For this, one must choose a group order larger than mnp 4 (for (rs T )Ax). Now the security is not anymore half the size of the group order but potentially half the size of the set from which t T is selected, that is at most the group order size minus that of np (for t T x). To be conservative we even propose, as an estimated security of the obtained protocol, to consider only half the size of α (that is the size of the group order minus that of mnp 4 ). In terms of efficiency, the improvement is shown in Table 1, last row. On the one hand, the key generation is now dominant and can be amortized only after about 10 matrix-vector multiplications. On the other hand, the verification time starts to be faster than the computation time. This is also shown in Figure 6 where the equivalent of the last row in Table 3 is shown for different matrix dimensions.

E Proven publicly delegatable protocol with negligible cryptographic operations

We first give the proof of Theorem 8, page 8, recalled hereafter, for the correctness, soundness and complexity of the protocol in Figures 4 and7. The flow of exchanges within our protocol is also illustrated in Figure 7. ) for j = 1, 2.

Finally, the last check is that these two parts, as well as the last one, which is e(g δ 1 ; g γ 2 ) v T •x = e(g 1 ; g 2 ) a(b+r)v T •x = e(T race(g δ•V X 1 ); g γ 2 ), are coherent with the definitions of ω and ζ above. Now, with a non-negligible probability , adversary A can pass the check for another y = y, by providing an adequate ζ . First, z T , s T 1 , s T 2 and C must be correct, as they are checked directly and independently by the Freivalds first four checks. Second, we have that e(ζ ; g 2 ) = e(g 1 ; g 2 ) abu T 0 y +abt T 0 x e(g δ 1 ; g γ 2 ) v T •x and therefore, we must also have e(ζ(ζ ) -1 ; g 2 ) = e(g 1 ; g 2 ) abu T 0 (y-y ) . As u 0 is a secret unknown to adversary A, for a random y the probability that u . The only other possibility is that adversary A was able to recover u T 0 . But that would directly implies that it has an advantage in the co-CDH: g η 1 = L η0 , g µ 2 = R µ0 .

In Table 4, we present more timings for the comparison between our protocol and, to our knowledge and according to Table 2, the best previously known from [START_REF] Elkhiyaoui | Efficient techniques for publicly verifiable delegation of computation[END_REF]. The associated speed-ups supports our claim of a Prover efficient protocol with a gain of two orders of magnitude.

  a randomized algorithm run by a Trustee which takes as input an element x in the domain of the outsourced function f . It returns σ x , an encoded version of the input x. • Compute(σ x , EK f ) → (σ y ): an algorithm run by the Prover to compute an encoded version σ y of the output y = f (x) given the encoded input σ x and the evaluation key EK x . • Verify(σ y , V K f ) → y or ⊥: given the encoded output σ y and the verification key V K f , the Verifier runs this algorithm to determine whether y = f (x) or not. If the verification passes it returns y otherwise it returns an error ⊥.

Figure 1 :

 1 Figure 1: Interactive protocol for Sparse-matrix vector multiplication verification under the co-CDH.

  operations in a group can be slightly faster than pairings. Moreover √ m 2 can be quite far off m. Therefore it might be interesting to use a non square vectorization b 1 ×b 2 , as long as b 1 b 2 ≥ m and b 1 + b 2 = Θ ( √ m) (and 0 padding if needed). Then we have U ∈ F b1×b2 p , µ ∈ F b1 p , η ∈ F b2 p , Y ∈ F b2×b1 p and z ∈ G b1 2 . The obtained protocol can compute e(g 1 ; g 2 ) u T •y with O ( √ m) cryptographic operations on the Verifier side and is given in Figure 2. Lemma 5. The protocol of Figure 2 for publicly delegation of a size m external group dot-product verification with rank-1 left hand side is sound, perfectly complete and requires the following number of operations where b 1 b 2 ≥ m and b 1 + b

1 Y

 1 group) operations; for z T v: O (b 1 ) cryptographic (group) operations; and for b1 i=1 e(z[i]; g µ[i] 2 ): O (b 1 ) cryptographic (pairings) operations. Then the preparation requires to compute g η 1 ∈ G b2 1 and g µ 2 ∈ G b1 1 , while the Prover needs to compute g η T for Y ∈ F b2×b1 p and b 1 b 2 = O (m).

  or not. If the test passes it returns T race(C) and if it fails it returns ⊥.

Figure 3 :Lemma 6 .

 36 Figure 3: Publicly delegatable protocol for the external dot-product.

2 1 Y, in O b 2 1 b 2 =

 12 ): O (b 1 ) cryptographic (pairings) operations; Then the preparation requires to compute w T ×U . This is O (b 1 b 2 = m) operations. Finally, the Prover needs to compute the matrix multiplication g U for U ∈ F b1×b2 p and Y ∈ F b2×b1 p O (mb 1 ). Therefore, one can take b 1 = O ( 3 √ m) and b 2 = O m 2/3 which gives only O m 2/3 cryptographic operations for the Verifier, and O m 4/3 cryptographic operations for the Prover. Now a dot-product can be cut in n k chunks of size k.

  dense matrix-vector product in a 256-bit finite field on a single core @3.4GHz Compute (y=Ax and overhead) server overhead KeyGen y=Ax Verify

Figure 5 :

 5 Figure 5: Protocol 4 performance.

Theorem 8 .

 8 Let A ∈ F m×n p whose matrix-vector products costs µ(A) arithmetic operations. Protocol 4 is sound under the co-CDH assumption. It is also perfectly complete and its number of performed operations is bounded as follows:

  T 0 (y -y ) = 0 is bounded by 1/|G 1 | and thus negligible. Thus adversary B can compute c ≡ u T 0 (y -y) -1 mod |G 1 | and (ζ/ζ ) c = g ab 1 . Therefore it breaks the co-CDH assumption with nonnegligible probability

Table 1 :

 1 Leading terms for the time and memory complexity bounds (exchange of A, x and y excluded).

	Scheme	[6]	[14]	[4]
	Mode	Public verif.	Public deleg.	Public deleg.
	Preparator (KeyGen)	2mn

Table 2 :

 2 Matrix-vector multiplication public verification over a 256-bit finite field with different protocols on a i7 @3.4GHz.

			1000×1000					2000×2000
		[12]	[6]	[14]	[4]	Fig. 4	[6]	[14]	[4]	Fig. 4
	KeyGen	141.68s 152.62s	-154.27s 0.80s	615.81s	-	612.72s 1.75s
	ProbGen	-	1.25s	2.28s	2.30s		-	2.13s	4.98s	4.56s	-
	Ax = y	20.14s	0.19s	0.19s	0.19s 0.19s	0.78s	0.78s	0.78s 0.78s
	Compute 188.60s 273.06s 433.88s 271.03s 2.26s 1097.96s 1715.46s 1079.71s 5.37s
	Verify	2.06s	26.62s	27.56s	0.33s		0.90s	52.60s	55.79s	0.62s	1.19s
			4000×4000					8000×8000
		[6]	[14]	[4]	Fig. 4		[6]		[14]	[4]	Fig. 4
	KeyGen	2433.10s	-2452.98s	4.89s	9800.42s	-	9839.26s 15.64s
	ProbGen	3.81s	13.29s	9.24s	-	7.41s	43.44s	18.46s	-
	Ax = y	3.28s	3.28s	3.28s	3.28s	13.30s	13.30s	13.30s 13.30s
	Compute 4360.43s 6815.40s 4329.46s 13.76s 17688.69s 27850.90s 17416.38s 37.00s
	Verify	103.14s	107.99s	1.20s	1.65s	211.07s	220.69s	2.37s	2.25s

Table 4 :

 4 Speed-up of our novel Protocol over a 256-bit finite field on a i7 @3.4GHz.

	Size	100	200	500	1000	2000	3000	4000
	[4]	2.77s 10.93s 67.93s 271.03s 1079.71s 2430.05s 4329.46s
	Fig. 7	0.17s	0.34s	0.98s	2.26s	5.37s	9.16s	13.76s
	Speed-up	17	32	69	120	201	265	315
	Size	5000	6000	7000	8000	9000	10000
	[4]	6790.15s 9780.24s 13309.61s 17416.38s 22002.51s 27175.12s
	Fig. 7	18.55s	24.03s	29.93s	37.00s	44.00s	51.97s
	Speed-up	366		407	445	471	500	523

* This work is partly funded by the OpenDreamKit Horizon 2020 European Research Infrastructures project (#676541).

A Fiore and Gennaro's protocol

For the sake of completeness, we present here the original protocol for matrix-vector verification in [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF], but with our rank-one update view. It stems from the fact that if s, t, ρ, τ are randomly generated vectors then the function g M [i,j] , where M = s • t T + ρ • τ T , is a pseudorandom function [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF]Theorem 3], provided that the Decision Linear assumption [6, Definition 3] holds (a generalization of the External Diffie-Hellman assumption for pairings).

• KeyGen: for A ∈ F m×n p , generate 3 multiplicative groups (G 1 , G 2 , G T ) of prime order p, with g 1 generating G 1 (resp. g 2 generating G 2 ), and a bilinear map e :

∈ G 1 , give W to the server and publish a = e(g α 1 , g 2 ) ∈ G T . Let x ∈ F n p be a query vector

B.2 Public verification in the random oracle model

Using Fiat-Shamir heuristic [START_REF] Fiat | How to prove yourself: Practical solutions to identification and signature problems[END_REF], the privately verifiable certificate of Section B.1 can be simulated noninteractively: uniformly sampled random values produced by the Verifier are replaced by cryptographic hashes (to prove security in the random oracle model) of the input and of previous messages in the protocol. Complexities are preserved, as producing cryptographically strong pseudo-random bits by a cryptographic hash function (e.g., like the extendable output functions of the SHA-3 family defined in [START_REF] Bertoni | Sponge-based pseudo-random number generators[END_REF][START_REF]standard: permutation-based hash and extendable-output functions[END_REF]), is linear in the size of both its input and output (with atomic operations often even faster than finite field ones):

• Preparator to Prover: A ∈ F m×n p .

• Verifier to Prover: x i ∈ F n p , for i = 1, . . . , k.

• Prover to Verifier: y i ∈ F m p , for i = 1, . . . , k.

• Verifier to Trustee: all the x i and y i .

• Trustee publishes and signs both u ∈ F m p and w ∈ F n p such that:

• Verifier public verification:

There is absolutely no overhead for the Prover; the cost for the Trustee is a single matrix-vector product for any k plus a cost linear in the input size; and the cost for the Verifier is O (nk). Using Fiat-Shamir heuristic this gives a possibility for an afterwards public verification (that is after the computations), but this not possible to test new vectors once u has been revealed.

C Proof of Theorem 4

We here give the proof of Theorem 4, page 4, recalled hereafter. Proof. For the correctness, we have that:

Then, by bilinearity, e(ζ i ; g 2 ) = e(g 1 ; g 2 ) hi+di = η i .

For the soundness, a malicious Prover can guess the correct output values, but this happens once in the number of elements of G T . Otherwise he could try to guess some matching h i and d i , but that happens less than one in the number of elements of F p . Finally, the Prover could produce directly ζ i . Suppose then it is possible to pass our verification scheme for some A, x and y = y = Ax. Then without loss of generality, we can suppose that the first coefficients of both vectors are different, y [1] = y [START_REF] Paulo | Pairing-friendly elliptic curves of prime order[END_REF] (via row permutations) and that y

Take a co-computational Diffie-Hellman problem (g c 1 , g d 2 ), where g cd 1 is unknown. Then denote by a = e(g c 1 ; g d 2 ) = e(g cd 1 ; g 2 ) and consider the vector z T = [a, e(1; 1), . . . , e(1;

Compute also the vector φ T = e(ω T ; g 2 )/χ T coefficient-wise. The latter correspond to φ T = e(g t T 1 ; g 2 ) for

. Now, if it is possible to break the scheme, then it is possible to compute ζ that will pass the verification for y as Ax, that is e(ζ ; g 2 ) = (z T y )(φ T x). Let h = u T y, d = t T x and h = u T y . Then e(ζ; g 2 ) = e(g h 1 ; g 2 )e(g d 1 ; g 2 ) and e(ζ ; V T requires O (n) group operations while g τi 1 , g ηi 1 , g ρi 2 , and

T and (g δ 1 )

T V . The complexity for the Verifier is

Finally for the soundness, assume that there is an adversary A that breaks the soundness of our protocol with non-negligible advantage for a matrix A ∈ F m×n p . In the following we will prove how an adversary B can use adversary A to break the co-CDH assumption with non-negligible advantage . Let assume that B was given a co-CDH sample (L = g a 1 , R = g b 2 ). First B simulates the soundness experiment to adversary A in the following manner: when A calls the oracle O KeyGen , adversary

p and a value r ∈ F p . We let u 0 be the vector representation of µ 0 • η T 0 and t 0 that of

T and (g δ 1 )

T * V = (g γ 1 )

T V and outputs:

1 ), (g

), (g γ 2 ), g γ 2 }. Thanks to the randomness and the decisional Diffie-Hellman assumption (DDH) in each group G i , as well as [START_REF] Fiore | Publicly verifiable delegation of large polynomials and matrix computations, with applications[END_REF]Theorem 3] for ω T = (L r ) v T , these public values are indistinguishable from randomly generated inputs. Further, we have

When adversary A calls the oracle O P robGen on input x, adversary B returns σ x = x. Therefore, if y = Ax and ζ = ω T x, then the verification will pass: indeed the first two checks will ensure that s

X when the third check ensures that z T = g η T

1

Y . This shows that:

T , (g δ 1 )

T V Figure 7: Exchanges in the proven publicly delegatable protocol with negligible cryptographic operations of Figure 4.


Iosif Baudouin Denis De Senneville
email: baudouin@imf.u-bordeaux2.fr

Alexandru Mendichovszky

M.Sc Sébastien Roujol

MD Isky Gordon

PhD Chrit Moonen

MD Nicolas Grenier

PhD Iosif Alexandru Mendichovszky

Improvement of MRI-functional measurement with automatic movement correction in native and transplanted kidneys

Keywords: MRI-Functional measurement, kidney perfusion, glomerular filtration rate, movement correction

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

F O R P E E R R E V I E W O N L Y F O R P E E R R E V I E W O N L Y F O R P E E R R E V I E W O N L Y

INTRODUCTION

Virtually all diseases of the kidney affect perfusion and glomerular filtration. Noninvasive and accurate measurement of both perfusion and glomerular filtration rate (GFR) could have a major impact in understanding renal physiopathology and for serial monitoring of the course of both acute and chronic kidney diseases. Dynamic contrast-enhanced (DCE) magnetic resonance imaging (MRI) is advocated to evaluate these functional parameters. However, publications in the literature show poor correlation when MRI-GFR has been compared with GFR measured by reference methods [START_REF] Buckley | Measurement of single kidney function using dynamic contrast-enhanced MRI: comparison of two models in human subjects[END_REF][START_REF] Lee | Renal function measurements from MR renography and a simplified multicompartmental model[END_REF]; this poor correlation precludes the use of DCE-MRI for GFR estimation in daily clinical practice. The inaccuracy is multi-factorial with unsolved problems regarding the ideal acquisition sequence; the dual MR effect (T1 and T2*) of contrast agents; the conversion of signal intensity into concentration; the pharmacokinetic models applied, as well as the difficulties in post-processing (segmentation and region of interest). DCE-MRI images are usually acquired during spontaneous breathing that will result in kidney movement.

Since respiratory-gated sequences (3) would lead to loss of temporal resolution, most groups studying DCE-MRI GFR have either repositioned images manually or ignored movement. However, movement causes artefacts in the pixel-based time-intensity analysis which will lead to inaccurate GFR quantification.

A compromise in choosing the acquisition parameters is required to achieve a sufficiently high signal-to-noise ratio (SNR). An ideal rapid isotropic 3D imaging of the moving kidneys is not achievable, and thus slices are usually oriented along the long axis of the kidney that will allow the best approximation of movement to be estimated. During DCE-MRI, voxel intensity changes may be due to one or more of the following factors: rapid, non-uniform movement of contrast within the different renal compartments; motion either due to respiration or physical

F O R P E E R R E V I E W O N L Y
movement, or low SNR due to the MRI acquisition sequence required for rapid acquisition.

Although movement can be addressed using image registration algorithms (4) applied in a postprocessing step, this is not an option with the kidney as the image amplitude and contrast change with time due to the transit of the contrast agent through the kidney following a bolus injection [START_REF] Yim | Registration of time series Gd-contrast enhanced magnetic resonance images for renography[END_REF].

The purpose of our study was to improve a 2D region tracking software for retrospective motion correction without sacrificing temporal resolution in quantitative renal DCE-MRI studies, code it into a radiologically useful software, and evaluate it using the Patlak-Rutland tracer kinetic model [START_REF] Patlak | Graphical evaluation of blood-to-brain transfer constants from multiple-time uptake data[END_REF][START_REF] Rutland | A single injection technique for subtraction of blood background in 131 Ihippuran renograms[END_REF][START_REF] Hackstein | Measurement of single-kidney glomerular filtration rate using a contrast-enhanced dynamic gradient-echo sequence and the Rutland-Patlak plot technique[END_REF][START_REF] Annet | Glomerular Filtration Rate: Assessment with dynamic contrast-enhanced MRI and a cortical compartment model in the rabbit kidney[END_REF]. This software was then applied to two different populations (healthy volunteers and renal transplant patients) in different institutions using different MRI scanners and different acquisition protocols. These protocols had different spatial resolutions and different levels of SNR. Motion during data acquisition may cause tracer kinetic model fitting errors in the post-processing step. Quantification of those model fitting errors may provide a quality criterion on motion correction and recent work demonstrated that it could even be used to drive the registration procedure [START_REF] Buonaccorsi | Comparison of the Performance of Tracer Kinetic Model-Driven Registration for Dynamic Contrast Enhanced MRI Using Different Models of Contrast Enhancement[END_REF]. In this study, quantification of motion correction accuracy is based on tracer kinetic models. As this quantification depends heavily on the chosen model, we selected to evaluate the improvement after movement correction using the Patlak-Rutland model (see Appendix 1). Improvements in GFR shown with this registration algorithm might also apply when using other tracer kinetic models.

All computations were performed in 2D on the middle slice of each kidney. [START_REF] Buckley | Measurement of single kidney function using dynamic contrast-enhanced MRI: comparison of two models in human subjects[END_REF] where is a criterion determining the accuracy of the registration. With image registration there is a compromise that has to be made between the complexity of the estimated spatial transformation T and the accuracy of the estimation strategy. Perturbations in image registration occur predominantly when the true motion field violates the brightness consistency model used for its approximation. In our case, because the contrast varies within the kidney with time, the algorithm has to be independent of the contrast changes within the kidney and also assumptions on the estimated kidney deformation have to be made. This work assumes that the kidney is a rigid body and its shape does not change during the MRI data acquisition. The new software takes into account that the two kidneys do not move equally with respiration and thus each kidney is dealt with separately during the correction process. Furthermore, the software has been developed to be virtually independent of operator interaction.

F O R P E E R R E V I E W O N L Y 2D motion

Proposed registration model

F O R P E E R R E V I E W O N L Y
The implemented registration model was proposed by Sun et al (11). This approach, based on kidney edge invariance, was found to be the most promising technique for kidney registration during bolus passage (see Appendix 2). In its original form, the algorithm was based only on the estimate of translation displacements with a pixel resolution which might be insufficient for low amplitude movement as in transplants. Our implementation allows the estimation of a rigid transformation (translation + rotation) with subpixel accuracy.

Original DCE MR images were loaded into a home based registration software (implemented in C++). Only the kidneys and aorta are required for functional renal analysis, so image masking was performed, eliminating un-necessary anatomical structures. A region of interest (ROI) encompassing the renal parenchyma, excluding the renal sinus, was manually drawn by the radiologist on an enhanced image, taken during the vascular phase, when enhancement of the renal cortex is maximal. This image is the reference image set for motion

F O R P E E R R E V I E W O N L Y
An exhaustive search was performed to determine parameters of T that maximises . A multi-step scheme was used to perform a subpixel registration while reducing enumeration of possible solutions. First, all possible values for the three parameters of T were explored with an accuracy of one pixel for translation and one degree of accuracy for rotation. Respiratory displacement of the kidney is mainly a cranio-caudal shift. Based on preliminary observations, the size of the search space was restrained to values of 31 pixels for cranio-caudal direction, 11 pixels for left-to-right direction and 5 degrees for in-plane rotations (these typical values were defined with respect to maximal kidney amplitude and image resolution in our data). Following this, an exhaustive search was performed with a step of 0.25 pixels for translation parameters to further enhance the accuracy of the registration to a sub-pixel level. The grey level intensity of pixels of the registered image was computed with a bilinear interpolation from the original brightness values. The ROI could then be propagated to all images of the time series.

If the motion had been accurately corrected, the kidney border would coincide with the contour of the original kidney ROI. An experienced operator reviewed the re-aligned data set and made manual small adjustements, based on a visual analysis, in case of small misalignements.

When the realignment was visually considered as impossible because of intra-scan or excessive out-of-plane movements, the image was rejected. The corrected images were exported back into their native format and loaded into a dedicated software package MIStar (Apollo Medical Imaging, Melbourne, Australia) for processing using the Patlak-Rutland model.

Patients

The movement correction algorithm was tested on DCE-MRI renal studies obtained in 2 different populations: 10 patients following kidney transplantation (age ranging from 24 to 63 years, mean 45) with variable renal function (creatinine clearance between 17 and 59.4 ml/min,

F O R P E E R R E V I E W O N L Y
mean 44.15 ml/min); and 10 healthy volunteers (age ranging from 23 to 36 years, mean 29.2).

Both studies received approval from the local ethics committees and all participants gave written informed consent for taking part in the studies.

MR acquisition protocol:

Following kidney transplantation, the patients underwent the MRI examination on a 1.5 T Philips system (ACS-NT, Philips Medical System, Best, The Netherlands) using a body phasedarray coil and a 3D SR-TFE pulse sequence, without fat saturation, with the following parameters: TR = 4.4 ms, TE = 2.5 ms, TI = 120 ms, FOV = 400×400 mm², FA = 10°, resulting in a temporal resolution of 1.5 seconds (per kidney volume). To ensure complete coverage of the kidney 5 slices (10 mm thick, no gap) were acquired for each dynamic volume, with an in-plane resolution of 128×50 pixels and a resulting voxel size of 3.2×8×10 mm 3 . During the functional scan a dose of 0.03 mmol (0.06 ml/kg) of Gd-DOTA (Dotarem®, Guerbet Group, Aulnay-sous-Bois, France) was injected, followed by a 20 ml flush of saline, both with a 2ml/sec injection rate, using an automatic injector (MedRad). This acquisition contained a total of 200 dynamic volumes over a total time of 300 seconds.

The healthy volunteers were scanned on a 1.5 T Siemens scanner (Avanto, Siemens Medical Solutions, Erlangen, Germany) with a dedicated abdominal TIM coil. The dynamic contrast-enhanced acquisition was performed using a gradient-echo 3D-FLASH pulse-sequence (VIBE) with the following parameters: TR = 1.63 ms, TE = 0.53 ms, flip angle = 17, strong fat saturation, PAT factor = 2 (GRAPPA), FOV = 400×325 mm², 18 slices covering the entire kidney, 7.5 mm slice thickness, no gap. The resulting voxel size was 3.1×3.1×7.5 mm 3 and each dynamic volume was acquired every 2.5 seconds. A 0.05 mmol (0.1 ml) /kg body weight dose of Gd-DTPA (Magnevist, Schering, Germany) was injected as a bolus at 2 ml/second injection rate

F O R P E E R R E V I E W O N L Y
using an automatic injector (Spectris). The contrast agent bolus was immediately followed by a 15 ml saline flush injected at the same speed. A total of 138 volumes were acquired in 345 seconds.

During the scans both patients and volunteers were asked to breathe normally and lie relaxed in the scanner. In both groups the slices were positioned in an oblique-coronal plane (along the long axis of the kidney) to minimise through-plane movement and ensure the presence of the aorta on at least one slice, needed for subsequent analysis using the Patlak-Rutland model (see Appendix 1).

Assessment of Movement Correction:

A manual ROI was generated for the cortex of each kidney on both sets of images. This or kidney [START_REF] Patlak | Graphical evaluation of blood-to-brain transfer constants from multiple-time uptake data[END_REF][START_REF] Rutland | A single injection technique for subtraction of blood background in 131 Ihippuran renograms[END_REF][START_REF] Hackstein | Measurement of single-kidney glomerular filtration rate using a contrast-enhanced dynamic gradient-echo sequence and the Rutland-Patlak plot technique[END_REF][START_REF] Annet | Glomerular Filtration Rate: Assessment with dynamic contrast-enhanced MRI and a cortical compartment model in the rabbit kidney[END_REF]. GFR was calculated using this model on a voxel-by-voxel basis providing two 2D functional GFR maps with GFR values (noted K p ) within the cortex, and a standard deviation map (noted (K p )). (K p ) relates the uncertainty on GFR computation and was thus used as a quality criterion in our study. For each kidney a mean of K p (noted p K ) and a mean of

F O R P E E R R E V I E W O N L Y
Kidney fMRI improvement with motion correction 10

(K p ) (noted ) ( p K
) over the cortical ROI were calculated using the time period of 60 -120 seconds post Gd contrast injection, considering time zero the first rise (more than 2 standard deviations) of the AIF signal from the baseline. These results were then compared between the uncorrected and movement corrected data.

As a number of questions relating to accurate GFR computation still remain unanswered (how to accurately evaluate kidney volumes? do transplant and native kidneys function identically as far as perfusion and filtration are concerned?), only relative variations of K p and ) ( p K were analysed.

Statistics

) ( p K
values from uncorrected and movement corrected data were compared using a Student t-test, considering as significant a p-value < 0.05.

RESULTS

Computation time on an Athlon 3.2 Ghz with 1.5 GB of RAM to correct a masked image series was 15 and 7 minutes for native and transplanted kidneys, respectively.

Visual assessment of the movement corrected data showed that manual repositioning was rarely required (in less than 1 % of all images in both groups) and few images were considered as too corrupted and rejected (less than 1 % of all images in both groups were rejected). Only a few seconds were thus necessary for this manual step.

Patients with renal transplants

F O R P E E R R E V I E W O N L Y
The typical amplitude of estimated renal displacement of transplanted kidneys didn't exceed one pixel, clearly indicating the importance of sub-pixel accuracy in the registration process.

Fluctuations in signal intensity time-curves, averaged over the cortical ROI were reduced in all cases following movement correction (illustrated in Figure 1). This is also reflected in the reduced dispersion of data on the Patlak-Rutland plot as shown in Figure 2. Figure 3 illustrates the 2D GFR and standard deviation maps before and after motion correction for the same patient. p K and

) ( p K
values for each patient before and after motion correction are presented in Figure 4. The average absolute variation of p K for the entire group was 6.4 % ± 4.8 % (max = 16.6 %). There was a significant reduction of

) ( p K
on motion corrected data sets (p = 0.003 Student t-test) compared to the non-corrected ones, with an average reduction of ) ( p K of 6.9 % ± 6.6 %. One patient (number 4) had a (K p ) of 21.4 % as shown on Figure 2.

Healthy volunteers with native kidneys

Right and left kidneys have been analyzed individually. Before motion correction, mean amplitude of estimated renal displacement of native kidneys was 6 pixels on the right side and 8 pixels on the left side. Figure 5 illustrates an example of the 2D GFR and standard deviation maps before and after motion correction. ) ( p K was 30.9 % ± 17.6 % (max = 60.8 % for patient number 4) for the right kidney, and 31.8 % ± 14 % (max = 55.3 % for patient number 1) for the left kidney.

DISCUSSION

The clinical use of DCE-MRI in renal studies has led to publications on quantification of the MR signal for both perfusion and GFR. None of the publications to date have applied any motion correction algorithm routinely to the data, yet this variable must be taken into account for any accurate quantification. Furthermore, there is no published evaluation of the accuracy of an automatic movement correction method applied to renal DCE-MRI. The registration method used in this study is based on the one proposed by Sun et al [START_REF] Sun | Integrated registration of dynamic renal perfusion MR images[END_REF]. To increase the accuracy, our algorithm estimates a rigid body transformation with subpixel accuracy. Changes in signal intensity due to bolus passage within the kidney do not affect the efficiency of our approach.

Visual assessment of the movement corrected data showed little movement in more than 98 % of tested data sets.

Accuracy of image registration algorithms can easily be quantified on synthetic or phantom data as the motion is fully controlled. However, actual organ displacement in-vivo is unknown. The approach whereby comparisons are made between corrected data using a postprocessing algorithm versus manual correction by several operators is questionable. Furthermore, using several operators as a reference method introduces important uncontrollable variables that

F O R P E E R R E V I E W O N L Y
GFR uncertainty in transplanted kidneys and 60.8 % in native kidneys). The mean of p K and ) ( p K variations for the transplant group are low (resp. 6.4 % and 6.9 %) compared to the standard deviation (resp. ± 4.8 % and ± 6.6 %). This relates to a wide dispersion of values around the mean (maximum resp. 16.6 % and 21.4 %) due to the different degrees of movement from one patient to another. Therefore, while a variation of less than 5 % was obtained in 4 patients for p K and in 5 patients for ) ( p K , a significant variation (greater than 10 %) was obtained on 2 patients for p K values and 3 patients for

) ( p K
values. In the native kidney group, p K variation > 5 % was obtained in 7 volunteers on the right kidney and in 9 volunteers on the left kidney; also a ) ( p K variation > 10 % was obtained in 9 volunteers on the right and in all volunteers on the left. These results suggest that motion correction is a necessary pre-requisite for quantification of DCE renal MRI in all cases.

Full analysis of the displacement of the kidneys requires a 3D image registration. The ideal image plane of renal DCE-MRI is along the long axis of the kidney which is the main axis of respiratory movement. The importance of a 3D correction algorithm can not be assessed from this study, but as the majority of the movement is along the long axis of the kidney, one may suggest that its added value might be small. As the transplant kidney lies at some distance from the diaphragm, movement correction in this clinical situation becomes less critical. This is reenforced by the results of this study where the improvement of GFR uncertainty is much smaller in the transplant group than in the normal volunteers. However, in transplanted kidneys number 3 and 9 some residual motion outside the imaging plane could explain why the motion correction didn't show improvement of GFR uncertainty (see Figure 4b).

F O R P E E R R E V I E W O N L Y

Comparison with other existing approaches for kidney registration

Althought navigator echoes can be used to estimate kidney displacement (13), motion information is in this case restricted to translational motion in the direction of the navigator.

Alternatively, kidney displacement can be estimated on MR images using image registration algorithms and three main classes of approaches can be distinguished in the literature (details for each approach are given in Appendix 2):

1) methods based on grey level intensity conservation, 2) Fourier based approaches, 3) methods based on geometrical characteristics invariance.

Methods based on grey levels intensity conservation and Fourier-based approaches are both very sensitive to the determination of the mask size in which the registration is computed.

An extension of this mask allows inclusion of surrounding anatomical structures to the registration process. Tissues closely surrounding the kidneys are not affected by bolus passage that may help the registration process, further efficiency of the registration process occurs by restraining computation to the desired part of the image. A delicate compromise must thus be found for the mask size in each patient study and each MR acquisition sequence. Critically, these two approaches are both sensitive to voxel intensity variation. This is a major drawback as the signal within the kidney changes significantly with time over the critical period required for analysis (due to the passage of the contrast agent through the kidney). Thus the change in voxel intensity may be due either to renal movement or to the normal passage of contrast within the kidney. This adds to the difficulty of absolute movement estimation (one identical image is always used as the reference image for the registration of the time series). Although these approaches to image registration can be sufficient in some cases (when contribution of contrast modification induced by the bolus is small compared to contribution of surrounding anatomical

F O R P E E R R E V I E W O N L Y
structures, as, for example, when no fat suppression is performed during the MR acquisition procedure), incremental estimation ( 14) is generally required (the displacement of the kidney is determined between subsequent images, the first image of each pair being the reference image).

For each image, a visual verification is required by the user and manual adjustment of the location of the mask with respect to the kidney is performed when needed. Automatic movement correction requires a more robust approach.

Contrary to the two previously mentioned approaches, methods based on geometrical characteristics invariance are fully stable with respect to the dynamic movement of the contrast agent through anatomical structures [START_REF] Yim | Registration of timeseries contrast enhanced magnetic resonance images for renography[END_REF]. This approach overcomes the problems related to voxel intensity variation during DCE-MRI acquisition and kidney motion estimation is not affected by the passage of contrast agents. Therefore, absolute motion estimation is feasible, allowing reduction in operator errors (especially for transplanted kidneys) and reduction in processing time required for movement correction by the radiologist. Limiting the mask size helps in restraining possible interferences generated by surrounding anatomical structures (induced, for example, by possible displacement outside the imaging plane) and increasing registration accuracy by limiting the computation to the desired part of the image.

No movement correction algorithm is perfect for DCE-MRI studies and the drawbacks of these techniques are their sensitivity to noise, change in kidney shape during acquisition or magnetic field homogeneity changes. Low SNR values in the renal parenchyma were only observed for dynamics acquired before bolus injection when using strong fat suppression in the volunteer MR protocol. Therefore, noise sensitivity does not represent a major limitation.

Changes in kidney shape during MR acquisition occur when motion outside the imaging plane and/or intra-scan motion artefacts are present. Such images are usually severely degraded and are removed from the analysis process. In this study, with 2 different MR acquisition protocols, this

F O R P E E R R E V I E W O N L Y
phenomenon occurred in less than 1 % of the total acquired images. Changes in magnetic field homogeneity at different time points in the respiratory cycle are difficult to correct and were not addressed in our study. Such changes may also lead to image distortion (in particular with long EPI acquisition train lengths).

In conclusion, the implemented registration method appeared to be efficient for movement correction for DCE-MRI in both native and transplanted kidneys. The correction method was independent of the MR acquisition parameters (SNR, spatial resolution, fat suppression) and operator intervention. Motion correction is a necessary pre-requisite to improve quantification of renal functional parameters. The results showed that the correction of native kidney displacements allowed a significant uncertainty reduction on the computed GFR with the Patlak-Rutland plot technique. However, the gain is modest in pelvic renal transplants due to low displacement amplitude. Other limiting factors such as the linearity of the signal with gadolinium concentration, the cortical volume measurement, and an accurate ROI positioning have to be solved. Then, a quantification of the error reduction in GFR measurement could be performed in future using a reference method and applied to a population with a large range of GFR values.

Although registration has been tested in 2D, all observations will be tested in the 3D case when technological progress in rapid MR acquisition sequence will allow sufficient spatial and temporal resolution.

F O R P E E R R E V I E W O N L Y APPENDIX 1: Patlak-Rutland model
The Patlak-Rutland plot technique describes a two compartment model with unilateral tracer flow from compartment 1 (the vascular space) into compartment 2 (nephron space) [START_REF] Patlak | Graphical evaluation of blood-to-brain transfer constants from multiple-time uptake data[END_REF][START_REF] Rutland | A single injection technique for subtraction of blood background in 131 Ihippuran renograms[END_REF][START_REF] Hackstein | Measurement of single-kidney glomerular filtration rate using a contrast-enhanced dynamic gradient-echo sequence and the Rutland-Patlak plot technique[END_REF][START_REF] Annet | Glomerular Filtration Rate: Assessment with dynamic contrast-enhanced MRI and a cortical compartment model in the rabbit kidney[END_REF]. The following assumptions are made:

-The interstitial space as a third space is neglected.

-Signal change is proportional to the concentration of gadolinium in a voxel.

-Prompt and complete gadolinium mixes inside the compartments.

-Hematocrit is constant in all renal vessels and the aorta.

-The gadolinium concentration in the aorta and the renal arteries is equal at any time.

The amount of gadolinium in the renal parenchyma R(t) can be expressed as the sum of gadolinium in the vascular space B(t) and the nephron Q(t) :

) ( ) ( ) ( t Q t B t R + = [3]
Two assumptions are then made:

-The amount of gadolinium the vascular space B(t) is proportional to the concentration of gadolinium in the aorta b(t).

-The amount of gadolinium filtered into the nephron is proportional to the integral of the gadolinium concentration curve of the aorta.

Those two assumptions can be mathematically expressed as follow:

= = t dx x b c t Q t b c t B 0 2 1 1 ) ( . ) ( ) ( . ) ( [4]
c 2 is equivalent to the gadolinium clearance from the vascular space into the nephron. The combination of equations [START_REF] Moricawa | Feasibility of simple respiratory triggering in mr-guided interventional procedures for liver tumors under general anesthesia[END_REF] and [START_REF] Maintz | A survey of medical image registration[END_REF] leads to:

F O R P E E R R E V I E W O N L Y + = t dx x b c t b c t R 0 2 1 ) ( . ) ( .
) ( [START_REF] Yim | Registration of time series Gd-contrast enhanced magnetic resonance images for renography[END_REF] When equation [START_REF] Yim | Registration of time series Gd-contrast enhanced magnetic resonance images for renography[END_REF] is divided by b(t), the Patlak-Rutland plot equation is obtained :

) ( ) ( . ) ( ) ( 0 2 1 t b dx x b c c t b t R t + = [6]
Equation [START_REF] Patlak | Graphical evaluation of blood-to-brain transfer constants from multiple-time uptake data[END_REF] is computed for each dynamics in the 60-120 seconds time frame. The obtained plot leads to a straight line with a slope (c 2 ) equal to the GFR (see figure 2). c 2 is then determined using a least square fit. Computation of c 2 is performed for each voxel of the renal parenchyma to obtain a 2D GFR map (noted K p ). Quality maps (noted (K p )) are then obtained by computing, for each voxel, the standard deviation of the differences between measured and fitted values. A mean of K p (noted p K ) and of (K p ) (noted ) ( p K ) over the cortex can then be calculated.

APPENDIX 2: Existing approaches for kidney registration

Three main classes of approaches for kidney tracking during bolus passage have been observed in the literature:

Methods based on grey levels intensity conservation

Those methods assume that M in equation [START_REF] Buckley | Measurement of single kidney function using dynamic contrast-enhanced MRI: comparison of two models in human subjects[END_REF] is a similarity criterion, computed on grey levels intensities, between the part of a reference image containing the kidney and the area it overlays in the image to register [START_REF] Giele | Movement correction of the kidney in dynamic MRI scans using FFT phase difference movement detection[END_REF][START_REF] Gupta | Fast method for correcting image misregistration due to organ motion in time-series MRI data[END_REF]. A number of similarity criteria have been proposed in the literature [START_REF] Maintz | A survey of medical image registration[END_REF]. High order criteria are robust to noise but, as a counterpart, are more sensitive to grey level intensities modification induced by the bolus. Inter-correlation coefficient is

F O R P E E R R E V I E W O N L Y
generally used as it offers the best compromise for kidney registration during perfusion [START_REF] Gupta | Fast method for correcting image misregistration due to organ motion in time-series MRI data[END_REF]. The objective is then to find parameters of the parametric spatial transformation that optimize this similarity criterion. For that purpose, exhaustive search or gradient driven approach are generally used.

Fourier based approaches

Those methods aim to express datasets in a new domain of representation and to exploit the new domain proprieties to estimate motion. For example, registration problems involving pure translation can be recovered by computing a phase difference in Fourier domain [START_REF] Giele | Movement correction of the kidney in dynamic MRI scans using FFT phase difference movement detection[END_REF].

Methods have been proposed to allow sub-pixel accuracy [START_REF] Shekarforoush | Subpixel image registration by estimating the polyphase decomposition of cross power spectrum[END_REF], and estimation of more complex displacements (a log-polar transformation is applied to the magnitude spectrum and the rotation and scale is recovered by computing phase difference in the log-polar space [START_REF] Reddy | An fft-based technique for translation, rotation, and scale-invariant image registration[END_REF][START_REF] Wolberg | Robust Image Registration Using Log-Polar Transform[END_REF]). Fourierbased approaches are very low time consuming and also more robust to estimate translation motion than methods based on grey levels intensity (see [START_REF] Giele | Movement correction of the kidney in dynamic MRI scans using FFT phase difference movement detection[END_REF]).

Methods based on geometrical characteristics invariance

Several template matching approaches based on geometrical characteristics invariance have been proposed in the literature. Gerig et al. have developed in 1991 a template matching between an estimated contour of the kidney and a manually drawn contour model [START_REF] Gerig | Semiautomated ROI analysis in dynamic MR studies. Part I: Image analysis tools for automatic correction of organ displacements[END_REF]. A linear transformation (rotation plus translation) was chosen to estimate the displacement between the two contours, and a bicubic interpolation was used for ensuring sub-pixel accuracy. The more recent and the most interesting method based on geometrical characteristics invariance has been proposed by Sun et al [START_REF] Sun | Integrated registration of dynamic renal perfusion MR images[END_REF]. The principle is the following: although the relative intensities between tissues vary with time, edges orientation along tissue boundary is always parallel across

F O R P E E R R E V I E W O N L Y
the image sequence. Note that, in its proposed form, the registration algorithm is combined with an auxiliary image segmentation step which is not evaluated in the current paper.

F O R P E E R R E V I E W O N L Y F O R P E E R R E V I E W O N L Y F O R P E E R R E V I E W O N L Y F O R P E E R R E V I E W O N L Y F O R P E E R R E V I E W O N L Y

  t. The objective was to relate the coordinate of each part of tissue in cur I with the corresponding tissue in ref I . A general definition for the problem of image registration can then be expressed as follow: registering

  ROI was different from the one used for movement correction and was drawn on the early enhanced image ref I . The arterial input function (AIF), required for the Patlak-Rutland analysis, was obtained from a manual ROI drawn on the aorta (volunteers) or iliac artery (transplant patients), just above the renal artery. The influence of inflow effects was minimised by the oblique-coronal positioning of the slices during data acquisition and no sign of inflow effects was observed during the Patlak analysis. For an objective in-vivo evaluation of the algorithm, we estimated single kidney GFR on the corrected and non-corrected data. Numerous publications have used the Patlak-Rutland model for estimation of the glomerular filtration of the cortex and

  volunteer before and after motion correction are presented in Figure 6. The average absolute variation of p K for the right kidney was 12.11 % ± 6.88 % (max = 25.6 %), and 11.6 % ± 6 % (max = 20.8 %) for the left kidney. A significant reduction of ) ( p K values was obtained on motion corrected data sets, with p = 0.003 (right kidney) and p = 0.002 left kidney (Student t-test). The average reduction of Kidney fMRI improvement with motion correction 12

  reduce reproducibility. This study has used an established tracer kinetic model (Patlak-Rutland model) which is independent of the user. The data following movement correction showed a reduced standard deviation with improvement of GFR uncertainty (up to 21.4 % reduction on

Figure 1 .

 1 Figure 1. Example of time intensity attenuation curves obtained on a kidney following a renal

Figure 2 .

 2 Figure 2. Example of Patlak-Rutland plots obtained on a transplanted kidney without (a) and

Figure 3 .

 3 Figure 3. Example of 2D GFR (a,b) and standard deviation (c,d) maps on transplanted kidney

Figure 4 .Figure 5 .

 45 Figure 4. Performances comparison between uncorrected and corrected data sets on the 10

Figure 6 .Figure 1 .Figure 2 .Figure 3 .Figure 4 .Figure 5 .Figure 6 .

 6123456 Figure 6. Performances comparison between uncorrected and corrected data sets on right (a,b)

d'Aquitaine, Philips Medical Systems, CDTU canceropôle network, Kidney Research UK.