A Message-Passing Approach to Phase Retrieval of Sparse Signals

Abstract

In phase retrieval, the goal is to recover a signal $\boldsymbol{x} \in \mathbb{C}^{N}$ from the magnitudes of linear measurements $\boldsymbol{Ax} \in \mathbb{C}^{M}$. While recent theory has established that  M ≈ 4 N intensity measurements are necessary and sufficient to recover generic $\boldsymbol{x}$, there is great interest in reducing the number of measurements through the exploitation of sparse $\boldsymbol{x}$, which is known as compressive phase retrieval. In this work, we detail a novel, probabilistic approach to compressive phase retrieval based on the generalized approximate message passing (GAMP) algorithm. We then present a numerical study of the proposed PR-GAMP algorithm, demonstrating its excellent phase-transition behavior, robustness to noise, and runtime. For example, to successfully recover  K-sparse signals, approximately $M \geq 2K\log _{2}(N/K)$ intensity measurements suffice when  K ≪  N and $\boldsymbol{A}$ has i.i.d Gaussian entries. When recovering a 6k-sparse 65k-pixel grayscale image from 32k randomly masked and blurred Fourier intensity measurements, PR-GAMP achieved 99% success rate with a median runtime of only 12. 6 seconds. Compared to the recently proposed CPRL, sparse-Fienup, and GESPAR algorithms, experiments show that PR-GAMP has a superior phase transition and orders-of-magnitude faster runtimes as the problem dimensions increase.

Appendices

Appendix A: Output Thresholding Rules

In this appendix, we derive expressions (10) and (12) that are used to compute the functions  g _out,  m and  g′ _out,  m defined in lines (D2) and (D3) of Table 1.

To facilitate the derivations in this appendix,¹³ we first rewrite  p  _Y |  Z ( y |  z) in a form different from (8). In particular, recalling that—under our AWGN assumption—the noisy transform outputs $u = z + w$ are conditionally distributed as $p(u\vert z) = \mathcal{N}(u;z,\nu ^{w})$, we first transform $u = ye^{j\theta }$ from rectangular to polar coordinates to obtain

$$\displaystyle{ p(y,\theta \vert z) = 1_{y\geq 0}1_{\theta \in [0,2\pi )}\,\mathcal{N}(ye^{j\theta };z,\nu ^{w})\,y, }$$

(21)

where  y is the Jacobian of the transformation, and then integrate out the unobserved phase $\theta$ to obtain

$$\displaystyle{ p_{Y \vert Z}(y\vert z) = 1_{y\geq 0}\,y\int _{0}^{2\pi }\mathcal{N}(ye^{j\theta };z,\nu ^{w})\,d\theta. }$$

(22)

We begin by deriving the integration constant

$$\displaystyle\begin{array}{rcl} & & C(y,\nu ^{w},\widehat{p},\nu ^{p}) \triangleq \int _{ \mathbb{C}}p_{Y \vert Z}(y\vert z)\,\mathcal{N}(z;\widehat{p},\nu ^{p})dz \\ & & = y\,1_{y\geq 0}\int _{0}^{2\pi }\int _{ \mathbb{C}}\mathcal{N}(ye^{j\theta };z,\nu ^{w})\mathcal{N}(z;\widehat{p},\nu ^{p})dzd\theta {}\end{array}$$

(23)

$$\displaystyle\begin{array}{rcl} & & = y\,1_{y\geq 0}\int _{0}^{2\pi }\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w} +\nu ^{p})d\theta,{}\end{array}$$

(24)

where we used the Gaussian-pdf multiplication rule¹⁴ in (24). Noting the similarity between (24) and (22), the equivalence between (22) and (8) implies that

$$\displaystyle\begin{array}{rcl} C(y,\nu ^{w},\widehat{p},\nu ^{p}) = \frac{2y} {\nu ^{w} +\nu ^{p}}\exp \bigg(-\frac{y^{2} + \vert \widehat{p}\vert ^{2}} {\nu ^{w} +\nu ^{p}} \bigg)I_{0}\bigg( \frac{2y\vert \widehat{p}\vert } {\nu ^{w} +\nu ^{p}}\bigg)\,1_{y\geq 0}.& &{}\end{array}$$

(25)

In the sequel, we make the practical assumption that  y > 0, allowing us to drop the indicator “1  _y ≥ 0” and invert  C.

Next, we derive the conditional mean

