---
abstract: '[ Image reconstruction from computed tomography (CT) measurement is a challenging statistical inverse problem since a high-dimensional conditional distribution needs to be estimated. Based on training data obtained from high-quality reconstructions, we aim to learn a conditional density of images from noisy low-dose CT measurements. To tackle this problem, we propose a hybrid conditional normalizing flow, which integrates the physical model by using the filtered back-projection as conditioner. We evaluate our approach on a low-dose CT benchmark and demonstrate superior performance in terms of structural similarity of our flow-based method compared to other deep learning based approaches. ]{}'
bibliography:
- 'references.bib'
---

Introduction {#sec:introduction}
============

Many important applications in medical imaging, such as computed tomography (CT) or magnetic resonance imaging (MRI), can be formulated as an inverse problem. The inverse problem consists in the reconstruction of an internal image of a patient based on radiological data. Many of these applications are ill-posed inverse problems, as small measurement errors can result in large errors in the reconstruction. In a classical way, an inverse problem is often formulated as follows: A forward operator $A: X \rightarrow Y$ maps the image $x^\dagger$ to (noisy) measurements $$\begin{aligned}
    y^\delta = Ax^\dagger + \mu,
    \label{eq:inverse_problem}\end{aligned}$$ where $\mu$ describes the noise. The research in inverse problems is focused in particular on developing algorithms for obtaining stable approximations of the true solution $x^\dagger$. In order to cover the uncertainties that occur especially with ill-posed problems, the theory of Bayesian inversion considers the posterior distribution $p(x|y^\delta)$ [@Dashti.2272013]. This posterior is the conditional density of the image $x$ conditioned on the measurements $y^\delta$.

The main task in statistical inverse problems is to approximate this high-dimensional conditional distribution. For high-dimensional, structured images, like they arise in CT, this is a challenging process. In the field of density estimation, conditional normalizing flows (NF) [@Winkler.11292019; @Ardizzone.742019] allow to learn expressive conditional densities by maximum likelihood training. Since the physical model is known in CT (Eq. \[eq:radon\_trafo\]), we propose a hybrid approach which integrates model-based reconstruction with conditional NFs.

![Ground truth samples from the LoDoPaB-CT dataset containing artifacts. These errors stem from the reconstruction technique that was used on the normal-dose measurements.[]{data-label="fig:ArtifactSamples"}](figures/lodopab_artifact_samples.pdf){width="\columnwidth"}

In many CT image reconstruction tasks the mean squared error (MSE) is used [@chen2017convnet_ct; @he2020iradon], which, however, has many known limitations [@Zhao.2017]. In the context of maximum likelihood estimation (MLE), the MSE loss arises from the assumption of i.i.d. standard Gaussian noise. However, this assumption is violated in CT training data since they are often obtained from reconstructions of high-dose or normal-dose measurements. E.g. the choice of the reconstruction algorithm can lead to artifacts, as shown in Figure \[fig:ArtifactSamples\]. This implies that the reconstruction error for individual pixels is no longer independent. We argue that these dependencies can be better captured by a flow-based model.

Our contributions are twofold: 1) We apply conditional normalizing flows to CT image reconstruction. 2) We propose a hybrid approach, which integrates the physical model by using the filtered back-projection as conditioner.

There are several deep learning approaches for implementing conditional generative models. Besides conditional generative adversarial networks (cGANs, [@Mirza.1162014]) and conditional variational autoencoders (cVAE, [@cVAE]) there are conditional normalizing flows [@Ardizzone.14.08.2018; @Winkler.11292019], which we will explore in this paper. These models can be understood as an extension of normalizing flow models [@Papamakarios.1252019]. Normalizing flows define a bijective mapping and have a tractable Jacobian determinant. Depending on the specific implementation, the inverse can be determined analytically (e.g. [@Dinh.27.05.2016], [@Kingma.792018]) or numerically (e.g. [@Behrmann.], [@Chen.20191207]). These models have in common that they define an invertible, differentiable transformation of the complex distribution $p(x)$ in image space into a base distribution $p(z)$. Using the change of variables formula, it is possible to evaluate the density of $p(x)$ using this transformation. This allows for training the model with the exact likelihood.

Related Work
------------

