System-informed magnetic resonance imaging (MRI) sequence design and control

Abstract

The present disclosure provides a nuclear magnetic resonance (NMR) imaging method, comprising at least one or more of: generating, via a first NMR simulation, an idealised target image assuming substantially idealised NMR scanner performance, and using an initial set of idealised NMR scanner control inputs (x); applying the idealised control inputs to an MR scanner system response model to obtain a set of predicted realised responses (B) representative of MR scanner operating parameters; generating, via a second NMR simulation, a predicted image using the set of predicted realised responses (B); and adapting the initial set of idealised NMR scanner control inputs (x) to give an optimised set of scanner control inputs in dependence on minimisation of a loss function that receives as inputs thereto at least the idealised target image and the predicted image and minimises errors therebetween.

Description

[0001] SYSTEM-INFORMED MAGNETIC RESONANCE IMAGING (MRI) SEQUENCE DESIGN AND CONTROL Technical Field

[0002] The present disclosure generally relates to the use and control of nuclear magnetic resonance (NMR) imaging systems and other magnetic resonance imaging (MRI) systems. More specifically, the present disclosure relates to the practical optimization of NMR and MRI systems and the need to ensure high accuracy imaging within the physical capabilities of NMR and MRI systems.

[0003] The present disclosure seeks to be applicable with all types of NMR and MRI scanners and seeks to combat negative effects on NMR and MRI scanner images such as eddy currents.

[0004] Background and Related Art

[0005] MRI scanners are engineered to give close to ideal performance such that acquisition methods can operate under this assumption. Eddy currents and amplifier non-linearities can significantly degrade image quality in MRI scans. Traditional compensation methods include pre-emphasis of input gradient currents but this creates a voltage overhead that is undesirable especially for lower performance hardware, which is becoming more widespread. Moreover, hardware imperfections often exist and can cause significant artefacts that degrade the quality of the acquired images. These imperfections are rarely considered when designing MR sequences and instead are overcome by applying separate mitigation strategies. This can lead to an inefficient use of resources, for example. We propose a new framework for sequence design: non-idealized system (NIS) optimization that accounts for hardware imperfections and learns novel strategies to overcome them by modifying gradient and RF properties as well as sampling properties.

[0006] Existing MRI systems (and NMR systems in general) treat a specification of the conditions to be realised during scanning as the primary specification, as shown in Figure 1. It is common to treat each sub-system (notably the main magnet, gradient, RF, data acquisition and patient handling systems) largely independently and as if they have idealised performance. This often includes a correction layer (e.g. pre-emphasis for the gradient system and non-linearity correction for the RF system) to pre-distort the drive signals so that the realised output is as close as possible to that specified for the task in hand. These corrections are often static, with settings determined in a completely separate calibration process, such as on a service day or at system set-up, but may also be updated during a pre-scan calibration step. Also, corrections may be applied to one sub-system to mitigate a non-ideal response of another subsystem. For example, shifting of the main clock frequency to correct for thermal drift in the main magnet, or for field shifts due to known eddy currents during gradient activity, or movements of the patient bed to keep a target region in the most favourable place in the scanner. Also, the existing state of the art favours templated timing structures with repeating loops whenever possible. This approach requires the scanner hardware, with its corrections, to approximate an ideal system that can achieve the specified scan conditions without persistent memory effects that cause responses at one moment to be influenced by the demand at a prior moment. Of course, it is widely recognised that this "ideal system with correction" is not perfectly achieved, but the response is usually to seek a further layer of correction.

[0007] In other words, different (MRI) sequences can generate different image contrasts that are used to identify specific pathology. The current design paradigm assumes a static system configuration that delivers ideal performance with fixed loop control systems that treat MRI sequences as repetitions of identical events. However, these sequences may not perform as expected due to hardware imperfections. For example, rapidly changing magnetic fields (required to generate MR images) can induce electrical (eddy) currents in nearby conducting media that degrade the quality of the acquired images1and impede diagnosis. Thus, a scenario arises where it would be beneficial to optimise the process of MRI or NMR sequences taking into account the non-ideal nature of the hardware associated with these systems. Therefore, allowing for much lower specification sub-systems that achieve satisfactory results by working in a unified framework and operate within the limits of the associated hardware, or achieving higher performance from existing specification subsystems. These optimisation systems may be implemented by utilising neural network or machine learning techniques.

[0008] MRzero is an existing proposed supervised self-learning pipeline that can discover MR sequences via knowledge of the Bloch equations only [1]. For sequence parameter optimization, this method requires a fully differentiable MR simulation process such that gradient descent can be used. It is this simulation process that we have modified to incorporate a realistic system model that can comprise any number of hardware imperfections. Loktyushin et al. performed simulations using an isochromat-based approach but here, we have adopted the more efficient phase distribution graph (PDG) method that is closely related to extended phase graphs (EPGs) but is made suitable for whole-sequence optimization with arbitrary gradients and timings and accounts for T2' relaxation [1].

[0009] Summary of the Invention Over the years, many research groups have developed techniques to separately overcome scanner hardware imperfections but this leads to an inefficient use of resources. Continuing with the example mentioned above, eddy currents are reduced by modifying the input current to the MR gradients (pre-emphasis) though this creates a voltage overhead. Individually addressing these imperfections can also lead to long and complicated acquisition protocols that impede clinical translation, so sacrifices are often required. The present disclosure represents a new total sequence design framework that starts with the actual system properties and can deliver maximum performance given all available resources. This way, we can directly account for hardware imperfections and learn novel acquisition strategies to overcome them. The resultant designed sequences will be inherently corrected for any imperfections we choose to include in the simulated MR signal model and so artefact-free images will be produced from a resource and timeefficient acquisition protocol.

[0010] The Non-idealized system (NIS) optimization approach is a new design framework that starts with actual system properties and delivers maximum performance given all available resources. It directly accounts for system imperfections and learns novel ways to overcome them.

[0011] The present disclosure largely considers a single but common hardware imperfection: eddy currents. However, this is only for exemplary reasons and the present disclosure seeks to optimize the imaging process in relation to a plurality of hardware imperfections. Rapidly changing magnetic fields are applied during an MR pulse sequence and induce currents, known as eddy current, in nearby conducting structures of the scanner. These create unwanted magnetic fields that have a detrimental effect on image quality [3]. Eddy currents accumulate as gradient waveforms ramp up and down but decay during constant periods. Since k-space locations are defined as the cumulative time integral of the area under a gradient waveform, linear eddy currents can cause a shift in actual k-space locations compared to the expected locations. Eddy current-induced imperfections have become more significant with the advent of high-performance gradient systems since eddy currents are larger for higher gradient amplitudes or faster slew rates. Pre-emphasis is the most typical correction method, whereby the waveform input to the gradient coil is intentionally distorted to counteract the subsequent eddy current perturbation. However, for some applications, lower performance scanners may not be able to achieve this voltage overhead due to tighter gradient constraints and therefore, more sophisticated solutions are necessary. In embodiments of the present invention, the required signal performance is specified using input parameters that enable forward modelling of a target image or other target signal. This is done using a forward model, that could be ideal system performance or other approximate model. Our approach is then to optimise the available control inputs (x), so that the actual system achieves as close as possible measured signals. That is, we seek to minimise the difference between the target image or target signal and the realised image or signal allowing free control of the inputs, subject to hardware constraints such as maximum voltage, current, temperature etc, or depletion of charge reservoir, available time for the examination etc. A feature of this approach is a unified control strategy with simultaneous and continuous adjustment of all control inputs, quite likely in a nonrepeating pattern.

[0012] The NIS optimization framework is presented with a common hardware imperfection (eddy currents) simulated in the forward MR signal model. Optimizations of gradient waveforms (and flip angles) are performed with different gradient (amplitude and slew rate) limits and a comparison is made to pre-emphasis.

