DIFFUSION-BASED VOICE CONVERSION WITH FAST MAXIMUM LIKELIHOOD SAMPLING SCHEME
Year: 2022

Abstract:
Voice conversion is a common speech synthesis task which can be solved in different ways depending on a particular real-world scenario. The most challenging one often referred to as one-shot many-to-many voice conversion consists in copying target voice from only one reference utterance in the most general case when both source and target speakers do not belong to the training dataset. We present a scalable high-quality solution based on diffusion probabilistic modeling and demonstrate its superior quality compared to state-of-the-art one-shot voice conversion approaches. Moreover, focusing on real-time applications, we investigate general principles which can make diffusion models faster while keeping synthesis quality at a high level. As a result, we develop a novel Stochastic Differential Equations solver suitable for various diffusion model types and generative tasks as shown through empirical studies and justify it by theoretical analysis.

heading: 1 INTRODUCTION
text: Voice conversion (VC) is the task of copying the target speaker’s voice while preserving the linguistic content of the utterance pronounced by the source speaker. Practical VC applications often require a model which is able to operate in one-shot mode (i.e. when only one reference utterance is provided to copy the target speaker’s voice) for any source and target speakers. Such models are usually referred to as one-shot many-to-many models (or sometimes zero-shot many-to-many models, or just any-to-any VC models). It is challenging to build such a model since it should be able to adapt to a new unseen voice having only one spoken utterance pronounced with it, so it was not until recently that successful one-shot VC solutions started to appear.
Conventional one-shot VC models are designed as autoencoders whose latent space ideally contains only the linguistic content of the encoded utterance while target voice identity information (usually taking shape of speaker embedding) is fed to the decoder as conditioning. Whereas in the pioneering AutoVC model (Qian et al., 2019) only speaker embedding from the pre-trained speaker verification network was used as conditioning, several other models improved on AutoVC enriching conditioning with phonetic features such as pitch and loudness (Qian et al., 2020; Nercessian, 2020), or training voice conversion and speaker embedding networks jointly (Chou & Lee, 2019). Also, several papers (Lin et al., 2021; Ishihara & Saito, 2020; Liu et al., 2021b) made use of attention mechanism to better fuse specific features of the reference utterance into the source utterance thus improving the decoder performance. Apart from providing the decoder with sufficiently rich information, one of the main problems autoencoder VC models face is to disentangle source speaker identity from speech content in the encoder. Some models (Qian et al., 2019; 2020; Nercessian, 2020) solve this problem by introducing an information bottleneck. Among other popular solutions of the disentanglement problem one can mention applying vector quantization technique to the content information (Wu et al., 2020; Wang et al., 2021), utilizing features of Variational AutoEncoders (Luong & Tran, 2021; Saito et al., 2018; Chou & Lee, 2019), introducing instance normalization layers (Chou & Lee, 2019; Chen et al., 2021b), and using Phonetic Posteriorgrams (PPGs) (Nercessian, 2020; Liu et al., 2021b).
The model we propose in this paper solves the disentanglement problem by employing the encoder predicting “average voice”: it is trained to transform mel features corresponding to each phoneme into mel features corresponding to this phoneme averaged across a large multi-speaker dataset. As for decoder, in our VC model, it is designed as a part of a Diffusion Probabilistic Model (DPM) since this class of generative models has shown very good results in speech-related tasks like raw waveform generation (Chen et al., 2021a; Kong et al., 2021) and mel feature generation (Popov et al., 2021; Jeong et al., 2021). However, this decoder choice poses a problem of slow inference because DPM forward pass scheme is iterative and to obtain high-quality results it is typically necessary to run it for hundreds of iterations (Ho et al., 2020; Nichol & Dhariwal, 2021). Addressing this issue, we develop a novel inference scheme that significantly reduces the number of iterations sufficient to produce samples of decent quality and does not require model re-training. Although several attempts have been recently made to reduce the number of DPM inference steps (Song et al., 2021a; San-Roman et al., 2021; Watson et al., 2021; Kong & Ping, 2021; Chen et al., 2021a), most of them apply to some particular types of DPMs. In contrast, our approach generalizes to all popular kinds of DPMs and has a strong connection with likelihood maximization.
This paper has the following structure: in Section 2 we present a one-shot many-to-many VC model and describe DPM it relies on; Section 3 introduces a novel DPM sampling scheme and establishes its connection with likelihood maximization; the experiments regarding voice conversion task as well as those demonstrating the benefits of the proposed sampling scheme are described in Section 4; we conclude in Section 5.

heading: 2 VOICE CONVERSION DIFFUSION MODEL
text: As with many other VC models, the one we propose belongs to the family of autoencoders. In fact, any conditional DPM with data-dependent prior (i.e. terminal distribution of forward diffusion) can be seen as such: forward diffusion gradually adding Gaussian noise to data can be regarded as encoder while reverse diffusion trying to remove this noise acts as a decoder. DPMs are trained to minimize the distance (expressed in different terms for different model types) between the trajectories of forward and reverse diffusion processes thus, speaking from the perspective of autoencoders, minimizing reconstruction error. Data-dependent priors have been proposed by Popov et al. (2021) and Lee et al. (2021), and we follow the former paper due to the flexibility of the continuous DPM framework used there. Our approach is summarized in Figure 1.

heading: 2.1 ENCODER
text: We choose average phoneme-level mel features as speaker-independent speech representation. To train the encoder to convert input mel-spectrograms into those of “average voice”, we take three steps: (i) first, we apply Montreal Forced Aligner (McAuliffe et al., 2017) to a large-scale multi-speaker LibriTTS dataset (Zen et al., 2019) to align speech frames with phonemes; (ii) next, we obtain average mel features for each particular phoneme by aggregating its mel features across the whole LibriTTS dataset; (iii) the encoder is then trained to minimize mean square error between output mel-spectrograms and ground truth “average voice” mel-spectrograms (i.e. input mel-spectrograms where each phoneme mel feature is replaced with the average one calculated on the previous step).
The encoder has exactly the same Transformer-based architecture used in Grad-TTS (Popov et al., 2021) except that its inputs are mel features rather than character or phoneme embeddings. Note that unlike Grad-TTS the encoder is trained separately from the decoder described in the next section.

heading: 2.2 DECODER
text: Whereas the encoder parameterizes the terminal distribution of the forward diffusion (i.e. the prior), the reverse diffusion is parameterized with the decoder. Following Song et al. (2021c) we use Itô calculus and define diffusions in terms of stochastic processes rather than discrete-time Markov chains.
The general DPM framework we utilize consists of forward and reverse diffusions given by the following Stochastic Differential Equations (SDEs):
dXt = 1
2 βt(X̄ −Xt)dt+
√ βtd −→ Wt , (1)
dX̂t =
( 1
2 (X̄ − X̂t)− sθ(X̂t, X̄, t)
) βtdt+ √ βtd ←− Wt , (2)
where t ∈ [0, 1], −→ W and ←− W are two independent Wiener processes in Rn, βt is non-negative function referred to as noise schedule, sθ is the score function with parameters θ and X̄ is n-dimensional vector. It can be shown (Popov et al., 2021) that the forward SDE (1) allows for explicit solution:
Law(Xt|X0) = N ( e− 1 2 ∫ t 0 βsdsX0 + ( 1− e− 12 ∫ t 0 βsds ) X̄, ( 1− e− ∫ t 0 βsds ) I ) , (3)
where I is n× n identity matrix. Thus, if noise follows linear schedule βt = β0 + t(β1 − β0) for β0 and β1 such that e− ∫ 1 0 βsds is close to zero, then Law (X1) is close to N (X̄, I) which is the prior in this DPM. The reverse diffusion (2) is trained by minimizing weighted L2 loss:
θ∗ = arg min θ L(θ) = arg min θ ∫ 1 0 λtEX0,Xt‖sθ(Xt, X̄, t)−∇log pt|0(Xt|X0)‖22dt, (4)
where pt|0(Xt|X0) is the probability density function (pdf) of the conditional distribution (3) and λt = 1− e− ∫ t 0 βsds. The distribution (3) is Gaussian, so we have
∇ log pt|0(Xt|X0) = − Xt −X0e− 1 2
∫ t 0 βsds − X̄(1− e− 12 ∫ t 0 βsds)
1− e− ∫ t 0 βsds
. (5)
At training, time variable t is sampled uniformly from [0, 1], noisy samplesXt are generated according to the formula (3) and the formula (5) is used to calculate loss function L on these samples. Note that Xt can be sampled without the necessity to calculate intermediate values {Xs}0<s<t which makes optimization task (4) time and memory efficient. A well-trained reverse diffusion (2) has trajectories that are close to those of the forward diffusion (1), so generating data with this DPM can be performed by sampling X̂1 from the prior N (X̄, I) and solving SDE (2) backwards in time.
The described above DPM was introduced by Popov et al. (2021) for text-to-speech task and we adapt it for our purposes. We put X̄ = ϕ(X0) where ϕ is the encoder, i.e. X̄ is the “average voice” mel-spectrogram which we want to transform into that of the target voice. We condition the decoder sθ = sθ(X̂t, X̄, gt(Y ), t) on some trainable function gt(Y ) to provide it with information about the target speaker (Y stands for forward trajectories of the target mel-spectrogram at inference and the ones of the training mel-spectrogram at training). This function is a neural network trained jointly with the decoder. We experimented with three input types for this network:
• d-only – the input is the speaker embedding extracted from the target mel-spectrogram Y0 with the pre-trained speaker verification network employed in (Jia et al., 2018);
• wodyn – in addition, the noisy target mel-spectrogram Yt is used as input; • whole – in addition, the whole dynamics of the target mel-spectrogram under forward
diffusion {Ys|s = 0.5/15, 1.5/15, .., 14.5/15} is used as input.
The decoder architecture is based on U-Net (Ronneberger et al., 2015) and is the same as in GradTTS but with four times more channels to better capture the whole range of human voices. The speaker conditioning network gt(Y ) is composed of 2D convolutions and MLPs and described in detail in Appendix H. Its output is 128-dimensional vector which is broadcast-concatenated to the concatenation of X̂t and X̄ as additional 128 channels.

