# Speak Like a Professional: Increasing Speech Intelligibility by Mimicking
Professional Announcer Voice with Voice Conversion

###### Abstract

In most of practical scenarios, the announcement system must deliver speech
messages in a noisy environment, in which the background noise cannot be
cancelled out. The local noise reduces speech intelligibility and increases
listening effort of the listener, hence hamper the effectiveness of
announcement system. There has been reported that voices of professional
announcers are clearer and more comprehensive than that of non-expert speakers
in noisy environment. This finding suggests that the speech intelligibility
might be related to the speaking style of professional announcer, which can be
adapted using voice conversion method. Motivated by this idea, this paper
proposes a speech intelligibility enhancement in noisy environment by applying
voice conversion method on non-professional voice. We discovered that the
professional announcers and non-professional speakers are clusterized into
different clusters on the speaker embedding plane. This implies that the
speech intelligibility can be controlled as an independent feature of speaker
individuality. To examine the advantage of converted voice in noisy
environment, we experimented using test words masked in pink noise at
different SNR levels. The results of objective and subjective evaluations
confirm that the speech intelligibility of converted voice is higher than that
of original voice in low SNR conditions.

Index Terms: Voice conversion, speech intelligibility, professional announcer,
speech-in-noise

##  1 Introduction

In most of the situation, public announcement systems must deliver speech
messages through adverse listening environments with various competing sounds
and reverberation. These negative conditions degrade intelligibility of the
announcement speech, and in some cases, damage the integrity of the intended
message. One solution to maintain the speech intelligibility is to increase
the intensity of playback volume to increase the signal-to-noise ratio (SNR).
However, this approach is only useful to some degrees due to the limitation of
playback equipment power and the comfort of listener. Consequently, various
approaches have been proposed to increase the speech intelligibility in noisy
environments without increasing the total power of speech. These include
modification of spectral properties [1, 2, 3], dynamic range compression [4,
5, 6], modification of speech modulation spectrum [7, 8], and time-scale
modification [9, 10].

It is known that the voice-related professions, such as professional
announcers, voice actor, and singer, can produce speech with the impression of
clearer and easier to hear than normal person [11, 12]. Moreover, recent
studies have shown that the speech from professional announcer can maintain
its intelligibility better than speech from non-expert person in very noisy
environment [13]. This phenomenon can be exploited to inspire speech
enhancement algorithms that seek to improve speech intelligibility in noise.
In this work, we aim to apply voice conversion technique to mimic a speaking
style of the professional announcer. Voice conversion refers to the process of
modifying voice personality without changing the linguistic information
conveyed in speech waveform. Furthermore, voice conversion can be applied to
control different attributes of voice style, i.e., gender and accent [14].
This study will clarify whether the announcer-adapted speech from voice
conversion model still inherits the noise-resistance property of natural
announcer speech. In addition, this study proposes a method to increase the
speech intelligibility without completely change the voice individuality,
which is useful in the situation where the identity of the speaker needs to be
preserved.

The structure of voice conversion model and the training procedure are
described in Section 2. Next, the detail for experimental settings is shown in
Section 3. Then the results of objective and subjective evaluations are
described in Section 4. Finally, we conclude our paper Section 5 with some
discussions for the future work.

##  2 StarGAN-v2 voice conversion model

Figure 1: Overview structure of StarGAN-v2 voice conversion model. The model
consists of a generator network to convert input mel-spectrogram, a
discriminator network for adversarial training, a style encoder to extract
speaker embedding, a pretrained F0 network for extracting F0 feature, a
pretrained speech recognition to extract linguistic features, and a pretrained
Parallel WaveGAN vocoder to generate waveform from mel-spectrogram.

In this paper, a voice conversion method based on StarGANv2 [15] model is
used. We follow the official implementation of StarGANv2-VC [16] which is
publicly available at 111https://github.com/yl4579/StarGANv2-VC. The overview
structure of voice conversion model is depicted in Fig. 1. The models consist
of 6 modules as follows:

  * •

Style encoder network consists of a stack of 4 residual network (resnet) [17]
layers an average pooling layer to extract the speaker embedding
$s_{\textrm{emb}}$ from mel-spectrogram. The dimension of speaker embedding
vector is set to 128.

  * •

F0 network is pretrained to predict the F0 value and voice/un-voice region
from the input mel-spectrogram. The F0 network is a stack of 2 convolutional
layers, 3 resnet layers and a recurrent bidirectional long-short term memory
layer. The output of the last resnet layer is used as the F0 features for
generator network.

  * •

Speech recognition (ASR) network is pretrained to predict the phoneme sequence
from mel-spectrogram. The ASR network is a joint CTC-attention VGG-BLSTM
network given by the Espnet Toolkit 222https://github.com/espnet/espnet [18].

  * •

Parallel WaveGAN [19] vocoder is pretrained to generate speech waveform from
input mel-spectrogram. The pretrained checkpoint is publicly available at
333https://github.com/kan-bayashi/ParallelWaveGAN.

  * •

The generator network (G) consists of two sub-modules: an encoder network and
a decoder network. The encoder takes the input mel-spectrogram and generate
the hidden feature vector. The concatenation of hidden feature vector and F0
feature vector is fed to the decoder network to generate converted mel-
spectrogram conditioned on the input speaker embedding.

  * •

Discriminator network consists of two sub-modules: a real-fake classifier
trained to discriminate real and converted mel-spectrogram, and a classifier
to predict the speaker identity of input mel-spectrogram.

###  2.1 Training procedure

The training data consists of utterances from 20 professional announcers from
ATR [20] dataset A-set and 20 non-expert speakers from ATR dataset C-set. All
the utterances are preprocessed by resampling to 24 kHz, removing leading and
trailing silence, and combining to 5-second chunks. There is total 22.234
utterances, in which 500 utterances are used for validation. The 80-band log-
mel spectrogram with band limited frequency range (0 to 8 kHz) is extracted
using short-time Fourier transform. The window length and frame shift are set
to 1024 and 256 respectively. We follow training strategy as described in [16]
with the same objective functions and hyper-parameters. The voice conversion
model is trained for 50 epochs with batch size of 48 using 2 Nvidia RTX3090
GPUs. The training process takes approximately 1 day to finish.

###  2.2 Visualization of speaker embedding

Figure 2: 2D visualization of first and second principal components of speaker
embedding. The red dots and eclipse shades denote the centroid and covariance
of each cluster.

The speaker embedding is a vector that encodes speaker individuality conveyed
in the input mel spectrogram into a compact vector. By analyzing the speaker
embedding using principal component analysis (PCA), we can factorize out the
dominant features of speaker individuality. Figure 2 plots the first and
second principal components of the speaker embedding after training. The first
component corresponds to the gender of the speakers, while the second
component corresponds to the voice type, which is non-expert voice or
professional announcer voice. This result suggests that the style of voice,
i.e., non-expert style or announcer style, can be controlled independently
from other voice attributes.

##  3 Experimental settings

###  3.1 Approach 1

From the study of Kobayashi et al. [13], it has been shown that speech from
professional announcers is perceptually clearer and easier to hear than that
of non-expert voice even in noisy environments. This phenomenon is mainly
originated from the differences of voice quality, clarity, projection, ring or
resonant. However, it is unclear that whether the announcer-adapted voice can
still possess this property of natural announcer voice. To investigate this
point, we transform speaker individuality of speaker M105 (non-expert in ATR
C-set) to that of target speaker MAU (professional announcer in ATR A-set)
using voice conversion model. The speaker embedding of MAU speaker is used to
synthesize the converted stimuli.

###  3.2 Approach 2

As can been clearly seen in Fig. 2, the second principal component of the
speaker embedding captures the difference between non-expert speaker and
professional announcer. This behavior is advantageous as the non-expert voice
can be converted to have the voice style of professional announcers only by
changing the second principal components of speaker embedding. Since the
professional announcers are more intelligible even in noisy environment, it is
expected that the second principal component can be used to increase the
intelligibility of non-expert voice. To clarify this point, we propose to
replace the second principal component of M105 speaker embeddings with the
average value calculated from the second principal component of all male
professional announcers. Then, the obtained speaker embedding is used to
synthesize converted stimuli. Different from Approach 1, the speaker
individuality of source speaker is partly changed as other principal
components of speaker embeddings are preserved.

###  3.3 Experiment speech stimuli

We select 520 Japanese words, each may contain 1 to 4 mora, from the ATR
Digital Voice Database A-set (ATR-A) and ATR Digital Voice Database C-set
(ATR-C) as the clean stimuli for target and source speaker. All speech
waveforms are resampled to 16 kHz sampling rate. There are 4 types of speech
stimuli in the experiments, which are denoted as follows:

  * •

Non-expert: Natural speech of non-professional speaker, which is collected
from speaker M105 in ATR-C set.

  * •

Announcer: Natural speech of professional announcer, which is collected from
speaker MAU in ATR-A set.

  * •

VC-1: converted speech from speaker M105 to speaker MAU (Approach 1) by voice
conversion model.

  * •

VC-2: converted speech by shifting the second principal components of speaker
embedding of speaker M105 (Approach 2) by voice conversion model.

To create the noisy stimuli, we mask the clean stimuli with pink noise at 5
different SNR levels: -9dB, -6dB, -3dB, 0dB, and $\infty$ (no noise). We
calculate the root-mean-square of the speech signal only in the voice region
and scale the noise signal to match with the desired SNR level. The voice
region is derived from the text transcription of ATR dataset. To avoid the
effect of different onset and offset timing between speech stimuli, duration
of each stimulus is adjusted to contain the same 200ms of leading noise and
200ms of trailing noise. In addition, speech stimuli are gated with two raised
cosine onset and offset windows of 40-ms to avoid overshoot distortion.

##  4 Evaluations

###  4.1 Objective evaluation

Figure 3: Average vowel space of non-expert, announcer, VC-1, and VC-2
stimuli, with respective areas: $2.64$, $3.60$, $2.91$ and $2.95$ ($\times
10^{5}$ $\textrm{Hz}^{2}$).