[0013] Without gradient constraints, NIS optimization returns some pre-emphasis behaviour and adapts waveforms to combat eddy currents with different time characteristics, giving negligible reconstruction errors. Pre-emphasis solutions for scanners with low amplitude limits require an increase in acquisition time whereas NIS optimization solutions automatically fit within any specified timings. Novel sequences are discovered for constrained scenarios for which eddy current compensation would ordinarily be unachievable. NIS optimization is equally successful for transient-state and steady-state gradient echo sequences; final images consistently have low error. Potential fine-tuning is possible via common machine learning techniques such as data augmentation.

[0014] Therefore, that is to say that the present disclosure seeks to optimize control inputs for the scanner system with the aim of producing the desired outcome (image or NMR signal) considering a full model of the system response. Essentially the non-ideal nature of the system is brought into the foreground and the system controls are explicitly optimized taking account of these in order to produce the target image / signal. This is in contrast with the existing approach, where hardware imperfections are corrected as an additional layer sitting between control inputs defined for image formation, and the instructions sent to hardware. Although the existing approach simplifies the process by decoupling hardware compensation from image formation, it does so at a cost in performance. A key benefit of the proposed method is that it allows a more diverse range of solutions, where conventional pre-compensation is not effective, since the focus is to produce the specified "target" image, not to hit predefined control inputs.

[0015] A first aspect of the present disclosure provides a nuclear magnetic resonance (NMR) imaging method, comprising at least one or more of: generating, via a first NMR simulation, an idealised target image assuming substantially idealised NMR scanner performance, and using an initial set of idealised NMR scanner control inputs (Xi); applying the idealised control inputs to an MR scanner system response model to obtain a set of predicted realised responses (B) representative of MR scanner operating parameters; generating, via a second NMR simulation, a predicted image using the set of predicted realised responses (B); and adapting the initial set of idealised NMR scanner control inputs (Xi) to give an optimised set of scanner control inputs (x0) in dependence on minimisation of a loss function that receives as inputs thereto at least the idealised target image and the predicted image and minimises errors therebetween.

[0016] This represents a step-change in MR sequence design. The framework inherently models and mitigates system imperfections, removing the need for separate correction strategies. In theory, any system factors can be included in the model. This "all-in-one" approach is in stark contrast to existing techniques that treat each scanner sub-system as already corrected to achieve ideal performance and designing only on that basis. By working with the drive signals directly and optimizing using imaging signal performance directly, scanner resources can be used as efficiently as possible. The present framework produces optimal sequences that simultaneously meet specified hardware constraints, account for modelled imperfections, and maintain sufficient reconstructed image quality.

[0017] In some aspects, the method may further comprise: operating an NMR scanner with the optimised set of scanner control inputs to generate an NMR image of a subject.

[0018] The above-mentioned aspects have the benefit of reduced hardware complexity and hardware strain. This in turn has the benefit of improved reach in relation to differing makes, models and price ranges of MRI scanners; as the present invention as currently disclosed can be utilised on many more MRI scanner systems then previous optimization systems that required large system overheads. Moreover, the present invention can be implemented on NMR or MRI scanner hardware that cannot be corrected by conventional methods. As such, this approach can allow for much lower specification sub-systems that achieve satisfactory results by working in a unified framework. In some aspects, the optimised set of scanner control inputs may have been optimised so as to be within corresponding system constraints of an NMR scanner to be operated using the optimised set of scanner control inputs (x0). Thus, the corresponding system constraints of the NMR scanner may comprise one or more selected from the group comprising: i) maximum voltage(s); ii) maximum current(s); iii) maximum temperature(s); iv) depletion of charge reservoir(s); v) available time to scan a subject; vi) maximum rate of change of voltage(s); vii) maximum rate of change of current(s); and / or viii) maximum power limit(s).

[0019] In some aspects, the initial set of idealised NMR scanner control inputs (Xi) and / or the optimised set of scanner control inputs (x0) may comprise one or more selected from but not limited to the group comprising: i) one or more demanded gradient waveforms (g(t)); ii) one or more radio frequency (RF) waveforms (rf); iii) one or more shim currents (Ishim); iv) one or more patient handling system controls

[0020] In some aspects, the predicted realised responses (B) may comprise one or more selected from but not limited to the group comprising: i) gradient fields, G; ii) RF fields, RF; iii) static fields, Bo; and iv) temperature response.

[0021] In some aspects, the MR scanner system response model may further comprise one or more physical models representing system imperfections or system constraints in the NMR scanner to be controlled with the scanner control inputs. Further, the MR scanner system response model may be: i) based on prior characterisation of the system, or ii) based on learned responses of resulting signals to changes in the control inputs, wherein the learned responses may be in the form of parametrised models or learned networks.

[0022] In some aspects, the NMR scanner may be recognised to be an imperfect NMR scanner system, wherein the NMR scanner system response model may be used to achieve control of the imperfect NMR scanner system, taking into account and using its imperfections to maximise possible performance.

[0023] In some aspects, after optimization the optimised set of scanner control inputs, x, may comprise a set of time resolved control inputs that may explicitly account for system imperfections to achieve an output that is as close as possible to the target output subject to hardware limits of the MR scanner.

[0024] In some aspects, the optimised control inputs may be arrived at by minimising a loss function. For example the loss function may be a composite loss function comprising: (i) an image loss term; (ii) a k-distance loss term; (iii) a voltage term; and (iv) a slew rate term. Further, the loss function may comprise the following:

[0025] Total Loss = Lj + Lk+ Lg+ Ls+ LRF[2] wherein: and wherein Vmaxand Smaxmay be voltage and slew rate limits respectively, w may be tunable regularization weights, IT may be the idealised target image, Ic may be the predicted image, kTmay be a target k-space, kc may be a predicted actual k-space, and gmaxand smaxma be voltage and slew rate limits respectively and w may be tunable regularization weights and g(t) may be defined initial gradient sequence waveforms.

[0026] In some aspects, the loss function may include at least one or more of the following: (i) an image (signal) loss that minimizes reconstruction error with respect to the target image or signal; (ii) a sampling loss that drives towards a desired sampling scheme; (iii) a hardware loss that penalizes solutions which exceed hardware limits of the MR scanner system.

[0027] In some aspects, by working with the set of scanner control inputs directly and optimizing using imaging signal performance directly, scanner resources may be used as efficiently as possible.

[0028] The above is advantageous as it reduces the required hardware performance, improves processing time whilst also ensuring high levels of NMR imaging quality and accuracy.

[0029] In some aspects, wherein sample points are located on a Cartesian grid or other specified sampling scheme there may be scope to shift or omit some sample points as part of the optimisation process in combination with optimising the control inputs to favour reconstructions of the NMR image of the subject. In some aspects, sample points are applied using a sampling scheme such that sample points can be shifted or removed to optimise reconstructions of the NMR image of the subject.

[0030] A second aspect of the present disclosure provides a system, comprising at least one or more of: a processor; one or more computer readable medium, storing instructions, which when executed by the processor, cause the processor to perform the method of any of the preceding aspects.

[0031] A third aspect of the present disclosure provides an NMR scanner having a set of scanner control inputs, the scanner control inputs being set to correspond to the values of the optimised set of scanner control inputs generated by the method of any of the preceding aspects.

[0032] Further features and advantages will be apparent from the appended claims.

[0033] Brief Description of the Drawings

[0034] Further features and advantages of the present invention will become apparent from the following description of an embodiment thereof, presented by way of example only, and with reference to the accompanying drawings, wherein like reference numerals refer to like parts, and wherein:

[0035] Figure 1 shows a flowchart of a classical gradient-only optimization system in accordance with the prior art;

[0036] Figure 2 shows an example flowchart of the gradient-only non-idealized system optimization in accordance with an embodiment of the present disclosure;

[0037] Figures 3A and 3B show the Eddy current responses, with two time constants simulated that represent the perturbation model used for forward MR signal simulation and the perturbed (realized) waveforms Goor B for both time constants and the unperturbed (demanded) equivalent waveforms g0or x, respectively, in accordance with an embodiment of the present disclosure;