heading: 2.3 RELATED VC MODELS
text: To the best of our knowledge, there exist two diffusion-based voice conversion models: VoiceGrad (Kameoka et al., 2020) and DiffSVC (Liu et al., 2021a). The one we propose differs from them in several important aspects. First, neither of the mentioned papers considers a one-shot many-tomany voice conversion scenario. Next, these models take no less than 100 reverse diffusion steps at inference while we pay special attention to reducing the number of iterations (see Section 3) achieving good quality with as few as 6 iterations. Furthermore, VoiceGrad performs voice conversion by running Langevin dynamics starting from the source mel-spectrogram, thus implicitly assuming that forward diffusion trajectories starting from the mel-spectrogram we want to synthesize are likely to pass through the neighborhood of the source mel-spectrogram on their way to Gaussian noise. Such an assumption allowing to have only one network instead of encoder-decoder architecture is too strong and hardly holds for real voices. Finally, DiffSVC performs singing voice conversion and relies on PPGs as speaker-independent speech representation.

heading: 3 MAXIMUM LIKELIHOOD SDE SOLVER
text: In this section, we develop a fixed-step first-order reverse SDE solver that maximizes the loglikelihood of sample paths of the forward diffusion. This solver differs from general-purpose Euler-Maruyama SDE solver (Kloeden & Platen, 1992) by infinitesimally small values which can however become significant when we sample from diffusion model using a few iterations.
Consider the following forward and reverse SDEs defined in Euclidean space Rn for t ∈ [0, 1]:
dXt = − 1
2 βtXtdt+
√ βtd −→ Wt (F ), dX̂t = ( −1
2 βtX̂t − βtsθ(X̂t, t)
) dt+ √ βtd ←− Wt (R), (6)
where −→ W is a forward Wiener process (i.e. its forward increments −→ Wt − −→ Ws are independent of −→ Ws for t > s) and ←− W is a backward Wiener process (i.e. backward increments ←− Ws − ←− Wt are independent of ←− Wt for s < t). Following Song et al. (2021c) we will call DPM (6) Variance Preserving (VP). For simplicity we will derive maximum likelihood solver for this particular type of diffusion models. The equation (1) underlying VC diffusion model described in Section 2 can be transformed into the equation (6-F) by a constant shift and we will call such diffusion models Mean Reverting Variance Preserving (MR-VP). VP model analysis carried out in this section can be easily extended (see Appendices D, E and F) to MR-VP model as well as to other common diffusion model types such as sub-VP and VE described by Song et al. (2021c).
The forward SDE (6-F) allows for explicit solution:
Law(Xt|Xs) = N (γs,tXs, (1− γ2s,t) I), γs,t = exp ( −1
2 ∫ t s βudu ) , (7)
for all 0 ≤ s < t ≤ 1. This formula is derived by means of Itô calculus in Appendix A. The reverse SDE (6-R) parameterized with a neural network sθ is trained to approximate gradient of the log-density of noisy data Xt:
θ∗ = arg min θ ∫ 1 0 λtEXt‖sθ(Xt, t)−∇log pt(Xt)‖22dt, (8)
where the expectation is taken with respect to noisy data distribution Law(Xt) with pdf pt(·) and λt is some positive weighting function. Note that certain Lipschitz constraints should be satisfied by coefficients of SDEs (6) to guarantee existence of strong solutions (Liptser & Shiryaev, 1978), and throughout this section we assume these conditions are satisfied as well as those from (Anderson, 1982) which guarantee that paths X̂ generated by the reverse SDE (6-R) for the optimal θ∗ equal forward SDE (6-F) paths X in distribution.
The generative procedure of a VP DPM consists in solving the reverse SDE (6-R) backwards in time starting from X̂1 ∼ N (0, I). Common Euler-Maruyama solver introduces discretization error (Kloeden & Platen, 1992) which may harm sample quality when the number of iterations is small. At the same time, it is possible to design unbiased (Henry-Labordère et al., 2017) or even exact (Beskos & Roberts, 2005) numerical solvers for some particular SDE types. The Theorem 1 shows that in the case of diffusion models we can make use of the forward diffusion (6-F) and propose a reverse SDE solver which is better than the general-purpose Euler-Maruyama one in terms of likelihood.
The solver proposed in the Theorem 1 is expressed in terms of the values defined as follows:
µs,t = γs,t 1− γ20,s 1− γ20,t , νs,t = γ0,s 1− γ2s,t 1− γ20,t , σ2s,t = (1− γ20,s)(1− γ2s,t) 1− γ20,t , (9)
κ∗t,h = νt−h,t(1− γ20,t)
γ0,tβth − 1, ω∗t,h = µt−h,t − 1 βth + 1 + κ∗t,h 1− γ20,t − 1 2 ,
(σ∗t,h) 2 = σ2t−h,t +
1 n ν2t−h,tEXt [Tr (Var (X0|Xt))] ,
(10)
where n is data dimensionality, Var (X0|Xt) is the covariance matrix of the conditional data distribution Law(X0|Xt) (so, Tr(Var (X0|Xt)) is the overall variance across all n dimensions) and the expectation EXt [·] is taken with respect to the unconditional noisy data distribution Law(Xt). Theorem 1. Consider a DPM characterized by SDEs (6) with reverse diffusion trained till optimality. Let N ∈ N be any natural number and h = 1/N . Consider the following class of fixed step size h reverse SDE solvers parameterized with triplets of real numbers {(κ̂t,h, ω̂t,h, σ̂t,h)|t = h, 2h, .., 1}:
X̂t−h = X̂t + βth
(( 1
2 + ω̂t,h
) X̂t + (1 + κ̂t,h)sθ∗(X̂t, t) ) + σ̂t,hξt, (11)
where θ∗ is given by (8), t = 1, 1− h, .., h and ξt are i.i.d. samples from N (0, I). Then:
(i) Log-likelihood of sample paths X = {Xkh}Nk=0 under generative model X̂ is maximized for κ̂t,h = κ∗t,h, ω̂t,h = ω ∗ t,h and σ̂t,h = σ ∗ t,h.
(ii) Assume that the SDE solver (11) starts from random variable X̂1 ∼ Law (X1). If X0 is a constant or a Gaussian random variable with diagonal isotropic covariance matrix (i.e. δ2 I for δ > 0), then generative model X̂ is exact for κ̂t,h = κ∗t,h, ω̂t,h = ω ∗ t,h and σ̂t,h = σ ∗ t,h.
The Theorem 1 provides an improved DPM sampling scheme which comes at no additional computational cost compared to standard methods (except for data-dependent term in σ∗ as discussed in Section 4.3) and requires neither model re-training nor extensive search on noise schedule space. The proof of this theorem is given in Appendix C. Note that it establishes optimality of the reverse SDE solver (11) with the parameters (10) in terms of likelihood of discrete paths X = {Xkh}Nk=0 while the optimality of continuous model (6-R) on continuous paths {Xt}t∈[0,1] is guaranteed for a model with parameters θ = θ∗ as shown in (Song et al., 2021c).
The class of reverse SDE solvers considered in the Theorem 1 is rather broad: it is the class of all fixed-step solvers whose increments at time t are linear combination of X̂t, sθ(X̂t, t) and Gaussian noise with zero mean and diagonal isotropic covariance matrix. As a particular case it includes Euler-Maruyama solver (κ̂t,h ≡ 0, ω̂t,h ≡ 0, σ̂t,h ≡ √ βth) and for fixed t and h → 0 we have κ∗t,h = ō(1), ω ∗ t,h = ō(1) and σ ∗ t,h = √ βth(1 + ō(1)) (the proof is given in Appendix B), so the optimal SDE solver significantly differs from general-purpose Euler-Maruyama solver only when N is rather small or t has the same order as h, i.e. on the final steps of DPM inference. Appendix G contains toy examples demonstrating the difference of the proposed optimal SDE solver and Euler-Maruyama one depending on step size.
The result (ii) from the Theorem 1 strengthens the result (i) for some particular data distributions, but it may seem useless since in practice data distribution is far from being constant or Gaussian. However, in case of generation with strong conditioning (e.g. mel-spectrogram inversion) the assumptions on the data distribution may become viable: in the limiting case when our model is conditioned on c = ψ(X0) for an injective function ψ, random variable X0|c becomes a constant ψ−1(c).