$$\displaystyle\begin{array}{rcl} \mathrm{E}_{Z\vert Y,P}\{Z\vert y,\widehat{p};\nu ^{p}\} = C(y,\nu ^{w},\widehat{p},\nu ^{p})^{-1}\int _{ \mathbb{C}}z\,p_{Y \vert Z}(y\vert z;\nu ^{w})\mathcal{N}(z;\widehat{p},\nu ^{p})dz.& &{}\end{array}$$

(26)

Plugging (22) into (26) and applying the Gaussian-pdf multiplication rule,

$$\displaystyle\begin{array}{rcl} & & \mathrm{E}_{Z\vert Y,P}\{Z\vert y,\widehat{p};\nu ^{p}\} \\ & & = C^{-1}y\int _{ 0}^{2\pi }\int _{ \mathbb{C}}z\,\mathcal{N}(z;ye^{j\theta },\nu ^{w})\mathcal{N}(z;\widehat{p},\nu ^{p})dzd\theta {}\end{array}$$

(27)

$$\displaystyle\begin{array}{rcl} & & = C^{-1}y\int _{ 0}^{2\pi }\int _{ \mathbb{C}}z\,\mathcal{N}\bigg(z; \frac{ye^{j\theta }/\nu ^{w} +\widehat{ p}/\nu ^{p}} {1/\nu ^{w} + 1/\nu ^{p}}, \frac{1} {1/\nu ^{w} + 1/\nu ^{p}}\bigg) \\ & & \quad \times \mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})dzd\theta {}\end{array}$$

(28)

$$\displaystyle\begin{array}{rcl} & & = C^{-1}y\int _{ 0}^{2\pi }\frac{ye^{j\theta }/\nu ^{w} +\widehat{ p}/\nu ^{p}} {1/\nu ^{w} + 1/\nu ^{p}} \mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta {}\end{array}$$

(29)

$$\displaystyle\begin{array}{rcl} & & = \frac{y/\nu ^{w}} {1/\nu ^{w} + 1/\nu ^{p}}C^{-1}y\int _{ 0}^{2\pi }e^{j\theta }\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta \\ & & \quad + \frac{\widehat{p}/\nu ^{p}} {1/\nu ^{w} + 1/\nu ^{p}}C^{-1}y\int _{ 0}^{2\pi }\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta {}\end{array}$$

(30)

$$\displaystyle\begin{array}{rcl} = \frac{y} {\nu ^{w}/\nu ^{p} + 1}C^{-1}y\int _{ 0}^{2\pi }e^{j\theta }\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta + \frac{\widehat{p}} {\nu ^{p}/\nu ^{w} + 1}.& &{}\end{array}$$

(31)

Expanding the $\mathcal{N}$ term, the integral in (31) becomes

$$\displaystyle\begin{array}{rcl} & & \int _{0}^{2\pi }e^{j\theta }\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta \\ & & = \frac{1} {\pi (\nu ^{w} +\nu ^{p})}\exp \bigg(-\frac{y^{2} + \vert \widehat{p}\vert ^{2}} {\nu ^{w} +\nu ^{p}} \bigg)\int _{0}^{2\pi }e^{j\theta }\exp \bigg( \frac{2y\vert \widehat{p}\vert } {\nu ^{w} +\nu ^{p}}\cos (\theta -\psi )\bigg)d\theta {}\end{array}$$

(32)

$$\displaystyle\begin{array}{rcl} & & = \frac{1} {\pi (\nu ^{w} +\nu ^{p})}\exp \bigg(-\frac{y^{2} + \vert \widehat{p}\vert ^{2}} {\nu ^{w} +\nu ^{p}} \bigg)e^{j\psi }\int _{ 0}^{2\pi }e^{j\theta '}\exp \bigg( \frac{2y\vert \widehat{p}\vert } {\nu ^{w} +\nu ^{p}}\cos (\theta ')\bigg)d\theta '{}\end{array}$$

(33)

$$\displaystyle\begin{array}{rcl} & & = \frac{2e^{j\psi }} {\nu ^{w} +\nu ^{p}}\exp \bigg(-\frac{y^{2} + \vert \widehat{p}\vert ^{2}} {\nu ^{w} +\nu ^{p}} \bigg)I_{1}\bigg( \frac{2y\vert \widehat{p}\vert } {\nu ^{w} +\nu ^{p}}\bigg),{}\end{array}$$

(34)