[0038] Figure 4 shows an example block diagram of the proposed framework in accordance with an embodiment of the present disclosure; Figures 5A-D show graphs of the target waveform and its perturbed counterparts for short, long and both eddy currents and present optimized demanded and equivalent realized waveforms for short, long and both eddy currents respectively, in accordance with an embodiment of the present disclosure;

[0039] Figures 6A-D show the NIS optimization results for an unconstrained scenario when simulating both eddy currents, in accordance with an embodiment of the present disclosure;

[0040] Figures 7A-D show comparison graphs between true pre-emphasis and NIS optimization; in accordance with an embodiment of the present disclosure;

[0041] Figures 8A-C show the NIS optimization results for amplitude limited cases (those in Figure 7 A) when modelling both eddy currents, in accordance with an embodiment of the present disclosure;

[0042] Figures 9A-C show optimization results for a fully constrained scenario with short eddy currents, long eddy currents and the combination of short and long eddy currents, in accordance with an embodiment of the present disclosure;

[0043] Figures 10A-D show optimization results for an unconstrained scenario when simulating both eddy currents but with different wkand seeks to show the effect of using a multiplicity of target images in order to avoid finding solutions that over fit to a specific target, in accordance with an embodiment of the present disclosure;

[0044] Figures 11A-F show NIS optimization results for an unconstrained steady-state scenario showing both optimization of RF flip angles as well as gradients, in accordance with an embodiment of the present disclosure;

[0045] Figures 12A-D show graphical representations of the reconstruction percentage error heat maps produced when the prediction of scanner performance and actual scanner performance are mismatched;

[0046] Figures 13A-B show graphs of the impact of enforcing zero gradients during the RF pulse period of each TR for unconstrained and fully constrained scenarios;

[0047] Figure 14 shows an example system of an MRI scanner that may be used in accordance with embodiments of the present disclosure; Figure 15 shows an example computer system block diagram that may be used for implementing the embodiment of the present disclosure;

[0048] Figure 16 shows a flow diagram of the framework in accordance with an embodiment of the present disclosure.

[0049] Detailed Description of the Preferred Embodiments

[0050] The present disclosure is a re-formulation of the MR sequence design problem; in which the approach focuses on the subsystem for generating magnetic field gradients used for spatial encoding, but it can be applied to all aspects of the MRI systems, NMR systems, magnetic resonance spectroscopy (MRS) systems and the associated expected end point of each is a globally optimized system that makes maximum use of all its available resources with each sub-system synergistically contributing to enhancing overall performance. Traditionally, designed gradient waveforms (constrained by specific scanner amplitude and slew rate limits) are input to a Bloch simulation from which image quality can be predicted. In most scanners, gradient systems are driven under voltage control, often with corrections as a hidden layer. Thus, the designer works with a specified gradient waveform but the gradient amplifier gives a different, hidden voltage waveform. To accommodate the hidden correction, the permitted gradient waveforms are often deliberately limited, reducing the apparent specification of the scanner.

[0051] Instead, the present disclosure directly formulates the sequence design problem in terms of voltages and / or currents and connects these to the actual realized gradients via a physical model. This model can contain any scanner imperfections that the designer wishes to account for; eddy currents are discussed in more detail as an example. The remaining workflow is unchanged compared to the classic optimization. The loss terms for the present disclosure can be chosen appropriately and tuned for the problem under consideration but we aim to maximize image reconstruction quality subject to actual scanner hardware properties and constraints.

[0052] Theory

[0053] Sequences are usually designed by considering a gradient waveform and respecting fixed amplitude and slew rate limits on the gradients themselves. In most scanners, gradient systems are driven under voltage control, often with corrections (e.g. pre-emphasis) as a hidden layer. Thus the designer works with a specified realized gradient waveform G(t) but gradient amplifier actually outputs a different demanded waveform g(t) . Instead we directly formulate the sequence design problem in terms of g(t) and connect this to G(t) via a physical model. The model could include many effects but here we focus on exponentially decaying eddy currents (Equation 1): anand Tnare amplitude and time constants for each eddy current term [3] and H(t) is a unit step function. Physical constraints apply to g(t) but G(t) defines the k-space trajectory for imaging, incorporated into forward NMR signal simulation in the optimization. The pipeline and a conventional method for comparison is shown in Figures 1, 2 and 3.

[0054] In Figure 1 an example of a classic optimization framework 1, in accordance with the prior art, is shown. The framework diagram shows the initial optimization which incorporates the Realized Gradients, G 10, the Bloch Simulation 11, the k-Space Data 12 and the Image 13. Post optimization, the framework then takes the Realized Gradients, G, 14 and feeds that through a Pre-emphasis 15 process before outputting the actual Demanded Gradients, g 16.

[0055] On the other hand, Figure 2 shows a Non-Idealized System (NIS) optimization framework 2 in accordance with embodiments of the present disclosure. As can be seen from a comparison of Figures 1 and 2; the NIS optimization framework 2 of the present disclosure does not perform a pre and post optimization procedure as was the case with the classic optimization framework 1. By working with the set of scanner control inputs directly and optimizing using imaging signal performance directly, scanner resources can be used as efficiently as possible. As such, the NIS optimization framework begins with the Demanded Gradients, g 20, and feeds this into the System Model 21. It is worth noting that the Demanded Gradients, g 20, would be considered as one of the system control inputs, x. From here the Realized Gradients, G, 22 are extracted. From the Realized Gradients, G 22, the output is processed in two ways. Firstly, the output is fed back to the Demanded Gradients, g, 20 and secondly, the Realized Gradients, G 22, are moved onto the Bloch Simulation 23. From the Bloch Simulation 23, the process continues to the k-Space Data 24, then the Image 25 and finally back to the Demanded Gradients, g 20. It is worth noting that the Realized Gradients, G 10, Bloch Simulation 11, k-Space Data 12, and Image Block 13 in Figures 1 perform the same function as the Realized Gradients, G 22, Bloch Simulation 23, k-Space Data 24, and Image Block 25 in Figure 2.

[0056] Following on from Figures 1 and 2, Figure 3A seeks to show Eddy current responses simulated in this work (time constants stated in legend) that represent the perturbation model used for forward MR signal simulation. Similarly, Figure 3B shows Perturbed (realized) waveforms (G_0) for both time constants listed and the unperturbed (demanded) equivalent waveforms (g_0).. Lxindicates the loss terms used during optimization in each case (defined in Equation 2 for NIS).

[0057] Figure 4 outlines an embodiment of the proposed architecture. Operation is initialized by defining what a desired "target" image or signal. This is flexibly achieved in the Specified sequence including sampling box in Figure 4, 40, and leads to an explicit generation of a target image or signal, 52. The target generation may be performed assuming idealized system performance or using a suitable approximation, a nominal system performance. Control inputs, as shown by the System control inputs box in Figure 4, are explicitly separated from realized fields, using a system response model. In the diagram, the control inputs x (including demanded gradient waveforms, g; RF waveforms, rf; and shim currents, ishim) are fed through the response model, which then produces predictions of the realized fields B (including gradient fields, G; RF fields, RF; static fields, Bo; and secondary factors such as temperature response which change performance): x System Response Model B

[0058] The realized fields B are then sufficient to perform a detailed simulation of the expected outcome, including the 'realised sampling' (S) which is a function of the realized gradients, and the NMR simulation requires knowledge of all field conditions. This can then be used to simulate the imaging data that would be achieved ('k-space data') which is reconstructed to give a predicted image or signal.

[0059] The system control inputs x are optimized by minimizing a loss function which includes the following terms:

[0060] (i) an image (signal) loss that minimizes reconstruction error with respect to the target image or signal

[0061] (ii) a sampling loss that drives towards a desired sampling scheme

[0062] (iii) a hardware loss that penalizes solutions which exceed hardware limits of the system