heading: 4 EXPERIMENTS
text: We trained two groups of models: Diff-VCTK models on VCTK (Yamagishi et al., 2019) dataset containing 109 speakers (9 speakers were held out for testing purposes) and Diff-LibriTTS models on LibriTTS (Zen et al., 2019) containing approximately 1100 speakers (10 speakers were held out). For every model both encoder and decoder were trained on the same dataset. Training hyperparameters, implementation and data processing details can be found in Appendix I. For mel-spectrogram inversion, we used the pre-trained universal HiFi-GAN vocoder (Kong et al., 2020) operating at 22.05kHz. All subjective human evaluation was carried out on Amazon Mechanical Turk (AMT) with Master assessors to ensure the reliability of the obtained Mean Opinion Scores (MOS). In all AMT tests we considered unseen-to-unseen conversion with 25 unseen (for both Diff-VCTK and Diff-LibriTTS) speakers: 9 VCTK speakers, 10 LibriTTS speakers and 6 internal speakers. For VCTK source speakers we also ensured that source phrases were unseen during training. We place other details of listening AMT tests in Appendix J. A small subset of speech samples used in them is available at our demo page https://diffvc-fast-ml-solver.github.io which we encourage to visit. The code will soon be published at https://github.com/huawei-noah/ Speech-Backbones.
As for sampling, we considered the following class of reverse SDE solvers:
X̂t−h = X̂t + βth
(( 1
2 + ω̂t,h
) (X̂t − X̄) + (1 + κ̂t,h)sθ(X̂t, X̄, gt(Y ), t) ) + σ̂t,hξt, (12)
where t = 1, 1 − h, .., h and ξt are i.i.d. samples from N (0, I). For κ̂t,h = κ∗t,h, ω̂t,h = ω∗t,h and σ̂t,h = σ∗t,h (where κ ∗ t,h, ω ∗ t,h and σ ∗ t,h are given by (10)) it becomes maximum likelihood reverse SDE solver for MR-VP DPM (1-2) as shown in Appendix D. In practice it is not trivial to estimate variance of the conditional distribution Law (X0|Xt), so we skipped this term in σ∗t,h
assuming this variance to be rather small because of strong conditioning on gt(Y ) and just used σ̂t,h = σt−h,t calling this sampling method ML-N (N = 1/h is the number of SDE solver steps). We also experimented with Euler-Maruyama solver EM-N (i.e. κ̂t,h = 0, ω̂t,h = 0, σ̂t,h = √ βth) and “probability flow sampling” from (Song et al., 2021c) which we denote by PF-N (κ̂t,h = −0.5, ω̂t,h = 0, σ̂t,h = 0).

heading: 4.1 SPEAKER CONDITIONING ANALYSIS
text: For each dataset we trained three models – one for each input type for the speaker conditioning network gt(Y ) (see Section 2.2). Although these input types had much influence neither on speaker similarity nor on speech naturalness, we did two experiments to choose the best models (one for each training dataset) in terms of speaker similarity for further comparison with baseline systems. We compared voice conversion results (produced by ML-30 sampling scheme) on 92 source-target pairs. AMT workers were asked which of three models (if any) sounded most similar to the target speaker and which of them (if any) sounded least similar. For Diff-VCTK and Diff-LibriTTS models each conversion pair was evaluated 4 and 5 times respectively. Table 1 demonstrates that for both Diff-VCTK and Diff-LibriTTS the best option is wodyn, i.e. to condition the decoder at time t on the speaker embedding together with the noisy target mel-spectrogram Yt. Conditioning on Yt allows making use of diffusion-specific information of how the noisy target sounds whereas embedding from the pre-trained speaker verification network contains information only about the clean target. Taking these results into consideration, we used Diff-VCTK-wodyn and Diff-LibriTTS-wodyn in the remaining experiments.

heading: 4.2 ANY-TO-ANY VOICE CONVERSION
text: We chose four recently proposed VC models capable of one-shot many-to-many synthesis as the baselines:
• AGAIN-VC (Chen et al., 2021b), an improved version of a conventional autoencoder AdaINVC solving the disentanglement problem by means of instance normalization;
• FragmentVC (Lin et al., 2021), an attention-based model relying on wav2vec 2.0 (Baevski et al., 2020) to obtain speech content from the source utterance;
As shown in (Kim et al., 2021), PPG-based VC models provide high voice conversion quality competitive even with that of the state-of-the-art VC models taking text transcription corresponding to the source utterance as input. Therefore, we can consider BNE-PPG-VC a state-of-the-art model in our setting.
Baseline voice conversion results were produced by the pre-trained VC models provided in official GitHub repositories. Since only BNE-PPG-VC has the model pre-trained on a large-scale dataset (namely, LibriTTS + VCTK), we did two subjective human evaluation tests: the first one comparing Diff-VCTK with AGAIN-VC, FragmentVC and VQMIVC trained on VCTK and the second one comparing Diff-LibriTTS with BNE-PPG-VC. The results of these tests are given in Tables 2 and 3 respectively. Speech naturalness and speaker similarity were assessed separately. AMT workers evaluated voice conversion quality on 350 source-target pairs on 5-point scale. In the first test, each pair was assessed 6 times on average both in speech naturalness and speaker similarity evaluation; as for the second one, each pair was assessed 8 and 9 times on average in speech naturalness and speaker similarity evaluation correspondingly. No less than 41 unique assessors took part in each test.
Table 2 demonstrates that our model performs significantly better than the baselines both in terms of naturalness and speaker similarity even when 6 reverse diffusion iterations are used. Despite working almost equally well on VCTK speakers, the best baseline VQMIVC shows poor performance on other speakers perhaps because of not being able to generalize to different domains with lower recording quality. Although Diff-VCTK performance also degrades on non-VCTK speakers, it achieves good speaker similarity of MOS 3.6 on VCTK ones when ML-30 sampling scheme is used and only slightly worse MOS 3.5 when 5x less iterations are used at inference.
Table 3 contains human evaluation results of Diff-LibriTTS for four sampling schemes: ML-30 with 30 reverse SDE solver steps and ML-6, EM-6 and PF-6 with 6 steps of reverse diffusion. The three schemes taking 6 steps achieved real-time factor (RTF) around 0.1 on GPU (i.e. inference was 10 times faster than real time) while the one taking 30 steps had RTF around 0.5. The proposed model Diff-LibriTTS-ML-30 and the baseline BNE-PPG-VC show the same performance on the VCTK test set in terms of speech naturalness the latter being slightly better in terms of speaker similarity which can perhaps be explained by the fact that BNE-PPG-VC was trained on the union of VCTK and LibriTTS whereas our model was trained only on LibriTTS. As for the whole test set containing unseen LibriTTS and internal speakers also, Diff-LibriTTS-ML-30 outperforms BNE-PPG-VC model achieving MOS 4.0 and 3.4 in terms of speech naturalness and speaker similarity respectively. Due to employing PPG extractor trained on a large-scale ASR dataset LibriSpeech (Panayotov et al., 2015), BNE-PPG-VC has fewer mispronunciation issues than our model but synthesized speech suffers from more sonic artifacts. This observation makes us believe that incorporating PPG features in the proposed diffusion VC framework is a promising direction for future research.
Table 3 also demonstrates the benefits of the proposed maximum likelihood sampling scheme over other sampling methods for a small number of inference steps: only ML-N scheme allows us to
use as few as N = 6 iterations with acceptable quality degradation of MOS 0.2 and 0.1 in terms of naturalness and speaker similarity respectively while two other competing methods lead to much more significant quality degradation.