Two objective metrics are used to evaluate the intelligibility of converted
speech: 1) average vowel space, and 2) extended short-time objective
intelligibility (eSTOI) [21].

####  4.1.1 Vowel space analysis

Several studies have reported that the expansion of vowel space corresponds to
an increase of speech intelligibility [22, 23]. Based on this information, we
compare the areas of the average vowel spaces derived from different types of
stimuli. The formant frequencies of 5 Japanese vowels (/a/, /e/, /i/, /o/, and
/u/) are extracted using Praat tool [24]. The locations of vowels in each
utterance are determined using the provided text transcription. Then the
average frequency of first and second formants of vowels are calculated across
all speech utterances. The vowel space is defined as the smallest polygon that
fits all the vowels. As can be seen from Fig. 3 the average vowel space of
professional announcer exhibits the largest area ($3.6\times 10^{5}$
$\textrm{Hz}^{2}$). This result aligns with the above assumption that the
speech intelligibility increases with a larger vowel space. Interestingly, the
converted voices from non-expert speaker (VC-1 and VC-2) show an expansion of
vowel space from the non-expert vowel space, from $2.64\times 10^{5}$
$\textrm{Hz}^{2}$ to $2.91\times 10^{5}$ $\textrm{Hz}^{2}$ and $2.95\times
10^{5}$ $\textrm{Hz}^{2}$ respectively. Moreover, the vowel space of VC-1 and
VC-2 appears to have similar shape to that of professional announcer. This
result indicates that the converted voice might have better intelligibility
than non-expert voice in noisy environment.

####  4.1.2 eSTOI measurements

To objectively measure the intelligibility of speech in noise, we calculate
the eSTOI of the speech stimuli at 4 SNR levels: -9dB, -6dB, -3dB and 0dB
using pySTOI python package444https://github.com/mpariente/pystoi. The clean
speech is used as the reference signal for eSTOI calculation. As can be seen
from Fig. 4, the announcer voice can resist to noisy environment better than
non-expert voice as expected. Moreover, VC-1 and VC-2 stimuli also show a
comparable performance to that of the announcer voice. This result confirms
the effectiveness of our proposed method.

Figure 4: Average eSTOI score of non-expert, announcer, VC-1, and VC-2 stimuli
across all utterances. The horizontal axis shows the SNR in dB. The eSTOI
score is in range [0, 1], with higher score indicates better intelligibility.
Figure 5: Average word correct answer rates and SE across participants. The
horizontal axis shows the SNR in dB and $\infty$ denotes that the speech
stimuli is presented without noise. Table 1: p-value of post-hoc pairwise
comparison between 4 types of stimuli using Tukey’s HSD test. Values that
reject the null hypothesis are in bold.

SNR | Non-expert - VC1 | Non-expert - VC2 | Non-expert - Announcer | Announcer - VC1 | Announcer - VC2 | VC1 - VC2
---|---|---|---|---|---|---
$-9$ dB | 0.001 | 0.001 | 0.003 | 0.77 | 0.77 | 0.90
$-6$ dB | 0.064 | 0.024 | 0.0021 | 0.47 | 0.71 | 0.9
$-3$ dB | 0.1 | 0.031 | 0.001 | 0.23 | 0.52 | 0.9
$0$ dB | 0.08 | 0.03 | 0.001 | 0.28 | 0.55 | 0.9
$\infty$ | 0.016 | 0.32 | 0.9 | 0.051 | 0.59 | 0.45

###  4.2 Subjective evaluations

We conducted listening tests to compare the intelligibility of 4 types of
speech stimuli. The listening test is conducted in sound-proof room to avoid
any interference. Speech stimuli are presented via a D/A converter (RME,
Fireface UCX), a headphone amplifier (STAX, SRM-1/MK2), and an electro-static
headphone (STAX, SR-404). The sound pressure level is fixed at 60dB LAeq,
measured on the clean speech stimuli by a sound level meter (Brüel&Kjær, Types
2250) via an artificial ear (Brüel&Kjær, Type 4153). Each participant listens
to a set of 300 different random words, which are equally distributed into 5
SNR levels and 4 types of speech stimuli. The stimulus is presented only once
for each trial and the order of presented speech stimuli is randomized for
each participant. The duration of the whole listening test is approximately 40
minutes, which is divided into 4 sections with 2 minutes break between each
section. There were 10 native Japanese participants, whose age ranging from 23
to 29 years old, joined our listening test.

Before each test, the listeners are provided with instruction and some sample
stimuli to get used to the sound level. Figure 5 reports the average word
correction rate across participants. One-way ANOVA test show that there were
statistical differences between 4 types of stimuli. A post-hoc pairwise
analysis using Tukey HSD test ($p<0.05$) was carried out to determine
statistical differences between pair of stimuli types in different SNR
conditions.

The results shown in Table 1 indicate that VC-2 and announcer stimuli are
significant different from non-expert stimuli in noisy conditions
($\textrm{SNR}\leq 0\textrm{dB}$). In addition, no statistical difference
between all 4 types of stimuli is found in clean condition
($\textrm{SNR}=\infty$). These results suggest that VC-2 is more effective
than VC-1 for enhancing speech intelligibility in noisy condition. The
possible reason for this difference might correspond to the increase amount of
distortion when speaker individuality is completely changed. However, further
analysis must be carried out to clarify this point.

##  5 Conclusions

The present work has proposed a speech intelligibility enhancement in noisy
environment using voice conversion technique. The results from objective
measurements and subjective evaluation confirm that adapting to announcer
voice can increase the intelligibility of non-expert speaker. By analyzing the
PCA of speaker embedding, it has been discovered that the announcer-speaking
style is an independent features of speaker individuality. By modifying the
second principal component of speaker embedding, we can manually control the
amount of announcer-speaking style, hence increasing the intelligibility of
speech in noisy environment without completely change the speaker
individuality. Statistical analysis shows that modifying the second principal
components yields the highest performance. Beside using announcer voice as the
target, the proposed method can be applied to mimic Lombard speech and Clear
speech, which are also speaking styles for increasing speech intelligibility
in noisy condition. For the future work, the proposed method can be extended
to generate converted speech adaptively to the noise condition in order to
further improve the speech intelligibility in noise.

##  6 Acknowledgements

This work was supported by the SCOPE Program of Ministry of Internal Affairs
and Communications (Grant number: 201605002).

## References

  * [1] W. Kleijn, J. B. Crespo, R. C. Hendriks, P. N. Petkov, B. Sauert, and P. Vary, “Optimizing speech intelligibility in a noisy environment: A unified view,” _IEEE Signal Processing Magazine_ , vol. 32, pp. 43–54, 2015.
  * [2] C. H. Taal and J. H. Jensen, “Sii-based speech preprocessing for intelligibility improvement in noise,” in _INTERSPEECH_ , 2013.
  * [3] C. H. Taal, R. C. Hendriks, and R. Heusdens, “Speech energy redistribution for intelligibility improvement in noise based on a perceptual distortion measure,” _Comput. Speech Lang._ , vol. 28, pp. 858–872, 2014.
  * [4] T.-C. Zorila and Y. Stylianou, “On spectral and time domain energy reallocation for speech-in-noise intelligibility enhancement,” in _INTERSPEECH_ , 2014.
  * [5] T.-C. Zorila, V. Kandia, and Y. Stylianou, “Speech-in-noise intelligibility improvement based on spectral shaping and dynamic range compression,” in _INTERSPEECH_ , 2012.
  * [6] C. Chermaz and S. King, “A sound engineering approach to near end listening enhancement,” in _INTERSPEECH_ , 2020.
  * [7] A. Amano-Kusumoto, T. Arai, K. Kinoshita, N. Hodoshima, and N. Vaughan, “Modulation enhancement of speech by a pre-processing algorithm for improving intelligibility in reverberant environments,” _Speech Commun._ , vol. 45, pp. 101–113, 2005.
  * [8] T. V. Ngo, R. Kubo, and M. Akagi, “Increasing speech intelligibility and naturalness in noise based on concepts of modulation spectrum and modulation transfer function,” _Speech Communication_ , 2021.
  * [9] Y. Tang and M. Cooke, “Subjective and objective evaluation of speech intelligibility enhancement under constant energy and duration constraints,” in _INTERSPEECH_ , 2011.
  * [10] V. Aubanel and M. Cooke, “Information-preserving temporal reallocation of speech in the presence of fluctuating maskers,” in _INTERSPEECH_ , 2013.
  * [11] H. Noh and D.-H. Lee, “How does speaking clearly influence acoustic measures? a speech clarity study using long-term average speech spectra in korean language,” _Clinical and Experimental Otorhinolaryngology_ , vol. 5, pp. 68 – 73, 2012.
  * [12] C. Kashimada, K. Ogita, T. Ishikawa, H. Hasegawa, and M. Ayama, “Effects of voice training on subjective evaluation of voice quality,” _The Journal of The Institute of Image Information and Television Engineers_ , vol. 63, pp. 1818–1823, 2009.
  * [13] M. Kobayashi and M. Akagi, “Intelligibility of announcer’s speech in noisy environments,” _IEICE Technical Report_ , vol. 119, pp. 95–99, 2020.
  * [14] T. V. Ho and M. Akagi, “Non-parallel voice conversion with controllable speaker individuality using variational autoencoder,” _2019 Asia-Pacific Signal and Information Processing Association Annual Summit and Conference (APSIPA ASC)_ , pp. 106–111, 2019.
  * [15] Y. Choi, Y. Uh, J. Yoo, and J.-W. Ha, “Stargan v2: Diverse image synthesis for multiple domains,” _2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR)_ , pp. 8185–8194, 2020.
  * [16] Y. A. Li, A. Zare, and N. Mesgarani, “Starganv2-vc: A diverse, unsupervised, non-parallel framework for natural-sounding voice conversion,” in _INTERSPEECH_ , 2021, pp. 1349–1353.
  * [17] K. He, X. Zhang, S. Ren, and J. Sun, “Deep residual learning for image recognition,” _2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR)_ , pp. 770–778, 2016.
  * [18] S. Watanabe, T. Hori, S. Karita, T. Hayashi, J. Nishitoba, Y. Unno, N. Yalta, J. Heymann, M. Wiesner, N. Chen, A. Renduchintala, and T. Ochiai, “Espnet: End-to-end speech processing toolkit,” _ArXiv_ , vol. abs/1804.00015, 2018\.
  * [19] R. Yamamoto, E. Song, and J.-M. Kim, “Parallel wavegan: A fast waveform generation model based on generative adversarial networks with multi-resolution spectrogram,” _ICASSP 2020 - 2020 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_ , pp. 6199–6203, 2020.
  * [20] A. T. R. I. International, “Digital voice database,” http://www.atr-p.com/.
  * [21] C. H. Taal, R. C. Hendriks, R. Heusdens, and J. R. Jensen, “A short-time objective intelligibility measure for time-frequency weighted noisy speech,” _2010 IEEE International Conference on Acoustics, Speech and Signal Processing_ , pp. 4214–4217, 2010.
  * [22] A. Amano-Kusumoto and J.-P. Hosom, “A review of research on speech intelligibility and correlations with acoustic features,” in _Center for Spoken Language Understanding, Oregon Health and Science University (Technical Report CSLU-011-001)_ , 2011.
  * [23] E. Godoy, M. Koutsogiannaki, and Y. Stylianou, “Approaching speech intelligibility enhancement with inspiration from lombard and clear speaking styles,” _Comput. Speech Lang._ , vol. 28, pp. 629–647, 2014.
  * [24] P. Boersma, “Praat, a system for doing phonetics by computer,” 2002.