[0063] The x may be initialized by the target image or signal specification, which could be in terms of conventional sequence parameters or by some other means that specifies the desired outcome. To assist optimization, it may be convenient to specify x as a discretized list of time resolved system input values or some other means to specify x over the full duration of system operation to achieve the desired output image or signal. After optimization, the system inputs, x, which will in general be a complete specification of time resolved control inputs that explicitly account for system imperfections to achieve an output that is as close as possible to the target output subject to hardware limits. In order to avoid over-fitting to one specific image, the method can be optimized using a range of potential image targets simultaneously or can be 'augmented' by changing the target image periodically during the optimization, as discussed in relation to Figure 10.

[0064] Further, Figure 4 shows an outline framework 4 of the proposed system in accordance with embodiments of the present disclosure. Fundamentally the system control inputs are fed through a model of the entire system response, to give a prediction of the true fields (and other conditions such as temperature, which affects performance). The outline framework 4 is first initialized by providing a required definition for the desired output target image or signal, this is achieved by the Specified sequence including sampling block 40. From the Specified sequence including sampling block 40 the idealised control inputs or control inputs are defined and the control inputs are explicitly separated from the realized fields. The control inputs may be the Demanded Gradients, g, or demanded RF, rf, as shown in the System control inputs in block 42. It is worth noting that the Demanded Gradients, g 20, would be considered as one of the system control inputs, x. Further control inputs may be applicable such as shim current, iShim, or temperatures. Block 41 then feeds into the System model or MR scanner system response model 43 which provides predictions of the realized fields B (including gradient fields, G; RF fields, RF; static fields, Bo; and secondary factors such as temperature response which change performance). The realized fields or realized responses B 44 may comprise gradient fields, G, RF fields, RF, static fields, Bo, temperature responses etc. The realized fields or realized responses B 44 are then sufficient to perform a detailed simulation of the expected outcome, including the Realised Sampling (S) 46 which is a function of the realized gradients, and the NMR simulation 47 which requires knowledge of all field conditions. From this, the simulation of the imaging data can be completed i.e., via the k-space data 48. This will be reconstructed 49 to give a predicted image, I 50. Moreover, the current framework uses the Specified Sequence including sampling block 40, in tandem, to produce an NMR simulation 51 or an idealised NMR simulation 51 resulting in a Target image, i 52 or an idealised NMR Target image, i 52. The Target image, i 52 is formed using an initial set of idealised NMR scanner control inputs, Xi. The reconstructed Predicted Image, I 50 is then compared with the Target image, i 52 and the system Hardware limits 54 to decipher the Loss terms 53 within the system. The loss terms include, but are not limited to, images losses, sampling degradation losses, k-distance losses, slew rate term, parameter losses, hardware losses etc. The loss terms 53 are then fed back into the system control inputs 41 for example the Demanded Gradients, g, and the Demanded RF, rf, to improve the accuracy of the system by adapting the initial set of idealised NMR scanner control inputs, Xi. Thus, the NMR scanner control inputs become optimised in dependence on minimsation of the loss function that receives as inputs thereto at least the idealised target image and the predicted image and minimises errors therebetween.

[0065] The NMR scanner can then utilize the optimized set of scanner control inputs to generate future NMR image of subjects. Using the continually optimized set of scanner control inputs will improve the accuracy of the generated NMR images. It is worth noting that the optimised set of scanner control inputs have been optimised so as to be within corresponding system constraints of an NMR scanner to be operated using the optimised set of scanner control inputs. This is such that the NMR scanner itself does not act as a constraint or limitation to the accuracy of the generated images. There are many possible NMR scanner constraints such as maximum voltage(s), maximum current(s), maximum temperature(s), depletion of charge reservoir(s), available time to scan a subject, maximum rate of change of voltage(s), maximum rate of change of current(s), maximum power limit(s) etc. In order to ensure the system acts within the parameters of the MRI scanners constraints the MR scanner system response model 43 comprises one or more physical models representing system imperfections or system constraints in the NMR scanner to be controlled with the optimised set of scanner control inputs. Accordingly, the MR scanner system response model 43 may be based on prior characterisation of the system or on learned responses of resulting signals to changes in the control inputs, wherein the learned responses are in the form of parametrised models or learned networks. As such, the NMR scanner can be described as imperfect, which differs to the frameworks often described within the prior art. Rather the present disclosure seeks to achieve control of the imperfect NMR scanner system, taking into account and using its imperfections to maximise possible performance of the realistic MRI scanner capabilities.

[0066] Methods

[0067] Starting from a base sequence, with defined demanded waveforms g(t), the sequence is first simulated assuming perfect system performance (i.e. G(t) = g(t)) to determine a 'target' k-space kTand via Bloch simulation, image IT. Subsequently, a true forward model is used to predict the actual k-space kcand image ic; crucially this image is always reconstructed via iFFT assuming k-space locations kT. Optimization of g(t) then uses a composite loss function (2A) including: (2B) an image term that minimizes reconstruction error with respect to the target image or signal; (2C) a k-location loss term; (2D) an amplitude term; (2E) a slew rate term; and (2F) a term that minimizes RF gradients since slice-selection is not simulated in the forward simulation method (a deliberate decision made to reduce the number of optimization variables and waveform complexity). The system also seeks to optimize a sampling loss that drives towards a desired sampling scheme and a hardware loss that penalizes solutions which exceed hardware limits of the MR scanner system..

