Optimization of high-channel count, switch matrices for multinuclear, high-field MRI

Modern magnetic resonance imaging systems are equipped with a large number of receive connectors in order to optimally support a large field-of-view and/or high acceleration in parallel imaging using high-channel count, phased array coils. Given that the MR system is equipped with a limited number of digitizing receivers and in order to support operation of multinuclear coil arrays, these connectors need to be flexibly routed to the receiver outside the RF shielded examination room. However, for a number of practical, economic and safety reasons, it is better to only route a subset of the connectors. This is usually accomplished with the use of switch matrices. These exist in a variety of topologies and differ in routing flexibility and technological implementation. A highly flexible implementation is a crossbar topology that allows to any one input to be routed to any one output and can use single PIN diodes as active elements. However, in this configuration, long open-ended transmission lines can potentially remain connected to the signal path leading to high transmission losses. Thus, especially for high-field systems compensation mechanisms are required to remove the effects of open-ended transmission line stubs. The selection of a limited number of lumped element reactance values to compensate for the for the effect of transmission line stubs in large-scale switch matrices capable of supporting multi-nuclear operation is non-trivial and is a combinatorial problem of high order. Here, we demonstrate the use of metaheuristic approaches to optimize the circuit design of these matrices that additionally carry out the optimization of distances between the parallel transmission lines. For a matrix with 128 inputs and 64 outputs a realization is proposed that displays a worst-case insertion loss of 3.8 dB.


