International audienceIn this work, we develop a computational framework for optimal design of experiment in parametric signal reconstruction. We apply this to the optimal design of one dimensional Q-space, Q-ball imaging and multiple Q-shell experimental design. We present how to construct sampling scheme leading to minimal condition number, and compare to state-of-the-art sampling methods. We show in particular a better noise performance of these scheme through Monte-Carlo simulations for the reconstruction of synthetic signal. This demonstrates the impact of this computational framework on acquisition in diffusion MRI

Caruyer, Emmanuel

Deriche, Rachid

INRIA a CCSD electronic archive server

HAL Id: hal-00747700https://hal.inria.fr/hal-00747700Submitted on 1 Nov 2012HAL is a multi-disciplinary open accessarchive for the deposit and dissemination of sci-entific research documents, whether they are pub-lished or not. The documents may come fromteaching and research institutions in France orabroad, or from public or private research centers.L’archive ouverte pluridisciplinaire HAL, estdestinée au dépôt et à la diffusion de documentsscientifiques de niveau recherche, publiés ou non,émanant des établissements d’enseignement et derecherche français ou étrangers, des laboratoirespublics ou privés.A Computational Framework for Experimental Designin Diffusion MRIEmmanuel Caruyer, Rachid DericheTo cite this version:Emmanuel Caruyer, Rachid Deriche. A Computational Framework for Experimental Design in Dif-fusion MRI. CDMRI - MICCAI Workshop on Computational Diffusion MRI, Oct 2012, Nice, France.￿hal-00747700￿A Computational Framework for ExperimentalDesign in Diffusion MRIEmmanuel Caruyer, Rachid DericheAthena Project-Team, Inria Sophia Antipolis - Méditerranée2004 Route des Lucioles - 06902 Sophia Antipolis Cedex - FranceAbstract. In this work, we develop a computational framework for opti-mal design of experiment in parametric signal reconstruction. We applythis to the optimal design of one dimensional q-space, q-ball imagingand multiple q-shell experimental design. We present how to constructsampling scheme leading to minimal condition number, and compare tostate-of-the-art sampling methods. We show in particular a better noiseperformance of these scheme through Monte-Carlo simulations for thereconstruction of synthetic signal. This demonstrates the impact of thiscomputational framework on acquisition in diffusion MRI.1 IntroductionDiffusion MRI investigates the properties of tissue microstructure from the anal-ysis of water molecules displacement. The diffusion characteristics, such as theensemble average propagator or the orientation distribution function, are relatedto the diffusion signal attenuation through continuous transforms. Since then,the first step in the processing pipeline is usually a parametric estimation of thediffusion signal, from a series of discrete measurements. The number of samplesin diffusion MRI is limited to keep the acquisition time compatible with in-vivomeasurements. Therefore, the choice of sampling points in the q-space is criticalfor a proper reconstruction and quantitative analysis of diffusion characteristics.The question of sampling efficiency has been widely studied for parametricestimation in diffusion MRI. Several approaches were proposed to uniformly ar-range points on the sphere in q-ball imaging (QBI), using an analogy betweenthe sampling directions and a pair of antipodal electric charges [1, 2], or geomet-ric constructions [3]. For the reconstruction of diffusion tensor MRI, the noiseperformance has been studied through the minimization of the condition number[4, 5]. In q-space MRI, several studies on multiple shell sampling [6–9] focused onthe efficiency of various sampling strategies, but they do not provide a methodto systematically improve the noise performance.In this work, we give a general method for optimal design of experiment inparametric signal reconstruction. We apply this to the optimal design on one di-mensional q-space experiment, q-ball imaging and multiple q-shell experimentaldesign. In the last section, we compare the proposed method to state-of-the-artsampling strategies.2 Theory2.1 Parametric estimation of the diffusion signalThe diffusion signal is approximated in a finite, orthonormal basis of functions∀q ∈ Ω, E(q) =R∑i=1cifi(q), (1)where Ω ⊂ R3. Depending on the application, we have Ω = R (1D diffusionsignal), Ω = S2 (QBI) or Ω = R3 (q-space imaging).Provided K measurements yk = E(qk) of the signal at wavevectors qk, thecoefficients ĉi are estimated by least squares. Put in matrix form, we writeĉ = (HTH)−1HTy. (2)H is the design matrix, and has entries Hki = fi(qk).2.2 Information matrix and optimal designWe present a general method to choose the sampling points qk to optimizethe noise performance. A useful index for noise performance and stability ofthe reconstruction is the condition number κ(M) = λmax(M)/λmin(M) of theinformation matrix M = 1/K HTH, where K is the number of measurements.The condition number is an upper bound to the error propagation from themeasurements to the coefficients estimates. The optimal value of κ(M) is 1, inwhich case the information matrix is proportional to the identity IR.The coefficients of the information matrix Mij can be interpreted as theapproximation of the continuous dot product 〈fi, fj〉Mij =1KK∑k=1fi(qk)fj(qk) ≈∫Ωfi(q)fj(q)dΩ(q) = δij (3)The basis is orthonormal, hence if this approximation is exact, M = IR, andthe associated condition number equals 1. This naturally introduces the notionof quadrature formula, and its generalization to higher dimension, called thecubature formula.Definition 1. A cubature formula for the integral I =∫Ωg(q)dΩ(q) is a col-lection of nodes qs and weights ωs such thatI =S∑s=1ωsg(qs) (4)If such a cubature formula exists for the integral in Eq. 3, then we place thesampling points at nodes qs, and the number of repetitions Ks at node qs isproportional to the weight ωs.3 MethodsIn this section, we derive cubature formulae for the simple harmonic oscillatorbasis [10] (1D diffusion signal), spherical harmonic basis [11] (QBI) and sphericalpolar Fourier basis [12] (q-space imaging).3.1 Optimal design in one dimensional q-space MRIThe simple harmonic oscillator basis for the reconstruction of real diffusion signalin one dimension [10] is given by Φi(q, u) = κi(u) exp(−2π2q2u2)L−1/2i (4π2q2u2),where u is a characteristic length, L−1/2i the generalized Laguerre polynomial ofdegree i and κi(u) a normalization constant.Put back into the general framework presented in Section 2.2, we have Ω = R,and the basis functions are fi = Φi. For the evaluation of the dot product inEq. 3, we use the substitution x = 4π2q2u2, so that〈Φi, Φj〉 = 2πuκi(u)κj(u)∫ ∞0L−1/2i (x)L−1/2j (x)x−1/2e−xdx. (5)When the basis is truncated to order N , the evaluation of Eq. 5 reduces tothe problem of Gauss-Laguerre quadrature [13]. The optimal samples are qs =√xs/2πu, with Ks repetitions, where Ks is proportional to xsexs/[L−1/2N (xs)]2and the nodes xs, s = 1 . . . N + 1 are the roots of L−1/2N+1 .3.2 Optimal design in q-ball imagingGeometric Electrostatic Optimal (L = 4)Fig. 1. Arrangements of K = 50 points on the unit sphere.The real, symmetric spherical harmonic basis {Ylm} truncated to order L hasdimension R = (L+1) ·(L+2)/2. Put back into the general framework presentedin Section 2.2, we have Ω = S2, and the basis functions are fi = Yi, wherei(l,m) = 1, 2, . . . , R for (l,m) = (0, 0), (2,−2), . . . , (L,L).This basis is equivalent to the basis of harmonic polynomial of degree L onS2, for which cubature formulae exist and are called spherical design.Definition 2. A spherical t-design [14] is a sequence of K points X = (uk),k = 1 . . .K on the unit sphere, such that the integral of any polynomial p(x, y, z)of degree at most t over the sphere is equal to the average value of the polynomialon X :1KK∑k=1p(ukx, uky, ukz) =∫S2p(ω)d2ω. (6)If the sampling directions uk form a spherical 2L-design, then the approximationin Eq. 3 is exact and κ(M) = 1. For the construction of a spherical 2L-designwith antipodal symmetry, we rely here on the equivalence criterion in [15]. Anexample of optimal direction set for order L = 4 is presented on Fig. 1, andcompared to a geometric [3] and an electrostatic [1, 2] arrangements of points.3.3 Optimal design in q-space imagingN = 1, L = 4,K = 66 N = 2, L = 4,K = 168Constant sampling Optimal sampling Constant sampling Optimal samplingFig. 2. 2-shell (left) and 3-shell (right) acquisition protocols: regular multiple q-shellwith constant number of points per shell, and optimal arrangement with minimal con-dition number for reconstruction in SPF basis.The truncated spherical polar Fourier (SPF) basis [12] is able to representthe diffusion signal in the whole q-space. To represent a continuous signal, whichverifies E(0) = 1, we have recently proposed a slightly modified version ofthe SPF basis [16]. We reconstruct the signal as E(q · u) = exp(−q2/2ζ) +∑nlm anlmCnlm(q · u). The basis functions are Cnlm(q · u) = Fn(q)Ylm(u), withFn(q) = χnq2ζexp(− q22ζ)L5/2n(q2ζ), (7)ζ is a scale factor, χn a normalization constant, and Ylm is the real sphericalharmonic function. When the radial and angular truncation orders are N andL, respectively, this basis has dimension R = N · (L + 1) · (L + 2)/2. Put backinto the general framework presented in Section 2.2, we have Ω = R3, andthe basis functions are fi = Ci, where i(n, l,m) = 1, 2, . . . , R for (n, l,m) =(0, 0, 0), (0, 2,−2), . . . , (N,L,L).For the construction of an optimal design for this basis, we build on thefindings of the previous two sections. We show that the radial part of the integralin Eq. 3 reduces to a Gauss-Laguerre quadrature problem, while the angular partreduces to a spherical design problem.Therefore we propose a design on N + 1 spheres in the q-space. The shell shas radius qs =√ζxs, where xs is the sth root of L5/2N+1. The number of pointsKs on shell s should be proportional to ωs = exp(−xs)/[xs(L5/2N (xs))2] Finally,the points on each sphere should form a spherical 2L-design. Example of pointssets generated with this method are depicted on Fig. 2. They are compared tomultiple shell sampling where the shell radii are evenly spaced, and the numberof points equal on each shell, as suggested in [7].4 Experiments and Results4.1 Quadrature formula in one dimensional q-space MRI0.0 0.5 1.0 1.5 2.0 2.5q ·d10-810-710-610-510-410-310-210-1100E(q)ground truthregular samplingLaguerre quadrature5 10 15 20 25 30K1.01.52.02.53.03.54.0(M)regularsamplingLaguerrequadratureFig. 3. Evaluation of Gauss-Laguerre quadrature for one dimensional q-space. (Left)An example of signal and its reconstruction. The blue squares and the green circlesrepresent the regular and quadrature samples respectively. The radii of the circles areproportional to the number of repeated acquisitions. (Right) Condition number of theinformation matrix.In this section, we show the feasibility of one dimensional q-space signalreconstruction from a set of measurements on a limited support size. We comparethe Gauss-Laguerre quadrature to a regular sampling on the range [0, qmax]We plot on Fig. 3 an example of reconstruction of a diffusion signal corre-sponding to the restricted diffusion between two parallel planes, separated bydistance d [10]. The truncation order in the basis was set to N = 5, and thecorresponding Gauss-Laguerre quadrature works on 6 nodes. The result for atotal of 20 acquisitions is visually identical to the reconstruction from a regularsampling. Besides, the associated information matrix is better conditioned forthe quadrature sampling. The reason why the condition number is not exactly1 in this case is that the quadrature weight ωs is approximated by the numberof repetitions at node qs, which is an integer.4.2 Evaluation of conventional schemes in q-ball imagingWe evaluate and report on Fig. 4 the noise performance of point sets gener-ated with electrostatic analogy [1, 2] and by geometrical construction [3], forthe reconstruction of SH coefficients of the diffusion signal. We compare thesesampling methods to the proposed, optimal point set based on spherical design.0 20 40 60 80 100 120 140K1.01.11.21.31.41.51.6(M)L=4geometricelectrostaticoptimal20 40 60 80 100 120 140K1.01.11.21.31.41.51.6(M)L=6geometricelectrostaticoptimalFig. 4. q-ball imaging: condition number of the information matrix corresponding tothe truncated SH up to order L = 4 (left) and L = 6 (right), of the electrostatic,geometric and optimal point sets. The geometric configurations are only provided forK ≥ 50, as this method is reported to be dedicated to large K by the author in [3]. Ourproposed, optimal design is based on spherical design, and therefore exists for K ≥ 24at order L = 4, and for K ≥ 46 at order L = 6 [15]. By construction, the conditionnumber associated is exactly 1.4.3 Multiple q-shell and three dimensional signal reconstruction40 80 120 160 200K12345678910(M)N=2constant samplingoptimal sampling40 80 120 160 200K5101520(M)N=3constant samplingoptimal samplingFig. 5. Condition number of multiple shell sampling, corresponding to the SPF basis.Using the results on Gauss-Laguerre quadrature and spherical designs, wegenerated optimal sampling schemes on multiple shells for the reconstruction inthe SPF basis. We compare this to the sampling strategy with shell radii evenlyspaced and constant number of points per shell proposed in [7]. The conditionnumber for the reconstruction in SPF basis is reported on Fig. 5.We also simulate both methods, for the sampling and reconstruction of asynthetic diffusion signal corresponding to a mixture of Gaussian. Visual recon-struction is reported on Fig. 6, and quantitative comparison on Fig. 7.-78 -39 0 39 78-78-3903978-78 -39 0 39 78-78-3903978-78 -39 0 39 78-78-3903978-78 -39 0 39 78-78-3903978-78 -39 0 39 78-78-3903978Ground truth signal Constant sampling Optimal samplingFig. 6. 2 shell (top) and 3 shell (bottom) sampling in action, for the reconstruction ofa synthetic diffusion signal corresponding to a mixture of Gaussian, simulating fibercrossing. The circles represent the sampling shells, and the line widths are proportionalto the number of points per shell.5 ConclusionIn this work, we develop a computational framework for optimal design of exper-iment in diffusion MRI. For the reconstruction of 1D, spherical and 3D signal, wepropose sampling scheme with minimal condition number. Monte-Carlo simula-tions confirm this result, as the signal to noise ratio of the parameters estimatedfrom optimal sampling scheme is improved with respect to conventional samplingscheme, for the same number of acquisitions.As a conclusion to this study, we claim that a sampling method is optimal fora reconstruction in a given basis and a given order. In addition to the technicaland physical limitations of the imaging system, the choice of the type of recon-struction must be taken into account when designing the acquisition protocol.40 60 80 100 120K0.180.200.220.240.260.280.30N=2constant samplingoptimal sampling40 60 80 100 120K0.180.200.220.240.26N=3constant samplingoptimal samplingFig. 7. Mean squared error for a synthetic diffusion signal corresponding to a mixtureof Gaussian, with Rician noise (SNR=25).References1. Jones, D., Horsfield, M., Simmons, A.: Optimal strategies for measuring diffusionin anisotropic systems by MRI. Magn Reson Med 42(39) (August 1999) 515 – 5252. Jansons, K.M., Alexander, D.C.: Persistent angular structure: new insights fomdiffusion magnetic resonance imaging data. Inverse Problems 19 (2003) 1031–10463. Koay, C.G.: A simple scheme for generating nearly uniform distribution of antipo-dally symmetric points on S2. J. Comput. Science 2(4) (2011) 377–3814. Papadakis, N.G., Xing, D., Huang, C.L., Hall, L.D., Carpenter, T.A.: A compar-ative study of acquisition schemes for DTI. J Magn Reson 137(1) (1999) 67–825. Skare, S., Hedehus, M., Moseley, M.E., Li, T.Q.: Condition number as a measure ofnoise performance of diffusion tensor data acquisition schemes with MRI. J MagnReson 147(2) (December 2000) 340–3526. Caruyer, E., Cheng, J., Lenglet, C., Sapiro, G., Jiang, T., Deriche, R.: Optimaldesign of multiple q-shells experiments for diffusion MRI. In: MICCAI Workshopon Computational Diffusion MRI - CDMRI’11, Toronto, Canada (September 2011)7. Assemlal, H.E., Tschumperlé, D., Brun, L.: Evaluation of q-space sampling strate-gies for the diffusion magnetic resonance imaging. In: MICCAI, London (2009)8. Portnoy, S., Ye, W., Entezari, A., Blackband, S.J., Vemuri, B.C.: A novel interlacedsampling scheme for multi-shell q-space MR microscopy. In: ISMRM. (2011)9. Daducci, A., McEwen, J., Van De Ville, D., Thiran, J.P., Wiaux, Y.: Harmonicanalysis of spherical sampling in diffusion MRI. In: ISMRM. (2011) 392910. Ozarslan, E., Koay, C., Basser, P.: Simple harmonic oscillator based estimationand reconstruction for 1D q-space MR. In: ISMRM. Volume 16. (2008) 3511. Descoteaux, M., Angelino, E., Fitzgibbons, S., Deriche, R.: Regularized, fast, androbust analytical q-ball imaging. Magn Reson Med 58(3) (2007) 497–51012. Assemlal, H., Tschumperlé, D., Brun, L.: Efficient and robust computation of pdffeatures from diffusion MR signal. Medical Image Analysis 13(5) (2009) 715–72913. Abramowitz, M., Stegun, I.A.: Handbook of mathematical functions. Dover (1970)14. Delsarte, P., Goethals, J., Seidel, J.: Spherical codes and designs. GeometriaeDedicata 6 (1977) 363–388 10.1007/BF03187604.15. Goethals, J.M., Seidel, J.J.: Cubature formulae, polytopes and spherical designs.In Coxeter, H., Davis, C., Grünbaum, B., Sherk, F., eds.: The Geometric vein: theCoxeter festschrift. Springer-Verlag, New York (1981) 203–21816. Anonymous: -. In: -. (-)

A Computational Framework for Experimental Design in Diffusion MRI

HAL-Rennes 1

https://hal.inria.fr/hal-00747700/document

A Computational Framework for Experimental Design in Diffusion MRI

Abstract

Similar works

Full text

Available Versions

INRIA a CCSD electronic archive server

HAL-Rennes 1