Optimized MLAA for quantitative non-TOF PET/MR of the brain

For quantitative tracer distribution in positron emission tomography, attenuation correction is essential. In a hybrid PET/CT system the CT images serve as a basis for generation of the attenuation map, but in PET/MR, the MR images do not have a similarly simple relationship with the attenuation map. Hence attenuation correction in PET/MR systems is more challenging. Typically either of two MR sequences are used: the Dixon or the ultra-short time echo (UTE) techniques. However these sequences have some well-known limitations. In this study, a reconstruction technique based on a modified and optimized non-TOF MLAA is proposed for PET/MR brain imaging. The idea is to tune the parameters of the MLTR applying some information from an attenuation image computed from the UTE sequences and a T1w MR image. In this MLTR algorithm, an αj parameter is introduced and optimized in order to drive the algorithm to a final attenuation map most consistent with the emission data. Because the non-TOF MLAA is used, a technique to reduce the cross-talk effect is proposed. In this study, the proposed algorithm is compared to the common reconstruction methods such as OSEM using a CT attenuation map, considered as the reference, and OSEM using the Dixon and UTE attenuation maps. To show the robustness and the reproducibility of the proposed algorithm, a set of 204 [18F]FDG patients, 35 [11C]PiB patients and 1 [18F]FET patient are used. The results show that by choosing an optimized value of αj in MLTR, the proposed algorithm improves the results compared to the standard MR-based attenuation correction methods (i.e. OSEM using the Dixon or the UTE attenuation maps), and the cross-talk and the scale problem are limited.

For quantitative tracer distribution in positron emission tomography, attenuation correction is essential. In a hybrid PET/CT system the CT images serve as a basis for generation of the attenuation map, but in PET/MR, the MR images do not have a similarly simple relationship with the attenuation map. Hence attenuation correction in PET/MR systems is more challenging. Typically either of two MR sequences are used: the Dixon or the ultra-short time echo (UTE) techniques. However these sequences have some well-known limitations. In this study, a reconstruction technique based on a modified and optimized non-TOF MLAA is proposed for PET/MR brain imaging. The idea is to tune the parameters of the MLTR applying some information from an attenuation image computed from the UTE sequences and a T1w MR image. In this MLTR algorithm, an j α parameter is introduced and optimized in order to drive the algorithm to a final attenuation map most consistent with the emission data. Because the non-TOF MLAA is used, a technique to reduce the cross-talk effect is proposed. In this study, the proposed algorithm is compared to the common reconstruction methods such as OSEM using a CT attenuation map, considered as the reference, and OSEM using the Dixon and UTE attenuation maps. To show the robustness and the reproducibility of the proposed algorithm, a set of 204 [ 18 F]FDG patients, 35 [ 11 C]PiB patients and 1 [ 18 F]FET patient are used. The results show that by choosing an optimized value of j α in MLTR, the proposed algorithm improves the results compared to the standard MR-based attenuation correction methods (i.e. OSEM using the Dixon or the UTE attenuation maps), and the cross-talk and the scale problem are limited.
Keywords: PET, MRI, attenuation correction, joint activity and attenuation estimation, maximum likelihood (Some figures may appear in colour only in the online journal)