where  ψ denotes the phase of $\widehat{p}$, and where the integral in (33) was resolved using the expression in [48, 9.6.19]. Plugging (34) into (31) gives

$$\displaystyle\begin{array}{rcl} \mathrm{E}_{Z\vert Y,P}\{Z\vert y,\widehat{p};\nu ^{p}\}& & = \frac{\widehat{p}} {\nu ^{p}/\nu ^{w} + 1} + \frac{ye^{j\psi }} {\nu ^{w}/\nu ^{p} + 1} \frac{I_{1}\big(\frac{2y\vert \widehat{p}\vert } {\nu ^{w}+\nu ^{p}} \big)} {I_{0}\big(\frac{2y\vert \widehat{p}\vert } {\nu ^{w}+\nu ^{p}} \big)},{}\end{array}$$

(35)

which agrees with (10).

Finally, we derive the conditional covariance

$$\displaystyle\begin{array}{rcl} \mathrm{var}_{Z\vert Y,P}\{Z\vert y,\widehat{p};\nu ^{p}\}& & = C(y,\nu ^{w},\widehat{p},\nu ^{p})^{-1}\int _{ \mathbb{C}}\vert z\vert ^{2}\,p_{ Y \vert Z}(y\vert z;\nu ^{w})\mathcal{N}(z;\widehat{p},\nu ^{p})dz \\ & & \quad -\vert \mathrm{E}_{Z\vert Y,P}\{Z\vert y,\widehat{p};\nu ^{p}\}\vert ^{2}. {}\end{array}$$

(36)

Focusing on the first term in (36), if we plug in (22) and apply the Gaussian-pdf multiplication rule, we get

$$\displaystyle\begin{array}{rcl} & & C(y,\nu ^{w},\widehat{p},\nu ^{p})^{-1}\int _{ \mathbb{C}}\vert z\vert ^{2}\,p_{ Y \vert Z}(y\vert z;\nu ^{w})\mathcal{N}(z;\widehat{p},\nu ^{p})dz \\ & & = C^{-1}y\int _{ 0}^{2\pi }\int _{ \mathbb{C}}\vert z\vert ^{2}\,\mathcal{N}\bigg(z; \frac{ye^{j\theta }/\nu ^{w} +\widehat{ p}/\nu ^{p}} {1/\nu ^{w} + 1/\nu ^{p}}, \frac{1} {1/\nu ^{w} + 1/\nu ^{p}}\bigg)dz \\ & & \quad \times \mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta {}\end{array}$$

(37)

$$\displaystyle\begin{array}{rcl} & & = C^{-1}y\int _{ 0}^{2\pi }\bigg(\big\vert \frac{ye^{j\theta }/\nu ^{w} +\widehat{ p}/\nu ^{p}} {1/\nu ^{w} + 1/\nu ^{p}} \big\vert ^{2} + \frac{1} {1/\nu ^{w} + 1/\nu ^{p}}\bigg)\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta {}\end{array}$$

(38)

$$\displaystyle\begin{array}{rcl} & & = C^{-1}y\int _{ 0}^{2\pi }\frac{\vert y\vert ^{2}/(\nu ^{w})^{2} + \vert \widehat{p}\vert ^{2}/(\nu ^{p})^{2} + 2y\vert \widehat{p}\vert /(\nu ^{w}\nu ^{p})\mathrm{Re}\{e^{j(\theta -\psi )}\}} {(1/\nu ^{w} + 1/\nu ^{p})^{2}} \\ & & \quad \times \mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta + \frac{1} {1/\nu ^{w} + 1/\nu ^{p}} {}\end{array}$$

(39)

$$\displaystyle\begin{array}{rcl} & & = \frac{\vert y\vert ^{2}/(\nu ^{w})^{2} + \vert \widehat{p}\vert ^{2}/(\nu ^{p})^{2}} {(1/\nu ^{w} + 1/\nu ^{p})^{2}} + \frac{1} {1/\nu ^{w} + 1/\nu ^{p}} \\ & & \quad + \frac{2y\vert \widehat{p}\vert /(\nu ^{w}\nu ^{p})} {(1/\nu ^{w} + 1/\nu ^{p})^{2}}C^{-1}y\mathrm{Re}\bigg\{e^{-j\psi }\int _{ 0}^{2\pi }e^{j\theta }\mathcal{N}(ye^{j\theta };\widehat{p},\nu ^{w}\! +\!\nu ^{p})d\theta \bigg\}{}\end{array}$$