heading: 4.3 MAXIMUM LIKELIHOOD SAMPLING
text: To show that the maximum likelihood sampling scheme proposed in Section 3 generalizes to different tasks and DPM types, we took the models trained by Song et al. (2021c) on CIFAR-10 image generation task and compared our method with other sampling schemes described in that paper in terms of Fréchet Inception Distance (FID).
The main difficulty in applying maximum likelihood SDE solver is estimating data-dependent term E[Tr (Var (X0|Xt))] in σ∗t,h. Although in the current experiments we just set this term to zero, we can think of two possible ways to estimate it: (i) approximate Var (X0|Xt) with Var (X̂0|X̂t = Xt): sample noisy data Xt, solve reverse SDE with sufficiently small step size starting from terminal condition X̂t = Xt several times, and calculate sample variance of the resulting solutions at initial points X̂0; (ii) use the formula (58) from Appendix C to calculate Var (X0|Xt) assuming that X0 is distributed normally with mean and variance equal to sample mean and sample variance computed on the training dataset. Experimenting with these techniques and exploring new ones seems to be an interesting future research direction.
Another important practical consideration is that the proposed scheme is proven to be optimal only for score matching networks trained till optimality. Therefore, in the experiments whose results are reported in Table 4 we apply maximum likelihood sampling scheme only when t ≤ τ while using standard Euler-Maruyama solver for t > τ for some hyperparameter τ ∈ [0, 1]. Such a modification relies on the assumption that score matching network is closer to being optimal for smaller noise.
Table 4 shows that despite likelihood and FID are two metrics that do not perfectly correlate (Song et al., 2021b), in most cases our maximum likelihood SDE solver performs best in terms of FID. Also, it is worth mentioning that although τ = 1 is always rather a good choice, tuning this hyperparameter can lead to even better performance. One can find randomly chosen generated images for various sampling methods in Figure 2.

heading: 5 CONCLUSION
text: In this paper, the novel one-shot many-to-many voice conversion model has been presented. Its encoder design and powerful diffusion-based decoder made it possible to achieve good results both in terms of speaker similarity and speech naturalness even on out-of-domain unseen speakers. Subjective human evaluation verified that the proposed model delivers scalable VC solution with competitive performance. Furthermore, aiming at fast synthesis, we have developed and theoretically justified the novel sampling scheme. The main idea behind it is to modify the general-purpose Euler-Maruyama SDE solver so as to maximize the likelihood of discrete sample paths of the forward diffusion. Due to the proposed sampling scheme, our VC model is capable of high-quality voice conversion with as few as 6 reverse diffusion steps. Moreover, experiments on the image generation task show that all known diffusion model types can benefit from the proposed SDE solver.

heading: A FORWARD VP SDE SOLUTION
text: Since function γ−10,tXt is linear in Xt, taking its differential does not require second order derivative term in Itô’s formula:
d(γ−10,tXt) = d ( e 1 2 ∫ t 0 βuduXt ) = e 1 2 ∫ t 0 βudu · 1
2 βtXtdt+ e
1 2 ∫ t 0 βudu · ( −1
2 βtXtdt+
√ βtd −→ Wt ) = √ βte 1 2 ∫ t 0 βudud −→ Wt.
(13)
Integrating this expression from s to t results in an Itô’s integral:
e 1 2 ∫ t 0 βuduXt − e 1 2 ∫ s 0 βuduXs = ∫ t s √ βτe 1 2 ∫ τ 0 βudud −→ Wτ , (14)
or
Xt = e − 12 ∫ t s βuduXs + ∫ t s √ βτe − 12 ∫ t τ βudud −→ Wτ . (15)
The integrand on the right-hand side is deterministic and belongs to L2[0, 1] (for practical noise schedule choices), so its Itô’s integral is a normal random variable, a martingale (meaning it has zero mean) and satisfies Itô’s isometry which allows to calculate its variance:
Var(Xt|Xs) = ∫ t s βτe − ∫ t τ βudu I dτ = ( 1− e− ∫ t s βudu ) I . (16)
Thus
Law(Xt|Xs) = N ( e− 1 2 ∫ t s βuduXs, ( 1− e− ∫ t s βudu ) I ) = N (γs,tXs, (1− γ2s,t) I) (17)

heading: B THE OPTIMAL COEFFICIENTS ASYMPTOTICS
text: First derive asymptotics for γ:
γt−h,t = e − 12 ∫ t t−h βudu = 1− 1
2 βth+ ō(h), (18)
γ20,t−h = e −
∫ t−h 0 βudu = e− ∫ t 0 βudue ∫ t t−h βudu = γ20,t(1 + βth) + ō(h), (19)
γ0,t−h = e − 12
∫ t−h 0 βudu = e− 1 2 ∫ t 0 βudue 1 2 ∫ t t−h βudu = γ0,t(1 + 1
2 βth) + ō(h), (20)
γ2t−h,t = e − ∫ t t−h βudu = 1− βth+ ō(h). (21)
Then find asymptotics for µ, ν and σ2:
µt−h,t =
( 1− 1
2 βth+ ō(h)
) 1− γ20,t − γ20,tβth+ ō(h)
1− γ20,t = 1− 1 2 βth− γ20,t 1− γ20,t βth+ō(h), (22)
νt−h,t = (γ0,t(1 + 1
2 βth) + ō(h))
βth+ ō(h)
1− γ20,t = γ0,t 1− γ20,t βth+ ō(h), (23)
σ2t−h,t = 1
1− γ20,t (βth+ ō(h))(1− γ20,t(1 + βth) + ō(h)) = βth+ ō(h). (24)
Finally we get asymptotics for κ∗, ω∗ and σ∗:
κ∗t,h = νt−h,t(1− γ20,t)
γ0,tβth − 1 = γ0,t−h(1− γ2t−h,t) γ0,tβth − 1
= (βth+ ō(h))((1 +
1 2βth)γ0,t + ō(h))
γ0,tβth − 1 = ō(1),
(25)
βthω ∗ t,h = µt−h,t − 1 + νt−h,t γ0,t − 1 2 βth = 1− 1 2 βth− γ20,t 1− γ20,t βth− 1− 1 2 βth+ 1 γ0,t ×
×
( γ0,t
1− γ20,t βth+ ō(h)
) + ō(h) = βth ( −1−
γ20,t 1− γ20,t + 1 1− γ20,t
) + ō(h) = ō(h),
(26)
(σ∗t,h) 2 = σ2t−h,t + ν 2 t−h,tEXt [Tr (Var (X0|Xt))] /n = βth+ ō(h)
+ γ20,t
(1− γ20,t)2 β2t h
2EXt [Tr (Var (X0|Xt))] /n = βth(1 + ō(1)). (27)

