Improving quantification in non-TOF 3D PET/MR by incorporating photon energy information

Hybrid PET/MR systems combine functional information obtained from positron emission tomography (PET) and anatomical information from magnetic resonance (MR) imaging. In spite of the advantages that such systems can offer, PET attenuation correction still represents one of the biggest challenges for imaging in the thorax. This is due to the fact that the MR signal is not directly correlated to gamma-photon attenuation. In current practice, pre-defined population-based attenuation values are used. However, this approach is prone to errors in tissues such as the lung where a variability of attenuation values can be found both within and between patients. A way to overcome this limitation is to exploit the fact that stand-alone PET emission data contain information on both the distribution of the radiotracer and photon attenuation. However, attempts to estimate the attenuation map from emission data only have shown limited success unless time-of-flight PET data is available. Several groups have investigated the possibility of using scattered data as an additional source of information to overcome re- construction ambiguities. This thesis presents work to extend the previous methods by using PET emission data acquired at multiple energy windows and incorporating prior information derived from MR. This thesis is organised as follows. We first cover both the literature and mathematical theory relevant to the framework. Then, we present and discuss results on the case of attenu- ation estimation from scattered data only, when the activity distribution is known. We then give an overview of several candidates for joint reconstruction, which reconstruct both the activity and attenuation from scattered and unscattered data. We present extensive results using simulated data and compare the proposed methods to state-of-the-art MLAA from a single energy window acquisition. We conclude with suggestions for future work to bring the proposed method into clinical practice

An Investigation of Stochastic Variance Reduction Algorithms for Relative Difference Penalised 3D PET Image Reconstruction

Penalised PET image reconstruction algorithms are often accelerated during early iterations with the use of subsets. However, these methods may exhibit limit cycle behaviour at later iterations due to variations between subsets. Desirable converged images can be achieved for a subclass of these algorithms via the implementation of a relaxed step size sequence, but the heuristic selection of parameters will impact the quality of the image sequence and algorithm convergence rates. In this work, we demonstrate the adaption and application of a class of stochastic variance reduction gradient algorithms for PET image reconstruction using the relative difference penalty and numerically compare convergence performance to BSREM. The two investigated algorithms are: SAGA and SVRG. These algorithms require the retention in memory of recently computed subset gradients, which are utilised in subsequent updates. We present several numerical studies based on Monte Carlo simulated data and a patient data set for fully 3D PET acquisitions. The impact of the number of subsets, different preconditioners and step size methods on the convergence of regions of interest values within the reconstructed images is explored. We observe that when using constant preconditioning, SAGA and SVRG demonstrate reduced variations in voxel values between subsequent updates and are less reliant on step size hyper-parameter selection than BSREM reconstructions. Furthermore, SAGA and SVRG can converge significantly faster to the penalised maximum likelihood solution than BSREM, particularly in low count data

An Investigation of Stochastic Variance Reduction Algorithms for Relative Difference Penalized 3D PET Image Reconstruction

Penalised PET image reconstruction algorithms are often accelerated during early iterations with the use of subsets. However, these methods may exhibit limit cycle behaviour at later iterations due to variations between subsets. Desirable converged images can be achieved for a subclass of these algorithms via the implementation of a relaxed step size sequence, but the heuristic selection of parameters will impact the quality of the image sequence and algorithm convergence rates. In this work, we demonstrate the adaption and application of a class of stochastic variance reduction gradient algorithms for PET image reconstruction using the relative difference penalty and numerically compare convergence performance to BSREM. The two investigated algorithms are: SAGA and SVRG. These algorithms require the retention in memory of recently computed subset gradients, which are utilised in subsequent updates. We present several numerical studies based on Monte Carlo simulated data and a patient data set for fully 3D PET acquisitions. The impact of the number of subsets, different preconditioners and step size methods on the convergence of regions of interest values within the reconstructed images is explored. We observe that when using constant preconditioning, SAGA and SVRG demonstrate reduced variations in voxel values between subsequent updates and are less reliant on step size hyper-parameter selection than BSREM reconstructions. Furthermore, SAGA and SVRG can converge significantly faster to the penalised maximum likelihood solution than BSREM, particularly in low count data

PET/MRI attenuation estimation in the lung:A review of past, present, and potential techniques