(40)

$$\displaystyle\begin{array}{rcl} & & = \frac{\vert y\vert ^{2}/(\nu ^{w})^{2} + \vert \widehat{p}\vert ^{2}/(\nu ^{p})^{2}} {(1/\nu ^{w} + 1/\nu ^{p})^{2}} + \frac{1} {1/\nu ^{w} + 1/\nu ^{p}} \\ & & \quad + \frac{2y\vert \widehat{p}\vert /(\nu ^{w}\nu ^{p})} {(1/\nu ^{w} + 1/\nu ^{p})^{2}}C^{-1}y \frac{2} {\nu ^{w} +\nu ^{p}}\exp \bigg(-\frac{y^{2} + \vert \widehat{p}\vert ^{2}} {\nu ^{w} +\nu ^{p}} \bigg)I_{1}\bigg( \frac{2y\vert \widehat{p}\vert } {\nu ^{w} +\nu ^{p}}\bigg){}\end{array}$$

(41)

$$\displaystyle\begin{array}{rcl} & & = \frac{\vert y\vert ^{2}/(\nu ^{w})^{2} + \vert \widehat{p}\vert ^{2}/(\nu ^{p})^{2}} {(1/\nu ^{w} + 1/\nu ^{p})^{2}} + \frac{1} {1/\nu ^{w} + 1/\nu ^{p}} + \frac{2y\vert \widehat{p}\vert /(\nu ^{w}\nu ^{p})} {(1/\nu ^{w} + 1/\nu ^{p})^{2}} \frac{I_{1}\big(\frac{2y\vert \widehat{p}\vert } {\nu ^{w}+\nu ^{p}} \big)} {I_{0}\big(\frac{2y\vert \widehat{p}\vert } {\nu ^{w}+\nu ^{p}} \big)},{}\end{array}$$

(42)

where (41) used (34) and (42) used (25). By plugging (42) back into (36), we obtain the expression given in (12).

Appendix B: EM Update for Noise Variance

Noting that

$$\displaystyle\begin{array}{rcl} \ln p(\boldsymbol{y},\boldsymbol{x};\nu ^{w})& & =\ln p(\boldsymbol{y}\vert \boldsymbol{x};\nu ^{w}) +\ln p(\boldsymbol{x};\nu ^{w}) {}\end{array}$$

(43)

$$\displaystyle\begin{array}{rcl} & & =\sum _{ m=1}^{M}\ln p_{ Y \vert Z}(y_{m}\vert \boldsymbol{a}_{m}^{\mathsf{H}}\boldsymbol{x};\nu ^{w}) + \text{const} {}\end{array}$$

(44)

$$\displaystyle\begin{array}{rcl}& & =\sum _{ m=1}^{M}\ln \Big(y_{ m}\int _{0}^{2\pi }\mathcal{N}(y_{ m}e^{j\theta _{m} };\boldsymbol{a}_{m}^{\mathsf{H}}\boldsymbol{x},\nu ^{w})d\theta _{ m}\Big) + \text{const},{}\end{array}$$

(45)

where (45) used the expression for  p  _Y |  Z from (22), we have

$$\displaystyle\begin{array}{rcl} & & \mathrm{E}\big\{\ln p(\boldsymbol{y},\boldsymbol{x};\nu ^{w})\big\vert \boldsymbol{y};\widehat{\nu ^{w}}[i]\big\} \\ & & \ =\int _{\mathbb{C}^{N}}p(\boldsymbol{x}\vert \boldsymbol{y};\widehat{\nu ^{w}}[i])\sum _{ m=1}^{M}\ln \Big(\int _{ 0}^{2\pi }\mathcal{N}(y_{ m}e^{j\theta _{m} };\boldsymbol{a}_{m}^{\mathsf{H}}\boldsymbol{x},\nu ^{w})d\theta _{ m}\Big)d\boldsymbol{x}.{}\end{array}$$

(46)

To circumvent the high-dimensional integral in (46), we use the same large system limit approximation used in the derivation of GAMP [27]: for sufficiently dense $\boldsymbol{A}$, as $N \rightarrow \infty$, the central limit theorem (CLT) suggests that $\boldsymbol{a}_{m}^{\mathsf{H}}\boldsymbol{x} = z_{m}$ will become Gaussian. In particular, when $\boldsymbol{x} \sim p(\boldsymbol{x}\vert \boldsymbol{y};\widehat{\nu ^{w}}[i])$, the CLT suggests that $\boldsymbol{a}_{m}^{\mathsf{H}}\boldsymbol{x} \sim \mathcal{N}(\widehat{z}_{m},\nu _{m}^{z})$, where