heading: C PROOF OF THE THEOREM 1
text: The key fact necessary to prove the Theorem 1 is established in the following Lemma 1. Let p0|t(·|x) be pdf of conditional distribution Law (X0|Xt = x). Then for any t ∈ [0, 1] and x ∈ Rn
sθ∗(x, t) = − 1
1− γ20,t
( x− γ0,tEp0|t(·|x)X0 ) . (28)
Proof of the Lemma 1. As mentioned in (Song et al., 2021c), an expression alternative to (8) can be derived for θ∗ under mild assumptions on the data density (Hyvärinen, 2005; Vincent, 2011):
θ∗ = arg min θ ∫ 1 0 λtEX0∼p0(·)EXt∼pt|0(·|X0)‖sθ(Xt, t)−∇log pt|0(Xt|X0)‖ 2 2dt, (29)
where Law (X0) is data distribution with pdf p0(·) and Law (Xt|X0 = x0) has pdf pt|0(·|x0). By Bayes formula we can rewrite this in terms of pdfs pt(·) and p0|t(·|xt) of distributions Law (Xt) and Law (X0|Xt = xt) correspondingly:
θ∗ = arg min θ ∫ 1 0 λtEXt∼pt(·)EX0∼p0|t(·|Xt)‖sθ(Xt, t)−∇log pt|0(Xt|X0)‖ 2 2dt. (30)
For any n-dimensional random variable ξ with finite second moment and deterministic vector a we have
E‖ξ − a‖22 = E‖ξ − Eξ + Eξ − a‖22 = E‖ξ − Eξ‖22 + 2〈E[ξ − Eξ],Eξ − a〉 + E‖Eξ − a‖22 = E‖ξ − Eξ‖22 + ‖Eξ − a‖22.
(31)
In our case ξ = ∇log pt|0(Xt|X0) and a = sθ(Xt, t), so E‖ξ − Eξ‖22 is independent of θ. Thus
θ∗ = arg min θ ∫ 1 0 λtEXt∼pt(·)‖sθ(Xt, t)− EX0∼p0|t(·|Xt) [ ∇log pt|0(Xt|X0) ] ‖22dt. (32)
Therefore, the optimal score estimation network sθ∗ can be expressed as
sθ∗(x, t) = Ep0|t(·|x) [ ∇log pt|0(x|X0) ] (33)
for all t ∈ [0, 1] and x ∈ supp {pt} = Rn. As proven in Appendix A, Law (Xt|X0) is Gaussian with mean vector γ0,tX0 and covariance matrix (1− γ20,t) I, so finally we obtain
sθ∗(x, t) = Ep0|t(·|x)
[ − 1
1− γ20,t (x− γ0,tX0)
] = − 1
1− γ20,t
( x− γ0,tEp0|t(·|x)X0 ) . (34)
Now let us prove the Theorem 1.
Proof of the Theorem 1. The sampling scheme (11) consists in adding Gaussian noise to a linear combination of X̂t and sθ∗(X̂t, t). Combining (11) and the Lemma 1 we get
X̂t−h = σ̂t,hξt + X̂t + βth
(( 1
2 + ω̂t,h
) X̂t + (1 + κ̂t,h)sθ∗(X̂t, t) ) = σ̂t,hξt
+ ( 1 + βth ( 1
2 + ω̂t,h
)) X̂t + βth(1 + κ̂t,h) ( − 1
1− γ20,t
( X̂t − γ0,tEp0|t(·|X̂t)X0
))
= σ̂t,hξt + ( 1 + βth ( 1
2 + ω̂t,h − 1 + κ̂t,h 1− γ20,t
)) X̂t + γ0,tβth(1 + κ̂t,h)
1− γ20,t Ep0|t(·|X̂t)X0,
(35)
where ξt are i.i.d. random variables from standard normal distribution N (0, I) for t = 1, 1− h, .., h. Thus, the distribution X̂t−h|X̂t is also Gaussian:
Law (X̂t−h|X̂t) = N ( µ̂t,h(κ̂t,h, ω̂t,h)X̂t + ν̂t,h(κ̂t,h)Ep0|t(·|X̂t)X0, σ̂ 2 t,h I ) , (36)
µ̂t,h(κ̂t,h, ω̂t,h) = 1 + βth
( 1
2 + ω̂t,h − 1 + κ̂t,h 1− γ20,t
) , (37)
ν̂t,h(κ̂t,h) = γ0,tβth(1 + κ̂t,h)
1− γ20,t , (38)
which leads to the following formula for the transition densities of the reverse diffusion:
p̂t−h|t(xt−h|xt) = 1√
2πσ̂nt,h exp
( − ‖xt−h − µ̂t,hxt − ν̂t,hEp0|t(·|xt)X0‖22
2σ̂2t,h
) . (39)
Moreover, comparing µ̂t,h and ν̂t,h with µt−h,t and νt−h,t defined in (9) we deduce that
ν̂t,h = νt−h,t ⇔ γ0,tβth(1 + κ̂t,h)
1− γ20,t = νt−h,t ⇔ κ̂t,h = κ∗t,h. (40)
If we also want µ̂t,h = µt−h,t to be satisfied, then we should have
1 + βth
( 1
2 + ω̂t,h − 1 + κ∗t,h 1− γ20,t
) = µt−h,t ⇔ ( µt−h,t − 1
βth − ω∗t,h + ω̂t,h
) βth+ 1 = µt−h,t,
(41)
i.e. ν̂t,h = νt−h,t and µ̂t,h = µt−h,t iff κ̂t,h = κ∗t,h and ω̂t,h = ω ∗ t,h for the parameters κ ∗ t,h and ω ∗ t,h defined in (10).
As for the corresponding densities of the forward process X , they are Gaussian when conditioned on the initial data point X0:
Law (Xt−h|Xt, X0) = N (µt−h,tXt + νt−h,tX0, σ2t−h,t I), (42)
where coefficients µt−h,t, νt−h,t and σt−h,t are defined in (9). This formula for Law (Xt−h|Xt, X0) follows from the general fact about Gaussian distributions appearing in many recent works on diffusion probabilistic modeling (Kingma et al., 2021): if Zt|Zs ∼ N (αt|sZs, σ2t|s I) and Zt|Z0 ∼ N (αt|0Z0, σ2t|0 I) for 0 < s < t, then
Law (Zs|Zt, Z0) = N
( σ2s|0
σ2t|0 αt|sZt +
σ2t|s σ2t|0 αs|0Z0, σ2s|0σ 2 t|s σ2t|0 I
) . (43)
This fact is a result of applying Bayes formula to normal distributions. In our case αt|s = γs,t and σ2t|s = 1− γ 2 s,t.
To get an expression for the densities pt−h|t(xt−h|xt) similar to (39), we need to integrate out the dependency on data X0 from Gaussian distribution Law (Xt−h|Xt, X0):
pt−h|t(xt−h|xt) = ∫ pt−h,0|t(xt−h, x0|xt)dx0 = ∫ pt−h|t,0(xt−h|xt, x0)p0|t(x0|xt)dx0
= EX0∼p0|t(·|xt)[pt−h|t,0(xt−h|xt, X0)], (44)
which implies the following formula:
pt−h|t(xt−h|xt) = 1√
2πσnt−h,t Ep0|t(·|xt)
[ exp ( −‖xt−h − µt−h,txt − νt−h,tX0‖ 2 2
2σ2t−h,t
)] . (45)
Note that in contrast with the transition densities (39) of the reverse process X̂ , the corresponding densities (45) of the forward process X are not normal in general.
Our goal is to find parameters κ̂, ω̂ and σ̂ that maximize log-likelihood of sample paths X under probability measure with transition densities p̂. Put tk = kh for k = 0, 1, .., N and write down this log-likelihood:
∫ p(x1, x1−h, .., x0) ( N−1∑ k=0 log p̂tk|tk+1(xtk |xtk+1) + log p̂1(x1) ) dx1dx1−h..dx0
= N−1∑ k=0 ∫ p(xtk , xtk+1) log p̂tk|tk+1(xtk |xtk+1)dxtk+1dxtk + ∫ p(x1) log p̂1(x1)dx1.
(46)
The last term does not depend on κ̂, ω̂ and σ̂, so we can ignore it. Let Rk be the k-th term in the sum above. Since we are free to have different coefficients κ̂t,h, ω̂t,h and σ̂t,h for different steps, we can maximize each Rk separately. Terms Rk can be expressed as
Rk = ∫ p(xtk , xtk+1) log p̂tk|tk+1(xtk |xtk+1)dxtk+1dxtk
= ∫ p(xtk+1)ptk|tk+1(xtk |xtk+1) log p̂tk|tk+1(xtk |xtk+1)dxtk+1dxtk
= EXtk+1 [∫ ptk|tk+1(xtk |Xtk+1) log p̂tk|tk+1(xtk |Xtk+1)dxtk ] .
(47)
From now on we will skip subscripts of µ, ν, σ, µ̂, ν̂, σ̂, κ̂, ω̂, κ∗ and ω∗ for brevity. Denote
Q(xtk , Xtk+1 , X0) = 1√
2πσn exp
( − ‖xtk − µXtk+1 − νX0‖22
2σ2
) log p̂tk|tk+1(xtk |Xtk+1). (48)
Using the formula (44) for the densities of X together with the explicit expression for the Gaussian density ptk|tk+1,0(xtk |Xtk+1 , X0) and applying Fubini’s theorem to change the order of integration, we rewrite Rk as
Rk = EXtk+1 [∫ ptk|tk+1(xtk |Xtk+1) log p̂tk|tk+1(xtk |Xtk+1)dxtk ] = EXtk+1 [∫ EX0∼p0|tk+1 (·|Xtk+1 ) [ ptk|tk+1,0(xtk |Xtk+1 , X0) log p̂tk|tk+1(xtk |Xtk+1) ] dxtk
] = EXtk+1 [∫ EX0∼p0|tk+1 (·|Xtk+1 )[Q(xtk , Xtk+1 , X0)]dxtk
] = EXtk+1EX0∼p0|tk+1 (·|Xtk+1 ) [∫ Q(xtk , Xtk+1 , X0)dxtk ] .
(49)
The formula (48) implies that the integral of Q(xtk , Xtk+1 , X0) with respect to xtk can be seen as expectation of log p̂tk|tk+1(ξ|Xtk+1) with respect to normal random variable ξ with mean µXtk+1 + νX0 and covariance matrix σ2 I. Plugging in the expression (39) into (48), we can calculate this integral:
Eξ [ − log √ 2π − n log σ̂ − ‖ξ − µ̂Xtk+1 − ν̂EX′0∼p0|tk+1 (·|Xtk+1 )X ′ 0‖22
2σ̂2
]
= − log √ 2π − n log σ̂ − Eξ‖ξ − µ̂Xtk+1 − ν̂Ep0|tk+1 (·|Xtk+1 )X ′ 0‖22
2σ̂2 .
(50)
Thus, terms Rk we want to maximize equal
Rk = − log √ 2π−n log σ̂−EXtk+1EX0∼p0|tk+1 (·|Xtk+1 ) Eξ‖ξ − µ̂Xtk+1 − ν̂Ep0|tk+1 (·|Xtk+1 )X ′ 0‖22
2σ̂2 (51)
Maximizing Rk with respect to (κ̂, ω̂, σ̂) is equivalent to minimizing EXtk+1Sk where Sk is given by
Sk = n log σ̂ + 1
2σ̂2 EX0∼p0|tk+1 (·|Xtk+1 )Eξ‖ξ − µ̂Xtk+1 − ν̂Ep0|tk+1 (·|Xtk+1 )X ′ 0‖22, (52)
where the expectation with respect to ξ ∼ N (µXtk+1 + νX0, σ2 I) can be calculated using the fact that for every vector â we can express Eξ‖ξ − â‖22 as
E‖ξ−Eξ+Eξ−â‖22 = E‖ξ−Eξ‖22+2〈E[ξ−Eξ],Eξ−â〉+E‖Eξ−â‖22 = nσ2+‖Eξ−â‖22. (53)
So, the outer expectation with respect to Law(X0|Xtk+1) in (52) can be simplified:
EX0∼p0|tk+1 (·|Xtk+1 ) [ nσ2 + ‖(µ− µ̂)Xtk+1 + νX0 − ν̂EX′0∼p0|tk+1 (·|Xtk+1 )X ′ 0‖22 ]
= nσ2 + EX0‖((µ− µ̂)Xtk+1 + νX0 − ν̂EX′0X ′ 0‖22 = nσ2 + (µ− µ̂)2‖Xtk+1‖22 + 2〈(µ− µ̂)Xtk+1 , (ν − ν̂)EX0X0〉+ EX0‖νX0 − ν̂EX′0X ′ 0‖22 = (µ− µ̂)2‖Xtk+1‖22 + 2〈(µ− µ̂)Xtk+1 , (ν − ν̂)EX0X0〉+ ν2EX0‖X0‖22 + ν̂2‖EX0X0‖22 + nσ2 − 2νν̂〈EX0X0,EX′0X ′ 0〉 = ‖(µ− µ̂)Xtk+1 + (ν − ν̂)EX0X0)‖22 + ν2EX0‖X0‖22 − ν2‖EX0X0‖22 + nσ2,
(54)
where all the expectations in the formula above are taken with respect to the conditional data distribution Law(X0|Xtk+1). So, the resulting expression for the terms Sk whose expectation with respect to Law (Xtk+1) we want to minimize is
Sk = n logσ̂ + 1
2σ̂2
( nσ2 + ‖(µ− µ̂)Xtk+1 + (ν − ν̂)E[X0|Xtk+1 ]‖22
+ ν2 ( E [ ‖X0‖22|Xtk+1 ] − ‖E[X0|Xtk+1 ]‖22 )) .
(55)
Now it is clear that κ∗tk+1,h and ω ∗ tk+1,h are optimal because µ̂tk+1,h(κ ∗ tk+1,h , ω∗tk+1,h) = µtk,tk+1 and ν̂tk+1,h(κ ∗ tk+1,h ) = νtk,tk+1 . For this choice of parameters we have
EXtk+1Sk = n log σ̂ + 1
2σ̂2
( nσ2 + ν2EXtk+1 [ E [ ‖X0‖22|Xtk+1 ] − ‖E[X0|Xtk+1 ]‖22 ]) . (56)
Note that E [ ‖X0‖22|Xtk+1 ] − ‖E[X0|Xtk+1 ]‖22 = Tr (Var (X0|Xtk+1)) is the overall variance of Law (X0|Xtk+1) along all n dimensions. Differentiating EXtk+1Sk with respect to σ̂ shows that the optimal σ∗tk+1,h should satisfy
n σ∗tk+1,h − 1 (σ∗tk+1,h) 3
( nσ2tk,tk+1 + ν 2 tk,tk+1 EXtk+1 [ Tr (Var (X0|Xtk+1)) ]) = 0, (57)
which is indeed satisfied by the parameters σ∗t,h defined in (10). Thus, the statement (i) is proven.
When it comes to proving that X̂ is exact, we have to show that Law (X̂tk) = Law (Xtk) for every k = 0, 1, .., N . By the assumption that Law (X̂1) = Law (X1) it is sufficient to prove that p̂tk|tk+1(xtk |xtk+1) ≡ ptk|tk+1(xtk |xtk+1) since the exactness will follow from this fact by mathematical induction. If X0 is a constant random variable, Law(X0|Xt) = Law (X0) also corresponds to the same constant, so Var (X0|Xt) = 0 meaning that σ∗t,h = σt−h,t, and the formulae (39) and (45) imply the desired result.
Let us now consider the second case when X0 ∼ N (µ̄, δ2 I). It is a matter of simple but lengthy computations to prove another property of Gaussian distributions similar to (43): if Z0 ∼ N (µ̄, δ2 I) and Zt|Z0 ∼ N (atZ0, b2t I), then Z0|Zt ∼ N (
b2t b2t+δ 2a2t µ̄ + δ 2at b2t+δ 2a2t Zt, δ2b2t b2t+δ 2a2t I)
and Zt ∼ N (µ̄at, (b2t + δ2a2t ) I). In our case at = γ0,t and b2t = 1− γ20,t, therefore
Law (X0|Xt) = N
( 1− γ20,t
1− γ20,t + δ2γ20,t µ̄+ δ2γ0,t 1− γ20,t + δ2γ20,t Xt, δ2(1− γ20,t) 1− γ20,t + δ2γ20,t I
) . (58)
So, Var (X0|Xt) does not depend on Xt and
(σ∗t,h) 2 = σ2t−h,t + ν2t−h,t n EXt [Tr (Var (X0|Xt))] = σ2t−h,t + ν2t−h,t δ2(1− γ20,t) 1− γ20,t + δ2γ20,t . (59)
Since Law (Xt|Xt−h), Law (Xt−h) and Law (Xt) are Gaussian, Bayes formula implies that Law (Xt−h|Xt) is Gaussian as well with the following mean and covariance matrix:
E[Xt−h|Xt] = γ0,t−h(1− γ2t−h,t) 1− γ20,t + δ2γ20,t µ̄+ γt−h,t(1− γ20,t−h + δ2γ20,t−h) 1− γ20,t + δ2γ20,t Xt, (60)
Var (Xt−h|Xt) = (1− γ2t−h,t)(1− γ20,t−h + δ2γ20,t−h)
1− γ20,t + δ2γ20,t I . (61)
The distribution Law (X̂t−h|X̂t) is also Gaussian by the formula (39), so to conclude the proof we just need to show that E[X̂t−h|X̂t = x] = E[Xt−h|Xt = x] and Var (X̂t−h|X̂t = x) = Var (Xt−h|Xt = x) for every x ∈ Rn for the optimal parameters (10). Recall that for κ∗t,h and ω∗t,h we have µ̂t,h(κ∗t,h, ω ∗ t,h) = µt−h,t and ν̂t,h(κ ∗ t,h) = νt−h,t. Utilizing the formulae (9), (39), (58) and the fact that γ0,t−h · γt−h,t = γ0,t (following from the definition of γ in (7)) we conclude that
E[X̂t−h|X̂t = x] = µ̂t,h(κ∗t,h, ω∗t,h)x+ ν̂t,h(κ∗t,h)Ep0|t(·|x)X0
= γt−h,t 1− γ20,t−h 1− γ20,t x+ γ0,t−h 1− γ2t−h,t 1− γ20,t
[ 1− γ20,t
1− γ20,t + δ2γ20,t µ̄+ δ2γ0,t 1− γ20,t + δ2γ20,t x
]
= γ0,t−h 1− γ2t−h,t
1− γ20,t + δ2γ20,t µ̄+ γt−h,t (1− γ20,t−h)(1− γ20,t) + δ2γ20,t(1− γ20,t−h) (1− γ20,t)(1− γ20,t + δ2γ20,t) x
+ γt−h,t δ2γ20,t−h(1− γ2t−h,t)
(1− γ20,t)(1− γ20,t + δ2γ20,t) x = γ0,t−h 1− γ2t−h,t 1− γ20,t + δ2γ20,t µ̄
+ γt−h,t (1− γ20,t−h)(1− γ20,t) + δ2γ20,t−h(1− γ20,t)
(1− γ20,t)(1− γ20,t + δ2γ20,t) x = γ0,t−h 1− γ2t−h,t 1− γ20,t + δ2γ20,t µ̄
+ γt−h,t 1− γ20,t−h + δ2γ20,t−h
1− γ20,t + δ2γ20,t x = E[Xt−h|Xt = x],
(62)
Var (X̂t−h|X̂t = x) = (σ∗t,h)2 I = ( σ2t−h,t + ν 2 t−h,t
δ2(1− γ20,t) 1− γ20,t + δ2γ20,t
) I
=
( (1− γ20,t−h)(1− γ2t−h,t)
1− γ20,t + γ20,t−h
δ2(1− γ2t−h,t)2
(1− γ20,t)(1− γ20,t + δ2γ20,t)
) I
= 1− γ2t−h,t (1− γ20,t)(1− γ20,t + δ2γ20,t) ( (1− γ20,t−h)(1− γ20,t) + δ2γ20,t(1− γ20,t−h) ) I
+ 1− γ2t−h,t (1− γ20,t)(1− γ20,t + δ2γ20,t) ( δ2γ20,t−h(1− γ2t−h,t) ) I
= 1− γ2t−h,t (1− γ20,t)(1− γ20,t + δ2γ20,t) ( (1− γ20,t−h)(1− γ20,t) + δ2γ20,t−h(1− γ20,t) ) I
= (1− γ2t−h,t)(1− γ20,t−h + δ2γ20,t−h)
1− γ20,t + δ2γ20,t I = Var (Xt−h|Xt = x).
(63)