# Data-iterative Optimization Score Model for Stable Ultra-Sparse-View CT
Reconstruction

Weiwen Wu  Yanyang Wang W. Wu, Y. Wang are with the Department of Biomedical
Engineering, Sun-Yat-sen University, Shenzhen Campus, Shenzhen, China. The
contribution of W. Wu and Y. Wang are equal. Wu is the corresponding author
(e-mail: wuweiw7@mail.sysu.edu.cn).

###### Abstract

Score-based generative models (SGMs) have gained prominence in sparse-view CT
reconstruction for their precise sampling of complex distributions. In SGM-
based reconstruction, data consistency in the score-based diffusion model
ensures close adherence of generated samples to observed data distribution,
crucial for improving image quality. Shortcomings in data consistency
characterization manifest in three aspects. Firstly, data from the
optimization process can lead to artifacts in reconstructed images. Secondly,
it often neglects that the generation model and original data constraints are
independently completed, fragmenting unity. Thirdly, it predominantly focuses
on constraining intermediate results in the inverse sampling process, rather
than ideal real images. Thus, we propose an iterative optimization data
scoring model. This paper introduces the data-iterative optimization score-
based model (DOSM), integrating innovative data consistency into the
Stochastic Differential Equation, a valuable constraint for ultra-sparse-view
CT reconstruction. The novelty of this data consistency element lies in its
sole reliance on original measurement data to confine generation outcomes,
effectively balancing measurement data and generative model constraints.
Additionally, we pioneer an inference strategy that traces back from current
iteration results to ideal truth, enhancing reconstruction stability. We
leverage conventional iteration techniques to optimize DOSM updates.
Quantitative and qualitative results from 23 views of numerical and clinical
cardiac datasets demonstrate DOSM’s superiority over other methods.
Remarkably, even with 10 views, our method achieves excellent performance.

Computed Tomography, image reconstruction, score-based generative model,
sparse-view, data consistency

##  1 Introduction

Computed Tomography (CT) has garnered extensive utilization in both medical
diagnosis due to its ability to provide practical and precise diagnostic
outcomes [1]. Sparse-view CT scanning emerges as a promising strategy to
decrease radiation dose, wherein only a fraction of projection data is
required for image reconstruction[2] [3]. However, the limited number of
measurement views leads to a reduced acquisition [4] of prior information
within imaged object, consequently yielding a deterioration in the image
quality [5]. Moreover, the diminished rank of the measurement matrix
introduces an augmented array of potential solutions[6], impeding convergence
and culminating in a heightened prevalence of uncertainties and inaccuracies
[7] [8].

Traditional reconstruction methods like filtered back-projection (FBP) yield
unsatisfactory results, marked by streaking artifacts and poor image quality
[2]. The advent of artificial intelligence has fostered the creation of
advanced deep-learning techniques for enhancing sparse-view CT reconstruction.
Notable examples of such techniques encompass FBP-ConvNet [9], Densenet
Deconvolution Network [10], residual encoder-decoder convolutional neural
network [11], generative adversarial network (CGAN) [12], and multi-domain
integrative swin transformer [13], among others. These deep-learning
approaches primarily hinge on supervised training with paired images to yield
superior outcomes.

Score-based generative models (SGMs) have risen to prominence due to their
exceptional capability to precisely sample from intricate distributions [14]
[15]. For example, denoising diffusion probabilistic models (DDPM)[14],
denoising diffusion restoration models (DDRM) [16], Denoising diffusion
implicit models (DDIM) [17], Stochastic Differential Equations (SDE) [15],
etc. Typical SGMs follow a two-stage process consisting of a forward stage
perturbing data to noise and a reverse stage converting noise back to data.
The reverse process is usually achieved by parameterized deep neural networks,
e.g., an SGM, which can be optimized with the training data. The power of SGM
seamlessly extends to various applications within the field of medical image
reconstruction [18] [19], particularly in scenarios such as sparse-view CT
reconstruction with remarkable achievements [20] [21]. Data consistency in a
score-based score-based model is to ensure that the generated samples or
images adhere closely to the observed data distribution, which plays an
important role for improving reconstructed image quality [22]. Data
consistency acts as a constraint during the sampling process to align the
generated samples with the actual data, improving the quality and realism of
the generated samples [23]. For a typical score-based generative model, the
data consistency is enforced [24] by iteratively adjusting the generated
samples based on the gradient of the log-likelihood of the observed data [15].
This process aims to reduce the discrepancy between the generated samples and
the real measurement data, resulting in better alignment and more accurate
modeling of the underlying data distribution [25]. By maintaining data
consistency [26], the model becomes more capable of generating samples that
are representative of the original data.

These existing SGMs for sparse-view CT reconstruction mainly focus on
performing the diffusion within sinogram-based and image-domain. In respect to
the sinogram-based SGM, the data consistency is formatted by utilizing the
generated full sinogram and the sampling location of the original sparse-view
data. For instance, a notable advancement involves the integration of a fully
unsupervised score-based generative model into the sinogram domain,
effectively enhancing sparse-view CT reconstruction capabilities [27].
Furthermore, a patch-based denoising diffusion probabilistic model was
tailored and developed [20]. In the image-domain context, the conventional
data consistency hinges on a sparse-sampling mask [18]. The data consistency
process involves employing partial projections derived from the reconstructed
image, wherein $\bm{y}_{t}=\bm{Ax}_{t}$, signifying the system matrix $\bm{A}$
and reconstructed image $\bm{x}_{t}$ at time step $t$. However, owing to the
inherent imperfections within the interim $\bm{y}_{t}$, inaccuracies manifest
in the estimated data
$\bm{y_{t}^{{}^{\prime}}}=mask*\bm{y}+(1-mask)*\bm{y}_{t}$ ($\bm{y}$ is the
sparse-view data), which in turn contributes to data consistency
discrepancies, subsequently introducing secondary artifacts in the
reconstruction process. Moreover, past investigations frequently updated the
data consistency and generative model as separate entities, disregarding the
intrinsic link between these modules [22][23]. These approaches not only
hinder the model’s ability to converge to a stable reconstruction outcome but
can also exacerbate the inherent instability within the reconstruction
process.

Figure 1: A comprehensive evaluation contrasting our proposed methodology with
contemporary state-of-the-art reconstruction techniques using 23 views’
measurements.

Overcoming the lack of precise data consistency constraints is a significant
challenge when striving for high-quality image reconstructions, especially in
the context of ultra-sparse-view CT reconstruction [5]. In this scenario, the
sparse-view sampling mask can introduce projection errors, leading to the
presence of erroneous artifacts [28]. To mitigate inaccuracies in projection,
a more effective approach involves utilizing measurement data solely for
establishing data consistency, defined as $\|\bm{y}-\bm{Ax}\|_{2}^{2}$. This
approach ensures the generation of reliable reconstruction outcomes.
Importantly, the focus of data consistency should be directed towards
$\bm{x}_{0}$ rather than $\bm{x}_{t}$ at time $t$. In simpler terms, when
$\bm{y}\neq\bm{Ax}_{t}$, it becomes insufficient to drive the generation
process towards satisfactory results. A viable solution is to enforce data
consistency as $\|\bm{y}-\bm{Ax}_{0}(\bm{x}_{t})\|_{2}^{2}$. However, a
significant challenge arises in establishing an explicit relationship between
$\bm{x}_{t}$ and $\bm{x}_{0}$. This is due to the fact that $p(\bm{x}_{0})$
follows a posterior distribution, which remains unknown during the inference
stage [22] [23]. Furthermore, it is imperative to devise a viable approach for
augmenting the coherence between the data consistency and generative update.
This strategy’s purpose should encompass the simultaneous constraint of both
data consistency and generative outcomes, effectively capturing the intricate
correlations inherent to these two modules. To fully harness the potential of
data consistency within the score-based model for tackling inverse problems,
pursuing a strategy that explores this avenue comprehensively is both valuable
and promising.

In this paper, we introduce an innovative data consistency component to the
SDE, which serves as a valuable constraint in the context of ultra-sparse-view
CT reconstruction. We term this framework the Data-iterative Optimization
Score-based Model (DOSM). The distinctive aspect of this new data consistency
element lies in its utilization of solely the original measurement data to
confine the outcomes of generation, effectively sidestepping any undesirable
structures. Furthermore, we establish a pioneering inference strategy that
traces back from $\bm{x}_{t}$ to $\bm{\hat{x}}_{0}$ using a solid theoretical
foundation at any time $t$. In contrast to prevailing state-of-the-art
techniques, our approach exhibits a noteworthy enhancement in both image
quality and quantitative outcomes (PSNR and SSIM), as indicated in Fig. 1. The
significant contributions made by this study can be succinctly outlined as
follows:

  * 1)

