Abstract
In MRI, at ultrahigh field, the use of parallel transmit radiofrequency (RF) arrays is very beneficial to better control spin excitation spatially. In that framework, the so-called “universal pulse” technique, proposed recently for head imaging at 7 tesla, gives access to “plug-and-play” nonadiabatic solutions exhibiting good robustness against intersubject variations in the resonant transmit fields. This new type of solution has been defined so far as the result of numerical pulse optimizations performed across a collection of RF field maps acquired on a small population sample (pulse design database). In this work, we investigate an alternative universal pulse design approach in the linear small tip angle regime whereby the database of RF field maps is first transformed into a second-order statistical model and which then exploits a statistical robust design formalism for the optimization of the RF and magnetic field gradient waveforms. Experimental validation with an eightfold transmit RF coil for 7 tesla brain imaging shows that this new approach brings some benefit in terms of computational efficiency. Hence, for a design database composed of 35 maps, the computation time initially of 50 min could be reduced down to 3 min. The proposed statistical approach thus enables integration of large databases, presumably necessary to ensure robust solutions. Finally, it provides means to compute flip angle statistics and, along with it, simple performance metrics for quality assurance (RF pulse performance) or guidance in the optimization of TX array architectures.
1. Introduction
In magnetic resonance imaging at high magnetic field, the apparition of standing wave effects [1–3] prevents uniform excitation of the spins across extended anatomical regions. For brain imaging, these effects are clearly visible at 3 tesla (T) and are often deemed unacceptable at field strength 4 T. Parallel transmission (pTX) is a promising technology to mitigate this problem as it provides degrees of freedom to “shim” the RF field distribution [4] (RF shimming). The transmit SENSE technique [5–7], enabled by the pTX technology, furthermore allows reducing the duration of multidimensional pulses [8] while preserving spatial definition. It thus offers one more solution to tackle the transmit field heterogeneity problem. The rationale here is to regard the problem of achieving a uniform magnetization profile as one particular problem that multidimensional RF pulses can solve, whereby the target magnetization profile is spatially uniform. The design of such pulses, here referred to as flip angle (FA) shimming, plays now an important role at UHF magnetic resonance to enable uniform spin excitations. Thus, the so-called fast-kz or spokes pulses [9, 10] and the kT-point pulses [11] are well-known slice-selective and nonselective FA shimming strategies taking advantage of the parallel transmission technology to produce highly uniform excitation profiles at a reasonable cost in terms of pulse duration as compared to RF shimming.
However, in addition to being nonuniform in space, the transmit RF field () also depends on the geometry and electric properties (dielectric permittivity and conductivity) of the sample and thus can vary significantly across individuals. Thus far, the usual technique to cope with this is to characterize the actual profiles using MRI-based measurements and to tailor RF shimming to this particular TX field distribution. This obviously can be time-consuming, in particular when the number of TX channels is high [12]. One alternative approach to handle the variability of the TX field distribution, called universal pulse (UP) [13], is to replace the subject-specific TX field distribution with a collection of such TX field distributions acquired on a database of subjects. The rationale here is to achieve, simultaneously across all subjects of the database, good excitation performances, and assume that this property generalizes to the population that the database is supposed to cover (e.g., the adult population). The validity of this approach has been demonstrated experimentally for various applications and hardware setups [14–16]. In those studies, a database of 10 to 20 individuals was deemed sufficient to yield good performances. However, the amount of data remained too limited to analyze quantitatively how pulse performance changes with the design database, and in particular with its size. Nevertheless, it is reasonable to assume that UP performance increases with the size of the design database; conversely, if a too small database is used to compute UPs, there is a risk of “data overfitting,” i.e., that the solution displays a higher FA shimming performance on the design database than its actual performance on the entire population.
To better control the actual performance of UP on the entire population, we propose in this work to model the transmit RF field as a stochastic process, whereby every new subject (giving rise to a new map) is a realization of this process. This representation appears very useful as it allows formulating the action of a FA shimming pulse in statistical terms. From the RF pulse design perspective, it also provides an adequate framework to apply the so-called statistical robust design theory, proposed by Alotto et al. [17], applied to the design of electromagnetic devices. For the latter application indeed, manufacturing errors of the components or uncontrollable environmental conditions such as temperature or humidity can lead to uncertainties in the response of the system. Applying classical optimization to such a problem can severely degrade the actual performance of the solution. A robust design approach mitigates this problem by accounting for the randomness of the response of the system in the optimization process using an appropriate statistical model. This idea can naturally be extended to the RF pulse engineering problem.
We thus present in this work a statistical representation of the TX field distribution in which, for every position , is assimilated to a Gaussian random process. We then derive a method to apply the aforementioned statistical robust design method to the design of FA shimming UPs, which uses the second-order statistics of instead of the collection of maps as input parameter. An experimental validation of the method is provided with data acquired with an eightfold transmit RF coil for 7 tesla head MRI. We demonstrate in particular (1) the utility of statistical robust design approach to boost computational efficiency, (2) the capacity of the method to enable the integration of a large database TX field maps, and (3) the utility of a statistical framework to compute flip angle statistics and to derive simple performance metrics for quality assurance.
2. Theory
2.1. Statistical Modelling of the Resonance TX Field
Let (tesla/volt) denote the spatially dependent ( referring to the position) resonant transmit RF magnetic field of the resonator of the -fold transmit array and Hz/T the gyromagnetic ratio of the proton. The so-called control field components (rad/s), , may then be arranged in a vector called the control field vector.
The control field vector being subject dependent is viewed as a stochastic variable where each new subject gives rise to a realization of . For every location, we assume that the mean value and the covariance matrix ( for the Hermitian conjugate) of the control field vector is well defined and we denote them, respectively, by and . In the following, the control field is assumed to have Gaussian statistics, i.e., .
Let be independent control field distributions (in practice TX field maps measured in different subjects). If is sufficiently large, this dataset can be exploited to estimate the mean and covariance matrix of the underlying random process by computing the sample mean and covariance matrix for this ensemble:
For the covariance matrix, we shall use the Ledoit–Wolf [18] correction of the sample covariance estimator when the condition ( in this study) is not satisfied.
Given a statistical model for , we now show that this information is sufficient to compute universal FA shimming pulses and to study the robustness of such pulses with respect to the intersubject variability of . We first propose a formalism to describe the action of FA shimming pulses. Although a generalization would be possible for large flip angles, this theory for now is limited to the small tip angle regime where there exists a linear relationship between the RF excitation and the spin’s tip angle [19].
2.1.1. Matrix Representation of the Spatial Distribution the Transmit RF Field
As the control field distribution is typically obtained by MR acquisition, it is common in practice that is not defined at each voxel of the underlying MR acquisition (existence of regions where the reconstruction of is not possible). We thus define by the set of voxel positions where a given realization of is known ( is subject and session dependent). In the following, we shall use the following matrix representation of the spatial distribution of the control field:where is the number of positions contained in . When there is no ambiguity about the choice of the region of interest, we shall simply use the notation .
2.2. Flip Angle Shimming Pulse
Let and be the RF and k-space trajectory waveforms of a given FA shimming pulse. Let be the duration of the RF and magnetic field gradient (MFG) waveforms and (rad/s) be the resonance offset at location due to the static field nonuniformity. Let us use below the notation . In the small tip angle approximation, the effective (time-averaged) control field created by this RF and magnetic field gradient (MFG) pulse pair is given by [8]where
Note that the flip angle created by the pulse is simply . As a result, the mean, , and variance, , of the effective control field are given by
In this model, the effect of off-resonances is taken into account but is not modelled as a random process . Equation (4) hence can be implemented as such if the static field offset is mapped in a subject-specific manner. Alternatively, one can assume , i.e., that the static field offset has a negligible influence on pulse performance, which is acceptable for short () pulses.
In the following, we present a method to compute universal FA shimming pulses which is entirely based on a second-order statistical model for the control field. Along with this statistically driven framework, we provide a method to evaluate the performance of a universal FA shimming pulse.
2.3. Statistical Robust Design of Universal Pulses
Let be a candidate FA shimming pulse to create a unitary effective control field . Let us denote bythe discrepancy between the unit target and , and by the probability that exceeds a certain error level . To optimize , we can exploit here a statistical robust design (SRD) approach [20]. Let us introduce a new design parameter , called SRD tolerance, and define the error level aswhere denotes the region of interest for the design of the pulse (any position outside being therefore ignored in this metric). The SRD problem may then be posed as the minimization problem:
Denoting now by the cumulative distribution function of , the condition gives
Hence, the problem described by equation (8) takes the equivalent form:where denotes the inverse function of . A useful simplification of this optimization problem is then obtained by replacing the maximum ( norm) by the norm for a certain choice of . This leads to the following new pulse optimization problem:
If, now, for a given pulse candidate , the following two conditions are satisfied for all positions in :then we can show (see Appendix A) that the pulse optimization problem (11) can be further simplified to take the form:where the dependence on the second-order statistics can be made explicit:
Retrospectively, the link with the database-driven objective function (simultaneous design across all subjects of the database of maps), denoted here by , can be explained as follows. The latter can be defined as the mean across the design database of the mean squares error of the effective control field:
Now, if we assimilate sample mean and statistical mean, we obtainwhich leads us to the same expression as the one obtained in equation (13) for . This shows the near equivalence (note however that the objective function defined in reference [15] is proportional to the average ( norm) of the FA-NRMSE rather than the norm and that the brain region for the evaluation of the FA-NRMSE was subject-specific rather than fixed as in equation (15)) of the database-driven and the statistical approaches when sample mean and statistical mean can be assimilated. This is acceptable in particular when the design database is large enough to be representative of the whole targeted population.
2.4. Derivation of a Lower Bound for the Coefficient of Variation of the Flip Angle
Let be a solution to the problem of minimizing the variance of the control field at position among the RF shim solutions that return on average, up to a phase term, a unitary control field:
In Appendix B, an explicit description of is obtained under the condition that
The expression for then readsand the variance of the value of the objective at the optimum is given by
Furthermore, it is shown in Appendix C that the locally optimal solution sets a performance limit for the coefficient of variation (CV) of the effective control field:which reads
Hence, under the condition , the lower bound for the CV of the effective control field, which we may also call ultimate FA-CV in the linear STA regime, and below denoted by , reads
The condition , on the other hand, has a clear physical interpretation: the mean control field magnitude dominates its standard deviation. The validity of this condition in practice will be tested subsequently using experimental data.
It is also of interest to define the ultimate FA-CV—in the sense of equation (23)—for the circularly polarized (CP) mode of the TX array. In this mode, a fixed relationship between the RF signals applied to the TX array ports is imposed in order to reproduce the behavior of a single port volume coil. The CP mode can be defined with the vector which describes the way the channels of the TX array interfere. The mean and variance of the CP mode control field, , is then given by and . Application of equation (23) on the CP mode control field statistics giveswhere
In what follows, we present an application of the proposed statistical model of the intersubject variability of the control field and their implication in terms of pulse design and flip angle shimming performance assessment.
3. Material and Methods
3.1. Control Field Second-Order Statistics
The data consisted of a collection of maps measured on different adult subjects (age years) which were recruited as part of previous studies, described in [14, 15, 21]. The MRI data were obtained on a parallel transmit-enabled whole-body Siemens 7 T system (Siemens Healthineers, Erlangen, Germany) equipped with an SC72 whole-body gradient insert (200 mT/m/ms and 100 mT/m maximum gradient slew rate and maximum gradient amplitude, respectively) and the 8Tx-32Rx Nova head coil (Nova Medical, Wilmington, MA, USA) for signal transmission and reception. This RF coil array possesses a CP mode for transmission in which the amplitude and phases at excitation ports 2 to 8 are imposed by the amplitude and phase applied to the first excitation port. This operation mode thus has only one degree of freedom (the complex amplitude applied to the first excitation port) and is designed to create the TX profile corresponding to a birdcage resonator. The parallel transmission system consisted of eight transmitters (1 kW peak power per channel) driven with the Siemens Tim Tx (Step 2.3) system. The protocol used to reconstruct the quantitative maps was acquired from a multislice interferometric turbo-FLASH acquisition [22, 23] (5 mm isotropic resolution, matrix size , TR = 20 s, TA = 4 min 40 s).
The control field second-order statistics of the TX array head coil were computed from the TX field maps of the first (chronology of data acquisition) 35 subjects of the database using the mean and covariance estimators, as described in the Theory section. For the estimation of the covariance matrix, the Ledoit–Wolf correction was used due to the limited number of subjects as compared to the dimension of the covariance matrix (64).
3.2. Coefficient of Variation of the Flip Angle
The lower bounds for the CV of the effective control field in the linear STA regime for the CP and pTX modes were computed from the estimated second-order statistics of the control field using equation (23). This enabled us in particular to compare the ultimate FA shimming performance when exploiting full pTX capability (independent RF signals applied to each port of the TX array) with the FA shimming performance of the CP mode of the TX array.
3.3. Statistical Robust Design of kT-Point Pulses
For that study, the spins were assumed to be at resonance everywhere, i.e., we assumed that the off-resonance effect on excitation was negligible. Computed FA shimming pulses were 7-kT-point [11] pulses (7 rectangular subpulses), wherein the duration of the RF subpulses and MFG blips were set to 140 μs and 50 μs, respectively. In all optimizations, the target flip angle, , was 10 deg. Furthermore, the RF peak power applied at each excitation port was limited to 500 W, the total RF energy (respectively, the RF energy per TX channel) of the pulse was assigned to not exceed 60 mJ (respectively, 12 mJ), and the -space displacement in , , and created by each MFG blip was imposed to not exceed 50 rad/m.
A first FA shimming pulse, referred to as SRD-UP, was designed from the second-order control field statistics of the Nova TX array head coil (estimated in Section 3.1) by solving numerically the problem (same as equation (13) but nonunitary target, duration of the pulse):
In the above expression, was defined as the set of positions that belonged to the brain region in at least half of the subjects of the database, i.e.,where denotes the mask of the brain for the -th subject of the database and is the indicator function of this region.
For comparison with SRD-UP, the solution (same parameterization and constraints as for the design of SRD-UP) of the database-driven UP objective was also computed. The database-driven UP objective was defined as follows:
This pulse is referred to as UP. Finally, a standard hard pulse where the transmit array is used in CP mode was defined. The amplitude of this pulse thus was set to return on average on the design database the target flip angle.
We then studied the spatial distribution of the mean and CV of the flip angle for UP, SRD-UP, and CP and compared it with the ultimate FA shimming performance . For convenience, this performance is attributed to a virtual pulse which would satisfy (see equations (4) and (19)) . In the following, this pulse is referred to as . Note that being a local solution, no RF and MFG waveforms can satisfy strictly this relationship unless the k-space trajectory defined by this pulse offers a complete coverage of the transmit k-space trajectory. This would lead in practice to a prohibitively long FA shimming pulse.
Using the five last subjects of the cohort (the 35 first constituting the design database), the flip angle normalized root mean squares error (FA-NRMSE) of UP, SRD-UP, and CP ( one of these pulses, , the FA created by pulse at location , ):was computed on each subject. For the FA-NRMSE computation, the FA maps, , were obtained by numerical integration of Bloch’s equations (whereby T1 and T2 relaxations were neglected) rather than using the linear small tip angle approximation (equation (3) For every subject of this test database, a subject-tailored kT-point pulse (minimization of the term of the sum in equation (28)), possessing the same parameterization as the above universal pulses (SRD-UP and UP), was also designed to minimize the FA-NRMSE. These solutions are referred to as , …, .
All pulse optimizations were performed using the second-order active-set algorithm, available in the Optimization Toolbox of MATLAB R2019a (The Mathworks, Natick, MA), which allows solving constrained optimization problems numerically. The active-set algorithm used the analytical expressions for the Jacobian of the objective function and the RF power, RF energy, and MFG constraints (Hoyos-Idrobo et al. [24]). The number of iterations was limited to . The algorithm however could stop before exceeding this number of iterations provided it reached a local minimum within default tolerance.
4. Results
4.1. Second-Order Field Statistics
In Figure 1, the second-order statistics of the control field for the Nova TX array are displayed for one axial slice through the brain. The mean control field thus gives rise to eight axial images, each one corresponding to one TX channel. We also provide the mean control field of the CP mode, . The covariance matrix is composed of the control field variance of every TX channel (eight diagonal images) as well as the covariance between channels (off-diagonal images).
(a)
(b)
4.2. Ultimate FA Shimming Performance
Figure 2 displays sagittal, coronal, and axial images of the lower bound for the CV of the effective control field (ultimate FA shimming performance) in CP mode and the pTX operation mode. Here, it can be seen that , i.e., that the pTX drive mode allows for a significant decrease in the CV of the effective control field as compared to the CP mode. In other words, for every voxel in the brain, there exists a pTX excitation mode that yields less intersubject variability than the CP mode. Quantitative comparisons of these performance boundaries are also provided in the form of cumulative histograms of and in Figure 3. Finally, histograms of the distributions of the ratio of CVs are provided in Figure 4. From this plot, it can be seen in particular that the CV of the effective control field can be reduced by a factor greater than two in some locations in the brain with pTX as compared to the CP mode.
4.3. FA Shimming Performance Simulation
The FA shimming performance analysis of , UP, SRD-UP, and for the standard hard pulse in the CP mode is shown in Figures 5 and 6. For UP and SRD-UP, in the brain (in the sense of (30)), the mean flip angle (statistical mean) is very close to the target (10) as it remains in the interval . By definition of (see equation (19)), for the virtual pulse , the mean flip angle is exactly equal to the target. As it can be expected, for the CP mode, the mean flip angle can deviate very strongly from the 10 target as it reaches for instance in the center . On the other hand, the median CV of the flip angle amounts to 12% for UP and CP and 10% for SRD-UP and is ultimately % for . In the “worst location” which was found to be the cerebellum (see Figure 5), the CV of the flip angle does not exceed 15% with but can exceed 20% with “real” pulses.
(a)
(b)
(a)
(b)
(c)
The result of the FA-NRMSE analysis made separately on the design and test subjects is presented in Figure 7. The results obtained on the test subjects are compliant with the FA-NRMSE statistics obtained on the design database in the sense that no significant offset (unpaired Wilcoxon rank-sum test) was found between the FA-NRMSE obtained in the design and test databases. On the test subjects, the performance was very comparable between , UP, and SRD-UP, with a median FA-NRMSE of 5 to 7%. The subject-tailored pulses returned FA-NRMSEs %.
(a)
(b)
4.4. Storage and Computation Time Efficiency
Denoting by the average number of voxels in the brain, the storage volume for the first- and second-order statistics of the control field and for a database of subjects are and , respectively. Thus, for (6 in our case), the statistical representation becomes more parsimonious. From the computation time point of view, the processing time (DELL Precision 7540, 32 Gb of RAM, 8 cores, Intel CORE i9 vPro 9th Gen) was 3270 s (669 iteration and 1851 evaluations of the objective function) and 160 s (489 iterations and 1471 evaluations of the objective function) for the computation of UP and SRD-UP, respectively. The computation time for one evaluation of the objective function (and the Jacobian) was 2.13 s and 0.26 s for UP and SRD-UP, respectively. As these workloads scaled linearly with the number of voxels ( in total for the 35 subjects of the database and for (see equation (27))), we found that the workloads per voxel were 5 μs for UP and 17.4 μs for the computation of UP and SRD-UP, respectively. As a result, but keeping in mind that these numbers may depend on the pulse parameterization (e.g., number of kT-points), the design of SRD-UP is computationally more efficient than the UP design when the size of the design database is greater than .
5. Discussion
In this work, we have shown that the statistical robust and the UP design methods are theoretically nearly equivalent when the size of the database is sufficiently large to embrace the variability of the resonant TX RF field across the population of interest. Despite its apparent complexity, the statistical robust design algorithm defined by equation (13) has a rather straightforward implementation as the general term of the summation across voxels depends explicitly on the mean vector and covariance matrix of the TX RF field (see equation (14)). As compared to the database-driven UP design methodology (Gras et al. [14]), an important difference is the involvement of the covariance matrix which has not equivalent in the original formulation. To our knowledge, the design of UPs using the SRD formalism, originally applied by Alotto and colleagues for the design of electromagnetic structures, has not been proposed so far for RF pulse design and could represent some advantages. First, as the SRD computation time does not depend on the size (number or subjects) of the database, this approach enables integration of potentially very large databases, which is a priori necessary to ensure optimal UP performance. Secondly, equipped with this model, one can likewise compute probabilities for obtaining a certain deviation from a target flip angle. Doing such estimations using a database of TX field maps would require unrealistically large number of subjects.
The theory presented in this work defines the performance of a FA shimming pulse as the coefficient of variation (std/mean) of the flip angle. Furthermore, in the linear STA regime, it is shown that this CV cannot be lower than a certain value which depends on the local second-order statistics of the control field. In the process of validating—and ultimately certifying—the use of pTX in a clinical setup, this performance measure may prove useful. One concrete application of this concept is the comparison of the CP mode ultimate performance () with the pTX performance (): for the Nova TX array, we found that the availability of eight independent TX channels to locally tune the drive mode so as to minimize the variance of the effective control allowed reducing the CV of the flip angle by a factor of , while achieving close to the target on average (see Figure 4). Interestingly, this theoretical analysis sheds some light onto the concept of UPs. Indeed, with statistical arguments, one can take advantage of the TX array architecture to mitigate both the spatial heterogeneity and the intersubject variability of the transmit RF field in the brain. Finally, the coefficient of variation metric for (see equation (23)), being characteristic of the TX array architecture, can potentially serve as a guide for coil-design purposes.
The present study was performed with the data collected on 40 subjects, from which a subgroup of 35 subjects was used to compute statistics. Given the high dimension of the covariance matrix of the control field for an eightfold TX array system, it was necessary to use the Ledoit–Wolf correction for the estimation of the covariance matrix of the control field (to obtain positive definite estimates). However, we can still expect a bias in the estimation of the covariance matrix. As a result, it is possible that the CV analysis reported here (lower bound and CV maps for the different kT-point pulses computed in this study) can deviate from the analysis that one would obtain with a larger database to estimate the control field second-order statistics.
The theory developed in this study is based on the assumption that the statistical distribution of the control field vector is Gaussian. This assumption is very convenient from the mathematical point of view as this property is inherited by the effective control field in the linear STA regime. With a few tens of subjects, this hypothesis has not yet been verified quantitatively, but given the complexity and the multiplicity of effects leading to the observed variability, the Gaussian model is a reasonable assumption a priori. Another hypothesis made in this study is the fact that is much greater than one. Here we notice that this is equivalent to the assumption that is much smaller than one. In Figure 2, was seen not to exceed 0.2, confirming that this hypothesis was approximately valid over the whole brain.
Extending this work to large flip angle pulse is one interesting prospect. However, in this case, the linear relationship between and exploited to derive a simple expression of the SRD objective does not hold. Here a pragmatic approach would be to replace exact calculation of the first and second moments of the flip angle (seen as a random variable) by numerical approximations thereof. A possible direction to perform such calculation is the use of the so-called unscented transform method [25] which seems here particularly adapted [20].
An important difference with the pulse design method reported in [13] is the deterministic nature of the off-resonance term in the current statistical design. In this study, we set for the design of UPs using the SRD formalism, i.e., we neglected off-resonance effects in the performance assessment. While acceptable for short excitation pulses, this can be a significant source of error in practice when the pulse duration becomes large [26]. Extending the statistical robust design approach to account for a random off-resonance term seems however feasible with some additional hypotheses. Hence, if it is assumed that is a centered Gaussian random variable independent of , it remains possible in equation (3) to compute the mean deviation of the effective control field with respect to the effective control field at resonance. This quantity takes the form:where denotes the effective control field at resonance, and where the expectation operator covers as well as statistics. A possible direction is to exploit this quantity as an additional objective to improve the robustness against off-resonances while taking into account the regional dependence of the variance of the resonance frequency offset.
In this work, for simplicity, the statistical robust pulse design implementation did not include specific absorption rate (SAR) constraints. However, as it was done in the past for the design of UPs, there would be no obstacle in adding explicit local SAR constraint to the statistical robust design procedure in the form of virtual observation points [27]. Interestingly also, the RF electric field distribution (-field) at the origin of the SAR is another imperfectly controlled physical quantity affecting RF safety management procedure, in particular at UHF and with transmit array architectures [28]. Like the field, the -field could advantageously be treated as well as a stochastic process, thereby allowing for a probabilistic model of the local 10 g SAR. Naturally, the dual problem consisting of minimizing the SAR under performance (flip angle error) constraints [29] could benefit as well from such a unified statistical model. One obstacle for the estimation of the -field statistics yet remains the absence of direct method to measure the -field in vivo and thus the necessity to resort to electromagnetic simulations on a large amount of head models to give access to them. Furthermore, due to the interweaved nature of the RF electric and magnetic fields, the values of the - and -field, seen here as random processes, are in fact necessarily correlated.
6. Conclusions
We have presented in this work a statistical robust approach to design calibration-free pTX universal pulses for brain imaging at ultrahigh field. Experimentally, this method has been shown to provide comparable FA-NRMSE results as the conventional database-driven UP design, but presents a significant advantage in terms of computational efficiency when the size of the design database exceeds 4. Furthermore, in contrast with the database-driven design of UP, the computation time of the statistical robust UP does not depend on the size of this database and consequently enables integration of a very large design database. This is an important criterion to ensure optimal UP performance. Furthermore, the theoretical lower bound for the coefficient of variation of the flip angle derived in this work, , appears as a valuable analytical performance limit to assess the FA shimming performance of a given pulse. Finally, this metric can be useful to characterize, and possibly also optimize, the performance of TX array architectures.
Appendix
A. Statistical Robust Pulse Design
In this section, we focus on the relaxed SRD problem (infinite norm replace by the norm) given by equation (11). Let us study the first- and second-order statistics of . From , we have: and Rician. Now, assuming that for all voxel positions, the mean effective control largely exceeds its standard deviation, i.e.,we can assume that owns a Gaussian distribution with the following mean and variance:
It follows from condition (A.1) that the distribution of is the so-called folded normal distribution and that the mean, , and the variance, , of are therefore given bywhere is defined as follows:and where denotes the cumulative density function of the standard normal distribution .
In Figure 8, we have plotted for three values of the SRD design variable the quantity:as a function of . The graph of as a function of is very flat in and, for and , we have (see Figure 8)
As a result, for a candidate pulse for which condition (A.6) is satisfied for all , then the objective in equation (11) simplifies to
Noting that (i) the leading term in equation (A.7) is independent of and (ii) thatthe pulse design problem thus simplifies to
We note that the condition given by equation (A.6) supposes that the candidate FA shimming pulse returns an effective control field that is “almost” uniformly unbiased ( everywhere). Therefore, it would not be realistic to assume that this condition is fulfilled for all . But, for the above simplification to be made, one only need to assume that this condition is fulfilled in the neighborhood of the solution of (A.9).
B. Lower Bound for the Coefficient of Variation of the Effective Control
Below we first solve the phase-constrained problem and then the phase-unconstrained problem. While the first problem can be solved analytically in general, the resolution of the second problem—the one that we focus on in this study—is conducted under the assumption that .
B.1. Phase-Constrained Problem
First, we look for the solution of the minimization problem:
Using the equality constraint, we have
The Lagrangian of the above problem reads
If is solution, thenand (true since has the Hermitian symmetry):
From equation (B.5), we obtain
From equations (B.6) and (B.5), we thus have
The solution to equation (B.1) thus reads
B.2. Phase-Unconstrained Problem
The phase-unconstrained problem (17) may be rewritten as follows:where . Assuming Gaussian statistics for (), then and so is Rician distributed with mean ( where and are the modified Bessel functions of the first kind with parameter and , respectively):
Now, if we suppose that(i.e., according to equation (18), ), then the distribution of reduces to where
Under this hypothesis, the problem (B.9) simplifies to
The Lagrangian of this optimization problem reads
The second-order optimality condition thus reads
Given now that is colinear to , a necessary condition to have is that is colinear with , and so, that there exists a scalar such that . From equation (B.16), we obtain (see equation (18)):
As a result, the solution to the phase-unconstrained problem is given (up to a uniform phase term across TX channels) by
Taking the above expression for in equation (B.11) (the condition under which the statistical distribution for can be assumed Gaussian), and observing that and , we obtain that the condition to be satisfied so that the distribution for is Gaussian is . Under this condition, we observe that the optimum of the phase-unconstrained problem converges to the optimum of the phase-constrained problem.
C. Locally Optimal FA Shimming Configuration
Given a FA shimming pulse , and a given position , let us define as follows:
This new pulse is simply obtained by scaling the RF waveform by the factor . thus satisfies
In other words, satisfies in the equality constraint of the local optimization problem given by equation (17). Thus we have by definition of (the solution of equation (17)):i.e.,and finally,
To complete the proof, we note that the function is a convex function and so that for any random complex random variable , . Applying this to yields
Data Availability
The datasets generated and/or analyzed during the current study and the pulse design source codes are available from the corresponding author on reasonable request.
Conflicts of Interest
The authors declare that they have no conflicts of interest.
Acknowledgments
This research received funding from the RETP program (Research Equipment and Technology Platforms) of the Leducq Foundation, the European Research Council under the European Union’s Seventh Framework Program (FP7/2013-2018), and ERC Grant Agreement #309674 and Proof of Concept #700812. The authors are grateful to Dr. B. Thirion for his valuable input regarding the estimation of the covariance matrix of the control field.