Introduction
For positron emission tomography (PET) or single photon emission tomography (SPECT), attenuation correction is an essential step in order to obtain a quantitative image of the tracer activity distribution. In PET, the attenuation sinogram can be obtained either directly with a radioactive transmission source or based on a computed tomography (CT) scan in hybrid PET/ CT. In the latter case, an additional step is required to extrapolate the CT attenuation image to the PET photon energy (511 keV) through a simple bilinear function (Kinahan et al 1998). In both cases, the attenuation sinogram is computed by forward projecting the attenuation image. For PET/MR imaging in general, the attenuation sinogram cannot be computed directly since the MR sequences provide no information about electron densities, leading to a lack of signal in the bone and air regions. Thus the MRI-based attenuation correction (MRAC) remains a challenging problem for a quantitative result.
Current MRAC algorithms can be classified into three categories: 1 -Template-and atlas-based methods, for which the general purpose is to compute a pseudo-CT image from a MR image of the patient, using a database of MR images and/ or a database of CT images (Kops and Herzog 2007, Hofmann et al 2008, Burgos et al 2014a, Johansson et al 2014. 2 -Segmentation-based methods, for which the goal is to derive an attenuation map from different MR sequences of the patient. Different methods have been proposed, for example, in Schlemmer et al (2008) only air and tissue are segmented and used for attenuation correction. For standard PET/MR systems, the common attenuation map is derived by the Dixon technique (Martinez-Möller et al 2009), that is dedicated for whole body and brain imaging, or the UTE technique (Keereman et al 2010) dedicated only for brain imaging. Recently a new promising technique named RESOLUTE (Ladefoged et al 2015) has been proposed based on the UTE technique. 3 -Emission-based methods, which combine the estimation of the emission map and the attenuation map in a simultaneous reconstruction algorithm framework. Censor et al (1979) were the first to propose a simultaneous calculation of the attenuation and activity coefficients. Manglos and Young (1993), Krol et al (1995Krol et al ( , 2001 and Nuyts et al (1999) improved this technique by using a maximum likehood approach. This algorithm is now commonly named MLAA (maximum likelihood reconstruction of attenuation and activity). It is often reported that MLAA in a non-TOF system generates some cross-talk in the reconstructed images between attenuation and emission data, and also the final tracer activity reconstruction is a scaled version of the true activity distribution (scale problem). In addition to these methods, Welch et al (1998) proposed to use the consistency conditions of the Radon transform, i.e. a set a mathematical conditions. Later, Bromiley et al (2001) extended this algorithm to include a 3-dimensional template. In a TOF system, the cross-talk problem is solved  but the constant factor  that scales the attenuation sinogram (the forward projection) to the correct absolute values is still unknown, which reflects the scaling problem. With TOF-PET data, an alternative is to estimate the activity image together with the sinogram of the attenuation factors: Rezaei et al (2014) and Defrise et al (2014) proposed an algorithm named MLACF (maximum likelihood attenuation correction factors) that jointly estimates the activity image and the attenuation sinogram (and therefore does not impose consistency of the attenuation correction factors). The scale problem in joint activity and attenuation could be solved by using the information from the MR images and applying an intensity prior on the estimated attenuation correction values in TOF-MLAA as in Mehranian and Zaidi (2015). Finally, the MLAA algorithm could be improved by applying some fixed external transmission sources as in Watson (2014) and Watson et al (2014), or by utilizing the inherent radioactive decay of the lutetium-176 in the LSO-crystals (Rothfuss et al 2014).
In PET/MR brain imaging, the vendor-supplied approach (Siemens) to reconstruct the emission image has been to use the Dixon sequence and to some extent the UTE technique for attenuation correction. Recently, Andersen et al (2014) showed that there is a strong spatial bias using the attenuation map based on the Dixon technique, because the lack of bones causes errors within or near the bone regions, and Dickson et al (2014) stated that the UTE technique still underestimated the activity. For non-TOF MLAA, the problem is severely underdetermined, so using the information from the MR-images to drive MLAA, particularly in the MLTR (maximum likelihood for transmission tomography) algorithm, is a desirable solution.
To achieve this, the local convergence of MLTR is manipulated. For this purpose, the MLTR algorithm proposed by Van Slambrouck and Nuyts (2014) is used, which provides a parameter j α for each voxel j, which affects the local convergence. In this study, a technique is proposed to limit the cross-talk within the brain region, and a way to compute an appropriate value for j α at each iteration is also proposed. For validation purposes, a large set of patients is used with different tracers to test the robustness and reproducibility of the method. The proposed algorithm is compared to OSEM (ordered subset) using a CT attenuation map (which is considered the reference), to OSEM using an attenuation map derived from the Dixon sequences and finally to OSEM using an attenuation map derived from the UTE sequences. The paper is organized as follows. In section 2, the proposed MLAA algorithm is described, the technique to mitigate the cross-talk is explained, and the j α driving is detailed. Then, in section 3, the information concerning the PET/MR system is given as well as the parameters of the reconstructed images. The method to study the cross-talk and determine the j α values is detailed. The patient characteristics are given, and finally the figures of merit analyzing the images are enumerated. The results of the proposed algorithm are shown in section 4 and compared to those of the other algorithms. Finally, in section 5, the results are discussed and some conclusions are drawn.

Materials and methods
In this section the different steps to generate an attenuation map and a corresponding PET image using simultaneous reconstruction of the activity and the attenuation for PET/MR brain imaging are explained. We initially describe the algorithms used in this study, and explain the methods used to incorporate both the T1-weighted MR image and the attenuation map derived from an ultrashort echo time sequences (UTE µ ) in the attenuation reconstruction of MLTR. In the two last sections, a technique to limit the cross-talk for MLAA in a non-TOF system is proposed and the initialisation of the estimated attenuation and activity is described.

MLAA
The study presented in this work is based on a maximum-likelihood activity and attenuation algorithm (MLAA) as proposed initially in Nuyts et al (1999). MLAA maximizes the Poisson log-likelihood, which is given by: with i the index of the projection lines, y i the measured sinogram, y ī the expected measurements computed from the emitter densities → λ and the linear attenuation → µ, both functions of position in the image. MLAA applies an alternated optimisation approach, using in every iteration MLEM to update → λ and MLTR to update → µ. In the following, the complete derivation of each algorithm is not given. A complete derivation of each algorithm is described in Nuyts et al (1999), Rezaei et al (2014), Van Slambrouck and Nuyts (2012) and (2014) for the estimation of the activity and the attenuation coefficients.
In this work, the Compton scatters and the random coincidences were taken into account, and both algorithms use three-dimensional (3D) PET data.
where ib η is the probability that a pair of photons emitted from the voxel b will be detected in the sinogram bin i (ignoring attenuation and detector sensitivity variations), l im is the intersection length of the projection line i with the voxel m, λ and μ are the estimated emitter density distribution and the estimated attenuation coefficients, respectively, Ā is the sinogram expectation of the attenuation factor (including the hardware attenuation factor of the MRI system included in the PET system, such as the RF coils, and the patient table), s represents the efficiencies, σ and ρ represent the scatter sinogram and the smoothed randoms sinogram respectively, and J is the total number of the voxels in the image.
Step (2) is also known as the forward-projection of the current PET image to the sinogram space.
In this study we implemented the common version of the attenuation and normalization weighted maximum likelihood expectation maximization with an ordinary Poisson model (OP-MLEM) algorithm (Shepp and Vardi 1982, Lange et al 1984, Politte and Snyder 1991 given by: with k the current iteration and I the total number of the pixel detectors.

MLTR algorithm.
To estimate the attenuation values, the maximum likelihood for transmission tomography (MLTR) algorithm is based on a gradient-ascent method including the scatter and the random terms. It is similar to the algorithm proposed and described in detail by Van Slambrouck and Nuyts (2014), and is given by: where l ij is the effective intersection length of the projection line i with voxel j. The factor j α is a design parameter, which can be set to any non-negative value and has an impact on the (position dependent) convergence of the algorithm. Equation (6) is slightly different from equation (15) in Van Slambrouck and Nuyts (2012) because an approximation was introduced to simplify the last factor in the denominator. This approximation assumes that y k ī ) ( ′ is equal to y i , as is expected close to convergence. As described in the next section, the scale factor j α can be chosen for each voxel j to improve and accelerate the convergence of the MLTR algorithm towards the correct attenuation value. It is important to note that the j α values are relative, since ĵ µ ∆ is weighted by the sum of all j α along the projection line in the denominator in the expression (6). Setting j α to 1 for all the voxels in the image is similar to the MLTR algorithm described in Nuyts et al (1998), and setting j α to ĵ µ produces the convex algorithm described by Lange and Fessler (1995). The value j α can be set above or below 1 to accelerate or slow down the updating.

Optimization.
In order to accelerate the MLAA algorithm, a block-iterative version of these algorithms was used as proposed by Hudson and Larkin (1994) for the emission estimation (OSEM). At the end of each subiteration, for the OSEM or the MLTR algorithm, a different subset was selected until all the subsets had been used.

Choosing the scale factor α j
In this study we proposed to drive the MLTR algorithm choosing an appropriate value of j α for each voxel of the image depending on the values of the T1w MR image and the attenuation map derived from the UTE µ (Keereman et al 2010). Both images were directly provided by the Siemens software at the end of the PET/MR acquisition.
The T1w image is related to the proton density and the spin relaxation mechanisms but, unlike a x-ray computed tomography (CT) image, it has no information about the electron density which is important for photon attenuation. In the case of brain imaging, the air and the cortical bone regions both have a very low signal intensity in the range (0-80) MR units (making it difficult to distinguish these regions from each other). The brain and soft-tissue regions have a high value (>100) with a good contrast.
All the details concerning the initializations of the first attenuation image used during the reconstruction are given section 2.4. In the MLTR algorithm (5), the update term ( ĵ µ ∆ ) is different for each voxel in the image and changes for each new subset. Three different cases should be distinguished: • ĵ µ ∆ = 0, the algorithm has already found an optimal attenuation value and it is not necessary to update further. • ĵ µ ∆ < 0, the algorithm decreased the attenuation value because the previous value of voxel j was too high according the likelihood. This case could be the air region in the head or around the patient, for instance.
In our experience, T1w MR images are typically of high quality, but do not provide enough information to separate air and bone. The UTE-μ image provides additional information useful for separating air and bone, but this image is not free of artifacts and may for example have regions with missing bone or excess bone. By combining image intensity information from both images to adjust j α , we attempt to complement the information from the T1w MR image with that of the UTE image, without incorporating the artifacts of the latter. In section 3, the optimization of the j α values X i (i = 0,...,3) to generate a correct attenuation value is given.

Cross-talk mitigation
The use of a simultaneous activity and attenuation reconstruction in a non-TOF system, like the Siemens mMR, generates a cross-talk effect in the images as shown by Nuyts et al (1999). This effect is possible because even close to the true solution some other solutions exist for a non-TOF PET. In Rezaei et al (2012), it was shown that the use of TOF information removes the cross-talk effect, but the final attenuation reconstruction is only determined up to a constant, that can be determined using a limited prior knowledge .
In order to limit the cross-talk and to mitigate the problem of the constant, our idea was to use the information included in the current estimate of the attenuation map from MLTR and the T1w image. Considering that the biggest region impacted by the cross-talk is the brain, a method to limit the effect in this specific region was developed. To do that, a kernel of N N × pixels was used, centered in each voxel of interest, checking for bone voxels in the current attenuation image and the signal intensity values in the T1w image. If only brain attenuation values were found in the kernel ( 0.11 µ < cm −1 ) and if the signal intensity in the T1w image was high (>70), the attenuation values for (center) the voxel j was confined to the range (0.096; 0.100) cm −1 . This range of values corresponds to the CT range brain attenuation. For the other voxels, MLTR updated the values normally. Because of this use of a priori information of expected brain tissue attenuation values, the scaling problem is also solved.

Map initialisation
The initial attenuation image is computed as close as possible to the true solution to avoid that MLTR algorithm converges to a local maximum of the likelihood. The initialisation consists of: (i) Creating a mask of the volume to reconstruct from the T1w image, as close as possible to the real support of the patient. To generate this mask a region-growing approach was used. Only one seed was used, placed at one of the corners of the T1w image. Then starting from this seed, all neighboring voxels inferior to a certain signal intensity value (80 in our case for every reconstruction) were set to 0 in the mask image, otherwise the values were set to 1. (ii) Using the previous mask, the estimated attenuation image values were set to 0.098 cm −1 (expected tissue attenuation value) inside the mask region, and 0.000 cm −1 outside. (iii) Considering there was no activity outside the patient, the mask previously computed was also used to initialize the estimated activity image. All the voxels inside the mask were set to 0.1 count per unit volume, and 0.0 outside. (iv) Finally, the first j α image was also initialized with this mask for the first subset in MLTR. All the voxels included in the mask were set to 1 and the voxels outside is 0.

Experiment
This section describes the scheme for evaluation of the proposed MLAA algorithm. The parameters to mitigate the cross-talk effect, and to optimize the j α values are also explained. The details about the patient characteristics and the set of figures of merit applied to evaluate the reconstruction results can be found at the end of this section.

Acquistion and image processing
The brain patient data sets were acquired in a fully-integrated PET/MR system (Biograph mMR, Siemens Healthcare) (Delso et al 2011). Each PET/MR acquisition was associated with a low dose CT from a PET/CT system (Biograph TruePoint or Biograph mCT(64), Siemens Healthcare) . In this study, the CT attenuation map and the PET reconstruction associated with it were considered as the reference for all the comparisons.
From the PET/MR acquisition, 2 attenuation maps were provided by Siemens. The first was based on the dixon-water-fat segmentation (DWFS) (Martinez-Möller et al 2009) that is the current industry standard method to correct for the attenuation in the PET/MR system. The dimension of this image was 126 192 128 × × with (∼2.60 2.60 3.12 × × ) mm 3 voxel size. The sequences were acquired with a repetition time (TR) of 3.6 ms and an echo time (TE) of 1.23/2.46 ms. The second attenuation map was derived from the UTE sequence (Keereman et al 2010) which, unlike DWFS, obtains a signal from the bone and therefore provides information that can be used for separating air and bone in the images. resampled to the Dixon attenuation map (Dixon µ ) to obtain the same voxel size. Before the reconstruction processing, Dixon µ , UTE µ and CT µ were resampled to the PET voxel size of (∼2.08 2.08 2.03 × × ) mm 3 . As explained previously, the estimated emission and the estimated attenuation suffer from a cross-talk effect. To evaluate the impact of the kernel technique to mitigate this effect, one [ 18 F]FDG brain patient data set was reconstructed with different kernel sizes and also without using the kernel technique. For this study the injected activity was 199 MBq and 10 min acquisition in 1 bed position. Kernel sizes used were: 3 3 × , 5 5 × , 7 7 × and 9 9 × . The emission image and the attenuation image were reconstructed with 5 iterations and 21 subsets (OSEM and MLTR). The dimension of the reconstructed volume was 150 150 127 × × and the voxel size (∼2.08 2.08 2.03 × × ) mm 3 . For this experiment 1 j α = was used everywhere. To quantify the impact of the kernel size, the sum of the squared differences was used, defined by: where j is the voxel of the image, X is the reconstructed attenuation or activity using MLAA and X corresponds to the CT µ in the case of the attenuation study, or to the activity computed by OSEM using a CT µ for the emission study. The region of interest (ROI) was defined as a large ellipsoid region within the brain.

Optimization of the α j parameters
As mentioned in the previous section, the goal of the method is to adjust the j α parameter in order to accelerate or reduce local convergence and by doing so, force MLAA towards a desirable solution. We had 4 parameters X i (i = 0,...,3) to optimize. To find the optimal j α values, systematic reconstructions of one patient were studied, assuming that the optimal j α values from that single patient could be applied to all the patients. To find the ideal value of j α , two sets of reconstructions were computed.
Finally, from the first set of reconstructions, the best combination of j α was defined as the one which minimized the difference between the attenuation corrected activity image and the reference activity image. 3 new values of j α were chosen, leading to 3 4 = 81 reconstructions. In this study, we limited the j α search values for X 0 1 ⩽ and for X i (i = 1,...,3) ⩾1. To evaluate the impact of the j α values on the reconstruction, we computed the sum of squared differences as defined in (7) in a ellipsoidal ROI. Only the estimated activity (not attenuation values) by MLAA reconstruction was analyzed by comparison to activity computed by OSEM using a CT µ .

Patient characteristics
To evaluate the robustness and the reproducibility of the method, 3 different clinical data sets were used. The sets consisted of 204 [ 18 F]FDG studies, 35 [ 11 C]PiB studies and one [ 18 F] FET study. The injected activities were 203 20 ± , 429 83 ± and 208 MBq respectively and the time between the tracer administration and the start of the scan was 665 184 ± , 1337 256 ± and 1200 s respectively. All the patients had a PET/CT (only the attenuation scan) examination before the PET/MR acquisition, in order to obtain the CT µ used as a reference, which is generally accepted in the literature. After each PET/MR examination, Dixon µ and UTE µ were accessible.

Evaluation
In order to quantify the accuracy of the proposed MLAA algorithm, its results were compared to those of: -OSEM using CT µ , used as a reference.
In all cases we reconstructed the data using 5 iterations and 21 subsets. The image dimension was 150 150 127 × × and the voxel size (∼2.08 2.08 2.03 × × ) mm 3 . We used a span 11 and a maximum ring difference of 60. For each patient, the normalization sinogram and the smoothed random sinogram were computed by the Siemens software (e7 tools). Since the single scatter simulation algorithm (Watson et al 1996, Watson 2000 from Siemens needs the attenuation map (as a background mask) to estimate the contribution of scattered photon pairs, the different attenuation maps used in the different methods could have a small impact on the scatter correction of the final PET image. In this study, we wanted to quantify the impact of the attenuation image only, and for this reason only the scatter sinogram computed from OSEM using CT µ was applied.
A set of different figures of merit were used to evaluate the accuracy of the MLAA method compared to the other approaches. All the PET reconstructed images and all the attenuation images were aligned to MNI space (Montreal Neurological Institute. The MNI defined a standard brain by using a large series of MRI scans on normal controls.) using the ICBM (International Consortium for Brain Mapping) 2009a template (Fonov et al 2009(Fonov et al , 2011. Then four figures of merits were deduced: • Normalized PET difference images (masked within the brain) of all the patients aligned to the MNI atlas. These values were computed as follows: where X denotes Dixon µ , UTE µ or MLAA map µ− . The resulting average images and the standard deviation images were obtained.
• From the aligned normalized PET images, we calculated the voxel-wise probability of having errors greater than ±5%. • The ROI analysis in the brain defined from the MNI atlas in 11 different regions. For each specific region in the brain, the percent difference was computed and plotted in a bar plot.

Cross-talk evaluation
In figure 1(a), the first row shows the impact of the kernel size in the MLAA reconstructed attenuation map. The image of the CT attenuation defined as the reference is also shown as the last image of the row. Figure 1(c) is a plot of the sum of the squared difference between the reconstructed emission using MLAA and the reconstructed emission using OSEM with the reference CT attenuation map. Figure 1(d) shows the plot of the sum of the squared differences between the reconstructed attenuation map from MLAA and the reference CT attenuation map. When the kernel technique is disabled, the cross-talk has the maximum impact. In figures 1(c) and (d) and visually in the images in figure 1(a), the kernel 3 3 × provides the largest reduction in cross-talk in the brain region. With a kernel size greater than 5, the effect of the kernel technique diminishes and the cross-talk effect reappears. For all the following reconstructions in this study, a 3 3 × kernel was used. within the selected domain is framed in a red box corresponding to X 0 = 0.001, X 1 = 50, X 2 = 1.5 and X 3 = 5. Because all the voxels within the mask were initially set to 0.098 cm −1 , meaning close to the solution for the brain attenuation, the parameter X 0 is small in order to slow down the conv ergence around this brain attenuation value. A relatively high value for X 1 is necessary to decrease the attenuation value and to find the air cavity regions. Both parameters X 2 and X 3 should be chosen not too high in order to find the optimal attenuation value for the bone values from the selected patient were reused for all the following reconstructions in this study assuming that these parameters would be nearly optimal for all the brain FDG PET data.
The optimisation procedure systematically selected the smallest value for X 0 , probably to minimize the cross-talk effect. We have only evaluated non-zero values for X 0 , because we did not want to stop the algorithm entirely from updating particular pixel values. Figure 4 compares the MLAA reconstructed attenuation map and emission image using a selected non-optimal set of values of j FDG ( ) α (X 0 = 1, X 1 = 5, X 2 = 5, X 3 = 500), the optimal values of j FDG ( ) α (defined above) and OSEM. In the case of the non-optimal j FDG ( ) α (first column figure 4), the MLAA algorithm failed to find the air cavity and the attenuation of the bone structure was too high leading to an overestimation of the reconstructed activity. In the case of the optimal j α (second column figure 4), the MLAA algorithm succeeded to find an air cavity similar to that of the CT µ . This further implies very similar PET images from reconstructions with MLAA and OSEM CT .   values, the second column shows a PET estimated and attenuation images using the optimal j FDG ( ) α values. Last column first row is the PET image computed from OSEM CT (last column, second row) that is the reference.

FDG patients evaluation
In figure 5, the normalized difference images images (according to equation (8)) for the FDG patients are illustrated in one central, transaxial slice. The result from the MLAA algorithm can be compared to OSEM DIXON and OSEM UTE . The results show that MLAA generated less error than the other methods, in particular close to and around the bone and close to the air cavity region at the frontal sinuses. The radial error, underestimation of activity values increasing from the center of the brain to the skull, which is large in the UTE μ-map and the Dixon μ-map, due to the missing bone as described in Andersen et al (2014), is very small when using the MLAA algorithm. Visually, a slight overestimation (∼2-3%) in MLAA within the brain region can be noted due to a small overestimation of the bone values in the estimated attenuation.
In figure 6, the standard deviation images, computed around the previous averaged results for the 204 FDG patients, are shown. The result documents that the MLAA algorithm seems to produce robust results with standard deviations at the same level as the other two methods in each region of the brain. Close to the air cavity region the variation is smaller than that of the OSEM UTE where the error from patient to patient could be very different.
In figure 7, the region of interest (ROI) analysis averaged on 204 FDG patients is shown. In the full brain region, the MLAA produces less error than the other reconstructions: the error is ∼−1.5% compared to ∼−6.5% for OSEM UTE and ∼−11.8% for OSEM DIXON . In all the regions, MLAA seems to obtain better results however, the errors in 2 regions (medulla and occipital) seem to be larger than other regions. With the results from the frontal cortex and the temporal regions, MLAA seems to be the more efficient algorithm to produce the correct activity near the air cavity region with an averaged error ∼−1%.
The probability map of errors superior to 5%, calculated from the 204 FDG patients, is shown in figure 8. Visual inspections shows that MLAA provides activity reconstructions   which are less likely to produce errors greater than 5% compared to the standard OSEM DIXON and OSEM UTE activity reconstructions. MLAA produces better results at air cavities than others, but there is some residual errors because of a slight overestimation of the bone values around the air cavity region. This result confirms that the MLAA algorithm has a higher probability to produce a more accurate emission map close to the bone and in the air cavity region compared to other methods.
In table 1, the average amount of brain having errors at different levels (5, 10 and 15%) is shown for the full brain region. For each case, the MLAA algorithm has a smaller affected amount of the full brain region than OSEM DIXON or OSEM UTE .

PiB and FET patient evaluation
In figure 9, the averaged percent difference images for the 35 PiB patients are shown. The results show that on average MLAA generates a uniform under-estimation (∼−0.5 − 1%) in the activity reconstructions. Note that the performance of MLAA is slightly superior to OSEM UTE (∼−3% averaged underestimation error in the support of the brain region). The  radial error is reduced with MLAA, and close to the bone and air cavity regions the MLAA algorithm seems to generate less error than the other techniques of reconstruction.
In figure

Discussion
In this work, we investigated and optimized the non-TOF MLAA algorithm by manipulating the step size of the MLTR attenuation update as described in Van Slambrouck and Nuyts (2014), in the context of PET/MR brain imaging. The algorithm relies only on the data provided by the Siemens mMR Biograph system, i.e. once established, it required no additional input from CT or external database (atlas). Currently, OSEM DIXON (Martinez-Möller et al 2009), including only water, fat, and soft-tissue values, is the basic, industry standard method for PET attenuation correction in a PET/MR system. For brain imaging, an attenuation map computed from the UTE sequences (Keereman et al 2010), including the bone values, is the more precise technique currently used in the clinic to correct for the attenuation. However, these techniques often yield systematic bias and fail locally in individual cases. Our goal was to propose a technique easy to implement, not based on an atlas technique, and efficient in computation time (∼10 min to reconstruct a brain with 5 iterations and 21 subsets, ∼2.0 mm voxel size and a dimension of 150 150 127 × × ). When standard MLAA is used with a non-TOF system, a cross-talk between attenuation and emission reconstructions seems inevitable (Nuyts et al 1999). In this work, the first step was to reduce or eliminate this cross-talk in the reconstructed images. A simple kernel technique was proposed and tested and a kernel size of 3 3 × was found to be optimal. Larger kernels were less effective, most likely because the probability of having a bone value in the kernel is higher, allowing MLAA to update normally and thereby generating some cross-talk. Another parameter having an impact on the cross-talk in this study was the update scaling parameter X 0 . Reducing its value in the soft-tissue and in the brain region made the algorithm slow down in this region limiting the cross-talk.
As explained in Rezaei et al (2012) and Defrise et al (2012), even with TOF information MLAA can provide an activity reconstruction which differs from the true activity by a global factor. In this study, we do not need to determine this constant since we impose a known attenuation in the brain region (where X 0 = 0.0005) and also by using a well chosen initial value in the first attenuation map ( 0.098 µ = cm −1 ) very close to the brain attenuation value. In this work, the j α parameter was studied in some detail. In order to manage properly each voxel of the attenuation image, the j α parameter was set to different values between the voxels as in Van Slambrouck and Nuyts (2014). For each of the tracers, the j α parameter values were optimized using a single patient study only, and therefore, we cannot claim that these values are optimal. Optimization on a group of patients was not considered because the procedure is very time consuming. However, similar results were obtained for the three tracers and good MLAA performance was observed on the entire set of patient studies. This indicates that the procedure is rather robust and yielded a useful set of parameter values. The optimization procedure produced a similar strategy for the three tracers: for pixels probably belonging to tissue, the updates are minimized, for pixels likely corresponding to air, the convergence is strongly accelerated, and for pixels suspected to correspond to bone, convergence is encouraged more if UTE μ is higher, which makes sense because a higher UTE μ values is evidence that the pixel corresponds to bone and not to air. For [ 18 F]FDG and [ 18 F] FET the j α values were very similar, the value of X 1 (air convergence) was higher for [11C]PiB than for the other tracers. Further studies may be performed to determine the impact, if any, of the statistics and the distribution of the tracer on the j α values. Based on the T1w MR image and the μ-map computed from the UTE sequences we encourage insertion of air and bone in the final μ-map computed by the MLAA algorithm. But we do so taking into account the sign of the MLTR update and by modifying the local convergence, which ensures that this insertion will not cause disagreement with the measurement.
With the optimized j α values, the proposed MLAA algorithm yielded fairly accurate results, and the performance seems to be rather independent of the tracer. Note that in some specific region of the brain, the MLAA algorithm performed markedly better than the UTE and Dixon based approaches. As it might be expected, this was particularly the case close to the air region cavity, because the MLTR algorithm is able to restore more precisely the shape of the air cavity and the bone structures found around the air cavity. Further studies are needed to determine in detail what caused the relatively high error in the medulla region. In the occipital region, an underestimation of ∼−5% was noted, due to some missing bones in this region. However, all the reconstruction methods had a high underestimation is this region, with MLAA being significantly better than OSEM DIXON and OSEM UTE .
In this study, the scatter correction was performed with the scatter sinogram computed from the CT μ-map to avoid the impact of the scatter in the results. Unlike OSEM DIXON or OSEM UTE , an initial attenuation map, as the first attenuation map used at the first subset of the proposed MLTR, should be applied to estimate the scatter producing a non-optimal scatter correction. In this work, the impact of the scatter was not studied, but in a recent study, Burgos et al (2014b) showed that using a μ-map like the UTE μ-map led to PET difference errors less than 1% in the brain region compared to OSEM using an ideal μ-map, which is the gold standard the CT μ-map. The uniform first μ-map used for MLAA had no bone, so it could be interesting to evaluate the error due to the scatter in the PET image. If the error is deemed too high, a new scatter correction could be done with the final μ-map computed from MLAA, and finish the reconstruction with the OSEM algorithm and the μ-map previously computed by MLAA.
In the proposed method, a few steps could be improved. Due to the uniform attenuation in the brain being fixed to the range (0.096; 0.100) cm −1 , a small overestimation of the PET value is obtained. Note that the mean attenuation of the ventricle is 0.095 cm −1 on CT. To improve this, the ventricle region could be segmented in the T1w MR image, then incorporated into the first attenuation map for MLTR, and j α could be set to 0 to keep the value of the ventricles constant. For the selection of j α , the raw values of the T1w MR images were used. Because the value could be different between patients (and more important, between systems and set-ups), a normalized value of the T1w MR image should be used instead. The impact in this homogenous study was probably limited.
Some additional information could be used to improve the PET image using a joint attenuation and emission reconstruction. Firstly, Watson (2014) and Watson et al (2014) proposed to use an external source around the patient to add some supplementary transmission data and improve the attenuation map computed by MLTR. Secondly (Berker et al 2014) proposed to use the information from the scattered coincidences and the photon energy measurement in order to improve the attenuation map computed by the MLAA algorithm. Finally, if the TOF is enabled in a PET/MR system, Rothfuss et al (2014) proposed to reconstruct an attenuation map using the inherent radioactive luthetium in the (LSO) crystal of the PET system as a transmission source. This information could be coupled with the MLAA algorithm to improve the attenuation and the PET image.

Conclusions
In this paper, a new approach using a non-TOF MLAA algorithm for PET/MR brain imaging is proposed, without any atlas information. The MLTR algorithm is guided by manipulating its local convergence, based on values from the UTE μ-map and the T1w MR image provided by the scanner. A parameter tuning method was proposed and applied for [ 18 F]FDG, [ 11 C]PiB and [ 18 F]FET tracers.
The algorithm is validated using more than 200 clinical patient data sets and reconstructed in 3D. It was found that the proposed MLAA algorithm had better results than OSEM DIXON or OSEM UTE for each tracer.