We ingeniously devise a novel data consistency term that harmoniously balances
the constraints from measurement data and the generative model. This seamless
integration facilitates the capture of intricate correlations and
significantly enhances their alignment. By fusing this innovative data
consistency approach with an SDE, we introduce an entirely new score-based
model for ultra-sparse-view CT reconstruction.

  * 2)

We develop an innovative and adaptable inference strategy that elegantly
traces from $\bm{x}_{t}$ back to its initial counterpart $\bm{\hat{x}}_{0}$.
Remarkably, as the number of iterations decreases, the accuracy of the
estimated $\bm{\hat{x}}_{0}$ improves significantly. The estimated
$\bm{\hat{x}}_{0}$ is introduced into the formulated data consistency for
formulating a novel unified reconstruction model.

  * 3)

We employ a suitable solving algorithm to optimize the established score-based
Model. To ensure solution stability, we incorporate the wisdom of the
conventional iteration reconstruction technique to optimize the update of data
consistency.

  * 4)

The efficacy and advantages of our proposed approach are meticulously examined
and evaluated. We validate our method using extensive sparse-view CT datasets,
encompassing both numerical simulations and real clinical cardiac data,
showcasing its superior reconstruction performance in comparison to state-of-
the-art alternatives.

##  2 Related Works

###  2.1  Stochastic Differential Equation Models

The fractional derivative-based generative model, known as SDE, stands as a
quintessential example within the realm of generative models. It has a
widespread utilization application in the domain of medical imaging for
addressing inverse problems. SDE gradually inject the noise onto the clean
data, shaping the training process of neural networks, which is termed as the
forward process. Subsequently, it delineates a perturbation process that
inverses the initial forward step, progressively denoising the noisy samples.
This iterative denoising process is referred to as the inverse process,
embodying the image reconstruction phase in image reconstruction. Notably,
within this framework, Y. Song et al[15] introduced a fractional score-based
Model to define the forward process. Throughout the entirety of the diffusion
process, data $\bm{x}$ is represented as $\bm{x}(t)=\bm{x}_{t}$, with
$t\in[0,1]$. Here, $\bm{x}_{0}$ signifies the distribution of the original
training data, while $\bm{x}_{T}$ approximates a spherical Gaussian
distribution and $T$ represents the number of diffusion step in forward
process. The formulation of the forward process is articulated as follows
[15]:

$d\bm{x}=f(\bm{x},t)dt+g(t)d\bm{w},$ (1)

where $f$ delineates the linear drift function, and $g(t)$ characterizes the
scalar diffusion coefficient. The vector $\bm{w}$ follows the standard
Brownian motion. The inverse process of Eq. (1) can be formulated as follows:

$\displaystyle d\bm{x}$
$\displaystyle=\left[{f(\bm{x},t)-g{{(t)}^{2}}{\nabla_{\bm{x_{t}}}}\log{p_{t}}(\bm{x})}\right]dt+g(t)d\bar{\bm{w}},$
(2)

where $dt$ represents an infinitesimal time step. In order to generate
reconstructed images through the inverse process of SDE, a time-dependent
score function ${\nabla{\bm{x_{t}}}}\log{p_{t}}(\bm{x})$ is required [18].
This scoring function estimator ${\bm{s_{\theta}}}$ can be obtained by means
of denoising scoring matching [29]:

$\displaystyle\begin{aligned}
\min_{\bm{\theta}}\mathbb{E}_{\bm{x}_{t}\mid\bm{x}_{0},\bm{x}_{0}}\left[\left\|\bm{s}_{\bm{\theta}}(\bm{x}_{t},t)-\nabla_{\bm{x}_{t}}\log
p\left(\bm{x}_{t}\mid\bm{x}_{0}\right)\right\|_{2}^{2}\right].\end{aligned}$
(3)

We estimate the score function $\theta^{*}$ through training a neural network.
The fractional scoring estimator ${\bm{s_{\theta}}}(\bm{x}_{t},t)$ can be
employed as a replacement for the scoring function in Eq. (2). Different
combinations of $f$ and $g$ can generate different types of SDE. For the
typical variation exploding (VE) -SDE, $f$ and $g$ are set as follows [15]:

$f=0,~{}~{}~{}g=\sqrt{\frac{{d\left[{{\sigma^{2}}(t)}\right]}}{{dt}}},$ (4)

where $\sigma{(t)}$ signifies the time-varying escalating scale function for
noise.

###  2.2 Score-based Sparse-View CT Reconstruction

Score-based Model shows good performance in sparse-view CT reconstruction [1].
Relying on the powerful generation ability of the score-based model, there are
some typical works including [27], [22] and etc. The sparse-view
reconstruction can be regarded as the classical inverse problem [18], which
aims at recovering an unknown signal from a set of observed signals.
Specifically, if $\bm{x}$ is an original signal and $\bm{y=Ax+\epsilon}$ is a
noisy observation given by a set of linear measurements. The sampling
procedure $\bm{A}$ is the radon transform. The sparse-view reconstruction is
equivalent to recovering the signal $\bm{x}$ from the measurement $\bm{y}$.
The $\bm{\varepsilon}$ denotes the measured noise. A common practice is to use
the least squares of linear models to solve image reconstruction [30]. Due to
the uncertainty of the inverse problem (where $\bm{A}^{T}\bm{A}$ is singular
or ill-conditioned, and the data is noisy), the standard practice is to
introduce regularization terms for reconstruction:

$\displaystyle\bm{x}=\underset{\bm{x}}{\operatorname{argmin}}\frac{1}{2}\|\bm{y}-\bm{A}\bm{x}\|_{2}^{2}+\lambda
R(\bm{x}),$ (5)

where $R(.)$ is a suitable regularization and it can be a sparse
transformation regularization.

###  2.3 Data Consistency Construction

The good reconstruction performance of the diffusion model not only depends on
the generation of the unconditional fractional function
${\nabla_{\bm{x_{t}}}}\log{p_{t}}(\bm{x})$, but also needs a suitable data
consistency strategy to guide the generation process. In solving CT imaging
problems, the data consistency strategy was first proposed by Y. Song et al
[18]. Specifically, the observed projection information, $y$, and the
sinusoidal graph generated at each step are completed.

$\displaystyle\bm{x}_{t}^{\prime}=\bm{A}^{T}[\bm{y}*mask+A\bm{x}_{t}*(1-mask)],$
(6)

where $\bm{x}_{t}^{\prime}$ is the result of adding measurement data
constraints. $mask$ represents to extract sparse-view data from the full-view
projection.

The other is a data consistency strategy based on posterior probability
sampling, which provides an easy-to-handle approximation for $p(\bm{y|x_{t}})$
in order to better solve the reconstructed inverse problem using the scoring
function [23].

$\displaystyle\nabla_{\bm{x}_{t}}\log
p\left(\bm{y}\mid\bm{x}_{t}\right)\simeq-\frac{1}{\sigma^{2}}\nabla_{\bm{x}_{t}}\left\|\bm{y}-{A}\left(\hat{\bm{x}}_{0}\left(\bm{x}_{t}\right)\right)\right\|_{2}^{2},$
(7)

where $\sigma$ is the measurement noise factor. $\hat{\bm{x}}_{0}$ is the
posterior mean calculated by the diffusion model and is a function of $x_{t}$.

To further improve the constraint generation process, Chung et al. proposed a
correction term for manifold constraints [22]. Hence, a constraint based on
manifold gradients is formulated to bound the gradient of the measurement term
within the confines of the data manifold.

$\displaystyle\nabla_{\bm{x}_{t}}\log
p\left(\bm{y}\mid\bm{x}_{t}\right)\simeq-\alpha\frac{\partial}{\partial\bm{x}_{t}}\left\|\bm{W}\left(\bm{y}-\bm{A}\hat{\bm{x}}_{0}\left(\bm{x}_{t}\right)\right)\right\|_{2}^{2},$
(8)

where $\alpha$ corresponds to the step length, and $\bm{W}$ denotes a
weighting factor. The authors in [22] additionally make use of the projected
estimate obtained from intermediate updates.

The first and third data consistency terms, due to inherent imperfections in
estimated data, introduce secondary artifacts during reconstruction. The
second data consistency term prioritizes derivative descent, making it
sensitive to perturbations and destabilizing the model. Additionally, these
data consistency strategies consider multiple task mappings but overlook CT-
specific physical priors. Simple FBP reconstruction limits the generative
model’s exploration of its potential. Importantly, these three data
consistency terms frequently update data and the generative model separately,
ignoring their intrinsic connection. This not only obstructs stable
convergence but can worsen inherent reconstruction instability.

##  3 Methods

###  3.1 Mathematical Model

Take into consideration a broad-ranging forward model for computed tomography
reconstruction, the ideal reconstruction solution of $\bm{{x_{0}}}$ should
satisfy:

$\bm{A}\bm{{x_{0}}}+\bm{\varepsilon}=\bm{y},$ (9)

Let’s remember that $\bm{A}$ embodies the scanner’s operation, referred to as
forward projection. Each individual row within $\bm{A}$ encapsulates the
coefficients of an equation that correlates with a singular ray. It describes
how the pixels are combined into ray sums. Here, $h_{ij}$ denotes an element
located at the $(i,j)-th$ position, which $i$ takes values from 1 to $I$, and
$j$ spans from 1 to $J$. $I$ corresponds to the total pixel count within the
reconstructed image, while $J$ is the number of X-ray pathways. In the score-
based model solving the CT reconstruction problem, the approximate fractional
function ${\bm{s_{\theta}}}$ obtained from the forward process training can
simulate the reverse denoising process in the reconstruction process [18].
${\bm{s_{\theta}}}$ is approximated as
${\nabla_{\bm{x_{t}}}}\log{p_{t}}(\bm{x})$. Since there is no explicit
correspondence [23] between $\bm{y}$ and $\bm{x}$, we cannot solve it directly
using the neural network of the score-based Model. A common practice [31] is
to solve the posterior probability according to Bayesian theory
$P(\bm{x|y})=p(\bm{x})p(\bm{y|x})/p(\bm{y})$. Specifically, the reconstructed
image $\bm{x}_{t}$ at the current time can be regarded as the inverse solution
through the observation data $\bm{y}$:

$\nabla_{\bm{x}_{t}}log\left(p\left(\bm{x}_{t}|\bm{y}\right)\right)=\nabla_{x_{t}}log\left(p\left(\bm{x}_{t}\right)\right)+\nabla_{\bm{x}_{t}}log\left(p\left(\bm{y}/\bm{x}_{t}\right)\right).$
(10)

Considering that $y$ and $x_{t}$, $x_{0}$ have no explicit correspondence,
they cannot be solved directly by Eq. (10). So we construct an approximate
relationship:

$\nabla_{\bm{x}_{t}}log\left(p\left(\bm{y}|\bm{x}_{t}\right)\right)\approx\nabla_{\bm{x}_{t}}log\left(p\left(\bm{y}|\hat{\bm{x}_{0}}\left(\bm{x}_{t}\right)\right)\right).$
(11)

Utilizing the diffusion model for addressing sparse-view CT reconstruction,
the mathematical formulation of the model can be expressed as follows:

$\min_{\bm{x}_{t}}M(\bm{x}_{t}),~{}~{}s.t.,\left\|\bm{A}\hat{\bm{x}_{0}}(\bm{x}_{t})-\bm{y}\right\|_{2}^{2}\leq{\tau}$
(12)

Here, the term $M\left(\bm{x}_{t}\right)$ signifies the solver configured
through Stochastic Differential Equations (SDE) at any time t. Furthermore,
Eq. (12) encompasses the solution for a condition-constrained diffusion model.
In essence, it can be transformed into an unconstrained optimization problem,
characterized by the subsequent expression:

$\bm{x}_{t-1}=\min_{\bm{x}_{t}}\frac{1}{2}\left\|\bm{A}\hat{\bm{x}_{0}}(\bm{x}_{t})-\bm{y}\right\|_{2}^{2}+\eta_{t}M(\bm{x}_{t})$
(13)

where $\eta_{t}$ embodies the equilibrium coefficient that finely tunes the
accentuation of the data consistency proportion at the temporal instant
denoted as $t$.

As introduced above, the posterior probability sampling of the score-based
model depends on $\bm{{x}}_{0}$, i.e., the underlying ground truth
reconstructed image. In this case, effective estimation of $\bm{{x}}_{0}$
becomes a critical step to conduct the posterior sampling. The estimation is
usually performed with the latest intermediate result of $\bm{x}_{t}$. For
example,
$\bm{\hat{x}}_{0}(\bm{x}_{t})=\bm{x}_{t}+\sigma_{t}^{2}\nabla_{\bm{x}_{t}}\log
p\left(\bm{x}_{t}\right)$ is applied in [32, 22] to achieve the denoised
result by computing the posterior expectation according to Tweedie’s formula
[33]. We argue that although this estimation enjoys simplicity, it’s sole
dependency on $\bm{x}_{t}$ may probably lead to an unstable and biased
estimation, resulting in a sub-optimal solution to the posterior sampling. To
obtain a more reliable estimation, this paper innovatively proposes an
effective method with $\bm{x}_{t}$s from multi-channels, as shown in Fig. 2.
Specifically, once the diffusion model is trained, the generation process will
be conducted in multiple, say $N$, channels starting from independently
sampled noise respectively. At each intermediate step, $\bm{{x}}_{0}$ can be
estimated by a weighted combination of the intermediate results from the $N$
channels as follows:

$\bm{\hat{x}}_{0}=\sum_{n=1}^{N}{w_{t}^{n}\bm{x}_{t}^{n}},$ (14)

where $\bm{x}_{t}^{n}$ indicates the intermediate result of the $n$-th channel
at the current $t$-th step of the reverse diffusion process and $w_{t}^{n}$
denotes the weight of $n^{th}$ channel of time $t$. While enjoying simplicity,
the estimation of $\bm{\hat{x}}_{0}$ in Eq. (14) admits several desirable
properties as follows.

Proposition: The starting pure noise images, $\bm{x}_{T}^{n}$ with
$n\in[1,\dots,N]$, can be considered as a set of denoised images from the same
ground truth $\bm{{x}}_{0}$ by gradually adding a sequence of independent
noise to $\bm{{x}}_{0}$, i.e.,
$\bm{x}_{t}=\bm{x}_{t-1}+\sqrt{\sigma_{t}^{2}-\sigma_{t-1}^{2}}\bm{z}_{t-1}$,
$\bm{z}_{t}\mathop{\sim}\limits^{i.i.d.}\mathcal{N}(\bm{0},\bm{I})$. Then
$\bm{\hat{x}}_{0}=\sum_{n=1}^{N}{w_{k}\bm{x}_{t}^{n}}$ is an unbiased estimate
of $\bm{{x}}_{0}$ when $N\to+\infty$ and $w_{t}^{n}=\frac{1}{N}$.

Remark: The variance expectation of the proposed estimation
$\bm{\hat{x}}_{0}=\frac{1}{N}\sum_{n=1}^{N}{\bm{x}_{t}^{n}}$ is expected to be
smaller than the variance expectation with single $\bm{x}_{t}^{n}$. Since
$\mathbb{E}\left(\text{Cov}\left(\bm{z}_{t_{1}}^{n_{1}},\bm{z}_{t_{2}}^{n_{2}}\right)\right)=0,\
\forall t_{1}\neq t_{2},\ n_{1}\neq
n_{2},\bm{z}_{t}^{n}\mathop{\sim}\limits^{i.i.d.}\mathcal{N}(\bm{0},\bm{I})$,
so
$\mathbb{E}\left(\text{Var}\left(\frac{1}{N}\sum_{n=1}^{N}{\bm{x}_{t}^{n}}\right)\right)=\frac{1}{N}\mathbb{E}\left(\text{Var}\left(\bm{x}_{t}^{n}\right)\right)$.
This means our proposed can improve the stability of reconstruction. Indeed,
the proposed estimation
$\bm{\hat{x}}_{0}=\frac{1}{N}\sum_{n=1}^{N}{\bm{x}_{t}^{n}}$ admit a reliable
solution and could effectively facilitate the enforcement of data consistency
term during the reconstruction process.

Considering the estimation of $\hat{x_{0}}$, we can establish the optimization
model as following:

$\bm{x}_{t-1}=\min_{\bm{x}_{t}}\frac{1}{2}\left\|\bm{A}\sum_{n=1}^{N}w_{t}^{n}\bm{x}_{t}^{n}-y\right\|_{2}^{2}+\eta_{t}M\left(\sum_{n=1}^{N}w_{t}^{n}\bm{x}_{t}^{n}\right)$
(15)

where $w_{t}^{n}$ represents the weight of $n^{th}$ noisy images from time $t$
and satisfy with $\sum_{n=1}^{N}w_{t}^{n}=1$.

###  3.2 Optimization Procedure

Eq. (15) demonstrates the VE-SDE prior is posed on the $n$ images rather than
a single image $\bm{x}_{t}^{n}$. As evident from our observations, Eq. (15)
unmistakably brings to light the inherent nexus between data consistency and
the diffusion generative model. A viable avenue for resolution emerges with
the introduction of an additional fidelity term aimed at solidifying the
alignment between these two pivotal terminologies. In this pursuit, we
introduce $u_{t}^{n}$ as a substitution for $x_{t}^{n}$, which is treated as a
bridge to connect $\hat{\bm{x}_{0}}$ and $M\left(\bm{x}_{t}\right)$. Actually,
we can simply the Eq. (15) into a series of separable optimization problems.
Consequently, the manifestation of Eq. (15) takes on a transformed
disposition, evolving into an intricately formulated optimization challenge
characterized by the subsequent expression:

$\displaystyle\left\\{\bm{x}_{t-1},\bm{u}_{t-1}\right\\}=$
$\displaystyle\min_{\left\\{\bm{x}_{t},\bm{u}_{t}\right\\}}\frac{1}{2}\left\|\bm{A}\sum_{n=1}^{N}w_{t}^{n}\bm{x}_{t}^{n}-\bm{y}\right\|_{2}^{2}$
(16)
$\displaystyle+\eta_{t}M\left(\sum_{n=1}^{N}w_{t}^{n}\bm{u}_{t}^{n}\right),$
$\displaystyle\text{ s.t. },\bm{x}_{t}^{n}=\bm{u}_{t}^{n},n=1,\ldots,N.$

Eq. (16) poses a problem constrained by the requirement
$\bm{x}_{t}^{n}=\bm{u}_{t}^{n}$, and under certain circumstances, it can be
transformed into an unconstrained optimization problem as below

$\displaystyle\left\\{\bm{x}_{t-1},\bm{u}_{t-1}\right\\}=$
$\displaystyle\min_{\left\\{\bm{x}_{t},\bm{u}_{t}\right\\}}\frac{1}{2}\left\|\bm{A}\sum_{n=1}^{N}w_{t}^{n}\bm{x}_{t}^{n}-\bm{y}\right\|_{2}^{2}$
(17)
$\displaystyle+\eta_{t}M\left(\sum_{n=1}^{N}w_{t}^{n}\bm{u}_{t}^{n}\right)+\frac{\beta_{t}}{2}\sum_{n=1}^{N}\left\|\bm{x}_{t}^{n}-\bm{u}_{t}^{n}\right\|_{2}^{2},$

where $\beta_{t}>0$ is a factor. It’s evident from Eq. (17) that the inclusion
of $\bm{u}_{t}$ serves a dual purpose. On one hand, it enforces the
reconstructed outcomes to align reasonably well with the measurements; on the
other hand, it imposes a constraint to confine the divergence between the
reconstruction outcomes and their generative counterparts within specific
bounds. Notably, the term
$\frac{\beta_{t}}{2}||\bm{x}_{t}-\bm{u}_{t}||_{2}^{2}$ constitutes a vital
component of data consistency. Consequently, the optimization of Eq. (16) can
be viewed as a process bifurcated into two distinct steps:

$\displaystyle\bm{x}_{t-1}=$
$\displaystyle\min_{\bm{x}_{t}}\frac{1}{2}\left\|\bm{A}\sum_{n=1}^{N}w_{t}^{n}\bm{x}_{t}^{n}-\bm{y}\right\|_{2}^{2}+\frac{\beta_{t}}{2}\sum_{n=1}^{N}\left\|\bm{x}_{t}^{n}-\bm{u}_{t}^{n}\right\|_{2}^{2}$
(18) $\displaystyle\bm{u}_{t-1}=$
$\displaystyle\min_{\bm{u}_{t}}\frac{1}{2}\sum_{n=1}^{N}\left\|\bm{x}_{t}^{n}-\bm{u}_{t}^{n}\right\|_{2}^{2}+\eta_{t}M\left(\sum_{n=1}^{N}w_{t}^{n}\bm{u}_{t}^{n}\right).$

Regarding Eq. (18), which introduces a new term to ensure data consistency.
Several iterative methods are available for achieving the best solution.
However, ensuring convergence is difficult while determining the appropriate
iterative step length presents a practical hindrance to optimization efforts.
The gradient-based descent technique is a popular solution avenue, commonly
embraced in recent diffusion models [22] [23]. However, these approaches often
demonstrate susceptibility to noise and perturbations, which is a common
occurrence in real-world scenarios. Thus, we have chosen to employ an advanced
optimization strategy for enhancing the performance of Eq. (18). Note that the
update of $\bm{x}_{t-1}$ comes from n images, we can treat
$\bm{x}_{t-1}=\sum_{n=1}^{N}{w_{t}^{n}\bm{x}_{t}^{n}}$ and obtain:

$\displaystyle\bm{x}_{t-1/2}=\min_{\bm{x}_{t}}\frac{1}{2}\left\|\bm{A}\sum_{n=1}^{N}w_{t}^{n}\bm{x}_{t}^{n}-y\right\|_{2}^{2}.$
(19)

By incorporating the wisdom of simultaneous iteration reconstruction technique
[34], in respect to the Eq. (LABEL:eqIII10) step, we can update it with:

$\displaystyle
x_{t-\frac{1}{2}}^{n}=x_{t}^{n}+DA^{T}E\left(y-A\sum_{n=1}^{N}w_{t}^{n}x_{t}^{n}\right),\quad
n=1,\ldots,N$ (20)

Here, $D$ and $E$ are diagonal matrices encompassing the reciprocal of the
column and row sums of the system matrix. Specifically, $d_{jj}$ corresponds
to one divided by the sum of elements in row $i$, i.e.,
$d_{jj}=1/{\sum_{i}h_{ij}}$, and similarly, $e_{jj}$ equals one divided by the
sum of elements in row $j$, i.e., $e_{jj}=1/{\sum_{j}h_{ij}}$. These matrices
serve to balance the impact of rays hitting individual pixels and pixels
intersected by each ray. The transpose of the matrix, denoted as $\bm{A}^{T}$,
facilitates the back-projection of projection images onto the reconstruction
area. It delineates the pixels influenced by a particular ray’s trajectory.
This mechanism enables the derivation of the updated $\bm{x}_{t-1}^{n}$ from
the previously updated $\bm{x}_{t-1/2}^{n}$, thus allowing for the subsequent
progression of the process.

$\displaystyle\bm{x}_{t-1}^{n}=\bm{x}_{t-1/2}^{n}+\beta_{t}(\bm{x}_{t}^{n}-\bm{u}_{t}^{n}),\quad
n=1,\ldots,N.$ (21)

If all the reconstruction procedure is finished, the estimated
$\hat{\bm{x}_{0}}=\sum_{n=1}^{N}w_{0}^{n}\bm{x}_{0}^{n}$. Regarding the second
step of Eq. (14), it outlines the diffusion model’s advancement process,
utilizing the VE-SDE approach for updates. Initially, the generation procedure
for the diffusion model is a generalized numerical solution derived from the
inverse SDE. In contrast to generating samples using the numerical SDE solver,
we opt for the Predictor-Corrector (PC) sampler [15]. This choice is driven by
its superior performance when dealing with VE-SDEs. In the context of PC
samplers, the predictor entails a numerical solver for the reverse-time SDE,
while the corrector encompasses any Markov Chain Monte Carlo (MCMC) technique
reliant solely on the scores. As previously mentioned, the SDE update
procedure is partitioned into predictor and corrector updates. Concerning the
predictor update, we can perform enhancements to refine the estimation of the
model parameters.

$\bm{u}_{t}^{n}=\bm{x}_{t}^{n}+\left(\sigma_{t}^{2}-\sigma_{t-1}^{2}\right)\bm{s}_{\theta^{*}}\left(\bm{x}_{t}^{n},t\right)+\sqrt{\sigma_{t}^{2}-\sigma_{t-1}^{2}}\bm{z},$
(22)

where the parameter $\sigma_{t}>0$ embodies a monotonically rising function in
relation to time $t$, as stipulated by [15]. $z~{}N\left(0,1\right)$ is random
noise that follows a Gaussian distribution, which is added by [15] and can
prevent the model from reaching the local optimization solution. Regarding the
corrector step, the update adheres to the subsequent formulation:

$\bm{u}_{t}^{n}=\bm{u}_{t}^{n}+\epsilon_{t-1}\bm{s}_{\theta^{*}}\left(\bm{x}_{t}^{n},t\right)+\sqrt{2\epsilon_{t-1}}\bm{z},$
(23)

where $\epsilon_{t-1}$ is the step size of at time $t-1$.

###  3.3 Overall Clarification

The DOSM algorithm’s overarching structure is depicted in Fig. 2, leveraging
the insights from the preceding analysis. We first train the VE-SDE model at
the training phase. For the test phase, with a set of input noise images
across $N$ channels, we harness both the data consistency principle and the
learned VE-SDE $S_{\theta}$ to derive $\bm{x}_{T-1}^{n}$ for $n=1,\ldots,N$.
This process sets the stage for computing $\hat{x}_{0}$ at the penultimate
time instance, specifically denoted as $\hat{x}_{0}^{T-1}$. Importantly, this
updated $\hat{x}_{0}^{T-1}$ becomes pivotal for enhancing data consistency at
the preceding time step of $T-2$. The iterative updating methodology results
in the acquisition of the ultimate estimation of $\hat{x}_{0}$.

In conventional data consistency approaches, the search direction remains
unconstrained, lacking a definitive objective for computing a dependable
gradient, which in turn contributes to the instability of the optimization
process. It is the reason why the other SGM-based has a relative poor
convergence, as the MCG and DPS in Fig. 1. This instability ultimately hampers
the accurate retrieval of the actual $\bm{x}_{0}$, as shown in Fig. 3(a). In
contrast, our proposed DOSM fully addresses the challenges in Fig. 3(b).
Firstly, it establishes a stable target, denoted as $\hat{x}_{0}$ and
validated by Proposition, right from the outset of the process. This grounded
objective provides a clear optimization direction, lending stability to the
overall procedure. Secondly, owing to the presence of $N$-channel noisy images
at time $t$, the DOSM inherently encompasses $N$ in distinct directions. This
amalgamation of directions not only enriches the search landscape but also
contributes significantly to fortifying the overall stability of the
reconstruction model. Collectively, these factors substantiate a better robust
convergence of our DOSM, as illustrated in Fig. 1.

Figure 2:  Overall architecture of DOSM network. (a) and (b) represent the training and testing of the network, respectively. The green line represents the diffusion of the backward SDE, the blue line represents the constraint of the data consistency policy, and the orange line represents the estimated $\bm{\hat{x}}_{0}$. Figure 3: Advantages of the DOSM reconstruction step (b) over the traditional diffusion model step (a). The equipotential surface is denoted as $M_{t}$, and the optimization process is centered around identifying the intersection solution between $M_{t}$ and the data consistency policy. The red and pink arrows indicate the directions of the diffusion model and the data consistency policy. The blue arrow indicates the direction in which the search optimization is advancing. Figure 4: Representative reconstructed results of simulated AAPM CT from 23 views by different methods. From top to bottom: reconstructed results, region of interest (ROI), and difference images. The columns from left to right are FBP, FISTA, DPS, MCG, Ours and Label. Figure 5: Similar to the illustration in Fig. 4 but from an alternate representative slice. Figure 6: Representative reconstructed results of clinical cardiae datasets from 23 views by different methods. From top to bottom: reconstructed results, region of interest, and different images. The columns from left to right: FBP, FISTA, DPS, MCG, Ours, and Label. Figure 7: Representative reconstructed results of simulated AAPM CT from 10 views by different methods. From top to bottom: reconstructed results, ROI, and difference images. Figure 8: Representative reconstructed results of clinical cardiae datasets from 10 views by different methods. From top to bottom: reconstructed results, ROI, and difference image. Figure 9: Visualized PSNR and SSIM quantitative results of AAPM CT datasets in 23 and 10 views. Table 1: The statistics results in terms of PSNR and SSIM from 23 and 10 views in clinical cardiac datasets | Clinical Cardiae Datasets
---|---
Method | 23 Views | 10 Views
| PSNR | SSIM | PSNR | SSIM
FBP | 20.30 | 0.374 | 15.16 | 0.209
FISTA |  22.29 | 0.618 | 21.84 | 0.508
CGAN | 22.57 | 0.623 | / | /
FBP-ConvNet | 24.31 | 0.721 | / | /
DPS | 31.27 | 0.763 | 28.25 | 0.657
MCG | 32.22 | 0.774 | 26.44 | 0.699
Ours |  33.93 | 0.881 | 32.11 | 0.874

##  4 Experiments

