Defining new sequences

Examples of sequence definitions can be found here. These are the actual sequence definitions used for the pre-implemented sequences in MRIBuilder. In this example, we will look at the file for DiffusionSpinEcho as an example.

We can see that there are several steps. First we define the sequence as a specific sub-type of Sequence:

const DiffusionSpinEcho = Sequence{:DiffusionSpinEcho}

Note the duplication of the sequence name (DiffusionSpinEcho). This name will have to be unique.

The next step is to define the sequence constructor:

function DiffusionSpinEcho(; scanner=DefaultScanner, parameters..., vars...)
    build_sequence(scanner) do
        (g1, g2) = dwi_gradients(...)
        seq = Sequence([
            :excitation => excitation_pulse(...),
            nothing,
            :gradient => g1,
            nothing,
            :refocus => refocus_pulse(...),
            nothing,
            g2,
            nothing,
            :readout => readout_event(...),
            nothing,
        ], name=:DiffusionSpinEcho, vars...)
        add_cost_function!(seq[6].duration + seq[7].duration)
        return seq
    end
end

The crucial part here are the individual parts used to build the sequence, defined as a vector:

[
    :excitation => excitation_pulse(...),
    nothing,
    :gradient => g1,
    nothing,
    :refocus => refocus_pulse(...),
    nothing,
    :gradient2 => g2,
    nothing,
    :readout => readout_event(...),
    nothing
]

We can see that this sequence is built in order by an excitation pulse, some unknown dead time (indicated by nothing), a gradient, more dead time, a refocus pulse, more dead time, another gradient, more dead time, and finally a readout.

Some of these components have been given specific names (e.g., :excitation => ...). This is optional, but can be useful to refer to individual components. There are helper functions available to create these components.

After creating the sequence object, we can now add secondary objectives to the cost function (using add_cost_function!). In this example, we have:

add_cost_function!(seq[6].duration + seq[7].duration)

If we check the order of the sequence component, we see that this minimises the sum of the duration of the second gradient and the wait block before this gradient. This cost function has been added to maximise the time between the second gradient and the readout (and hence minimise the effect of eddy currents on the readout). Note that this is a secondary cost function that will only take effect if it does not interfere with any user-defined constraints and cost functions (see sequence optimisation). Some secondary cost functions will be automatically defined for you within the individual components (e.g., a trapezoidal gradient has a secondary cost function to maximise the slew rate). There is even a tertiary cost function, which minimises the total sequence duration.

The next step is to define summary variables that the user can constrain when setting up a specific instance of this sequence:

@defvar begin
    diffusion_time(ge::DiffusionSpinEcho) = start_time(ge, :gradient2) - start_time(ge, :gradient)
    echo_time(ge::DiffusionSpinEcho) = 2 * (variables.effective_time(ge, :refocus) - variables.effective_time(ge, :excitation))
end

For this sequence, we define the variables.diffusion_time as the time between the start of the first and second gradient pulse, and the variables.echo_time as twice the time between the refocus and excitation pulses. These variables need to be defined within a @defvar block.

In addition to these sequence-specific summary variables, there are also a lot of variables already pre-defined on individual components, such as the variables.slice_thickness of the RF pulse or the variables.gradient_strength of the gradient pulses. To access these summary variables on a sequence-level, we need to tell MRIBuilder for which RF pulses/gradients we are interested in computing these variables:

get_pulse(seq::DiffusionSpinEcho) = (excitation=seq[:excitation], refocus=seq[:refocus])
get_gradient(seq::DiffusionSpinEcho) = (gradient=seq[:gradient], gradient2=seq[:gradient2])
get_readout(seq::DiffusionSpinEcho) = seq.readout

Note that we can indicate we are interested in multiple RF pulses/gradients by supplying them as a named tuple ((excitation=..., refocus=...)).

Setting this allows us to get RF pulse or gradient-specific properties by calling the sequence, for example:

using MRIBuilder
sequence = DiffusionSpinEcho(bval=1., TE=:min, slice_thickness=2.)
variables.flip_angle(sequence)
(excitation = 89.99999999999999, refocus = 180.0000000000496)

Here we can see that we get the variables.flip_angle for each of the two RF pulses defined using get_pulse above.

The final component to defining summary variables is to define one or more default coherence pathways using get_pathway:

get_pathway(seq::DiffusionSpinEcho) = Pathway(seq, [90, 180], 1, group=:diffusion)