heading: D REVERSE MR-VP SDE SOLVER
text: MR-VP DPM is characterized by the following forward and reverse diffusions:
dXt = 1
2 βt(X̄ −Xt)dt+
√ βtd −→ Wt , (64)
dX̂t =
( 1
2 βt(X̄ − X̂t)− βtsθ(X̂t, X̄, t)
) dt+ √ βtd ←− Wt. (65)
Using the same method as in Appendix A, we can show that for s < t
Law (Xt|Xs) = N (γs,tXs + (1− γs,t)X̄, (1− γ2s,t) I), γs,t = e− 1 2
∫ t s βudu. (66)
With the following notation:
µs,t = γs,t 1− γ20,s 1− γ20,t , νs,t = γ0,s 1− γ2s,t 1− γ20,t , σ2s,t = (1− γ20,s)(1− γ2s,t) 1− γ20,t (67)
we can write down the parameters of Gaussian distribution Xs|Xt, X0:
E[Xs|Xt, X0] = X̄ + µs,t(Xt − X̄) + νs,t(X0 − X̄), Var (Xs|Xt, X0) = σ2s,t I . (68)
The Lemma 1 for MR-VP DPMs takes the following shape:
sθ∗(x, X̄, t) = − 1
1− γ20,t
( x− (1− γ0,t)X̄ − γ0,tEp0|t(·|x)X0 ) = − 1
1− γ20,t
( (x− X̄)− γ0,t ( Ep0|t(·|x)X0 − X̄ )) .
(69)
The class of reverse SDE solvers we consider is
X̂t−h = X̂t + βth
(( 1
2 + ω̂t,h
) (X̂t − X̄) + (1 + κ̂t,h)sθ(X̂t, X̄, t) ) + σ̂t,hξt, (70)
where t = 1, 1 − h, .., h and ξt are i.i.d. samples from N (0, I). Repeating the argument of the Theorem 1 leads to the following optimal (in terms of likelihood of the forward diffusion sample paths) parameters:
κ∗t,h = νt−h,t(1− γ20,t)
γ0,tβth − 1, ω∗t,h = µt−h,t − 1 βth + 1 + κ∗t,h 1− γ20,t − 1 2 ,
(σ∗t,h) 2 = σ2t−h,t +
1 n ν2t−h,tEXt [Tr (Var (X0|Xt))] ,
(71)
which are actually the same as the optimal parameters (10) for VP DPM. It is of no surprise since MR-VP DPM and VP-DPM differ only by a constant shift.