Deep learning methods have been successfully applied to many ill-posed inverse problems such as CT [@Arridge.2019]. In particular, end-to-end learned methods have been used. Those methods can be classified in three main groups: post-processing [@chen2017convnet_ct], fully-learned [@he2020iradon] and learned iterative algorithms [@adler2017]. These end-to-end methods have in common that they learn a parameterized operator ${T_\theta:Y \rightarrow X}$ by optimizing the parameters $\theta$ using training data $\{ (y_i^\delta, x_i^\dagger) \}_{i=1}^N$. Usually, this is done by minimizing the MSE between the ground truth data $x_i^\dagger$ and the reconstruction $T_\theta(y_i^\delta)$ as $$\begin{aligned}
    \hat{\theta} \in \arg \min_{\theta \in \Theta} \frac{1}{N} \sum_{i=1}^N \Vert T_\theta(y_i^\delta) - x_i^\dagger \Vert^2.\end{aligned}$$

Recently, Deep Image Priors were used for CT, achieving promising results in the low-data regime [@baguer2020computed]. Similar to our approach is the work of [@Adler.11142018], who employed a Wasserstein GAN to draw samples from the conditional distribution. However, in this approach it is not possible to evaluate the likelihood of the generated samples. [@Ardizzone.14.08.2018] have used invertible neural networks to approximate the conditional distribution and to analyze inverse problems. In a subsequent paper this concept was extended to conditional invertible neural networks (cINNs) which yielded good performance in the field of conditional image generation [@Ardizzone.742019].