Introduction
Modern magnetic resonance imaging (MRI) systems rely heavily on parallel imaging techniques using multi-channel, phased array coils. Currently, they are routinely employed to improve image quality, to reduce the scan time or to obtain more information (e.g. in multinuclear or multi-parametric measurements) within a given time frame. Moreover, reducing the scan time brings about further benefits in terms of increasing patient compliance (shorter scan times are more comfortable for the patient) and enables an increased workload of expensive MRI machinery (patient throughput can be increased). In fact, several applications would not be feasible at all without parallel imaging techniques. These include, but are not limited to, imaging the beating heart [1,2], fast functional/diffusion imaging [3][4][5][6][7] and correcting unwanted artifacts, e.g. aliasing and ghosting [8,9]. With the increase in the number of elements in phased array coils [10][11][12], more and more parallel receive channels are required per MRI system. Unfortunately, this has a significant negative impact in terms of machine cost due to the requirement for additional receiver units and a technologically complex patient table [13], on patient comfort and operator workflow arising from the bulky cable bundles on the local RF coils [13], and on the large amounts of data that require handling prior to final image reconstruction. Thus, the number of receive antennas (coils)-or at least the number of physical connection points in the patient table [14]present in an MRI examination frequently exceeds the number of receivers available in the system. The problem is aggravated when large fields-of-view are to be covered and multiple coil arrays are used during a scan, e.g. in spine imaging and moving table acquisitions.
An early commercial solution to the problem was the incorporation of mode matrices in the receiver coil array [15]. These make use of the fact that the majority of the receive signal's power is contained in a small number of RF modes. The concept can be compared to the use of circular polarized receive coils which employ two quadrature coil elements and combines these into a single receive signal. Thus, high-channel-count coil systems can be connected to fewer receiver channels, while avoiding significant penalties in the receive signal-to-noise ratio (SNR). However, the mode matrix approach significantly reduces the available acceleration factor for a given number of receive elements as it actually reduces the physical number of receive elements [16]. This factor is particularly significant in high-field systems, which are increasing in availability, as they potentially allow higher speed-up factors [17]. Thus, most modern MRI scanners have a large number of physical coil connectors distributed on the patient table that need to be flexibility connected to the available receivers with the established routing depending on the coils selected for the desired imaging task e.g. when full body coverage is desired [14] or multiple different type of local coils are employed in a single imaging session [13].
Signal routing from the physical coil connectors commonly located on the patient table to the receivers located in the technical room outside the Faraday shield in clinical MRI systems is usually achieved by using a module termed either the "switch matrix" (a term originally used in telecommunications and computer communication networks (e.g. [18] or "matrix switch" (e.g. [19]. The switch matrix is commonly located close to the magnet to avoid long coaxial cable bundles but physically behind the initial low noise preamplifier in order to maintain the highest possible SNR. A large number of implementations are feasible and can generally be categorized as crossbar switches or multiplex matrix switches [20]. Both of these, as well as hybrid topologies, have been used in MRI applications [21,22]. Crossbar switch matrix implementations can either be achieved with the use of integrated circuits configured as single-pole double-throw (SPDT) switches, e.g. based on GaAs, GaN [23], or MEMS [24] switches, or by employing simpler single-pole single-throw (SPST) configurations, e.g. by using PIN-diodes. While the former implementation seems to be technologically more robust, as it allows to disconnect the remaining transmission line length from the active signal path and therefore does not contain remaining physically connected open-ended transmission lines, it requires the receive signal to pass through a high number of active switches. This potentially leads to high signal attenuation. As an example, a well-designed GaAs IC with suitably high linear power handling capability for MRI receive signals easily displays an insertion loss of 0.35 dB (e.g. SKY13351-378LF-Skyworks Solutions, Inc., USA). The receive signal has to cross this device (M-1) x N times in the worst case with M being the number of input lines and N the number of output lines.
Recently, PIN-diode controlled switch matrices have been shown to be feasible even in high-field MRI applications by either compensating open transmission line impedances with suitable lumped element terminations in one direction [25] or in both directions [26]. In addition, a method to reduce the PCB footprint for large size matrices using a combination of switch types has been presented [27].
To date, switch matrix design has mostly been concerned with routing the signals of the proton Larmor frequency as the number of available high receive channel arrays operating at non-proton frequencies at clinical field strength is limited to research coil implementations, e.g. [28,29]. However, the routing of X-nuclei (non-proton, resonating at different frequencies) signals might still be required if the coil connectors on the patient table do not provide dedicated plugs for the X-nucleus' receive coil arrays.
A number of switch matrix implementations with variable degrees of flexibility in routing configuration have been described in the literature: sparse or cascaded matrix designs are described in [13,30,31] and a full switch matrix is described in [32]. While they are routinely employed in MRI scanners operating at clinical field strength of 1.5 T and 3 T, in some commercial 9.4 T and 7 T systems, the switch matrix was removed and the signal was routed directly without the possibility to switch it to different receiver units. However, with the increased availability of high-field scanners, the benefits of X-nuclei imaging is likely to attract increased interest [33,34]. As multichannel, X-nucleus coil arrays also provide improved SNR compared to volume coils, as is also the case in proton imaging, there is compelling motivation to also use coil arrays in these SNR-starved applications. In fact, multichannel X-nucleus arrays have been under active development for a number of years, and one can compare exemplary implementations in [35,36] or [37] for a review on the topic. In this context, designing a switch matrix suitable for operation with protons and the most common X-nuclei at high-field is desirable. It would, for example, enable the combined operation of a CP mode proton coil with a 32-channel sodium array in a system equipped with 32 receive channels or the connection of more advanced double-tuned arrays as described in [38]. Note, that switch matrix implementations operating at the proton frequency as well as at one or more X-nucleus frequency has been previously described in [39].
In this manuscript we present a novel design methodology for a low-loss, high channelcount switch matrix based on SPST switches employing PIN diodes [26] capable of operating at, for example, 400 MHz (proton -1 H), 376 MHz (Fluorine -19 F), 162 MHz (Phosphorus -31 P), 106 MHz (Sodium -23 Na) and 54 MHz (Oxygen -17 O), these being the resonance frequencies of the said nuclei on a 9.4 T MR system. The method is based on multi-parameter optimization strategies and uses analytical transmission line formulas to derive optimum transmission line spacing and compensation elements on both horizontal and vertical switch matrix transmission lines. Fig 1 shows the topology of a switch matrix employing a single PIN-diode at each junction. Note that single PIN-diode switches can make a connection between an input and an output, but cannot break it-they are acting in a SPST manner. Therefore, the transmission line stubs remain electrically connected. Thus, attenuation of the signal occurs, with the level of attenuation depending on the impedance transformation of the open transmission line to the location where the PIN-diode connects upper and lower transmission lines. The remaining vertical and horizontal stub lines will contribute to the signal attenuation across the switch matrix as described in [26].

Methods
In order to demonstrate the implementation challenges, the attenuation arising from the unterminated transmission lines for a low channel-count (5 inputs x 5 outputs) switch matrix has been measured on a prototype board [26] for the connections shown in color in Fig 2 for a number of nuclei frequently used in high field MRI. Even for the small-sized matrix, the task of determining suitable compensation elements is tedious and evaluation of optimal compensation element values quickly becomes infeasible. Thus, for the sake of simplicity, Fig 2 reports results obtained for a single compensation element of 91 pF. In this prototype "switching" from unterminated to terminated stub lines was accomplished by soldering the termination impedance on the PCB. Only a limited number of switch matrix configurations with different stub line lengths were evaluated experimentally in order to demonstrate the negative influence of unterminated transmission lines on the IL per se. In a product implementation connection of the designated termination resistance will be implemented via SPST PIN diode switches.
The worst-case IL in the small laboratory sample matrix evaluated in Fig 2 is 24 dB in the uncompensated case, while it is only 1.94 dB for a connection from input #1 to output channel #5 on the proton frequency when compensation elements are employed. To investigate the influence of IL on SNR and in order to establish the requirements for the circuit, optimization experiments were carried out on a 9.4 T animal scanner [40] by inserting a variable attenuator at the location of the proposed switch matrix. Single slice gradient echo (TE = 10 ms, TR = 100 ms, FA = 15) images were acquired using a linear-polarized birdcage resonator and a 63 ml sample (per 1000 g distilled water 5 g NaCl and 1.25 g NiSo 4 x 6 H 2 O) to measure SNR as a function of attenuation as described in [41], Method 2. The results of the measurement are given in Fig 3.

Fig 2. Signal routes (colored) evaluated at different frequencies (left) and insertion loss measured on the bench for the different nuclei considered
here. If the insertion loss was above 2 dB, a single value compensation (91 pF) was attached at one or both transmission lines and the insertion loss remeasured for all possible combination of compensation elements. All values are in dB; all input lines are designated as "Input" and all output lines as "Channel" with consecutive numbering. Transmission line stubs are labeled with alphabetical capital letters A to J. For the sake of clarity, only a limited number of matrix configurations have been evaluated experimentally, e.g. the switch was configured with fixed connections from input #1 to output #5, input #2 to output #4 and so forth to demonstrate the effect of transmission line stubs of different length for different resonant frequencies. In the table on the right, the value for the insertion loss shown in black indicates values measured for the unterminated case (that is, without the 91 pF attached to the transmission line stub). Red values are the ones representing IL for the terminated case (where the 91 pF are connected to the transmission line stubs). For the terminated case, the location of the termination capacitance is given underneath the measured IL values.
https://doi.org/10.1371/journal.pone.0237494.g002 The influence of signal attenuation on image quality for multiple nuclei was evaluated experimentally on the same NR system using an interchangeable set of quadrature birdcage coils with integrated transmit/receive switches and preamplifiers [42] and the prototype switch matrix shown in Fig 1 placed behind the coil connectors. Imaging of a phantom (50 ml cylindrical tube filled with doped water and 150 mM NaF (Sigma-Aldrich, German)) was carried out using the FLASH sequence at the 1 H, 19  As can be seen from Figs 3 and 4, small values of attenuation already have a negative influence on the SNR of the image. Therefore, the switch matrix should be optimized for as low IL as possible. Consequently, for larger switch matrix sizes, which are likely to be encountered in the modern MRI systems, an algorithm is required that determines optimum compensation element values in an automated way. This optimization problem can use a single cost function (C) and minimizes the maximum insertion loss (IL) encountered for all switch configurations given that the optimum stub termination is selected for a given signal routing scheme. For a M x N switch matrix this is expressed as with M and N being the number of input and output lines, respectively, and f k being the operating frequencies. Eq (1) is similar to using a weighted Tchebycheff model with weighting coefficients set to equal one as described in [43]. Due to the large number of variables equaling possible terminations for each stub line, metaheuristic approaches [44] are suitable for optimization of the switch matrix circuit. Metaheuristics can be classified as single solution approaches which act on a single candidate solution, e.g. simulated annealing [45] and variant local search methods, or population-based approaches [46]. The latter includes natureinspired algorithms, such as evolutionary algorithms, particle swarm optimization [47], firefly algorithms [48] and artificial bee colonies [49]. While metaheuristic population-based approaches are not guaranteed to find the globally optimum solution, a recent comparison has shown that they are competitive with deterministic Lipschitz algorithms [50] in many respects. This, along with their relative ease of use, are the motivation for applying them to the switch matrix optimization for multi-nuclear, high-field MRI. IL can be computed as the sum of the mismatch loss and transmission losses over the transmission lines used in the switch matrix, e.g. compare [51]. Mismatch loss is computed by first transforming the compensation impedance of the two stub lines to the active switch junction. This uses the impedance transformation formula for transmission lines of characteristic impedance Z 0 terminated with the impedance Z L along the stub line length l with γ being the complex propagation constant of the transmission line [52].At this point, the parallel impedance of the stub lines and the 50 Ohm load from the receiver are computed and further transformed to the switch matrix input. This enables the impedance, and hence the mismatch loss at this point, to be calculated. The microstrip transmission lines used in the calculations have been designed to display a characteristic impedance of 50 Ohms when built on an FR4 substrate with a thickness of 1.5 mm (relative dielectric constant ε r = 4.6, copper trace thickness t = 0.035 mm, substrate thickness h = 1.5 mm, microstrip trace width w = 2.725 mm, loss tangent of substrate tan δ = 0.022, surface roughness-rms deviation of the conductor surface from a plane-Rough = 0.055, and relative conductivity with respect to copper Rho = 1). The transmission line parameters were calculated online (http://mcalc.sourceforge.net) using a method which employs the formulae derived by Hammerstad and Jensen [53] and is conveniently summarized in e.g. [54].

Implementation of circuit optimization
Circuit optimization was carried out with Python [55] using the platypus package [56]. The cost function is implemented in two subroutines-one computing the optimum termination and returning minimal insertion loss for a single switch junction, while the outer subroutine computes the maximum insertion loss over all junctions. The cost function is parameterized so as to allow the maximal insertion loss for a number of operating frequencies to be found, as well as having the choice to place more than one compensation element per transmission line stub. It should be noted that the range of termination elements allowed was limited from 300 nH inductive to 200 pF capacitive, which covers available lumped, high-frequency elements with suitable self-resonances and quality factors. This range was also mapped onto a real valued space, between +1 and -1 for equidistant gridding during program execution to ensure that standardized initialization of the search space remained feasible. The final goal of the study was to design a multi-nuclear switch matrix with 128x64 channels. The dimensions were selected according to the maximum number of receiver channels available in commercial 7 T systems and the maximum number of receive elements in an experimental coil array described so far [57].I Initial investigations of algorithm performance were carried out on differently sized matrices for the sake of simplicity.
In an initial evaluation, several multi-objective optimization algorithms based on evolution strategies-a variant of evolutionary computing which is inspired by acting upon populations underlying variation and selection in generational loops [58]-were compared in order to evaluate their performance for the switch matrix optimization task. This was done for a 64x32 sized matrix operating at the proton frequency of 400 MHz and allowing two compensation elements for each stub. For this investigation, the transmission lines were assumed to be lossless. The strategies investigated were: a non-dominated sorting genetic algorithm II (NSGA-II) [59] and its more recent version NSGA-III [60], a covariance matrix adaptation evolution strategy (CMA-ES) [61], a generalized differential evolution (GDE3) [62], an indicator-based evolutionary algorithm (IBEA) [63], a multi-objective evolutionary algorithm based on decomposition (MOEA/D) [64], a multi-objective (OMOPSO) [65] and speed-constrained multi-objective particle swarm optimizer (SMPSO) [66], a strength Pareto evolutionary algorithm (SPEA2) [67] and an epsilon multi-objective evolutionary algorithm (ε-MOEA) [68]. In addition to using different optimization methods, the strategies employed also differ in their approaches to evaluating the Pareto optimality of the population. For example, the NSGA derivatives employ dominance depth, while SPEA uses dominance count/rank; a different approach is to use performance measures for the selection step, e.g. in IBEA [43]. For an indepth discussion of the properties of the different algorithms, the reader is referred to standard textbooks, e.g. [43]. A list of active projects implementing multi-objective optimization in the Python programming language is given in [69]. All algorithms were executed with 10 random seeds and ran for a maximum of 2000 iterations.
Using the optimization algorithm selected during the initial evaluation, a multi-nuclear switch matrix with 16 inputs and 8 outputs was designed for operation at the 1 H, 19 F, 31 P, 23 Na, and 17 O Larmor frequencies of a 9.4 T system. This time, attenuation on the microstrip lines was accounted for by using a complex propagation constant for each frequency. In addition, a variable, but unique, spacing between neighboring transmission lines in the range 10 mm to 60 mm was considered during the optimization to allow the circuit geometry to be optimized in combination with the compensation elements. The lower limit was imposed so as to avoid closely spaced lines with high coupling between neighboring lines. Compensation elements were allowed from a maximum of 300 nH inductive to 200 pF capacitive, as above. This more realistic design example was evaluated with the three best-performing algorithms found during the initial comparison.
Finally, a large sized multi-nuclear matrix with 128 input lines and 64 output lines was optimized using the same constraints as given above. This time, only the ε-MOEA algorithm was investigated as it provided the best results for the two test cases investigated beforehand. The number of algorithm iterations was increased until no further improvement in the overall cost could be obtained. Table 1 shows the highest insertion loss for any path across a 64x32 sized matrix obtained from optimizing the compensation elements with the different algorithms. It can be seen that NSGA-II, SPEA and ε-MOEA performed best in this task, while other algorithms converged to significantly worse solutions, despite being executed for a number of seed values.

Results
The three algorithms that converged to an insertion loss below around 1.2 dB were further evaluated on a more realistic multi-nuclear switch design task. This task was based on a 16 x 8 switch matrix and the insertion losses for each nucleus and each switch configuration are shown in S1 to S3 Figs. The overall results for the optimization with a maximum number of algorithm executions of 4000 times are given in Table 2.
The results for the full-sized switch matrix (128 x 64 ports) designed with ε-MOEA are given in Fig 5. Table 3

Discussion
We have shown that SNRs of both, proton and X-nucleus signals can be significantly affected when using a switch matrix to route receive signals from the coil connectors on the patient table to the receivers of an MR system. Thus, suitable circuit topologies that reduce signal attenuation in the switch matrix as much as possible are required. In this manuscript we demonstrate that circuit optimization of a large, multi-nuclear switch matrix is feasible even when using a limited number of compensation elements to compensate for the influence of open-ended transmission lines of varying length. Optimization, although using a single cost function only, was carried out with a variety of multi-objective optimization strategies which employ either swarm like features or genetic principles and are thus suitable candidates for metaheuristic, non-convex problems. The non-convexity of the optimization problem at hand arises from the magnitude operation in the IL formulation. These algorithms have been For better visibility, the switch matrix was rotated so that the output ports are shown on the left-hand side (abscissa) ranging from 1 (top) to 64 (bottom). The inputs are arranged along the ordinate and increase from right to left. Each color patch gives the maximum IL for the connection between the corresponding input and output port. All data are simulation results calculated using the transmission line formulas from Hammerstad et al. [53].
https://doi.org/10.1371/journal.pone.0237494.g005 employed successfully in other fields of MRI, e.g. for B 1 shimming [70] using different single cost functions and a particle swarm approach and in coil circuit optimization [71] with three optimization targets and employing CMA-ES. In the presented optimization task, the optimization was able to generate a circuit topology that caused an addition IL of only 0.6 dB above that of the attenuation of the longest transmission line. With a maximum overall IL of below 4 dB, the degradation in image SNR is negligible, as we demonstrated with attenuation measurements on a 9.4 T scanner.
In our investigation, we found ε-MOEA to perform best out of the tested candidate algorithms. While this is a lesser-known evolution strategy [72] and was shown to be less performant in an biobjective testbed [73] it consistently outperformed the other algorithms implemented in the platypus collection. The good performance of the ε-MOEA algorithm might be due to the fact that, in contrast to generational algorithms which evolve the entire population at every iteration, ε-dominance archives ensure convergence and diversity throughout search for a proper selection of the ε parameter [74]. However, further investigation is required to back up this hypothesis.
From initial tests using multiple cost functions-single cost functions for all switch matrix junctions-it quickly became evident that one combined target function performed superior. This is in line with [70], where the formulation of multiple optimization goals was found to be feasible but was not chosen in the final implementation. The number of algorithm executions for the optimization to converge of around 16,000 is in the range reported in [50]. While the metaheuristic approaches employed are not guaranteed to converge to the global optimum, a comparison with the IL on the longest transmission line arising from attenuation only lends reasonable credit only to the results obtained with the proposed optimization strategy.
It should be noted, therefore, that for the investigations presented here, no experimental validation has been carried out for the full-scale matrix, since the intention was solely to investigate a feasible solution of the given optimization problem. However, from the investigations carried out in prior work [26,27], we feel that agreement between simulation and measurement has been sufficiently demonstrated. Further validation is given in S4 Fig which compares IL for a 6x3 matrix optimized with the proposed algorithm with those obtained from a circuit simulator. The impact on of insufficient transmission line compensation on MRI image quality has been investigated both previously and in this work and allows a compromise between the envisaged switch matrix implementation and the associated image degradation to be found.
While this work focuses on the classical topology of a PIN-diode based switch matrix, which results in a rather large PCB footprint-e.g. the 128 x 64 switch matrix designed in this work measures approximately 1250 mm x 660 mm-the workflow can be used in combination with alternative matrix implementations. One possibility could, for example, be to use the switch matrix designed in this work to optimize the sub-matrices presented in [27] for stacked matrix designs or to reduce the flexibility of the switch matrix by using a multiplexer circuit that combines a number of input lines prior to routing through the switch matrix, e.g. as suggested in [21,25]. It is certainly also feasible to change the number of allowed compensation elements, e.g. using a single termination only as described previously or to increase the number of elements to three to decrease the maximum IL further. Both plots are normalized with respect to their maximum IL in color coding. While the pattern looks similar the circuit simulator predicts slightly larger insertion losses in all cases, which is probably due to using slightly different formulae for calculating transmission line properties. (PNG)