heading: E REVERSE SUB-VP SDE SOLVER
text: Sub-VP DPM is characterized by the following forward and reverse diffusions:
dXt = − 1
2 βtXtdt+
√ βt(1− e−2 ∫ t 0 βudu)d −→ Wt , (72)
dX̂t =
( −1
2 βtX̂t − βt
( 1− e−2 ∫ t 0 βudu ) sθ(X̂t, t) ) dt+ √ βt(1− e−2 ∫ t 0 βudu)d ←− Wt. (73)
Using the same method as in Appendix A, we can show that for s < t
Law (Xt|Xs) = N (γs,tXs, ( 1 + γ40,t − γ2s,t(1 + γ40,s) ) I), γs,t = e − 12 ∫ t s βudu. (74)
Note that for s = 0 this expression simplifies to
Law (Xt|X0) = N (γ0,tX0, (1− γ20,t)2 I). (75)
With the following notation:
µs,t =γs,t ( 1− γ20,s 1− γ20,t )2 , νs,t = γ0,s 1 + γ40,t − γ2s,t(1 + γ40,s) (1− γ20,t)2 ,
σ2s,t = (1− γ20,s)2(1 + γ40,t − γ2s,t(1 + γ40,s))
(1− γ20,t)2
(76)
we can write down the parameters of Gaussian distribution Xs|Xt, X0:
E[Xs|Xt, X0] = µs,tXt + νs,tX0, Var (Xs|Xt, X0) = σ2s,t I . (77)
The Lemma 1 for sub-VP DPMs takes the following shape:
sθ∗(x, t) = − 1 (1− γ20,t)2 ( x− γ0,tEp0|t(·|x)X0 ) . (78)
The class of reverse SDE solvers we consider is
X̂t−h = X̂t + βth
(( 1
2 + ω̂t,h
) X̂t + (1 + κ̂t,h) ( 1− e−2 ∫ t 0 βudu ) sθ(X̂t, t) ) + σ̂t,hξt, (79)
where t = 1, 1 − h, .., h and ξt are i.i.d. samples from N (0, I). Repeating the argument of the Theorem 1 leads to the following optimal (in terms of likelihood of the forward diffusion sample paths) parameters:
κ∗t,h = νt−h,t(1− γ20,t) γ0,tβth(1 + γ20,t) − 1, ω∗t,h = µt−h,t − 1 βth + (1 + κ∗t,h)(1 + γ 2 0,t) 1− γ20,t − 1 2 ,
(σ∗t,h) 2 = σ2t−h,t +
1 n ν2t−h,tEXt [Tr (Var (X0|Xt))] .
(80)

heading: F REVERSE VE SDE SOLVER
text: VE DPM is characterized by the following forward and reverse diffusions:
dXt = √ (σ2t ) ′ d −→ Wt , (81)
dX̂t = − ( σ2t )′ sθ(X̂t, t)dt+ √ (σ2t ) ′ d ←− Wt. (82)
Since for s < t
Xt = Xs + ∫ t s √ (σ2u) ′ d −→ Wu , (83)
similar argument as in Appendix A allows showing that
Law (Xt|Xs) = N ( Xs, I · ∫ t s ( σ2u )′ du ) = N (Xs, (σ2t − σ2s) I). (84)
With the following notation:
µs,t = σ2s − σ20 σ2t − σ20 , νs,t = σ2t − σ2s σ2t − σ20 , σ2s,t = (σ2t − σ2s)(σ2s − σ20) σ2t − σ20 , (85)
we can write down the parameters of Gaussian distribution Xs|Xt, X0:
E[Xs|Xt, X0] = µs,tXt + νs,tX0, Var (Xs|Xt, X0) = σ2s,t I . (86)
The Lemma 1 for VE DPMs takes the following shape:
sθ∗(x, t) = − 1
σ2t − σ20
( x− Ep0|t(·|x)X0 ) . (87)
Repeating the argument of the Theorem 1 leads to the following optimal (in terms of likelihood of the forward diffusion sample paths) reverse SDE solver:
X̂t−h = X̂t + (σ 2 t − σ2t−h)sθ(X̂t, t) + σ∗t,hξt, (88)
where
(σ∗t,h) 2 = (σ2t − σ2t−h)(σ2t−h − σ20) σ2t − σ20 + 1 n ν2t−h,tEXt [Tr (Var (X0|Xt))] , (89)
t = 1, 1− h, .., h and ξt are i.i.d. samples from N (0, I).