Here the coherence Pathway sets out what a specific set of spins might experience during the sequence. In this case the sequence experiences two RF pulses and is excited by the first pulse and flipped by the second ([90, 180]). It is then observed during the first readout (1). For such a Pathway we can compute:

By defining a default pathway for the sequence, the user can now put constraints on any or all of these variables.

Component helper functions

MRIBuilder.Parts.HelperFunctions.dwi_gradientsMethod
dwi_gradients(; type, optimise=false, refocus=true, orientation=[1, 0, 0], group=:diffusion, variables...)

Returns two diffusion-weighting gradients that are guaranteed to cancel each other out.

Parameters

  • type: A symbol describing the type of gradient. One of:
    • :trapezoid: Pulsed trapezoidal gradient profile. See Trapezoid for the relevant variables.
    • :instant: instantaneous gradient (e.g., to test short-pulse approximations). See InstantGradient for the relevant variables.
  • optimise: Whether to optimise this readout event in isolation from the rest of the sequence. Use this with caution. It can speed up the optimisation (and for very complicated sequences make it more robust), however the resulting parameters might not represent the optimal solution of any external constraints (which are ignored if the readout is optimised in isolation).
  • scanner: Used for testing. Do not set this parameter at this level (instead set it for the total sequence using build_sequence).
source
MRIBuilder.Parts.HelperFunctions.excitation_pulseMethod
excitation_pulse(; parameters...)

Create an excitation RF pulse.

By default there is no slice-selective gradient. To enable slice selection min_slice_thickness has to be set to a number or to :min. If min_slice_thickness is not set or is set to :min, then either bandwidth or duration should be set, otherwise the optimisation might be unconstrained (ignore this for shape=:instant).

Parameters

For an InstantPulse (i.e., shape=:instant), only the flip_angle, phase, and group will be used. All other parameters are ignored.

  • optimise: set to true to optimise this RF pulse separately from the embedding sequence.

Pulse parameters

  • shape: The shape of the RF pulse. One of :sinc (for SincPulse), :constant/:hard (for ConstantPulse), or :instant (for InstantPulse). Default is :sinc for slice-selective pulses or :instant otherwise.
  • flip_angle: size of the flip due to the RF pulse in degrees (default: 90).
  • phase: angle of the RF pulse in the x-y plane in degrees (default: 0).
  • frequency: frequency of the RF pulse relative to the Larmor frequency in kHz (default: 0).
  • bandwidth: width of the RF pulse in Fourier space in kHz (default: free variable).
  • duration: duration of the RF pulse in ms (default: free variable).
  • Nzeros: sets the number of zero crossings for a SincPulse (default: 3). Can be set to a tuple of two numbers to set a different number of zero crossings before and after the pulse maximum.
  • group: name of the group of which the RF pulse is part. This is used to add transformations after the sequence is optimised.

Slice selection

  • slice_thickness: minimum slice thickness that should be possible without adjusting the sequence timings in um (not mm!) (default: no slice selection). Can be set to :min to indicate that this should be minimised given the scanner constraints and user values for bandwidth or duration.
  • rephase: set to false to disable the spin rephasing after the RF pulse.
  • rotate_grad: name of the parameter with which the slice selection gradient will be rotated after sequence optimisation (default: :FOV).
  • scanner: overrides the global_scanner for this part of the sequence. Recommended to set only if not part of a larger sequence.
source
MRIBuilder.Parts.HelperFunctions.gradient_spoilerMethod
gradient_spoiler(; optimise=false, orientation=[0, 0, 1], rotate=:FOV, scale=:spoiler, spoiler=1., duration=:min, variables...)

Returns two DWI gradients that are guaranteed to cancel each other out.

Parameters

  • orientation: Orientation of the gradient (default: slice-select direction).
  • rotate: in which coordinate system is the orientation defined (default: :FOV).
  • scale: variable controlling how the gradients should be scaled after optimsation (default: :spoiler).
  • optimise: Whether to optimise this readout event in isolation from the rest of the sequence. Use this with caution. It can speed up the optimisation (and for very complicated sequences make it more robust), however the resulting parameters might not represent the optimal solution of any external constraints (which are ignored if the readout is optimised in isolation).
  • scanner: Used for testing. Do not set this parameter at this level (instead set it for the total sequence using build_sequence).

Variables

  • variables.spoiler: maximum spoiler scale (before applying any reductions due to scale).
  • Any other parameters expected by Trapezoid.