$$\displaystyle\begin{array}{rcl} \widehat{z}_{m}[i]& \triangleq & \sum _{n=1}^{N}a_{ mn}\widehat{x}_{n}[i],{}\end{array}$$

(47)

$$\displaystyle\begin{array}{rcl}\nu _{m}^{z}[i]& \triangleq & \sum _{ n=1}^{N}\vert a_{ mn}\vert ^{2}\nu _{ n}^{x}[i],{}\end{array}$$

(48)

such that $\widehat{x}_{n}[i]$ and  ν  _n  ^x[ i] are the mean and variance of the marginal posterior pdf $p(x_{n}\vert \boldsymbol{y};\widehat{\nu ^{w}}[i])$. In this case,

$$\displaystyle\begin{array}{rcl} & & \mathrm{E}\big\{\ln p(\boldsymbol{y},\boldsymbol{x};\nu ^{w})\big\vert \boldsymbol{y};\widehat{\nu ^{w}}[i]\big\} \\ & & =\sum _{ m=1}^{M}\int _{ \mathbb{C}}\mathcal{N}(z_{m};\widehat{z}_{m}[i],\nu _{m}^{z}[i])\,\ln \int _{ 0}^{2\pi }\mathcal{N}(y_{ m}e^{j\theta _{m} };z_{m},\nu ^{w})d\theta _{ m}dz_{m}.{}\end{array}$$

(49)

From (14), we see that any solution $\widehat{\nu ^{w}}[i\! +\! 1]> 0$ is necessarily a value of  ν  ^w that zeros the derivative of the expected log-pdf. Thus, using the expected log-pdf expression from (49),