heading: G TOY EXAMPLES
text: In this section we consider toy examples where data distribution X0 is represented by a single point (corresponding to the case (ii) of the Theorem 1) and by two points (corresponding to more general case (i)). In the first case the point is unit vector i = (1, 1, .., 1) of dimensionality 100, in the second one two points i and −2i have the same probability. We compare performance of two solvers, Euler-Maruyama and the proposed Maximum Likelihood, depending on the number N ∈ {1, 2, 5, 10, 100, 1000} of solver steps. The output of the perfectly trained score matching network sθ∗ is computed analytically and Gaussian noise with variance ε ∈ {0.0, 0.1, 0.5} is added to approximate the realistic case when the network sθ we use is not trained till optimality. We considered VP diffusion model (6) with β0 = 0.05 and β1 = 20.0.
The results of the comparison are given in Table 5 and can be summarized in the following:
• for both methods, larger N means better quality;
• for both methods, more accurate score matching networks (smaller ε) means better quality;
• for large number of steps, both methods perform the same;
• it takes less number of steps for the proposed Maximum Likelihood solver to converge with a good accuracy to data distribution than it does for Euler-Maruyama solver;

heading: N = 1 conv / div div / div div / div div / div div / div div / div
text: Also, in the second example where X0 ∈ {i,−2i} the Maximum Likelihood SDE solver reconstructs the probabilities of these two points better than Euler-Maruyama which tends to output “i-samples” (which are closer to the origin) more frequently than “−2i-samples”. E.g. for ε = 0.0 and N = 10 the frequency of “i-samples” generated by Euler-Maruyama scheme is 54% while this frequency for Maximum Likelihood scheme is 50% (500k independent runs were used to calculate these frequencies).

heading: H SPEAKER CONDITIONING NETWORK
text: The function x · tanh(softplus(x)) is used as a non-linearity in the speaker conditioning network gt(Y ). First, time embedding te is obtained by the following procedure: time t ∈ [0, 1] is encoded with positional encoding (Song et al., 2021c), then resulting 256-dimensional vector t′ is passed through the first linear module with 1024 units, then a non-linearity is applied to it and then it is passed through the second linear module with 256 units. Next, noisy mel-spectrogram Yt for wodyn input type or Yt concatenated with {Ys|s = 0.5/15, 1.5/15, .., 14.5/15} for whole is passed through 6 blocks consisting of 2D convolutional layers each followed by instance normalization and Gated Linear Unit. The number of input and output channels of these convolutions is (1, 64), (32, 64), (32, 128), (64, 128), (64, 256), (128, 256) for wodyn input type and the same but with 16 input channels in the first convolution for whole input type. After the 2nd and 4th blocks MLP1(te) and MLP2(te) are broadcast-added where MLP1 (MLP2) are composed of a nonlinearity followed by a linear module with 32 (64) units. After the last 6th block the result is passed through the final convolution with 128 output channels and average pooling along both time and frequency axes is applied resulting in 128-dimensional vector. All convolutions except for the final one have (kernel, stride, zero padding) = (3, 1, 1) while for the final one the corresponding parameters are (1, 0, 0). Denote the result of such processing of Y by c for wodyn and whole input types.
Clean target mel-spectrogram Y0 is used to obtain 256-dimensional speaker embedding d with the pre-trained speaker verification network (Jia et al., 2018) which is not trained. Vectors d, c and t′ are concatenated (except for d-only input type where we concatenate only d and t′), passed through a linear module with 512 units followed by a non-linearity and a linear module with 128 units. The resulting 128-dimensional vector is the output of the speaker conditioning network gt(Y ).

heading: I TRAINING HYPERPARAMETERS AND OTHER DETAILS
text: Encoders and decoders were trained with batch sizes 128 and 32 and Adam optimizer with initial learning rates 0.0005 and 0.0001 correspondingly. Encoders and decoders in VCTK models were trained for 500 and 200 epochs respectively; as for LibriTTS models, they were trained for 300 and
110 epochs. The datasets were downsampled to 22.05kHz which was the operating rate of our VC models. VCTK recordings were preprocessed by removing silence in the beginning and in the end of utterances. To fit GPU memory, decoders were trained on random speech segments of approximately 1.5 seconds rather than on the whole utterances. Training segments for reconstruction and the ones used as input to the speaker conditioning network gt(Y ) were different random segments extracted from the same training utterances. Noise schedule parameters β0 and β1 were set to 0.05 and 20.0.
Our VC models operated on mel-spectrograms with 80 mel features and sampling rate 22.05kHz. Short-Time Fourier Transform was used to calculate spectra with 1024 frequency bins. Hann window of length 1024 was applied with hop size 256.
For Diff-LibriTTS models we used simple spectral subtraction algorithm in mel domain with spectral floor parameter β = 0.02 as post-processing to reduce background noise sometimes produced by these models. Noise spectrum was estimated on speech fragments automatically detected as the ones corresponding to silence in source mel-spectrogram.

heading: J DETAILS OF AMT TESTS
text: For fair comparison with the baselines all the recordings were downsampled to 16kHz; we also normalized their loudness. In speech naturalness tests workers chosen by geographic criterion were asked to assess the overall quality of the synthesized speech, i.e to estimate how clean and natural (human-sounding) it was. Five-point Likert scale was used: 1 - “Bad”, 2 - “Poor”, 3 - “Fair”, 4 - “Good”, 5 - “Excellent”. Assessors were asked to wear headphones and work in a quiet environment. As for speaker similarity tests, workers were asked to assess how similar synthesized samples sounded to target speech samples in terms of speaker similarity. Assessors were asked not to pay attention to the overall quality of the synthesized speech (e.g. background noise or incorrect pronunciation). Five-point scale was used: 1 - “Different: absolutely sure”, 2 - “Different: moderately sure”, 3 - “Cannot decide more same or more different”, 4 - “Same: moderately sure”, 5 - “Same: absolutely sure”.

Review summaries:
This paper proposes a method to perform voice conversion through recent methods in Diffusion Based Probability Modeling (DPM) through Scochastic Differential Equations (SDE), building upon recent work Grad-TTS and Glow-TTS. The architecture is an encoder-decoder setup, trained separately, described in upcoming paragraphs. Viewed at a high level, the encoder maps input to 'noise', and the decoder inverts noise to output, in keeping with the SDE formalism. To transform average features to output voice, speaker info is conditioned on the average voice features before feeding to decoder.   The encoder learns to map input voice features to an average mel spectrogram. In order to do this, the Montreal Forced Aligner is used to align mel frames to phonemes. The average output voice (for encoder training) is obtained by aligning and averaging out corresponding phonemes in the dataset.   The main novelty in this paper is claimed to be in the decoder setup solving the reverse SDE. The authors derive a modified SDE that does MLE inference on the original reverse SDE. Through this approach, the claims are that the likelihoods are better, and its estimation is 'faster' (ostensibly, owing to getting around step size limitations).   Evaluations are carried out for VCTK and LibriTTS, and quite well done, and the samples provided show that the approach works. Comparisons are made against several other approaches, include subjective evaluations and FID scores, and look for similarity and naturalness metrics. They also compare against different flavors of SDE setups (Variance Exploding, Preserving, etc.).   

The paper proposes a diffusion probabilistic model-based voice conversion method for one-shot voice conversion scenario.  The proposed method can generate high-quality converted speech compared to state-of-the-art approaches.  Furthermore, to improve the real-time factor, a novel stochastic differential equations solver is proposed which makes the diffusion model faster.  The proposed solver is also suitable for other generative tasks. 

This paper tackles the problem of one-shot many-to-many voice conversion (with both unseen source and target speakers) using a diffusion model on mel-spectrograms. The authors propose a dedicated architecture able to generalize to unseen speakers naturally (without relying on Phoetic Posteriorgrams (PPGs) as in previous works) by conditioning the diffusion model on "average-phoneme" spectrograms together with the target speech. They also address the important problem of the errors arising from using the Euler-Maruyama scheme with large discretization steps by proposing a novel SDE solver, allowing to perform voice conversion with as few as 6 diffusion steps. They demonstrate the applicability of this solver on other modalities (CIFAR-10) using on unconditional model. Extensive and rigorous experiments are conducted on the voice conversion task and an accompanying web page is provided with numerous and convincing voice conversions.