Spherical-Radial Transformation Methods

These methods use a transformation to a spherical-radial (SR) coordinate system. Let ${\bf y}= r{\bf z}$, with $\vert\vert{\bf z}\vert\vert _2 = 1$, so that ${\bf y}^t{\bf y}= r^2$ and $d{\bf y}= r^{m-1}d{\bf z}$. Deák (1980-90) used this transformation as the basis for several methods for MVN problems, and the methods described in this section can be considered as generalizations of Deák's methods. After the SR transformation, the MVT problem becomes

$$\begin{eqnarray*} {\bf T}({\bf a},{\bf b},\mbox{\boldmath$\Sigma$}, \nu) &=& \fr... ...bf z}\vert\vert _2 = 1}F({\bf a},{\bf b},C,\nu,{\bf z})d{\bf z}, \end{eqnarray*}$$

where $F({\bf a},{\bf b},C,\nu,{\bf z})$ can be written in the form,

$$F({\bf a},{\bf b},C,\nu,{\bf z}) = \frac{2\Gamma(\frac{\nu+m... ... z})} \frac{r^{m-1}}{(1+\frac{r^2}{\nu})^{\frac{\nu+m}{2}}}dr.$$

We assume that the $F$ values can be quickly computed using standard statistical software. If we let ${\bf v}= C{\bf z}$, the limits for the $r$-variable integration are given by

$$\rho_l({\bf z}) = \max\{0,\max\limits_{v_i>0}\{a_i/v_i\}, \... ...imits_{v_i>0}\{b_i/v_i\}, \min\limits_{v_i<0}\{a_i/v_i\} \}\}.$$

For a given radial direction, ${\bf z}$, the $\rho$ limits are the distances from the origin to the two points where the vector with direction ${\bf z}$ intersects the boundary of the integration region.

A Monte-Carlo algorithm for the MVT problem uses

$${\bf T}({\bf a}, {\bf b}, \mbox{\boldmath $\Sigma$}, \nu) \a... ...{1}{N}\sum\limits_{k=1}^N F({\bf a},{\bf b},C,\nu,{\bf z}_k),$$

where the ${\bf z}_k$ points are uniformly random from $U_m$, the surface of the $m$-sphere. Sets of these points can easily be generated from sets of points uniform on $C_{m-1}$, the unit hypercube $[0,1]^{m-1}$, using the transformation (see Fang and Wang, 1994), from a point ${\bf w}\in C_{m-1}$ to a point ${\bf z}\in U_m$, defined by

$$\begin{eqnarray*} z_{m-2i+2}({\bf w}) & = & \sin(2\pi w_{m-2i+1})\sqrt{1-w_{m-2i... ...{m-2i}^\frac{2}{m-2i}} \prod_{k=1}^{i-1} w_{m-2k}^\frac{1}{m-2k} \end{eqnarray*}$$

for $i = 1, 2, \ldots l$, where $l=\lfloor\frac{m}{2}\rfloor - 1$, and ending with

$$z_2({\bf w}) = \sin(2\pi w_1)\prod_{k=1}^l w_{m-2k}^\frac{1}... ...\bf w}) = \cos(2\pi w_1)\prod_{k=1}^l w_{m-2k}^\frac{1}{m-2k},$$

when $m$ is even, or ending with

$$z_3({\bf w}) = (2w_1 - 1 )\prod_{k=1}^l w_{m-2k}^\frac{1}{m-2k},$$

$$z_2({\bf w})=2\sin(2\pi w_2)\sqrt{w_1(1-w_1)}\prod_{k=1}^lw_... ...\pi w_2)\sqrt{w_1(1-w_1)}\prod_{k=1}^lw_{m-2k}^\frac{1}{m-2k},$$

when $m$ is odd. This transformation has a constant Jacobian, so a Monte-Carlo algorithm for the MVT problem based on uniform $C_{m-1}$ points uses

$${\bf T}({\bf a}, {\bf b}, \mbox{\boldmath$\Sigma$}, \nu) \ap... ...m\limits_{k=1}^N F({\bf a},{\bf b},C,\nu,{\bf z}({\bf w}_k)),$$ (9)

