Compressed-sensing reconstruction

Description

This example describes how to perform a compressed-sensingreconstruction of a CS-2 accelerated acquisition.

Loading Package

using LazyArtifacts # loading data
using SEQ_BRUKER_a_MP2RAGE_CS_360
using CairoMakie # plotting

In addition we load the package internally used to perform the reconstruction

using SEQ_BRUKER_a_MP2RAGE_CS_360.MRIReco
using SEQ_BRUKER_a_MP2RAGE_CS_360.MRIReco.RegularizedLeastSquares

Download the datasets

if you run the literate example offline change the following line by : `MP2artifacts = artifact"MP2RAGEdata"

datadir = Main.MP2_artifacts
@info "The test data is located at $datadir."
[ Info: The test data is located at /home/runner/.julia/artifacts/04cd4c29bb9e2aeb5384fbc70a9af0e1a37ca369.

If you want to perform your own reconstruction, you can change the following line in order to point to another a bruker dataset

path_bruker = joinpath(datadir, "MP2RAGE_CS2")
"/home/runner/.julia/artifacts/04cd4c29bb9e2aeb5384fbc70a9af0e1a37ca369/MP2RAGE_CS2"

Compressed-sensing reconstruction

In order to use an advanced reconstruction we will pass some parameters that will be used by the reconstruction package MRIReco.jl

We have to create a parameter dictionnary that will be used. If you need more information about it take a look at MRIReco.jl

CS = Dict{Symbol,Any}()
CS[:sparseTrafo] = "Wavelet" #sparse trafo
CS[:reg] = L1Regularization(100.)       # regularization
CS[:solver] = FISTA    # solver
CS[:iterations] = 30

d = reconstruction_MP2RAGE(path_bruker; mean_NR=true,paramsCS = CS)
Dict{Any, Any} with 7 entries:
  "LUT"            => Float32[93.0 0.5; 94.0 0.5; … ; 5287.0 -0.5; 5288.0 -0.5]
  "params_reco"    => Dict{Symbol, Any}(:iterations=>30, :solver=>FISTA, :reg=>…
  "params_MP2RAGE" => ParamsMP2RAGE(800.0, 2284.1, 7.0, 5000.0, 128, 7.0, 7.0)
  "im_reco"        => ComplexF32[20.6719+0.527998im 10.2103+14.8792im … 68.1409…
  "T1map"          => Float32[2686.0 2301.0 … 1311.0 1511.0; 2327.0 1187.0 … 10…
  "params_prot"    => BrukerFile: /home/runner/.julia/artifacts/04cd4c29bb9e2ae…
  "MP2RAGE"        => Float32[-0.374886 -0.295793 … 0.0963395 -0.00548892; -0.3…

for comparison purposes let's perform the undersampled reconstruction (without the paramCS keyword)

d_under = reconstruction_MP2RAGE(path_bruker; mean_NR=true)
Dict{Any, Any} with 7 entries:
  "LUT"            => Float32[93.0 0.5; 94.0 0.5; … ; 5287.0 -0.5; 5288.0 -0.5]
  "params_reco"    => Dict{Any, Any}(:iterations=>1, :reconSize=>(128, 128, 96)…
  "params_MP2RAGE" => ParamsMP2RAGE(800.0, 2284.1, 7.0, 5000.0, 128, 7.0, 7.0)
  "im_reco"        => ComplexF32[61.5519+126.668im 14.7966+115.13im … 53.1283-5…
  "T1map"          => Float32[2463.0 2481.0 … 1006.0 1777.0; 1924.0 1701.0 … 12…
  "params_prot"    => BrukerFile: /home/runner/.julia/artifacts/04cd4c29bb9e2ae…
  "MP2RAGE"        => Float32[-0.333003 -0.336819 … 0.260251 -0.124936; -0.1811…

We can check the results

begin
  f = Figure(size=(500,400))
  ax=Axis(f[1,1],title="TI₁ undersampled")
  h=heatmap!(ax,abs.(d_under["im_reco"][:,:,60,1,1]),colormap=:grays)

  ax=Axis(f[1,2],title="TI₁ CS")
  h=heatmap!(ax,abs.(d["im_reco"][:,:,60,1,1]),colormap=:grays)


  ax=Axis(f[2,1],title="T₁ map undersampled")
  h=heatmap!(ax,d_under["T1map"][:,:,60,1],colorrange = (500,2000))

  ax=Axis(f[2,2],title="T₁ map CS")
  h=heatmap!(ax,d["T1map"][:,:,60,1],colorrange = (500,2000))

  for ax in f.content   # hide decoration befor adding colorbar
    hidedecorations!(ax)
  end

  Colorbar(f[2,3],h,label = "T₁ [ms]", flip_vertical_label=true)
  f
end
Example block output

Write results in BIDS format

Results can be written following most of the qBIDS format recommandation

subject_name = "sub_01_cs"
dir_path = "" # directory path where the files will be create
write_bids_MP2RAGE(d,subject_name,dir_path)

This page was generated using Literate.jl.