[0068] Total Loss = L[ + Lk+ Lg+ Ls+ LRF[2A]

[0069] IT is the idealised target image, Ic is the predicted image, kTis a target k-space, kc is a predicted actual k-space, and gmaxand smaxare voltage and slew rate limits respectively and w are tunable regularization weights and g(t) are defined initial gradient sequence waveforms.

[0070] Optimizations used the Adam [4] algorithm with an initial learning rate of 0.001 for all sequence parameters. Throughout, a spoiled gradient echo (GRE) sequence with flip angle (FA) = 5° and TR = 10ms was used as a starting point. A numerical brain phantom with feasible Bo, Bi, Ti, T2, T2', PD and isotropic diffusion coefficient values was used for signal simulation using PDGs. Time was discretized into steps At = 0.1ms, with the first 2ms of each TR reserved for RF pulses (no gradients) and an acquisition window of 3.2ms (0.5ms after the RF pulse). Images were simulated on a 32x32 grid, and all TR periods were simulated to allow for longer-term eddy currents whose effect persists over multiple periods. 32 TRs is insufficient to enter a steady-state, so a separate sequence comprising 52 dummy cycles and FA = 15° is also considered. 52 cycles are chosen as a compromise between achieving a reasonable steady-state and reducing computation time. ITin this scenario is generated using 300 dummies (i.e. a definite steady-state) and RF flip angles are additionally optimized such that lccan obtain the correct contrast.

[0071] Unconstrained (gmax= 50mT, Smax= 500mT / m / ms) optimizations are first run using a transient-state sequence and different eddy current properties: (i) mono-exponential perturbation with = 0.5ms (a short eddy current); (ii) mono-exponential perturbation = 50ms (a long eddy current); and (iii) bi-exponential perturbation with both short and long eddy currents (TT= 0.5ms and T2= 50ms). Results are analyzed by examining optimal demanded (gop) and realized (Gop) gradient waveforms, images (lc), k-space locations (kc) and loss term values. MR systems with lower performance hardware are represented by assigning different amplitude limits. These amplitude limited optimizations are run with gmax= {11, 15, 20, 25, 30}mT / m but with slew rate unconstrained. The optimal results are separately analyzed and compared to corresponding pre-emphasis solutions in terms of image reconstruction errors (NRMSE with respect to IT), effective partial Fourier factors and sequence timings.

[0072] NIS optimization is able to discover a solution that gives low reconstruction NRMSE even in scenarios where conventional pre-emphasis would fail. To demonstrate this, fully constrained optimizations are run for the three different sets of eddy current properties - (i), (ii) and (iii) defined above - but with demanded gradient amplitude limited to gmax= 12.5mT / m and its slew rate to Smax= 125mT / m / ms. gop, Gop, final lcand final kcare plotted.

[0073] The influence of Lkis explored by running uncontrolled optimizations with different wk. To improve the universality of the discovered solutions at lower wk, data augmentation techniques are incorporated in the form of random rotations of ITover a range of ±90°. 20 rotated versions of ITare generated and switched between every five iterations. Reconstruction NRMSE, optimal gradient waveforms (gop, Gop) and final kcare compared. Optimization of an unconstrained steady-state GRE sequence (described above) is run with variation of RF flip angles and g(t) during both dummy and acquisition periods. Results are analyzed as for the transient-state case.

[0074] The stability of NIS optimization solutions is investigated by introducing errors in eddy current parameters compared to those used during optimization. This is done by running the optimal demanded waveforms (gop) through the forward MR signal model once more but using different combinations of 'misestimated' a and T (each varied by ±25% from their respective ground-truth values). This analysis is performed separately for short and long eddy currents, and for the unconstrained and fully constrained scenarios. Heat maps of NRMSE between the reconstructed images in each misestimation case and lTare plotted. Lastly the repeatability of unconstrained and fully constrained optimization solutions and influence of LRFin each case, are explored by comparing Gop.

[0075] Note that unconstrained optimizations are run without restarts and to 20000 iterations to ensure convergence. Though restarts of the Adam algorithm hindered convergence for the unconstrained scenarios, they were added for the amplitude limited and fully constrained scenarios because a lower total loss became accessible. The restart strategy had 8 resets of the algorithm at iterations 500, 1500, 3000, 5000, 7500, 10500, 14000 and 18000. Results

[0076] Optimizations completed without gradient constraints reveal that the NIS optimization framework provides solutions that resemble pre-emphasis, especially during the readout gradients. Interestingly other gradient lobes are not pre-emphasized as this is found to be unnecessary to obtain a good image reconstruction. These solutions are preferential to standard pre-emphasis because they demand lower gradient amplitudes for the spoiler lobe, for instance. Time responses in Figure 5 correspond to scenarios when a short eddy current (T = 0.5ms), a long eddy current (T = 50ms) and their combination are simulated. Demanded (unperturbed) and realized (perturbed) waveforms are plotted separately for the different time constant combinations and can be compared to the target waveforms (subscript 0 in Figure 5A). Within each gradient lobe, a sharper response is apparent when T is shorter whilst the solution for a longer T simply increases the amplitude of each lobe versus the target waveform. Thus, Figure 5A is the target x-waveform and its perturbed counterparts for short, long and both eddy currents. Figure 5B, 5C and 5D present optimized demanded and equivalent realized waveforms for short, long and both eddy currents respectively.

[0077] Figure 6 considers the most scanner-realistic scenario (i.e. the combination of short and long eddy currents) and summarizes results from a NIS optimization without gradient constraints. Following 20000 iterations, the optimization has almost perfectly recovered the target k-space trajectory from the perturbed equivalent and gives negligible reconstruction error. Convergence of the weighted loss terms (w, = 1, wk= 20e-6 and wRF= le4) is plotted for the three terms invoked for this scenario. The optimized demanded waveforms show pre-emphasis features to mitigate short-term and long-term eddy current effects; the latter are countered with slight non-zero negative lobes after the spoiler and during the RF pulse that increase during the acquisition (i.e. with TR number). Thus, Figure 6A relates to the k-space locations for the target image, at iteration zero (perturbed) and after convergence (optimized) (from left to right). Figure 6B shows the target image on the left alongside difference images between this and the image reconstructions obtained at iteration zero and at convergence. Voxel values are normalized to the maximum target image value. Figure 6C shows optimized demanded and realized waveforms for the first five TRs. Finally, Figure 6D shows the loss terms values during optimization; Ls= Le= 0 for all unconstrained scenarios.

[0078] Further to the above, the target image (IT) is generated using the solid gray waveform shown in Figure 7A. The prewinder gradient has an amplitude of -25.0mT / m, the readout gradient is 7.34mT / m, and the spoiler is 17.6mT / m. By default, the target sequence has TE = 3.1ms and TR = 10ms. Pre-emphasis solutions usually require gradient amplitudes to increase with respect to the initial waveforms but these may be unobtainable for a scanner that has reduced hardware performance with tighter gradient constraints. This can be overcome either by: (i) increasing TE to ensure gradient lobe area is sufficient or (ii) using partial Fourier whilst keeping TE constant. These methods are outlined in Figure 7; as gradient amplitude limit decreases, the required TE to maintain pre-emphasis increases. NIS optimization inherently yields a solution that fits into any specified timing structure (TE = 3.1ms throughout) but has an associated increase in image NRMSE, whilst pre-emphasis theoretically gives a perfect reconstruction.

[0079] Figure 7D shows the partial Fourier factor required to maintain pre-emphasis whilst keeping TE constant. This is compared to the effective partial Fourier factor apparent when using NIS optimization. As the gradient amplitude limit decreases, the first few readout samples might be misplaced. The corresponding k-space lines that are misplaced by more than 0.6Akxwith respect to those in kTare identified and removed prior to image reconstruction. The image is then zero-filled to enable calculation of NRMSE with respect to IT. The number of removed samples is plotted as an "effective partial Fourier factor" and is consistently below that of pre-emphasis. Variable TE pre-emphasis solutions for {11, 15, 20}mT / m limited cases are also shown. Specifically, TE = 4.16ms (TR = 12.84ms), TE = 3.73ms (TR = 11.47ms) and TE = 3.44ms (TR = 10.51ms) for the three scenarios respectively. To summarise, Figure 7A shows the gradient changes necessary to achieve a true pre-emphasis solution at different gradient amplitude limits. This calculation has been made to an accuracy of lOps. Figure 7B shows the equivalent changes in TE required for a wider range of limits. Figure 7C shows that NIS optimization does not require a change in acquisition time but NRMSE does vary with the amplitude limit. Figure 7D shows the partial Fourier factors required to maintain pre-emphasis at a constant TE versus the effective factor generated by NIS optimizations.

[0080] Gradient constraints can make a pre-emphasis solution unobtainable within a reasonable acquisition time; image quality will deteriorate since eddy currents are not fully mitigated. Alternatively, NIS optimization consistently yields solutions that counter eddy currents within hardware limits and with low reconstruction error. Figure 8 shows optimized sequences obtained for the amplitude limited cases that Figure 7 introduces. As preemphasis is forbidden by the chosen amplitude limits, more novel waveforms are obtained with a large negative spoiler gradient between the readout gradient and RF pulse of the next repetition. For smaller limits, readout samples are discarded from the reconstructions as described above. Though unconventional, these solutions give low NRMSE with respect to ITafter the required zero-filling: 2.62%, 2.08% and 1.68% for 11, 15 and 20 mT / m- limited. Further, Figure 8 shows the NIS optimization results for amplitude limited cases (those in Figure 7A) when modelling both eddy currents. One sample is removed for the Figure 8A, i.e., the 20mT / m-limited case, three for Figure 8B i.e., 15mT / m-limited and four for Figure 8C i.e., llmT-limited. Difference maps between the target and optimized reconstructions are shown for each case on the left-hand side, where all images are normalized by the maximum value in the target. Optimized gradients and (truncated) k- space locations are alongside.

[0081] Figure 9 displays results of NIS optimizations when amplitude and slew rate constraints (12.5mT / m and 125mT / m / ms) are simultaneously active. Novel gradient waveforms are obtained for each scenario that mitigate eddy currents and maintain low reconstruction error (2.05%, 2.25% and 2.28% for short, long and both eddy currents respectively). Misplaced k-space locations were again systematically identified and respective readout samples removed from the reconstructions (three for short and long, and four for both); this explains the truncated k-space location plots in Figures 5 and 6 where sampling points do not fully extend to kx= -15 as is the case for kT. Reconstruction NRMSEs when not discarding the samples are 2.53%, 2.04% and 2.86% respectively. To summarise, Figure 9 shows the optimization results for a fully constrained scenario with Figure 9A being related to short eddy currents, Figure 9B being related to long eddy currents and Figure 9C being related to the combination of short and long eddy currents.

[0082] Equation 2 above, includes a k-space location term (Lk) that encourages sample points to lie on the regular Cartesian grid used to generate the target image. However, it is worth noting that many other sampling schemes could be used such as spiral, radial, variable density, 2D, 3D etc. Figure 10 investigates the importance of this loss term and shows that the inclusion of known machine learning techniques, such as data augmentation, during the optimization partially mitigates its requirement. Target image rotation reduces the potential for over-fitting to any one target as indicated by the smoother optimal gradient waveforms. Thus, Figure 10 shows the optimization results for an unconstrained scenario when simulating both eddy currents but with different wk: 2xl0'5for Figure 10B and 2xl0'7for Figures 10C and 10D, whilst data augmentation ('aug') is incorporated in Figure 10 D. For each case, optimal demanded and realized waveforms for three TRs and respective k-space rows are displayed. Image 1 (denoted by dots in Figure 10A) is the one used during non-augmented optimizations. Data augmentation produces lower NRMSE across all target image rotations.

[0083] Figure 11 shows NIS optimization results when simulating a GRE acquisition but with acquisition in a steady-state. Without constraints, optimal demanded gradient waveforms show similar pre-emphasis behaviour to those in a transient-state optimization (Figure 6). Since the number of dummy cycles used to generate the target image and those included in the optimization do not match, flip angle is additionally optimized to tune image contrast. Strong variation occurs during dummy TRs, but FA is unchanged from its initial value of 15° during the acquisition. Moreover, Figure 11 shows the NIS optimization results for an unconstrained steady-state scenario where both eddy currents are simulated. Figure 11A shows k-space locations for the target image, at iteration zero (perturbed) and after convergence (optimized). The normalized target image is shown in Figure 11B alongside difference images between this and the normalized initial and final reconstructions respectively. Optimal demanded waveforms and their realized equivalents are shown in Figures 11C and 11D for the first five dummy and acquisition period TRs. Figure HE shows the loss terms values during optimization. Figure HF shows the optimized flip angle scheme, color-coded to match Figures 11C and HD.

[0084] Realistically, the eddy current parameters modelled during the optimization process are unlikely to perfectly match those encountered during an experiment on a real-life scanner. For example, after measurement of these properties some deviation may occur due to scanner hardware updates. Figure 12 shows the impact of incorrectly estimating eddy current parameters (time constant r and amplitude a) on image reconstruction quality. Errors are generally less pronounced when misestimating the shorter eddy current parameters with this component being more sensitive to changes in T than a. The opposite is true for the longer eddy current term, which appears to be mostly insensitive to T misestimation. Unsurprisingly the fully constrained solution is more susceptible to inaccuracies in the MR signal model, but over a reasonable range of parameter errors, NIS solutions are robust. Therefore, Figure 12 shows the reconstruction percentage error heat maps produced when the prediction of scanner performance and actual scanner performance are mismatched. NIS optimization results for unconstrained (Figures 12A and 12B) and fully constrained (Figures 12C and 12D) scenarios using both eddy currents are taken as a basis. Separate maps are plotted for short and long eddy currents, with corresponding ground-truth time constants stated in each subplot title. Errors are calculated between the target image and reconstructions obtained when misestimating eddy current parameters according to the plotted x- and y-axis combinations.

[0085] Since slice or slab-selection is not simulated in the current MR signal simulation method, the loss term stated in Equation 2F is necessary. Figure 13 shows the impact of this term for unconstrained and fully constrained scenarios, with unperturbed demanded waveforms displayed. The unconstrained case is minorly impacted by this loss term; gradient waveforms are maximally different (4.7mT / m) for the first TR prewinder lobe. More different optimal waveforms are obtained in the fully constrained scenario as a new solution is found. However, this is not unexpected since re-runs of fully constrained NIS optimizations produce solutions with almost identical Iopand kopbut very different Gop, implying that indistinguishable local minima exist for this particular scenario. Thus, Figure 13 shows the impact of enforcing zero gradients during the RF pulse period of each TR for unconstrained in Figure 13A and fully constrained in Figure 13B. Realized waveforms are plotted for an optimization that uses LRF(non-zero wRF) and one that does not (zero wRF); waveforms from a repeat non-zero wRFoptimization are also shown. RF pulses occur during the darker coloured boxes.

[0086] Discussion

[0087] The present disclosure introduces a new sequence design framework that models and inherently accounts for MR scanner imperfections. We consider eddy current effects in this proof-of-concept study and show how our method can optimize gradient waveforms to yield a non-artefacted image even under severe hardware constraints, chosen here to be gradient amplitude and slew rate limits. Crucially, NIS optimization provides novel acquisition schemes that successfully counter eddy currents even in scenarios when traditional pre-emphasis would not be possible.

[0088] Our framework adapts to the time constant of the eddy current being modelled during optimization and shows some expected pre-emphasis behaviour when no gradient constraints are applied. The prewinder and readout gradients are somewhat preemphasized in Figure 5 but spoiler amplitude is not increased because its area is already sufficient to achieve spoiling and changes are not required to counter the simulated eddy currents. The time response of these solutions is however different for different T; gradient amplitudes are increased and decay away sharply within the duration of a single gradient lobe for the short scenario, whereas Gpongand G^'1® are simply offset from one another. To combat a bi-exponential eddy current perturbation, these strategies are combined (see Figures 5-6). In addition, a small negative lobe between the spoiler and subsequent RF increases amplitude per repetition to counter the long eddy current that spans multiple TRs. For this unconstrained case, convergence is quickly achieved with LRFplateauing before 3000 iterations and little image quality improvement seen beyond iteration 7500; further improvements are only driven by the loss term Lk(that aligns sampling points to the target Cartesian grid). The target image is almost perfectly recovered with residual errors of 0.17%.

[0089] However, given the proliferation of lower cost and thus more widely available MRI scanners, the higher gradient amplitudes associated to pre-emphasis may no longer be achievable. These systems can span from commercially available low-field scanners such as the Free. Max (Siemens Healthcare, Erlangen, Germany) to table-top scanners that inevitably have reduced gradient performance compared to clinical systems. In this instance, complete prewinding or spoiling could be unachievable unless gradient duration is increased to obtain the necessary gradient moment, as shown in Figures 7A and 7B. For example, a gradient system limited to 15mT / m would require TE of our chosen sequence to be increased from 3.1ms to 3.73ms and thus TR from 10ms to 11.47ms. For the low resolution 2D GRE acquisition used here, total acquisition time increases by a not insignificant 14.7%.

[0090] On the other hand, NIS optimization yields a solution that automatically fits within the original sequence timings (i.e. TE = 3.1ms and TR = 10ms is maintained) without significantly impacting image reconstruction quality. In particular, as gradient amplitude is more tightly constrained, NRMSE increases but still only reaches 2.9% with a stringent limit of llmT / m (Figure 7C). Reconstruction errors are negligible for limits beyond 20mT / m after which the effective partial Fourier factor for NIS optimization is zero (Figure 7D). Even before this point however, NIS optimization outperforms a fixed TE preemphasis implementation across a range of amplitude limits. Though the extent of k-space is reduced for the lowest amplitude limits, which slightly impacts image contrast, image NRMSEs < 3% are reported for each amplitude limited scenario. Some imperfections remain in Figure 8 but since these are located at the edges of k-space, their contribution is minimal. To avoid a forbidden pre-emphasis solution, NIS optimization discovers more innovative gradient waveforms that simultaneously achieve a low reconstruction error, satisfy hardware constraints and mitigate eddy currents. In particular, a large negative spoiler is discovered for each scenario and prewinder duration increases as their amplitude is more constrained.

[0091] Simultaneously constraining gradient amplitude and slew rate to values (12.5mT / m and 125mT / m / ms) that are significantly below the maximum values of the target sequence (25mT / m and 323mT / m / ms), ordinarily would mean that eddy current compensation is unachievable. However, NIS optimization is able to obtain a solution that removes any artefacts associated to the simulated eddy current terms regardless of time constant by discovering novel demanded waveforms. Following removal of mislocated readout samples, NRMSEs < 2.3% for each scenario in Figure 9. Small and more slowly ramping spoiler gradients are chosen for the T = 0.5ms case, and besides some gradient blips prior to each RF pulse, waveform shape is fairly regular and unchanged from one TR to the next. Larger variation across repetitions is seen when the T = 50ms eddy current is simulated since its impact extends across several TRs. Optimal waveforms for the "long" case feature flat-top lobes with irregular spoiler amplitudes. The solution for the "both" case is not a linear combination of its component solutions; few inter-TR patterns exist. Long-term eddy current compensation primarily involves modifying the spoiler lobe and subsequent gradient during the 2.3ms prior to the next RF pulse; prewinders and readouts are often unchanged.

[0092] To obtain a consistent NIS optimization solution, Lkis included with a default weighting of wk= 2e-5. This enforces sampling locations to lie on a regular Cartesian grid such that almost perfect alignment to those used to generate ITis obtained for an unconstrained scenario. However, it is worth noting that many other sampling schemes can be used in conjunction with the present invention; this could be include sampling schemes such as spiral, radial, variable density, 2D, 3D etc. NRMSEs for this scenario span 0.1-0.2% for the range of rotated target images in Figure 10A. However, since this requires prior knowledge of kT, making an optimization independent of Lkwould make our framework more adaptable. Reducing Lkproduces an optimal solution with more chaotic k-space locations, whereby readout amplitudes vary more rapidly. This is demonstrated in Figure 10C, where wkhas been reduced by two orders of magnitude compared to in Figure 10B; prewinder and spoiler portions of the optimal demanded waveforms are unchanged. Despite this, a low reconstruction NRMSE is obtained because the framework is over-fitting to the specific ITused for optimization, indicated by the pink dot. To reduce over-fitting, data augmentation can be added into the optimization; a simple rotation of ITevery five iterations has been used in Figure 10D, though this method could be developed to include other transformations.

[0093] Without changing wkdata augmentation suppresses fluctuations of the readout gradients and produces more regular k-space locations (especially for outer kylines) such that optimal waveforms are more similar to those in the "high wk" case. Though reconstruction NRMSE for this scenario is poorer for the first target image (that used for "low wk" optimization) compared to the "low wk" case, NRMSE across all other images in Figure 10A (rotations used for "low wk+ aug" optimization) are reduced. This solution will be more universal and produce a reasonable reconstruction for any of the tested brain phantom orientations, whereas the "low wk" result will only produce a non-artefacted image for the one orientation the optimization was trained on. Consequently, these results suggest data augmentation reduces both over-fitting and the dependence of NIS optimizations on the Lkterm.

[0094] Other sequences and their parameters (beyond demanded gradient waveforms) can be tuned using NIS optimization. To demonstrate this, our GRE sequence is driven to a steady-state using a simple dummy cycle approach to obtain the expected Ti-weighting. Optimization of an unconstrained scenario yields pre-emphasis trends again (Figure 11D), as expected. Dummy cycle gradients (Figure 11C) are not modified during optimization, though the corresponding RF pulse flip angles are varied to achieve the correct image contrast (Figure 11F). Here the number of dummy cycles used to generate IT(300) differs from the number used during optimization (52); the latter was chosen as a trade-off between entering a sufficient steady-state and the number of optimization variables. The observed flip angle variation is reminiscent of other literature strategies that can drive sequences to a steady-state. However, the acquisition flip angle scheme does not stray from the starting value of 15°. Despite significantly more degrees-of-freedom versus equivalent transient-state optimizations (and thus computational time), good convergence for all loss terms is achieved in Figure HE (long before the optimization was terminated).

[0095] On this topic, computational expense was a limiting factor for all optimizations. A low resolution image (matrix size of 32x32) is chosen to reduce simulation time; an optimization of 3200 gradient timepoints for each axis (6400 total; corresponding to the transient-state scenarios considered here) takes 160 minutes on a 20(40)xIntel (R) Xeon (R) Silver 4210 2.20GHz CPU, 251 GB RAM, 32GB NVIDIA Tesla V100 GPU (NVIDIA, Santa Clara, CA, USA). The steady-state optimization takes approximately 7s per iteration versus 0.5s for the transient-state scenarios.

[0096] Given that real-life MRI systems exhibit instabilities, it is important to characterize the consistency of the NIS optimization solutions when the eddy current parameters are mismatched from those used during optimization. Figures 12A-B demonstrate that the discovered pre-emphasis-like solution for the unconstrained scenario performs well for a range of (<Z, T) inaccuracies. NRMSEs remain below 3% even when the short-term eddy current amplitude is misestimated by ±25%. Greater dependence is seen for r on shortterm eddy current NRMSE, with a maximum 6.8% error obtained for simultaneous -25% a and +25% T misestimation. NRMSEs are generally larger for the longer-term eddy current, and so this should be prioritized when establishing the system model to be used during NIS optimization, a dominates in Figure 12B with NRMSE mostly invariant of r. Similar trends are seen for the fully constrained scenarios but errors are generally higher. Short-term eddy currents are less significant again with maximum NRMSEs of 6% for 10% misestimations of either parameter in Figure 12C. Long-term eddy current NRMSEs can reach the same magnitude as the eddy current parameter misestimation amount, highlighting a need for accurate characterization when using more challenging hardware configurations. The complexities of fully constrained optimizations is further demonstrated in Figure 13. Quite different solutions are obtained depending on whether wRFis active or not, suggesting that the NIS optimization seeks to utilize the RF period for eddy current compensation. Figure 13A shows that re-runs of an unconstrained scenario yield almost identical results, suggesting a globally optimal solution is found. However, re-runs of fully constrained scenarios (Figure 13B) give different results during the spoiler; prewinder / readout regions are mostly unchanged. This indicates that the post-readout period in each TR is the most susceptible to changes during the optimization, and that some local minima exist when gradient constraints are active. Each solution is equally valid though since corresponding k-space locations and reconstruction NRMSEs are indistinguishable. Despite this search-space degeneracy, NIS optimization consistently provides artefact-free image reconstructions. Future work will integrate control of RF and gradient systems more completely with more accurate simulation methods, obviating the need for some elements of the loss function such as the RF loss term (equation 2F), which actually only operates solely on the gradient system.

[0097] In a further embodiment of the present disclosure, the framework may involve 'a-priori' optimization for a range of desired outcomes (output images or signals) given known system performance, which can itself be characterized by performing detailed off-line measurements of system response. There is no requirement for the system response to have a particular form (e.g. linear system), only that it should behave predictably. The system model may be pre-trained based on extensive prior data or be continuously adapting to include data from the latest system operations. All prior data acquisitions contain information about relationships between control inputs and actual signal outputs.

[0098] Figure 14 shows an example an MRI scanner 14, the patient 142 is located in a whole body MRI system 141 equipped with N radiofrequency transmit antennas or coils and gradient coils 143. The portion of the device external to the patient is monitored by M dedicated sensors 145. The sensors are placed at several different positions on the external part of the device. They measure the amplitude and phase of any radiofrequency currents induced on the interventional device by the transmit field. These measurements are relayed to the computer 147. The computer 147 controls the overall operation of the MRI system such as the control of the transmitter control system 144 which subsequently control the respective radiofrequency coils and gradient coils 143.

[0099] Figure 15 is a block diagram of a typical general purpose computer system 147, as seen in Figure 14, the computer system 147 can form the processing platform for framework, as described, according to the present disclosure. The computer system 12 comprises a central processing unit (CPU), random access memory (RAM), and input / output ports (I / O) into which data can be received and output therefrom as is well known in the art. Additionally included is the radio frequency coils and gradient coils required for the function of the MRI scanner, a visual display unit and a network access point for sending or receiving data including the receiving parameter settings, sensor data, historic data, training data etc. Although, many further appliances could be used if necessary.

[0100] The computer system 147 also includes some non-volatile storage 150, such as a hard disk drive, solid-state drive, or NVMe drive. Stored on the non-volatile storage 150 is a number of executable computer programs together with data and data structures required for their operation or training. Overall control of the system 147 is undertaken through the execution of the Control Program 1521, the System Model or MR Scanner System Response Model Program 1522, the Forward Model 1523 and the Loss Function model 1524. The other data contained in the non-volatile storage 150 is the required training data 1525.

[0101] The computer-implemented system may applicable with a plethora of devices i.e., devices that are capable of processing data and controlling external systems. The devices should be capable of at least receiving, processing, analysing, and storing of data that is extracted from sensors relating to the MRI scanner. Therefore, devices that may be applicable for use in this system are not limited to, but at least include:

[0102] • desktop computers

[0103] • laptop computers

[0104] • mobile telephones

[0105] • tablet computers

[0106] Figure 16 shows a flow diagram 16 of the framework in accordance with the present disclosure. In particular the flow diagram relates to the generating and adapting of the set of NMR scanner control inputs, x, based on the outcomes of the loss function. The flow diagram starts at S160. This step will include any required initialisation of the computer system i.e., the processor, associated memory etc. In sl61, an idealised target image is generated, via a first NMR simulation, assuming substantially idealised NMR scanner performance, and using an initial set of NMR scanner control inputs (Xi). Then, in sl62 the initial set of control inputs (Xi) are applied to an MR scanner system response model to obtain a set of predicted realised responses (B) representative of MR scanner operating parameters. Next, in sl63, a predicted image using the set of predicted realised responses, B (as obtained in sl62), are generated via a second NMR simulation. At this point, the flowchart goes in two directions. Firstly, the flowchart feedback via S165 which adapts the initial set of NMR scanner control inputs (Xi) in sl61 to give an optimised set of scanner control inputs (x0) in dependence on minimisation of a loss function that receives as inputs thereto at least the idealised target image and the predicted image and minimises errors therebetween. On the other hand, S163 feeds into sl64 where the NMR scanner is operated with the set of scanner control inputs to generate an NMR image of a subject in the NME scanner.

[0107] In a further embodiment of the present disclosure, the framework may use a degree of 'online optimization' to allow the user to change parameters of the image sequence (for example, resolution and field of view, or contrast properties). Fast optimization methods would then be employed, taking predefined solutions as starting points.

[0108] It is worth noting that optimization methods are not constrained to be classical only - any method for optimizing the control inputs (for example using neural networks to approximate the end-to-end performance) can be deployed.

[0109] It will be understood that the above list is non-exhaustive, and that the method and system described herein is applicable to many technical problem domains to which machine learning models may be applied.

[0110] Various modifications, whether by addition, substitution, or deletion will be apparent to the intended reader to provide further embodiments of the present disclosure, any and all of which are intended to be encompassed by the appended claims.

[0111] References

[0112] 1. Loktyushin A, Herz K, Dang N, et al. MRzero - Automated discovery of MRI sequences using supervised learning. Magn. Reson. Med. 2021;86:709-724 doi: 10.1002 / mrm.28727.

[0113] 2. Endres J, Dang N, Glang F, Loktyushin A, Weinmuller S, Zaiss M. Phase distribution graphs for differentiable and efficient simulations of arbitrary MRI sequences. In: Proc. Inti. Soc. Mag. Reson. Med. 30. ; 2022.

[0114] 3. Bernstein MA, King KF, Xiaohong JZ. Handbook of MRI Pulse Sequences. Academic Press; 2004. doi : https: / / doi.org / 10.1016 / B978-0-12-092861-3.X5000-6.

[0115] 4. Kingma DP, Ba JL. Adam: a method for stochastic optimization. arXiv 2014; 1412.6980: 1-15.

Claims

Claims1. A magnetic resonance (MR) imaging method, comprising: generating, via a first MR simulation, an idealised target image assuming substantially idealised MR scanner performance, and using an initial set of MR scanner control inputs (Xi); applying the initial set of MR scanner control inputs to an MR scanner system response model to obtain a set of predicted realised responses (B) representative of MR scanner operating parameters; generating, via a second MR simulation, a predicted image using the set of predicted realised responses (B); and adapting the initial set of MR scanner control inputs (Xi) to give an optimised set of scanner control inputs (x0) in dependence on minimisation of a loss function that receives as inputs thereto at least the idealised target image and the predicted image and minimises errors there between.

2. A method according to claim 1, and further comprising: operating an MR scanner with the optimised set of scanner control inputs (x0) to generate an MR image of a subject.

3. A method according to any of the preceding claims, wherein the optimised set of scanner control inputs (x0) have been optimised so as to be within corresponding system constraints of an MR scanner to be operated using the optimised set of scanner control inputs.

4. A method according to claim 3, wherein the corresponding system constraints of the MR scanner comprise one or more selected from the group comprising: i) one or more maximum voltages; ii) one or more maximum currents; iii) one or more maximum temperatures; iv) depletion of one or more charge reservoirs; v) one or more available times to scan a subject; vi) one or more maximum rates of change of voltage; vii) one or more maximum rates of change of current; and / or viii) one or more maximum power limits.