source
MRIBuilder.Parts.HelperFunctions.interpret_image_sizeMethod
interpret_image_size(fov, resolution, voxel_size, slice_thickness)

Combines the user-provided information of the image size to produce a tuple with:

  • slice_thickness: if not set explicitly, will be set to the third element of voxel_size
  • resolution_z: number of voxels in the slice-select direction.
  • keywords parameters for the readout_event functions with the fov, resolution, and voxel_size in the x- and y- dimensions.
source
MRIBuilder.Parts.HelperFunctions.readout_eventMethod
readout_event(; type, optimise=false, variables...)

Adds a readout event to the sequence.

Parameters

  • type: A symbol describing the type of readout. It will default to :epi if a resolution has been set and :instant otherwise. Can be set to one of the following:
    • :epi: EPI readout. See EPIReadout for the relevant variables.
    • :instant: single isolated readout event SingleReadout (e.g., for NMR). Does not expect any variables.
  • optimise: Whether to optimise this readout event in isolation from the rest of the sequence. Use this with caution. It can speed up the optimisation (and for very complicated sequences make it more robust), however the resulting parameters might not represent the optimal solution of any external constraints (which are ignored if the readout is optimised in isolation).
  • scanner: Used for testing. Do not set this parameter at this level (instead set it for the total sequence using build_sequence).
source
MRIBuilder.Parts.HelperFunctions.refocus_pulseMethod
refocus_pulse(; parameters...)

Create an excitation RF pulse.

By default there is no slice-selective gradient. To enable slice selection slice_thickness has to be set to a number or to :min. If slice_thickness is not set or is set to :min, then either bandwidth or duration should be set, otherwise the optimisation might be unconstrained (ignore this for shape=:instant).

Parameters

  • optimise: set to true to optimise this RF pulse separately from the embedding sequence.

Pulse parameters

For an InstantPulse (i.e., shape=:instant), only the flip_angle, phase, and group will be used. All other parameters are ignored.

  • shape: The shape of the RF pulse. One of :sinc (for SincPulse), :constant/:hard (for ConstantPulse), or :instant (for InstantPulse). Default is :sinc for slice-selective pulses or :instant otherwise.
  • flip_angle: size of the flip due to the RF pulse in degrees (default: 180).
  • phase: angle of the RF pulse in the x-y plane in degrees (default: 0).
  • frequency: frequency of the RF pulse relative to the Larmor frequency in kHz (default: 0).
  • bandwidth: width of the RF pulse in Fourier space in kHz (default: free variable).
  • duration: duration of the RF pulse in ms (default: free variable).
  • Nzeros: sets the number of zero crossings for a SincPulse (default: 3). Can be set to a tuple of two numbers to set a different number of zero crossings before and after the pulse maximum.
  • group: name of the group of which the RF pulse is part. This is used to add transformations after the sequence is optimised.

Slice selection and spoilers

  • slice_thickness: minimum slice thickness that should be possible without adjusting the sequence timings in um (not mm!) (default: no slice selection). Can be set to :min to indicate that this should be minimised given the scanner constraints and user values for bandwidth or duration.
  • spoiler: set to the spatial scale on which the spins should be dephased in mm. For rotating spoilers, this does include the contribution from the slice select gradient as well.
  • rotate_grad: name of the parameter with which the slice selection and spoiler gradient will be rotated after sequence optimisation (default: :FOV).
  • scanner: overrides the global_scanner for this part of the sequence. Recommended to set only if not part of a larger sequence.
  • orientation: vector with orientation of slice select gradient and pulses (defaults: z-direction).
source
MRIBuilder.Parts.HelperFunctions.saturation_pulseMethod
saturation_pulse(; pulse=(), binomial_order=1, variables...)

Creates a saturation pulse consisting of repeating instances of pulse_type.

Parameters

  • pulse: Parameters/variables used to define the basic RF pulse used in the saturation pulse. These can be:
  • binomial_order: How many repeats of pulse_type there should be in a single binomial pulse.

Variables

  • nrepeat: number of repeats of the binomial pulse. This should typically be set by the user.
  • interblock_delay: time between the start of each pulse. Should be short compared with the T1, but long enough to limit SAR.
  • duration: total duration of the saturation pulse.
source

Optimisation helper functions

MRIBuilder.BuildSequences.build_sequenceMethod

Wrapper to build a sequence.

Use as