with all $w_{i,k} \sim Uniform(0,1)$. This method can be also used for integration regions that are not hyper-rectangles. All that is needed for other regions is an efficient method for computing the intersection points for the boundary of the integration region for each radial direction. Lohr (1990) has described this type of generalization for MVN problems.

The use of various types of antithetic variates, as described by Deák (1990) for improving the convergence of SR MVN methods, can improve the convergence of this type of method. Let $Z$ be an $m \times m$ uniformly random (with Haar measure, see Stewart, 1980) orthogonal matrix with columns $\{{\bf z}_j\}$, and define

$$S_n(Z) = \frac{1}{2^n{m \choose n}} \sum_{{\bf s}}\sum_{1 \l... ...\bf b},C,\nu,\frac{\sum_{l=1}^n s_l {\bf z}_{j_l}}{\sqrt{n}}),$$

where ${\bf s}= (s_1, s_2, ..., s_n) = (\pm 1, ..., \pm 1)$ and the outer sum is taken over the $2^n$ possible sign combinations for the components of ${\bf s}$. The sample points used by $S_n(Z)$ are very evenly spread over the surface of the unit $m$-sphere. For MVN problems, Deák found that the larger values for the parameter $n$ (which must satisfy $n \leq m$) can provide values for $S_n$ with significantly smaller variances. But the larger the $n$ value, the higher the computational cost for $S_n$, so these two features of the $S_n$ sums must be balanced, for practical computations. Deák recommended values of $n= 1, 2$ or $3$ for typical computations. MVT estimates based on $S_n$ are obtained using

$${\bf T}({\bf a}, {\bf b}, \mbox{\boldmath$\Sigma$}, \nu) \approx \frac{1}{M}\sum\limits_{k=1}^M S_n(Z_k).$$ (10)

One problem with the SR transformation algorithms is that $F({\bf a},{\bf b},C,\nu,{\bf z})$ as a function of ${\bf z}$, although continuous, is not very smooth, because of sharp corners of the integration region defined by ${\bf a}\leq rC{\bf z}\leq {\bf b}$. This problem resulted in approximations with large variation and slower convergence for SR algorithms for MVN problems (see Genz, 1992). For some combinations of ${\bf a}$ and ${\bf b}$, (e.g. if ${\bf0}< {\bf a}< {\bf b}$ ) many $F$ values will be zero, and this can cause further reductions in the efficiency of the SR algorithms.

SR methods for the NCMVT problem can easily be constructed by combining a numerical integration method for the $s$ variable (in equation (3)) with an SR method for the inner MVN integral. Using the SR transformation for the inner integral we have

$$\begin{eqnarray*} {\bf T}({\bf a}, {\bf b}, \mbox{\boldmath$\Sigma$}, \nu, \mbox... ... a},{\bf b},C,\nu,\mbox{\boldmath$\delta$},{\bf z},s)d{\bf z}ds, \end{eqnarray*}$$

with $G$ appropriately defined. A simple Monte-Carlo algorithm for the NCMVT problem uses

$${\bf T}({\bf a}, {\bf b}, \mbox{\boldmath $\Sigma$}, \nu, \m... ...}\sum\limits_{k=1}^N G({\bf a},{\bf b},C,\nu,{\bf z}_k, s_k),$$

where $\{{\bf z}_k\}$ are uniformly random from the surface of $m$-sphere and $s_k\sim \chi_\nu$.

Second formulations of the MVT and NCMVT problems, using the SR coordinate system, reverse the order of the radial and spherical integrations. The result for the MVT problem is

$$\begin{eqnarray*} {\bf T}({\bf a},{\bf b},\mbox{\boldmath$\Sigma$}, \nu) &=& \f... ...rac{r^2}{\nu})^{\frac{\nu+m}{2}}} H({\bf a},{\bf b},C,\nu,r)dr, \end{eqnarray*}$$

where $H$ is appropriately defined. Algorithms for the general MVT and NCMVT problems using these formulations can suffer from convergence problems similar to those for the first formulations. However, Somerville (1997-99) has found the MVT second formulation to be useful for some confidence interval computation applications.