5. A method according to any of the preceding claims, wherein the initial set of MR scanner control inputs (Xi) and / or the optimised set of scanner control inputs (x0) comprise one or more selected from the group comprising: i) one or more demanded gradient waveforms (g(t)); ii) one or more radio frequency (RF) waveforms (rf); iii) one or more shim currents (Ishim); and iv) one or more patient handling system controls.

6. A method according to any of the preceding claims, wherein the predicted realised responses (B) comprise one or more selected from the group comprising: i) gradient fields, G; ii) RF fields, RF; iii) static fields, Bo; and iv) temperature response.

7. A method according to any of the preceding claims, wherein the MR scanner system response model comprises one or more physical models representing system imperfections or system constraints in the MR scanner to be controlled with the optimised set of scanner control inputs (x0).

8. A method according to claim 7, wherein the MR scanner system response model is: i) based on prior characterisation of the system, or ii) based on learned responses of resulting signals to changes in the control inputs, wherein the learned responses are in the form of parametrised models or learned networks.

9. A method according to any of the preceding claims, wherein the MR scanner is recognised to be an imperfect MR scanner system, wherein the MR scanner system response model is used to achieve control of the imperfect MR scanner system, taking into account and using its imperfections to maximise possible performance.

10. A method according to any of the preceding claims, wherein after optimization the optimised set of scanner control inputs (x0) comprise a set of time resolved control inputs that explicitly account for system imperfections to achieve an output that is as close as possible to the target output subject to hardware limits of the MR scanner.