build_sequence(scanner;) do
    ...
end

Within the code block you can create one or more sequences, e.g.

seq = Sequence(
    SincPulse(flip_angle=90, phase=0, duration=2., bandwidth=:max)
    nothing.,
    SingleReadout
)

You can also add any arbitrary constraints or objectives using one of:

  • set_simple_constraints!
  • apply_simple_constraint!
  • add_cost_function!

As soon as the code block ends the sequence is optimised (if optimise=true) and returned.

Parameters

  • scanner: Set to a Scanner to limit the gradient strength and slew rate. When this call to build_sequence is embedded in another, this parameter can be set to nothing to indicate that the same scanner should be used.
  • optimise: Whether to optimise and fix the sequence as soon as it is returned. This defaults to true if a scanner is provided and false if no scanner is provided.
  • n_attempts: How many times to restart the optimiser (default: 20)? Decrease if you want to quickly check feasibility of a simple sequence. Increase if convergence fails for a complex sequence.
  • kwargs...: Other keywords are passed on as attributes to the optimiser_constructor (e.g., set print_level=3 to make the Ipopt optimiser quieter).
source

Pathways API

MRIBuilder.Pathways.PathwayType
Pathway(sequence::Sequence, pulse_effects::Vector{:Symbol/Number}, readout_index=1; group=nothing)

Describes how a specific spin/isochromat might experience the sequence.

Only a single pathway through the RF pulses is considered, so that at every point in time the spins are in one of the following four states:

  • +longitudinal: initial relaxed state
  • +transverse: excited state. During this time gradients will affect the variables.area_under_curve (or variables.qval) and variables.bval.
  • -longitudinal: inverse state
  • -transverse: inverse excited state. During this time all gradients will have the inverse effect compared with +transverse.

The RF pulses cause mappings between these different states as described below.

Parameters

  • sequence: MRI Sequence to be considered.
  • pulse_effects: How each RF pulse affects the spins. This can be one of the following:
    • :skip/:ignore/0: This RF pulse leaves the spins unaffected.
    • :refocus/:invert/180: Flips the sign of the spin state (i.e., +longitudinal <-> -longitudinal, +transverse <-> -transverse)
    • :excite/90: Takes spin state one step along the following sequence +longitudinal -> +transverse -> -longitudinal -> -transverse -> +longitudinal
    • :neg_excite/270/-90: Inverse step compared with :excite.
  • readout_index: After encountering the number of pulses as defined in pulse_effects, continue the Pathway until the readout given by index is reached. If set to 0 the Pathway is terminated immediately after the last RF pulse.
  • group: which gradient grouping to consider for the net_dephasing and bmat. If not set, all gradients will be considered (using their current alignment).

Attributes

Over the pathway the following values are computed. Each can be accessed by calling the appropriate function:

Timings

  • variables.duration_state: The total amount of time spent in a specific state in this pathway (+longitudinal, +transverse, -longitudinal, or -transverse)
  • variables.duration_transverse: The total amount of time the spins spent in the transverse plane in ms. This can be used to quantify the expected effect of T2-decay.
  • variables.duration_dephase: The total amount of time the spins spent in the +transverse relative to -transverse state in ms. The absolute value of this can be used to quantify the expected effect of T2'-decay.

Effect of gradients

The area under curve, q-values, and b-values are computed separately for each group of gradients (depending on the group keyword set during construction).

source
MRIBuilder.Pathways.interpret_pulse_effectsMethod
interpret_pulse_effects(number_or_symbol)

Interpret the various numbers and symbols that can be passed on to a Pathway.

The result will be one of:

  • :ignore (if input is 0, :ignore, or :skip).
  • :excite (if input is 90 or :excite).
  • :refocus (if input is 180, :refocus, or :excite).
  • :neg_excite (if input is -90, 270, or :negexcite).
source
MRIBuilder.Pathways.update_gradient_tracker_till_time!Method
update_gradient_tracker_till_time!(walker::PathwayWalker, key, new_time)
update_gradient_tracker_till_time!(tracker::GradientTracker, new_time)

Update the bmat for any time passed since the last update (assuming there will no gradients during that period).

The bmat is updated with the outer produce of qvec with itself multiplied by the time since the last update.

When called with the first signature the tracker will be created from scratch if a tracker with that key does not exist.

source
MRIBuilder.Pathways.update_walker_gradient!Method
update_walker_gradient!(gradient_block::GradientWaveform, walker::PathwayWalker, gradient_start_time::VariableType; overlapping_pulses=[], overlapping_readouts=[])