Background on Computed Tomography {#sec:Background}
=================================

Computed tomography allows for a non-invasive acquisition of the inside of the human body, which makes it one of the most important tools in modern medical imaging [@buzug2008computed_tomography]. CT is a primary example of an inverse problem. The determination of the interior distribution cannot be achieved directly. It has to be inferred from the measured attenuation of X-rays sent through the body.

Current research focuses on reconstruction methods for low-dose CT measurements to reduce the health risk from radiation [@shan2019nn_vs_vendor_ct; @baguer2020computed]. One strategy to reduce the dose is to measure at fewer angles. This can result in undersampled measurements and therefore in the existence of ambiguous solutions to the inverse problem. Another option is to lower the intensity of the X-ray. This leads to increased Poisson noise on the measurements and adds to the instability of the inversion. In this paper, we test our reconstruction model for the lower intensity case.

The basic principle of a CT machine with a parallel beam geometry can be described by the 2D Radon transform [$A: X \rightarrow Y$]{} [@radon1986radon_trafo]: $$A x(s,\varphi) = \int_{\mathbb{R}} x\left( s \begin{bmatrix} \cos(\varphi) \\ \sin(\varphi) \end{bmatrix} + t \begin{bmatrix} -\sin(\varphi) \\ \phantom{-}\cos(\varphi) \end{bmatrix} \right) \, \mathrm{d}t.
    \label{eq:radon_trafo}$$ It is an integration along a line, which is parameterized by the distance $s \in \mathbb{R}$ and angle $\varphi \in [0,\pi]$ (cf. Figure \[fig:parallel\_beam\]). For a fixed pair $(s,\varphi)$ this results in the log ratio of initial and final intensity at the detector for a single X-ray beam (Beer-Lambert’s law). The whole measurement, called *sinogram*, is the collection of the transforms for all pairs $(s,\varphi)$. The task in CT is to invert this process to get a reconstruction of the body. The inversion of the Radon transform is an ill-posed problem since the operator is linear and compact [@natterer2001math_tomography]. The consequences is an instable inverse mapping, which amplifies even small measurement noise.

A common inversion model is the filtered back-projection (FBP) [@Shepp1974FBP]. The reconstruction for position $(i,j)$ is calculated by a convolution over $s$ and an integration along $\varphi$ as $$\begin{aligned}
    x(i,j) = \int_0^{\pi} y(s, \varphi) \star h(s) \vert_{s = i \cos(\varphi) + j \sin(\varphi)} \, \mathrm{d}\varphi.\end{aligned}$$ In general, $h$ is chosen as a high-pass filter such as the Ram-Lak filter [@Ramachandran2236]. In reality, we can only measure for a finite number of angles and distances. In this discrete setting the FBP only works well for a high number of measurement angles. Otherwise severe streak artifacts can appear in the reconstruction.

Normalizing Flow
----------------

Normalizing flows are a special class of probabilistic models, which transform a simple base distribution via a number of invertible transformations into a desired target density [@Papamakarios.1252019]. They have been used for variational inference [@Rezende.5212015], density estimation [@Dinh.27.05.2016] and generative modeling [@Kingma.792018; @Chen.20191207]. Recently the field of application was extended to inverse problems [@Ardizzone.14.08.2018]. Normalizing flows offer the possibility of exact likelihood training.

Normalizing flows can be implemented using deep invertible neural networks. Research on these invertible networks has developed a number of different building blocks. [@Behrmann.] define a class of invertible transformations based on the architecture of residual networks [@He.10.12.2015]. We focus on models based on affine coupling layers [@Dinh.30.10.2014; @Dinh.27.05.2016]. In particular, we focus on the Glow architecture [@Kingma.792018]. This model extends the coupling layer approach by two components, ActNorm and invertible 1x1 Convolutions. These components are arranged in a multi-scale architecture [@Dinh.27.05.2016]: Every $L$ layers an invertible squeeze operation is applied to reduce the spatial dimension, and a part of the channels is split off to reduce the computational cost.

Methods {#sec:methods}
=======

Problem Setting {#sec:probSetting}
---------------

To estimate conditional densities, data pairs from measurements $y^\delta$ and ground truth images $x^\dagger$ are required. In computed tomography (CT) it is not possible to obtain actual ground truth data, because no picture can be taken of the interior of the human body. Instead of using ground truth images we use reconstructions based on high-dose measurements $y^{\delta_1} = A x^\dagger + \mu_1$, i.e. $x^{\delta_1} = T_\text{FBP}(y^{\delta_1})$. Because $x^{\delta_1}$ is the output of an reconstruction algorithms, it can contain artifacts and differ from the actual image $x^\dagger$, see Figure \[fig:ArtifactSamples\] for an example. In the next step we simulate low-dose CT measurements using this reconstruction as $y^{\delta_2} = A x^{\delta_1} + \mu_2$, where $\mathrm{Var}[\mu_2] \geq \mathrm{Var}[\mu_1]$, since low-dose measurements are more prone to measurement noise. The training set then consists of data pairs $\{y^{\delta_2}, x^{\delta_1} \}$. An example of such a dataset is LoDoPaB-CT [@leuschner2019lodopabct], which we use to benchmark our proposed conditional flow.

Normalizing Flow with FBP Conditioning
--------------------------------------

From a statistical point-of-view, an inverse problem can also be seen as a generating process $x \sim p(x|y)$ [@Dashti.2272013; @Arridge.2019]. The task in such a statistical inverse problem is to estimate this conditional distribution. We are using conditional normalizing flows (NF) [@Winkler.11292019] to approximate the target density $p(x|y)$. The conditional NF is composed of a series of invertible transformations $F= f_K \circ \dots \circ f_1$. Here, every individual transformation is parameterized by $\theta$ and receives a conditional input $y$: $f_i = f_{\theta_i}(\cdot, y)$. This transformation defines a transport map, which converts the initial density into a simple, easy-to-sample density $p_Z$. This model defines a conditional density $q(x|y,\theta)$ and using the change-of-variables formula the conditional density can be calculated: $$\begin{aligned}
    q(x|y;\theta) = p_{Z}(F_\theta(x;y)) \left\vert \det\left(\frac{\partial F_\theta(x;y)}{\partial x}\right) \right\vert.\end{aligned}$$ We denote the Jacobian for one data point $x_i, y_i$ with $J_i = \frac{\partial F_\theta(x_i;y_i)}{\partial x}$. Instead of directly using the measurements $y_i$ as conditional inputs, we propose to employ a reconstruction, e.g. the filtered back-projection $\hat{x}_i = T_{\text{FBP}}(y)$.

Assume a dataset $\{(x_i, y_i)\}_{i=1}^N$ of measurements $y_i$ and reconstructions $x_i$. To approximate the target density $p(x|y)$ a conditional NF $q(x|y, \theta)$ can be trained using the negative log-likelihood as a loss function. Using a standard normal distribution, i.e. $p_Z \sim \mathcal{N}(0, I)$, this amounts to minimizing $$\begin{aligned}
    \mathcal{L}(\theta) = \frac{1}{N} \sum_{i=1}^N \frac{\Vert F_\theta(x_i;T_{\text{FBP}}(y_i))\Vert_2^2}{2}  - \log \left\vert \det J_i \right\vert .\end{aligned}$$

Conditional coupling layers
---------------------------

We are using the conditional coupling layer from [@Ardizzone.742019] to construct a conditional invertible neural network (cINN), which is an extension of the affine coupling layer from [@Dinh.27.05.2016]. We propose to integrate the model-driven approach of inverse problems by using the filtered back-projection $\hat{x} = T_\text{FBP}(y)$ as conditional input instead of the raw sinogramm measurements $y$. The input $u=[u_1, u_2]$ to an affine coupling layer is split into two parts and both parts are transformed individually: $$\begin{aligned}
    v_1 &= u_1 \odot \exp(s_1(u_2, \hat{x})) + t_1(u_2, \hat{x}) \\ \newline
    v_2 &= u_2 \odot \exp(s_2(v_1, \hat{x})) + t_2(v_1, \hat{x}).\end{aligned}$$ The transformations $s_1, s_2, t_1, t_2$ do not need to be invertible and are modelled as convolutional neural networks (CNNs). The inverse of an affine coupling layer is: $$\begin{aligned}
    u_1 &= (v_1 - t_1(u_2,\hat{x})) \odot \exp(-s_1(u_2, \hat{x})) \\ \newline
    u_2 &= (v_2 - t_2(v_1,\hat{x})) \odot \exp(-s_2(v_1,\hat{x})).\end{aligned}$$ The log-determinant of the Jacobian for one affine coupling layer can be calculated as the sum over $s_i$, i.e. $\sum_j s_1(u_2, \hat{x})_j + \sum_j s_2(v_1, \hat{x})_j$. A deep invertible network can be built as a sequence of multiple such layers, with a permutation of the dimensions after each layer.

The conditional input $\hat{x}$ is added as an extra input to each transformation in the coupling layer. In practice, an additional conditioning network $H$ is added, so instead of $\hat{x}$ the output $H(\hat{x})$ is used. This conditioning network $H$ is under no architectural constraints and can contain all usual elements (i.e. BatchNorm, pooling layer, etc.) of a CNN.

![Reconstruction and standard deviation of cINN. $1000$ Samples were used for the reconstruction. The top row shows the ground truth image and the corresponding FBP.[]{data-label="fig:Reconstruction"}](figures/ConditionedMean.pdf)

Results {#sec:results}
=======

Sampling from the model is a two-step process: First, a sample $z$ is drawn from the base density $p_Z$. Second, this sample is transformed by the inverse to obtain an image. By repeatedly sampling $z_j$ for a fixed input $y^\delta$ we thus estimate the conditional mean as $$\begin{aligned}
    \hat{x} = \mathbb{E}_{z}[F^{-1}(z,T_{\text{FBP}}(y^\delta))] \approx \frac{1}{n} \sum_{j=1}^n F^{-1}(z_j, T_{\text{FBP}}(y^\delta)).\end{aligned}$$ We evaluate our model on the LoDoPaB-CT dataset [@leuschner2019lodopabct]. For this dataset we are in the case of oversampling, so we expect a uni-modal distribution. This enables the choice of the conditional mean as the reconstruction method. If a highly multi-modal distribution is expected, the conditional mean is not the optimal choice. To measure the error between reconstruction and ground truth, the PSNR and the SSIM [@wang2004ssim] are evaluated. Both are common quality metrics for the evaluation of CT and MRI reconstructions [@joemai2017ssim_in_ct; @adler2017; @zbontar2018dataset_fastmri]. The PSNR is a pixel-wise metric which is defined via the MSE. The SSIM is a structural metric, which compares local patterns of pixels and is not calculated on a per pixel basis.

Implementation {#subsec:implementation}
--------------

We follow the multi-scale architecture design of RealNVP [@Dinh.27.05.2016]. After each block, consisting of 6 coupling layers, downsampling is performed. The downsampling is done using the Haar downsampling from [@Ardizzone.742019] and the variant used in [@Jacobsen.2018]. The dimensions have to be permuted after each coupling layer. This is done using the invertible 1x1 convolutions from [@Kingma.792018]. The model is implemented using the library FrEIA[^1]. A conditioning network was used to further extract features from the filtered back-projection. This conditioning network was trained together with the full flow-based model. Details on the implementation can be found in the supplementary material.

LoDoPaB-CT Dataset {#subsec:lodopab}
------------------

We evaluate our method on the low-dose parallel beam (LoDoPaB) CT dataset [@leuschner2019lodopabct], which contains over two-dimensional CT images and corresponding simulated low photon count measurements. The ground truth images $x^{\delta_1}$ are human chest CT reconstructions from the LIDC/IDRI database [@armato2011lidc_idri], cropped to $\num{362}\times\num{362}$ pixels. Projections are computed using parallel beam geometry with angles and beams. Poisson noise is applied to model a low photon count ($\mu_2$ in Section \[sec:probSetting\]).

Evaluation on LoDoPaB-CT
------------------------

We have evaluated our model on the LoDoPaB-CT dataset. First we examined the dependence of PSNR and SSIM on the number of samples for the conditional mean. The results are shown in Figure \[fig:SSIMPSNR\] (appendix). Both PSNR and SSIM increase with a higher number of samples. This allows for a trade-off between quality of reconstruction and time. For our evaluation we have chosen a conditional mean with 1000 samples. Table \[tab:results\_lodopab\] shows the scores on the test dataset. The comparison includes several classical and deep learning approaches.

The Filtered Backprojection (FBP), TV regularization and Deep Image Prior (DIP) + TV [@baguer2020computed] work without additional training. The FBP can also be combined with a U-Net, which acts as post-processing network [@jin2017cnn_imaging]. The approach to fully learn the inversion process is followed by the iRadonMap [@he2020iradon]. Learned Primal Dual [@adler2017] and Learned Gradient Descent [@adler2017iterative_nn] are learned iterative schemes.

In terms of PSNR our model is comparable to other state-of-the-art deep learning approaches, despite not explicitly trained to minimize the MSE between the prediction and the ground truth images. Regarding the SSIM our model outperforms all other approaches. This further underlines our hypothesis that using the more flexible flow objective enable to incorporate structural properties.

Conclusion {#sec:conclusion}
==========

We have investigated how flow-based models can be applied as a conditional density estimator for the reconstruction of low-dose CT images. Using this generative approach, we were able to obtain high-quality reconstructions that outperformed all other deep learning approaches in terms of structural similarity. So far only coupling-based INNs were used, but future work should explore other architectures such as i-ResNets [@Behrmann.] for this conditional density estimation task. Furthermore, our hybrid approach that integrates the physical model into the conditioning could enable the use of more advanced reconstruction algorithms. Thus, conditional flows are a promising avenue for statistical model-based inverse problems such as CT reconstruction.

Acknowledgements {#sec:acknowledgements}
================

Johannes Leuschner and Maximilian Schmidt acknowledge the support by the Deutsche Forschungsgemeinschaft (DFG) within the framework of GRK 2224/1 “$\pi^3$: Parameter Identification – Analysis, Algorithms, Applications”.

Appendix
========

Model Architecture
------------------

The model was trained for stochastic gradient steps of batchsize with the Adam-Optimizer [@kingma2014adam] using a weight decay of . The last layer in the subnetworks of each coupling layer is initialized with zero. This initializes the model as a whole with the identity.

A multiscale architecture was used for implementation of the cINN. The model includes 6 resolution levels. After each level a part of the channels is split off and passed on to the output. After each resolution level downsampling is performed. Downsampling was performed using the iRevNet variant [@Jacobsen.2018] as well as the Haar Downsampling by [@Ardizzone.742019]. The input size of the CT-images is $1 \times 352 \times 352$. The full cINN was build as follows.

  cINN                                        Output size
  ------------------------------------------- ---------------------------
  iRevNet-Downsampling                        $4 \times 176 \times 176$
  level 1 conditional section                 $4 \times 176 \times 176$
  iRevNet-Downsampling                        $16 \times 88 \times 88$
  Split: $8 \times 88 \times 88$ to output    $8 \times 88 \times 88$
  level 2 conditional section                 $8 \times 88 \times 88$
  iRevNet-Downsampling                        $32 \times 44 \times 44$
  Split: $16 \times 44 \times 44$ to output   $16 \times 44 \times 44$
  level 3 conditional section                 $16 \times 44 \times 44$
  iRevNet-Downsampling                        $64 \times 22 \times 22$
  Split: $32 \times 22 \times 22$ to output   $32 \times 22 \times 22$
  level 4 conditional section                 $32 \times 22 \times 22$
  Haar-Downsampling                           $128 \times 11 \times 11$
  Split: $96 \times 11 \times 11$ to output   $32 \times 11 \times 11$
  level 5 conditional section                 $32 \times 11 \times 11$
  Split: $28 \times 11 \times 11$ to output   $4 \times 11 \times 11$
  level 6 Dense-conditional section           $484$

A conditioning network was used to extract features from the conditional input. Similar to the cINN, this network consists of 6 resolution levels. The output from the resolution level of the conditioning network is used as the conditioning input for the respective resolution level in the cINN. If not specified otherwise, a kernel size of $k=3$ is used. In addition, batch normalization (BN) is applied.

  Conv2d: $1 \rightarrow 3$, stride=2 + LeakyReLU
  -----------------------------------------------------------
  Conv2d: $32 \rightarrow 64$ + LeakyReLU
  Conv2d: $64 \rightarrow 128$ + BN + LeakyReLU
  Conv2d: $128 \rightarrow 64$ + BN + LeakyReLU
  Conv2d: $64 \rightarrow 32$ + BN + LeakyReLU
  Conv2d: $32 \rightarrow 4$ + BN ($\rightarrow$ level $1$)

  LeakyRelu
  -----------------------------------------------------------
  Conv2d: $4 \rightarrow 32$, stride=2 + BN + LeakyReLU
  Conv2d: $32 \rightarrow 32$ (k=1) + LeakyReLU
  Conv2d: $32 \rightarrow 32$ + BN + LeakyReLU
  Conv2d: $32 \rightarrow 8$ + BN ($\rightarrow$ level $2$)

  LeakyRelu
  ------------------------------------------------------------
  Conv2d: $8 \rightarrow 32 $ (k=1) + LeakyReLU
  Conv2d: $32 \rightarrow 64$, stride = 2+ LeakyReLU
  Conv2d: $64 \rightarrow 16$ + BN ($\rightarrow$ level $3$)

  LeakyRelu
  ------------------------------------------------------------
  Conv2d: $16 \rightarrow 64 $ (k=1) + LeakyReLU
  Conv2d: $64 \rightarrow 64$, stride = 2+ LeakyReLU
  Conv2d: $64 \rightarrow 32$ + BN ($\rightarrow$ level $4$)

  LeakyRelu
  -------------------------------------------------------------------
  Conv2d: $32\rightarrow 96 $ (k=1) + LeakyReLU
  Conv2d: $96 \rightarrow 128$, stride = 2+ LeakyReLU
  Conv2d: $128 \rightarrow 32$ (k=1) + BN ($\rightarrow$ level $5$)

  LeakyRelu
  ----------------------------------------------------------
  Conv2d: $32\rightarrow 64 $, stride=2 + LeakyReLU
  Conv2d: $64 \rightarrow 256$, stride = 2+ LeakyReLU
  Average Pooling + Flatten + BN ($\rightarrow$ level $6$)

To implement the subnetworks in the coupling layers a CNN variant and a fully connected variant were used. The input channels are denoted by $c_{in}$ and the output channels by $c_{out}$.

  CNN-subnetwork (k=1) or (k=3)
  ----------------------------------------------
  Conv2d: $c_{in}\rightarrow 64 $, + LeakyReLU
  Conv2d: $64 \rightarrow 92$ + LeakyReLU
  Conv2d: $92 \rightarrow c_{out}$

  Dense-subnetwork
  ----------------------------------------------
  Dense: $c_{in}\rightarrow 512 $, + LeakyReLU
  Dense: $512 \rightarrow 512$ + LeakyReLU
  Dense: $512 \rightarrow c_{out}$

Using this two variants of subnetworks the conditional sections are implemented as follows.

[|l|l|]{}\
Coupling (CNN-subnet k=1) &\
Glow $1 \times 1$ convolution &\
Coupling (CNN-subnet k=3) &\
Glow $1 \times 1$ convolution &\

[|l|l|]{}\
Random permutation &\
Dense-subnetwork &\

After each downsampling a small unconditioned subnetwork is used:

  CNN-subnetwork (without conditional input)
  ---------------------------------------------------
  Conv2d: $c_{in}\rightarrow 64$ (k=1), + LeakyReLU
  Conv2d: $64 \rightarrow 64$ (k=1) + LeakyReLU
  Conv2d: $64 \rightarrow c_{out}$ (k=1)

The downsampling section is built as follows:

  Downsample section (Haar or iRevNet)
  -----------------------------------------
  Haar or iRevNet downsampling
  Glow $1 \times 1$ convolution
  Coupling (unconditional CNN-subnetwork)
  Glow $1 \times 1$ convolution
  Coupling (unconditional CNN-subnetwork)

Evaluation of the Conditional Mean
----------------------------------

We have used the conditional mean as a reconstruction for the CT image. Figure \[fig:SSIMPSNR\] shows the performance in relation to the number of samples used.

![ PSNR and SSIM for the validation set of the LoDoPaB-CT dataset. The PSNR is colored in red and the SSIM is colored in blue.[]{data-label="fig:SSIMPSNR"}](figures/PSNR_SSIM_bySample.pdf){width="\columnwidth"}

Additional Examples
-------------------

![image](figures/Reconstruction.pdf) \[fig:Reconstruction2\]

[^1]: https://github.com/VLL-HD/FrEIA