In this study, we conduct our experimental investigations by harnessing the
Python programming language within the PyTorch framework. All experiments are
implemented on a high-performance computing system equipped with an NVIDIA RTX
A6000 48GB graphics processing unit. The training procedure of VE-SDE follows
with the guidelines recommended by Song et al[18]. Here, the Adam optimization
algorithm is used and the learning rate is set to $2\times 10^{-4}$. In
configuring the noise variance, we establish fixed values of
$\sigma_{\text{min}}=0.01$ and $\sigma_{\text{max}}=378$. The number of
iterations within the sampling process is set as $2000$. Furthermore, a series
of comparison methods are chosen, including FISTA [35], CGAN [12], FBP-ConvNet
[9], DPS [23] and MCG [22]. To quantitatively assess and contrast the efficacy
of the reconstruction outcomes, we employ two metrics: the Peak Signal-to-
Noise Ratio (PSNR) and the Structural Similarity (SSIM). Elevated values of
PSNR and SSIM correspond to heightened levels of reconstruction quality.
Furthermore, our code and datasets are provided for review with the link.

###  4.1 Data Description

We first train and test our model using a simulated dataset coming from the
AAPM 2016 CT Low Dose Challenge [36]. For the simulation study, the clinical
data from nine patients’ are used for training (4,000 slices), and the rest of
the patient’s data is used for testing (400 slices). To generate the 23 views
and 10 views data, an equiangular fan-beam projection geometry is employed. We
systematically generated data from both 23 and 10 views to facilitate the
sparse-view CT reconstruction process. Regarding the scanning geometry
configuration, the distances from the rotation center to the source and
detector are 1500mm and 500 mm. The detector width is 41.3 cm with 720
detector elements.

The real cardiac clinical datasets used in [37] is utilized in this
experiment. The curved cylindrical detector array encompassed a total of 880
individual units, while the complete scan is composed of 2200 views. The
field-of-view (FOV) diameter covers an expansive 49.8 × 49.8 $cm^{2}$, and the
image matrix featured dimensions of 512 × 512 pixels. Notably, the distance
from the X-ray source to the center and detector are set as 53.85 $cm$ and
103.68 $cm$. To validate the universality and exceptional performance of our
approach, the model trained using simulation data is seamlessly applied to the
testing of clinical cardiac data.

###  4.2 Experimental results

####  4.2.1 Simulate Dataset Results

We first compare the reconstruction results with 23 views on the simulated
clinical datasets, as shown in Figs. 4 and 5. Upon scrutiny, it becomes
evident that the FBP results exhibit an excessive proliferation of artifacts
that obscure image intricacies, rendering it challenging to extract
informative content beyond the skeletal structure. FISTA, as a classical
iterative reconstruction technique, manages to enhance image quality compared
to FBP, yet elusive finer details still persist. FBP-ConvNet, an exemplar of
supervised deep learning approaches, brings about substantial improvements in
image quality, accentuating large structures and edges. However, a notable
drawback arises in the form of inaccurately predicted reconstruction details,
a consequence of the profound information loss due to the ultra-sparse-view
measurements.

Incorporating data consistency into the generation process, the score-based
model significantly enhances the efficacy of image reconstruction. In the
sparse-view CT reconstruction landscape, DPS [23] and MCG [22] stand as
prevailing benchmarks. In comparison to traditional and typical deep learning-
based methods, the score-based methodology consistently attains superior
results, yielding heightened image quality. A closer examination within the
region of interest (ROIs) reveals that the image edges delineated by yellow
circles in the DPS result become blurred. Similarly, the MCG outcomes exhibit
conspicuous artifacts and noise, particularly within the extracted ROI of Fig.
4. This artifact’s presence emanates from the intrinsic integration of
estimated projection within MCG, causing an overly constrained data
consistency enforcement during iterations. Additionally, image edges denoted
by yellow arrows are also subjected to blurring.

In contradistinction to these contenders, our proposed approach showcases the
most remarkable image quality, accompanied by a discernible preservation of
image features and details. This assertion is fortified by the discernible
clarity of the yellow circles and arrows within the ROIs. The difference map
from our method demonstrates the closest proximity to zero, unequivocally
signifying the supremacy of our approach over other baseline methods.
Moreover, to demonstrate the out-performance of our method, the profiles along
the yellow location with all reconstruction algorithms are shown in Fig. 10.
One can see that the profile reconstructed by our method is closer to the
ground truth.

Likewise, we proceed to assess the efficacy of our methodology involving only
10 views for image reconstruction, as presented in Fig. 7. The DPS and MCG
baselines have a pronounced presence of noise within the reconstruction
outcomes under such ultra-sparse-view condition, which is further confirmed by
the extracted ROI. That demonstrates the data consistency strategies within
these two methods struggle to adequately guide the generation process of the
score-based model under ultra-sparse-view settings. Within the ROI featured in
Fig. 7, the textural intricacies captured within the yellow circles, as
reconstructed by our method, significantly outperform those achieved by the
alternative approaches. The difference map can also clearly demonstrate the
gap between the reconstructed image and the ground truth obtained by our
method is the smallest. The reconstruction’s quantitative results from both 23
and 10 views in Fig. 9 consistently demonstrate excellent reconstruction
performance with the best PSNRs and SSIMs of our method.

####  4.2.2 Cardiac Clinical Dataset Results

To further evaluate the efficacy and superior performance of our approach, we
perform experimental validation on real clinical cardiac CT data. We employ
the trained VE-SDE in above simulated clinical datasets for this real-world
study. The reconstruction results from 23 and 10 views are showcased in Figs.
6 and 8, respectively. It becomes apparent that the FBP results are marred by
an excessive proliferation of artifacts, making it challenging to distinguish
beneficial structures and details from these artifacts. While FISTA improves
image quality in comparison to FBP, it falls short in capturing finer details.
Conversely, FBP-ConvNet excels in enhancing image quality by highlighting
prominent structures and edges. However, inaccuracies arise due to the extreme
paucity of view measurements. This point is illustrated by the yellow circles
in Fig. 6.

Consistently outperforming traditional deep learning-based reconstructions,
the score-based methodology consistently delivers elevated image quality. Upon
closer examination of the specific ROI in Fig. 6, it’s evident that DPS
results exhibit blurring at edges denoted by yellow circles. Furthermore, it
introduces artifacts marked by yellow arrows. Similarly, the outcomes of MCG
display conspicuous artifacts and noise, particularly within the extracted
ROIs of Figs. 6 and 8. While MCG introduces details and features, they are
erroneous structures and artifacts that ultimately degrade image quality.
Additionally, the image edges indicated by yellow arrows in Fig. 6 are
inaccurate in the MCG results. Notably, the image details and edges introduced
by MCG in Fig. 8 are also inaccurate. The discrepancy map in Figs. 6 and 8
vividly illustrate that the disparity between the reconstructed image and the
ground truth obtained by our method is minimal.

Moreover, profile analyses along the yellow location with all reconstruction
algorithms in Fig. 10 reveal that the profile reconstructed by our method
closely approximates the ground truth. The quantitative reconstruction results
from both 23 and 10 views, presented in Tables I, consistently underscore the
exceptional reconstruction performance of our method, evidenced by superior
PSNRs and SSIMs.

Figure 10: The intensity profiles along the specified yellow line in the
reconstructed images. (a) represent 23 views of simulated dataset results. (b)
represent 23 views of clinical cardiae datasets results.

###  4.3 Ablation experiment

In this section, ablation studies are performed to probe the effectiveness of
the various modules of the proposed method. These experiments are carried out
on the 23 views reconstruction task of simulated clinical data. We randomly
selected 100 test data for ablation experiments.

Figure 11: Comparison of the DOSM reconstruction performance under different
experimental factors. We perform factor optimization using PSNR and SSIM. (a),
(b), and (c) respectively represent the selection of parameters N for
estimating $\bm{\hat{x}}_{0}$, the number of iterations, and the parameters of
$\beta$.

####  4.3.1 Selection of $N$

The impact of $N$ on estimating $\bm{\hat{x}}_{0}$ is crucial. As the number
of N increases, the constructed $\bm{\hat{x}}_{0}$ estimator can achieve
better results. Figure 11(a) demonstrates that larger values of $N$ yield
better quantitative outcomes, while an N of 5 maintains stable performance.
However, higher $N$ increases computational costs. Our study finds a balance
between performance and efficiency by choosing $N$ as 5.

####  4.3.2 Iteration number setting within data consistency

To provide effective guidance for the data consistency strategy during the
reconstruction process, the process underwent multiple iterations to ensure
the optimal alignment with data consistency requirements. The quantitative
evaluation results of our approach are shown in Fig. 11(b), illustrating the
impact of different iteration settings. The results indicate that as the
iteration number increases, the reconstruction outcomes are enhanced. It’s
worth noting that the performance is stable with the data consistency settings
when the iteration is set as 20 iterations. The selection of 20 iterations
within our data consistency strategy achieves a balance between reconstruction
performance and computational efficiency.

####  4.3.3 Selection of $\beta$

As shown in Eq. (17), the new data consistency strategy introduces an
additional term to harmonize the impact of the current noisy data and the
reconstruction outcomes generated through VE-SDE. Here, all $\beta_{t}$ are
set as the same $\beta$. We delved further into the influence of $\beta$, as
illustrated in Fig. 11(c). Evidently, an optimal factor can enhance the
quantitative outcomes. In our investigation, this factor is designated as
$\beta=0.1$.

Figure 12: 4 views of the reconstruction. results

##  5 Discussions and Conclusion

While our proposed methods have demonstrated superior performance compared to
other comparative techniques, it is imperative to address certain subtleties
that demand attention. First, we acknowledge that our method surpasses
classical score-based models in terms of computational time cost. This
divergence emerges from the essential need to estimate $\bm{x}_{0}$ for
establishing the data consistency. It’s worth noting that this estimation
process can be parallelized, effectively mitigating the computational overhead
in practical implementation. Our proposed method has exhibited strong
performance across both 23 and 10 views. However, a pertinent question arises:
Can this prowess extend to even fewer views, such as 4 views. As illustrated
in Fig. 12, our approach continues to successfully reconstruct a good image
even with only 4 views, while refining the image details remains a challenge.
The pursuit of more advanced data consistency strategies and sophisticated
score-based generative models emerges as an intriguing avenue for 4 views of
CT reconstruction. In this study, we underscored the significance of achieving
a higher-quality estimation of $\bm{x}_{0}$ holds significance. The potential
to utilize a prior image as an improved $\bm{x}_{0}$ for guidance in
reconstruction holds promise.