Update the walker's qvec and bmat based on the given gradient_block.

The following steps will be taken:

  • Do nothing if walker.is_transverse is false
  • increase the appropriate walker.bmat by the outer product of walker.qvec with itself multiplied by the time since the last gradient
  • update the appropriate walker.qvec and walker.bmat based on the gradient waveform. This will require appropriate qvec/bmat functions to be defined for the gradient building block.
  • update walker.last_gradient_time to the time at the end of the gradient.

This requires variables.bmat_gradient and variables.qvec to be implemented for the GradientWaveform.

source
MRIBuilder.Pathways.update_walker_pulse!Method
update_walker_pulse!(walker::PathwayWalker, pulse_effects::Vector, pulse_time)

Apply the first element of pulse_effects to the walker at the given pulse_time.

The following steps will be taken if the first pulse_effect is not :ignore

  • if walker.is_transverse is true before the pulse, increase the walker.bmat by the outer product of walker.qvec with itself multiplied by the time since the last gradient
  • update walker.duration_states with time since last pulse.
  • update walker.last_pulse_time
  • update walker.is_transverse, and walker.is_positive based on the first value in pulse_effects.
  • if walker.is_positive changed in the previous step than the walker.qvec needs to be flipped.
  • remove the first element from pulse_effects.

If the first element is :ignore the only effect is that the first element is removed from pulse_effects.

source
MRIBuilder.Pathways.update_walker_till_time!Function
update_walker_till_time!(walker::PathwayWalker, new_time[, gradient_group])

Updates all parts of a PathwayWalker up to the given time.

This updates the walker.duration_states and the bmat for each gradient tracker. If gradient_group are provided, then only the gradient tracker matching that group will be updated. If that gradient tracker does not exist, it will be created.

This function is used to get the walker up to date till the start of a gradient, pulse, or final readout.

source
MRIBuilder.Pathways.walk_pathway!Method
walk_pathway!(bb::Sequence/BuildingBlock, walker::PathwayWalker, pulse_effects::Vector, nreadout::Ref{Int}, start_time)

Computes the effect of a specific ContainerBlock (starting at start_time) on the PathwayWalker.

For individual pulses and gradients, the following behaviour is implemented:

For overlapping gradients/pulses, one should first encounter the part of the gradient before the variables.effective_time of the pulse, then apply the pulse, and then the rest of the gradient. This effective time can be passed on to update_walker_gradient! to allow part of the gradient waveform to be applied.

The function should return true if the Pathway has reached its end (i.e., the final readout) and false otherwise.

source
MRIBuilder.Variables.variables.area_under_curveFunction
area_under_curve(pathway::Pathway)

Return net displacement in k-space (i.e., spoiling) experienced by the spins following a specific Pathway.

Only gradients active while the spins are in the transverse plane are considered.

Returns a NamedTuple with the area_under_curve for all gradient groups.

source
MRIBuilder.Variables.variables.bmatFunction
bmat(pathway::Pathway)

Return 3x3 diffusion-weighted matrix experienced by the spins following a specific Pathway in rad^2 ms/um^2.

Only gradients active while the spins are in the transverse plane are considered.

Returns a NamedTuple with the bmat for all gradient groups.

source
MRIBuilder.Variables.variables.bvalFunction
bval(pathway::Pathway)

Return size of diffusion-weighting experienced by the spins following a specific Pathway in rad^2 ms/um^2.

Only gradients active while the spins are in the transverse plane will contribute to the diffusion weighting.

Returns a NamedTuple with the bval for all gradient groups.

source
MRIBuilder.Variables.variables.duration_stateFunction
duration_state(pathway::Pathway, transverse::Bool, positive::Bool)

Returns how long the Pathway spent in a specific state.

The requested state can be set using transverse and positive as follows:

  • transverse=false, positive=true: +longitudinal
  • transverse=true, positive=true: +transverse
  • transverse=false, positive=false: -longitudinal
  • transverse=true, positive=false: -transverse
source
MRIBuilder.Variables.variables.net_dephasingFunction
net_dephasing(pathway::Pathway)

Return net displacement vector in k-space/q-space experienced by the spins following a specific Pathway in kHz/um.

Only gradients active while the spins are in the transverse plane are considered.

Returns a NamedTuple with the qvec for all gradient groups.

source