Positron emission tomography/magnetic resonance imaging (PET/MRI) potentially offers several advantages over positron emission tomography/computed tomography (PET/CT), for example, no CT radiation dose and soft tissue images from MR acquired at the same time as the PET. However, obtaining accurate linear attenuation correction (LAC) factors for the lung remains difficult in PET/MRI. LACs depend on electron density and in the lung, these vary significantly both within an individual and from person to person. Current commercial practice is to use a single-valued population-based lung LAC, and better estimation is needed to improve quantification. Given the under-appreciation of lung attenuation estimation as an issue, the inaccuracy of PET quantification due to the use of single-valued lung LACs, the unique challenges of lung estimation, and the emerging status of PET/MRI scanners in lung disease, a review is timely. This paper highlights past and present methods, categorising them into segmentation, atlas/mapping, and emission-based schemes. Potential strategies for future developments are also presented

STIR Software for Tomographic Image Reconstruction

<h1>Summary of changes in STIR release 5.2.0</h1> <p>This version is 100% backwards compatible with STIR 5.0 as far as usage goes. However, there are changes in the output of scatter estimation and ECAT8 normalisation, see below for more information.</p> <h2>Overall summary</h2> <p>Of course, there is also the usual code-cleanup and improvements to the documentation. See also <a href="https://github.com/UCL/STIR/milestone/6">the 5.2 milestone on GitHub</a>.</p> <p>Overall code management and assistance by Kris Thielemans (UCL and ASC). Other main contributors were Daniel Deidda (NPL) and Markus Jehl (Positrigo).</p> <h2>Patch release info</h2> <ul> <li>5.2.0 released 30/10/2023</li> </ul> <h2>Summary for end users (also to be read by developers)</h2> <h3>Bug fixes</h3> <ul> <li>Scatter estimation was setting initial activity image to 1 at set-up, effectively ignoring the initial image, aside from geometric info.</li> <li>Setting SPECTUB resolution model with STIR python or SIRF divided slope by 10 in error. The problem did not occur when set using parameter file</li> </ul> <h3>Changed functionality</h3> <ul> <li>The ECAT8 normalisation (used for the Siemens mMR) code now takes the 4th component <em>axial effects</em> into account. These normalisation factors are therefore different (even up to ~10%). This gives improved axial uniformity in the images. The use of the axial effects can be switched off by adding setting <code>use_axial_effects_factors:=0</code> to the parameter file (see an example in <code>examples/Siemens-mMR/correct_projdata_no_axial_effects.par</code>), or the class member of the same name. In addition, the Siemens normalisation header is now read (using a new class <code>InterfileNormHeaderSiemens</code>) such that hard-coded variables for the Siemens mMR have been removed. Further testing of this functionality is still required however. <a href="https://github.com/UCL/STIR/pull/1182/">PR #1182</a>.</li> <li>Interfile header parsing now correctly identifies keywords that contain a colon by checking for <code>:=</code>.</li> <li>The <code>set_up()</code> method of the ray-tracing projection matrix now skips further processing if it was already called with data of the same characteristics. This will means that any cached data will be re-used, potentially leading to a speed-up when re-using it from Python. <a href="https://github.com/UCL/STIR/pull/1281/">PR #1281</a>.</li> </ul> <h3>New functionality</h3> <ul> <li><p>The <code>Discretised Shape3D</code> shape/ROI has now an extra value <code>label index</code>. For ROIs, this allows using a single volume with multiple ROIs encoded as labels, such as output by ITKSnap and many others. When used as a shape in <code>generate_image</code>, it could be used to extract a single ROI from such a label image. <a href="https://github.com/UCL/STIR/pull/1196/">PR #1196</a>.</p> </li> <li><p>Global variables in SPECTUB have been substituted by class members, such that multiple SPECTUB projectors can be used. <a href="https://github.com/UCL/STIR/pull/1169/">PR #1169</a>.</p> </li> <li><p>Global variables in PinholeSPECTUB have been substituted by class members, such that multiple PinholeSPECTUB projectors can be used. <a href="https://github.com/UCL/STIR/pull/1212/">PR #1212</a>.</p> </li> <li><p>Scatter estimation is now smoothed in axial direction for BlocksOnCylindrical scanners. <a href="https://github.com/UCL/STIR/pull/1172/">PR #1172</a>.</p> </li> <li><p><code>InverseSSRB</code> now works for BlocksOnCylindrical after a rewrite. <a href="https://github.com/UCL/STIR/pull/1172/">PR #1172</a>. /</p> </li> <li><p>Parallelised function <code>set_fan_data_add_gaps_help</code> across segments to reduce computation time. <a href="https://github.com/UCL/STIR/pull/1168/">PR #1168</a>.</p> </li> <li><p>New utility <code>SPECT_dicom_to_interfile</code> which reads a DICOM file with SPECT projection data and extracts the data and writes one or more Interfile 3.3 headers (still somewhat preliminary). <a href="https://github.com/UCL/STIR/pull/1182/">PR #1182</a>.</p> </li> <li><p>The new <code>stir_timings</code> utility is mostly useful for developers, but you could use it to optimise the number of OpenMP threads to use for your data. <a href="https://github.com/UCL/STIR/pull/1237/">PR #1237</a>.</p> </li> <li><p>New classes <code>SegmentIndices</code>, <code>ViewgramIndices</code> and <code>SinogramIndices</code>, used by <code>ProjData</code> related classes, as opposed to having to specify all the elements directly, e.g. in C++</p> <pre><code> auto sinogram = proj_data.get_sinogram(sinogram_indices);</code></pre> <p>This makes these functions more future proof, in particular for TOF. The older functions are now deprecated. Note that as <code>Bin</code> is now derived from <code>ViewgramIndices</code>, instantations of <code>Bin</code> can now be used to specify the indices as well in most places. There is still more work to do here, mostly related to the symmetries. <a href="https://github.com/UCL/STIR/pull/1273/">PR #1273</a>.</p> </li> </ul> <h4>Python (and MATLAB)</h4> <ul> <li>Examples use <code>stir.ProjData.read_from_file</code> as opposed to <code>stir.ProjData_read_from_file</code>. The former is supported since SWIG 3.0, and the <a href="https://swig.org/Doc4.1/Python.html#Python_nn20">default from SWIG 4.1</a>.</li> <li>Addition of <code>DetectionPosition</code> and <code>DetectionPositionPair</code>.</li> <li><code>bin.time_frame_num</code> is now no longer a function in Python, but acts like a variable (as the other <code>Bin</code> members).</li> <li>Addition of <code>RadionuclideDB</code></li> </ul> <h3>New examples</h3> <ul> <li><code>examples/python/construct_projdata_demo.py</code> illustrates constructing a <code>ProjDataInMemory</code></li> </ul> <h3>Changed functionality</h3> <ul> <li>Scatter estimation was resetting the activity image to 1 before each iteration. This led to cases where the reconstructed image (and therefore the scatter estimation) did not converge, especially when using a small number of sub-iterations. Now, the reconstructed image is continuouslu updated between scatter iterations by default. This should also allow users to use less sub-iterations, therefore saving some time for the scatter estimation. The old behaviour can be re-enabled by setting <code>restart_reconstruction_every_scatter_iteration</code> to true either via a parameter file or via the <code>set_restart_reconstruction_every_scatter_iteration()</code> function. <a href="https://github.com/UCL/STIR/pull/1160/">PR #1160</a>.</li> <li>energy resolution functions and keywords have now more documentation. <code>Scanner::check_consistency</code> also checks if the energy resolution is less than 20 (as it is FWHM/reference_energy). <a href="https://github.com/UCL/STIR/pull/1149/">PR #1149</a>.</li> <li>Errors now throw <code>std::runtime_error</code> instead of <code>std::string</code>. <a href="https://github.com/UCL/STIR/pull/1131/">PR #1131</a>.</li> <li>The parameter <code>use_view_offset</code> was removed from the <code>interpolate_projdata</code> functions. View-offset is now always taken into account. <a href="https://github.com/UCL/STIR/pull/1172/">PR #1172</a>.</li> <li>The info, warning and error calls are thread safe now (which makes them slower), and the logging output in <code>distributable.cxx</code> was changed from verbosity 2 (which is the STIR default) to verbosity 3. This is to reduce the default output during iterative reconstructions. <a href="https://github.com/UCL/STIR/pull/1243/">PR #1243</a>.</li> <li>The <code>Succeeded</code> class has a new method <code>bool succeeded()</code> enabling more concise code (avoiding the need for comparing with <code>Succeeded::yes</code> which is especially verbose in Python).</li> <li>The example files for the Siemens mMR now use lower min/max thresholds for the (single) scatter scale. This gives better results, see <a href="https://github.com/UCL/STIR/issues/1163/">Issue #1163</a>. <a href="https://github.com/UCL/STIR/pull/1279/">PR #1279</a>.</li> </ul> <h3>Deprecated functionality and upcoming changes to required tool versions</h3> <ul> <li>The following functions (previously used for upsampling the scatter estimate) have been made obsolete or replaced, and will be removed in STIR version 6.0.0: <code>interpolate_axial_position</code>, <code>extend_sinogram_in_views</code> and <code>extend_segment_in_views</code></li> <li>Constructors/functions in <code>ProjData</code> related classes that explicitly use <code>axial_pos_num</code>, <code>view_num</code> etc in their arguments are now deprecated, and should be replaced by their respective versions that use <code>SegmentIndices</code>, <code>ViewgramIndices</code> or <code>SinogramIndices</code>. The former will not be compatible with TOF information that will be introduced in version 6.0.0.</li> <li>Use of the AVW library to read Analyze files will be removed in 6.0, as this has not been checked in more than 15 years. Use ITK instead.</li> <li>GE VOLPET and IE support will be removed in 6.0, as we have no files to test this, and it's obsolete anyway.</li> <li>STIR version 6.0.0 will require C++ 14 (currently we require C++ 11, but already support C++ 20) and CMake 3.14.</li> </ul> <h3>Build system and dependencies</h3> <ul> <li>CMake 3.12 is now required on Windows.</li> <li>We now use CMake's <a href="https://gitlab.kitware.com/cmake/community/-/wikis/doc/tutorials/Object-Library">OBJECT library feature</a> for the registries. This avoids re-compilation of the registries for every executable and therefore speeds-up building time. Use of STIR in an external project is not affected as long as the recommended practice was followed. This is now documented in the User's Guide. <a href="https://github.com/UCL/STIR/pull/1141/">PR #1141</a>.</li> <li>The <code>error</code> and <code>warning</code> functions are now no longer included from <code>common.h</code> and need to be included manually when used (as was already the case for <code>#include "stir/info.h"</code>). <a href="https://github.com/UCL/STIR/pull/1192/">PR #1192</a>.</li> <li>add .h and .i as dependencies for SWIG generated wrappers to make sure they get rebuild. (Currently adding all .h files, which is too much, but CMake needs a fix before we can do this properly). <a href="https://github.com/UCL/STIR/pull/1218/">PR #1218</a>.</li> </ul> <h3>Changes for developers</h3> <ul> <li>moved all functionality in <code>CListEventCylindricalScannerWithDiscreteDetectors</code> to template class <code>CListEventScannerWithDiscreteDetectors</code> (templated in <code>ProjDataInfoT</code>). This enables re-use for generic/blocksoncylindrical scanners. <a href="https://github.com/UCL/STIR/pull/1222/">PR #1222</a>.</li> <li>rewritten <code>ProjDataInMemory</code> to avoid streams, causing a speed-up of some operations, and removing a limit of total size of 2GB. <a href="https://github.com/UCL/STIR/pull/1260/">PR #1260</a>.</li> </ul> <h3>Known problems</h3> <ul> <li>See <a href="https://github.com/UCL/STIR/labels/bug">our issue tracker</a>.</li> </ul> <h3>Minor (?) bug fixes</h3> <ul> <li>Small change in scatter simulation to how non-arccorrected bins are computed. Added a check in the construction of non-arccorrected projdata that the number of tangential bins is not larger than the maximum non-arccorrected bins. <a href="https://github.com/UCL/STIR/pull/1152/">PR #1152</a>.</li> <li><code>extend_segment_in_views</code> does not handle view offsets correctly and does not work for BlocksOnCylindrical scanners <a href="https://github.com/UCL/STIR/issues/1177/">issue #1177</a>. A new function <code>extend_segment</code> was added that works for Cylindrical and BlocksOnCylindrical and allows extension in tangential and axial direction as well. <a href="https://github.com/UCL/STIR/pull/1172/">PR #1172</a>.</li> <li><code>sample_function_on_regular_grid</code> did not handle offset correctly in all places <a href="https://github.com/UCL/STIR/issues/1178/">issue #1178</a>. <a href="https://github.com/UCL/STIR/pull/1172/">PR #1172</a>.</li> <li>Ray tracing projection for BlocksOnCylindrical scanner geometries contained a bug where some bins were swapped across oblique segments <a href="https://github.com/UCL/STIR/issues/1223/">issue #1223</a>. This sometimes lead to large artifacts in reconstructions. <a href="https://github.com/UCL/STIR/pull/1231/">PR #1231</a>.</li> </ul> <h3>Documentation changes</h3> <ul> <li>Updated the STIR developers guide, which was quite out-of-date w.r.t. C++ features etc.</li> </ul> <h3>recon_test_pack changes</h3> <ul> <li>Updated headers of most images and projection data to avoid warnings.</li> </ul> <h3>Other changes to tests</h3> <ul> <li><code>test_Scanner.cxx</code> tests for energy resolution, <a href="https://github.com/UCL/STIR/pull/1149/">PR #1149</a>.</li> <li>New file <code>test_interpolate_projdata</code>, <a href="https://github.com/UCL/STIR/pull/1141/">PR #1141</a>.</li> </ul>If you use this software, please cite it using the metadata from this file