The proposed general framework exhibits remarkable reconstruction prowess for
ultra-sparse-view CT reconstruction. Key attributes can be succinctly outlined
as follows. Firstly, leveraging insights from conventional iterative
reconstruction methods enhances reconstruction stability and elevates image
quality. Moreover, an augmented data consistency paradigm is developed,
seamlessly integrating $\bm{x}_{0}$ estimation. This two-fold approach
significantly bolsters model robustness, where improved $p(\bm{x}_{0})$
directly correlates with heightened image quality. Lastly, we introduce a
discrepancy term to align the original data-driven prior with the SGM
counterpart. Through extensive experiments on diverse datasets, including
simulations and real clinical scenarios, our approach consistently
demonstrates superior performance and unwavering stability.

##  6 Acknowledgments

We would like to express our sincere thanks to Jianjia Zhang and Zirong Li
from the School of Biomedical Engineering of Sun Yat-sen University for
polishing the expression of research.

## References

  * [1] A. Kazerouni, E. K. Aghdam, M. Heidari, R. Azad, M. Fayyaz, I. Hacihaliloglu, and D. Merhof, “Diffusion models for medical image analysis: A comprehensive survey,” _arXiv preprint arXiv:2211.07804_ , 2022.
  * [2] M. J. Willemink and P. B. Noël, “The evolution of image reconstruction for ct—from filtered back projection to artificial intelligence,” _European radiology_ , vol. 29, pp. 2185–2195, 2019.
  * [3] Y. Han and J. C. Ye, “Framing u-net via deep convolutional framelets: Application to sparse-view ct,” _IEEE transactions on medical imaging_ , vol. 37, no. 6, pp. 1418–1429, 2018.
  * [4] J. Bian, J. H. Siewerdsen, X. Han, E. Y. Sidky, J. L. Prince, C. A. Pelizzari, and X. Pan, “Evaluation of sparse-view reconstruction from flat-panel-detector cone-beam ct,” _Physics in Medicine & Biology_, vol. 55, no. 22, p. 6575, 2010.
  * [5] W. Wu, D. Hu, C. Niu, H. Yu, V. Vardhanabhuti, and G. Wang, “Drone: dual-domain residual-based optimization network for sparse-view CT reconstruction,” _IEEE Transactions on Medical Imaging_ , vol. 40, no. 11, pp. 3002–3014, 2021.
  * [6] S. Niu, G. Yu, J. Ma, and J. Wang, “Nonlocal low-rank and sparse matrix decomposition for spectral ct reconstruction,” _Inverse problems_ , vol. 34, no. 2, p. 024003, 2018.
  * [7] A. Ben-Tal and A. Nemirovski, “Robust optimization–methodology and applications,” _Mathematical programming_ , vol. 92, pp. 453–480, 2002.
  * [8] A. Raj, Y. Bresler, and B. Li, “Improving robustness of deep-learning-based image reconstruction,” in _International Conference on Machine Learning_. PMLR, 2020, pp. 7932–7942.
  * [9] K. H. Jin, M. T. McCann, E. Froustey, and M. Unser, “Deep convolutional neural network for inverse problems in imaging,” _IEEE Transactions on Image Processing_ , vol. 26, no. 9, pp. 4509–4522, 2017.
  * [10] Z. Zhang, X. Liang, X. Dong, Y. Xie, and G. Cao, “A sparse-view ct reconstruction method based on combination of densenet and deconvolution,” _IEEE transactions on medical imaging_ , vol. 37, no. 6, pp. 1407–1417, 2018\.
  * [11] H. Chen, Y. Zhang, M. K. Kalra, F. Lin, Y. Chen, P. Liao, J. Zhou, and G. Wang, “Low-dose ct with a residual encoder-decoder convolutional neural network,” _IEEE Transactions on Medical Imaging_ , vol. 36, no. 12, pp. 2524–2535, 2017\.
  * [12] M. U. Ghani and W. C. Karl, “Deep learning-based sinogram completion for low-dose ct,” in _2018 IEEE 13th Image, Video, and Multidimensional Signal Processing Workshop (IVMSP)_. IEEE, 2018, pp. 1–5.
  * [13] J. Pan, H. Zhang, W. Wu, Z. Gao, and W. Wu, “Multi-domain integrative swin transformer network for sparse-view tomographic reconstruction,” _Patterns_ , vol. 3, no. 6, 2022.
  * [14] J. Ho, A. Jain, and P. Abbeel, “Denoising diffusion probabilistic models,” in _Advances in Neural Information Processing Systems_ , H. Larochelle, M. Ranzato, R. Hadsell, M. Balcan, and H. Lin, Eds., vol. 33. Curran Associates, Inc., 2020, pp. 6840–6851.
  * [15] Y. Song, J. Sohl-Dickstein, D. P. Kingma, A. Kumar, S. Ermon, and B. Poole, “Score-based generative modeling through stochastic differential equations,” in _International Conference on Learning Representations_ , 2021\.
  * [16] B. Kawar, M. Elad, S. Ermon, and J. Song, “Denoising diffusion restoration models,” _Advances in Neural Information Processing Systems_ , vol. 35, pp. 23 593–23 606, 2022.
  * [17] J. Song, C. Meng, and S. Ermon, “Denoising diffusion implicit models,” _arXiv preprint arXiv:2010.02502_ , 2020.
  * [18] Y. Song, L. Shen, L. Xing, and S. Ermon, “Solving inverse problems in medical imaging with score-based generative models,” in _International Conference on Learning Representations_ , 2022.
  * [19] H. Li, G. Ditzler, J. Roveda, and A. Li, “Descod-ecg: Deep score-based diffusion model for ecg baseline wander and noise removal,” _IEEE Journal of Biomedical and Health Informatics_ , pp. 1–11, 2023.
  * [20] W. Xia, W. Cong, and G. Wang, “Patch-based denoising diffusion probabilistic model for sparse-view ct reconstruction,” _arXiv preprint arXiv:2211.10388_ , 2022.
  * [21] Z. Li, Y. Wang, J. Zhang, W. Wu, and H. Yu, “Two-and-a-half order score-based model for solving 3d ill-posed inverse problems,” _arXiv preprint arXiv:2308.08511_ , 2023.
  * [22] H. Chung, B. Sim, D. Ryu, and J. C. Ye, “Improving diffusion models for inverse problems using manifold constraints,” _Advances in Neural Information Processing Systems_ , vol. 35, pp. 25 683–25 696, 2022.
  * [23] H. Chung, J. Kim, M. T. Mccann, M. L. Klasky, and J. C. Ye, “Diffusion posterior sampling for general noisy inverse problems,” in _The Eleventh International Conference on Learning Representations_ , 2022.
  * [24] T. Karras, M. Aittala, T. Aila, and S. Laine, “Elucidating the design space of diffusion-based generative models,” _Advances in Neural Information Processing Systems_ , vol. 35, pp. 26 565–26 577, 2022.
  * [25] A. Kazerouni, E. K. Aghdam, M. Heidari, R. Azad, M. Fayyaz, I. Hacihaliloglu, and D. Merhof, “Diffusion models in medical imaging: A comprehensive survey,” _Medical Image Analysis_ , p. 102846, 2023.
  * [26] H. Chung, B. Sim, and J. C. Ye, “Come-closer-diffuse-faster: Accelerating conditional diffusion models for inverse problems through stochastic contraction,” in _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , 2022, pp. 12 413–12 422.
  * [27] B. Guan, C. Yang, L. Zhang, S. Niu, M. Zhang, Y. Wang, W. Wu, and Q. Liu, “Generative modeling in sinogram domain for sparse-view CT reconstruction,” _arXiv preprint arXiv:2211.13926_ , 2022.
  * [28] E. Y. Sidky, C.-M. Kao, and X. Pan, “Accurate image reconstruction from few-views and limited-angle data in divergent-beam ct,” _Journal of X-ray Science and Technology_ , vol. 14, no. 2, pp. 119–139, 2006.
  * [29] Y. Song and S. Ermon, “Generative modeling by estimating gradients of the data distribution,” in _Advances in Neural Information Processing Systems_ , H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, Eds., vol. 32. Curran Associates, Inc., 2019.
  * [30] M. V. Zibetti, C. Lin, and G. T. Herman, “Total variation superiorized conjugate gradient method for image reconstruction,” _Inverse Problems_ , vol. 34, no. 3, p. 034001, 2018.
  * [31] H. Chung, D. Ryu, M. T. McCann, M. L. Klasky, and J. C. Ye, “Solving 3d inverse problems using pre-trained 2d diffusion models,” in _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , 2023, pp. 22 542–22 551.
  * [32] S. Lee, H. Chung, M. Park, J. Park, W.-S. Ryu, and J. C. Ye, “Improving 3d imaging with pre-trained perpendicular 2d diffusion models,” _arXiv preprint arXiv:2303.08440_ , 2023.
  * [33] B. Efron, “Tweedie’s formula and selection bias,” _Journal of the American Statistical Association_ , vol. 106, no. 496, pp. 1602–1614, 2011.
  * [34] J. Trampert and J.-J. Leveque, “Simultaneous iterative reconstruction technique: Physical interpretation based on the generalized least squares solution,” _Journal of Geophysical Research: Solid Earth_ , vol. 95, no. B8, pp. 12 553–12 559, 1990.
  * [35] A. Beck and M. Teboulle, “A fast iterative shrinkage-thresholding algorithm for linear inverse problems,” _SIAM journal on Imaging Sciences_ , vol. 2, no. 1, pp. 183–202, 2009.
  * [36] https://ctcicblog.mayo.edu/2016-low-dose-ct-grand-challenge/.
  * [37] Q. Xu, H. Yu, X. Mou, L. Zhang, J. Hsieh, and G. Wang, “Low-dose x-ray ct reconstruction via dictionary learning,” _IEEE transactions on medical imaging_ , vol. 31, no. 9, pp. 1682–1697, 2012.