$$\displaystyle\begin{array}{rcl} 0& & =\sum _{ m=1}^{M}\int _{ \mathbb{C}}\mathcal{N}(z_{m};\widehat{z}_{m}[i],\nu _{m}^{z}[i])\frac{\int _{0}^{2\pi } \frac{\partial } {\partial \nu ^{w}}\mathcal{N}(y_{m}e^{j\theta _{m}};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])d\theta _{m}} {\int _{0}^{2\pi }\mathcal{N}(y_{m}e^{j\theta _{m}'};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])d\theta _{m}'} dz_{m}.{}\end{array}$$

(50)

Plugging the derivative expression (see [39])

$$\displaystyle\begin{array}{rcl} & & \frac{\partial } {\partial \nu ^{w}}\mathcal{N}(y_{m}e^{j\theta _{m} };z_{m},\widehat{\nu ^{w}}[i\! +\! 1]) \\ & & = \frac{\mathcal{N}(y_{m}e^{j\theta _{m}};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])} {\widehat{\nu ^{w}}[i\! +\! 1]^{2}} \big(\vert y_{m}e^{j\theta _{m} } - z_{m}\vert ^{2} -\widehat{\nu ^{w}}[i\! +\! 1]\big),{}\end{array}$$

(51)

into (50) and multiplying both sides by $\widehat{\nu ^{w}}[i\! +\! 1]^{2}$, we find

$$\displaystyle\begin{array}{rcl} \widehat{\nu ^{w}}[i\! +\! 1]& & = \frac{1} {M}\sum _{m=1}^{M}\int _{ \mathbb{C}}\mathcal{N}(z_{m};\widehat{z}_{m}[i],\nu _{m}^{z}[i]) \\ & & \quad \times \frac{\int _{0}^{2\pi }\vert y_{m}e^{j\theta _{m}} - z_{m}\vert ^{2}\mathcal{N}(y_{m}e^{j\theta _{m}};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])d\theta _{m}} {\int _{0}^{2\pi }\mathcal{N}(y_{m}e^{j\theta _{m}'};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])d\theta _{m}'} dz_{m}{}\end{array}$$

(52)

$$\displaystyle\begin{array}{rcl} & & = \frac{1} {M}\sum _{m=1}^{M}\int _{ \mathbb{C}}\mathcal{N}(z_{m};\widehat{z}_{m}[i],\nu _{m}^{z}[i]) \\ & & \quad \times \int _{0}^{2\pi }\vert y_{ m}e^{j\theta _{m} } - z_{m}\vert ^{2}p(\theta _{ m};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])d\theta _{ m}dz_{m}{}\end{array}$$

(53)

with the newly defined pdf

$$\displaystyle\begin{array}{rcl} p(\theta _{m};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])& \triangleq & \frac{\mathcal{N}(y_{m}e^{j\theta _{m}};z_{ m},\widehat{\nu ^{w}}[i\! +\! 1])} {\int _{0}^{2\pi }\mathcal{N}(y_{m}e^{j\theta _{m}'};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])d\theta _{m}'}{}\end{array}$$

(54)

$$\displaystyle\begin{array}{rcl} & \propto & \exp \Big(-\frac{\vert z_{m} - y_{m}e^{j\theta _{m}}\vert ^{2}} {\widehat{\nu ^{w}}[i\! +\! 1]} \Big) {}\end{array}$$

(55)

$$\displaystyle\begin{array}{rcl}& \propto & \exp \big(\kappa _{m}\cos (\theta _{m} -\phi _{m})\big)\text{ for }\kappa _{m} \triangleq \frac{2\vert z_{m}\vert y_{m}} {\widehat{\nu ^{w}}[i\! +\! 1]},{}\end{array}$$

(56)

where  ϕ  _m is the phase of  z  _m (recall (5)). The proportionality (56) identifies this pdf as a von Mises distribution [49], which can be stated in normalized form as

$$\displaystyle\begin{array}{rcl} p(\theta _{m};z_{m},\widehat{\nu ^{w}}[i\! +\! 1])& & = \frac{\exp (\kappa _{m}\cos (\theta _{m} -\phi _{m}))} {2\pi I_{0}(\kappa _{m})}.{}\end{array}$$

(57)

Expanding the quadratic in (53) and plugging in (57), we get

$$\displaystyle\begin{array}{rcl} \widehat{\nu ^{w}}[i\! +\! 1]& & = \frac{1} {M}\sum _{m=1}^{M}\int _{ \mathbb{C}}\mathcal{N}(z_{m};\widehat{z}_{m}[i],\nu _{m}^{z}[i])\bigg(y_{ m}^{2} + \vert z_{ m}\vert ^{2} \\ & & \quad - 2y_{m}\vert z_{m}\vert \int _{0}^{2\pi }\cos (\theta _{ m} -\phi _{m})\,\frac{\exp (\kappa _{m}\cos (\theta _{m} -\phi _{m}))} {2\pi I_{0}(\kappa _{m})} d\theta _{m}\bigg)dz_{m} {}\end{array}$$

(58)

$$\displaystyle\begin{array}{rcl} & & = \frac{1} {M}\sum _{m=1}^{M}\int _{ \mathbb{C}}\mathcal{N}(z_{m};\widehat{z}_{m}[i],\nu _{m}^{z}[i]) \\ & & \quad \times \bigg (y_{m}^{2} + \vert z_{ m}\vert ^{2} - 2y_{ m}\vert z_{m}\vert R_{0}\bigg(\frac{2\vert z_{m}\vert y_{m}} {\widehat{\nu ^{w}}[i\! +\! 1]} \bigg)\bigg)dz_{m},{}\end{array}$$

(59)

where  R ₀(⋅ ) is the modified Bessel function ratio defined in (13) and (59) follows from [48, 9.6.19]. To proceed further, we make use of the expansion $R_{0}(\kappa ) = 1 -\frac{1} {2\kappa } - \frac{1} {8\kappa ^{2}} - \frac{1} {8\kappa ^{3}} + o(\kappa ^{-3})$ from [50, Lemma 5] to justify the high-SNR approximation

$$\displaystyle\begin{array}{rcl} R_{0}(\kappa ) \approx 1 -\frac{1} {2\kappa },& &{}\end{array}$$

(60)

which, when applied to (59), yields

$$\displaystyle\begin{array}{rcl} \widehat{\nu ^{w}}[i\! +\! 1] \approx \frac{2} {M}\sum _{m=1}^{M}\int _{ \mathbb{C}}\big(y_{m} -\vert z_{m}\vert \big)^{2}\mathcal{N}(z_{ m};\widehat{z}_{m}[i],\nu _{m}^{z}[i])dz_{ m}.& &{}\end{array}$$

(61)

Although (61) can be reduced to an expression that involves the mean of a Rician distribution, our empirical experience suggests that it suffices to assume  ν  _m  ^z[ i] ≈ 0 in (61), after which we obtain the much simpler expression given in (15).