11. A method according to any of the preceding claims, wherein the loss function is a composite loss function comprising: (i) an image loss term; (ii) a k-distance loss term; (iii) a voltage term; and (iv) a slew rate term.

12. A method according to claim 11, wherein the loss function comprisesTotal Losswherein:and wherein Vmaxand Smaxare voltage and slew rate limits respectively, w are tunable regularization weights, IT is the idealised target image, Ic is the predicted image, kTis a target k-space, kc is a predicted actual k-space, and gmaxand smaxare voltage and slew rate limits respectively and w are tunable regularization weights and g(t) are defined initial gradient sequence waveforms.

13. A method according to any of the preceding claims, wherein the loss function includes at least one or more of the following:(i) an image (signal) loss that minimizes reconstruction error with respect to the target image or signal(ii) a sampling loss that drives towards a desired sampling scheme(iii) a hardware loss that penalizes solutions which exceed hardware limits of the MR scanner system.

14. A method according to any of the preceding claims, wherein by working with the set of scanner control inputs directly and optimizing using imaging signal performance directly, scanner resources can be used as efficiently as possible.

15. A method according to any of the preceding claims, wherein sample points are applied using a sampling scheme such that sample points can be shifted or removed to optimise reconstructions of the MR image of the subject.

16. A system, comprising: a processor; one or more computer readable medium, storing instructions, which when executed by the processor, cause the processor to perform the method of any of the preceding claims.

17. A MR scanner having a set of scanner control inputs, the scanner control inputs being set to correspond to the values of the optimised set of scanner control inputs generated by the method of any of claims 1 to 15 or the system of claims 16.

Images