Equiprobable Approximation to Bivariate Lognormal Returns

This section starts by defining a convenient notation to represent two rate-of-return shocks that are distributed according to a multivariate lognormal that allows for nonzero covariances. It next presents a computationally simple method for constructing a numerical approximation to that joint distribution. The final section uses both the numerical approximation, and more accurate (but enormously slower) standard numerical integration tools to assess the accuracy of the Campbell-Vicera analytical approximation to the solution to the optimal portfolio choice problem described in Portfolio-Multi-CRRA.

1Statistical Theory¶

Consider a set of two normally distributed risks

which are statistically independent ($\ShkMeanOneLog_{1,t+1} \perp \ShkMeanOneLog_{2,t+1}$) even though the scale of shock 1’s standard deviation is determined by the proportionality factor $\scale$ multiplied by the value of shock 2’s standard deviation, $\sigAll$. (This permits us to increase or decrease the size of both risks by changing $\sigAll$ and to adjust the relative size of risk 1 vs risk 2 by adjusting $\scale$.)

The $\ShkMeanOneLog$ are interpreted as the logarithms of level variables $\ShkMeanOne$, and the means of the log variables have been chosen such that the expectations of the levels are independent of the size of the risk (cf. Math Facts). That is, defining $\ShkMeanOne_{i} = \Ex_{t}[e^{\ShkMeanOneLog_{i,t+1}}]$ for $i \in \{1,2\}$:

From the first shock we can construct a log rate-of-return variable that can be represented equivalently in either of two ways:

where $\ShkLogZeroLog_{1,t+1}~\sim~\mathcal{N}(0.,(\scale \sigAll)^{2})$; the notation for $\ShkLogZeroLog$ is motivated by the fact that the addition of the extra term cancels the nonzero mean of the original $\ShkMeanOneLog$. Then

where note to avoid confusion that $\riskyAlt_{1} \neq \Ex_{t}[{{\risky}}_{1,t+1}]$ while $\log \Ex_{t}[\Risky_{1,t+1}] \neq \risky_{1}=\Ex_{t}[\risky_{1,t+1}]$.

Using $\Ex_{t}[\ShkMeanOneLog_{2,t+1}] = -0.5 \sigAll^{2}$ and $\Ex_{t}[(\omega/\scale)\ShkMeanOneLog_{1,t+1}] = -0.5 ((\omega/\scale)\sigAll)^{2}$, we can by analogy define a second return

for some constants $\omega$ and $\zeta$.

Since $(\omega/\scale)\ShkMeanOneLog_{1,t+1}$ is the only component of ${\risky}_{2,t+1}$ that is correlated with ${\risky}_{1,t+1}$,

Thus, the parameter $\omega$ controls the covariance between the risky returns. If we set $\omega = 0$ then $\risky_{1,t+1} \perp \risky_{2,t+1}$ (the returns are independent).

Next we want to find the value of $\zeta$ such that the expected level of the return is unaffected by $\sigAll$ (so that we will be able to explore independently the distinct effects of the components of each shock and their covariance):

regardless of the values of $\scale$ and $\sigAll$. Using (5), we therefore need:

Using standard facts about lognormals (cf. Math Facts), and for convenience defining $\hat{\omega}= (\sigAll/\scale) \omega$, we have

which means that we can rewrite (5) directly as

Hence, from the independent mean-one lognormally distributed shocks $\ShkMeanOne_{1}$ and $\ShkMeanOne_{2}$, we have constructed a pair of jointly lognormally distributed shocks whose covariance is controlled by the parameter $\omega$, whose relative and absolute variances are controlled by the parameters $\scale$ and $\sigAll$, and whose means are $\risky_{1}$ and $\risky_{2}$.

To sum up, the process can be described in either of two ways:

with covariance matrix

where a final useful result that follows from (12) and (6) is

2Computational Theory¶

To reduce clutter, define $\tilde{\Risky}_{i}=\Risky_{i,t+1}$ and interpret $\Ex$ as $\Ex_{t}[]$, so that we can write the expectation of some function $\hFunc$ that depends on the realization of the return shocks as:

where $\FFunc(\tilde{\Risky}_{1},\tilde{\Risky}_{2})$ is the cumulative distribution function for a multivariate lognormal with the covariance matrix defined by (12).^[1]

Standard numerical computation software can compute this double integral, but at such a slow speed as to be unusable for many purposes.

Computation of the expectation can be massively speeded up by advance construction of a numerical approximation to $\FFunc(\tilde{\Risky}_{1},\tilde{\Risky}_{2})$.

Such approximations often take the approach of replacing the distribution function with a discretized approximation to it; appropriate weights $w[i,j]$ are attached to each of a finite set of points indexed by $i$ and $j$, and the approximation to the integral is given by:

where various methods are used for constructing the weights $w[i,j]$ and the nodes (corresponding to the $\{i,j\}$ pairs). The matrices $\hat{\Risky}_{1}$ and $\hat{\Risky}_{2}$ contain the conditional mean values of $\tilde{\Risky}_{1}$ and $\tilde{\Risky}_{2}$ associated with each of the regions.

Perhaps the most popular such method is Gauss-Hermite interpolation (see Judd (1998) for an exposition, or Kopecky & Suen (2010) for a recent candidate for a better choice). Here, we will pursue a particularly simple and intuitive alternative: Equiprobable discretization. In this method, $m=n$ and boundaries on the joint CDF are determined in such a way as to divide up the total probability mass into submasses of equal size (each of which therefore has a mass of $n^{-2}$). This is conceptually easier if we represent the underlying shocks as statistically independent, as with $\ShkMeanOne_{1}$ and $\ShkMeanOne_{2}$ above; in that case, each submass is a square region in the $\ShkMeanOne_{1}$ and $\ShkMeanOne_{2}$ grid. We then compute the average value of $\ShkMeanOne_{1}$ and $\ShkMeanOne_{2}$ conditional on their being located in each of the subdivisions of the range of the CDF. Since $\ShkMeanOne_{1}$ and $\ShkMeanOne_{2}$ are IID, the representation of the approximating summation is even simpler than (15):

where $\Risky_{i}(\ShkMeanOne_{1},\ShkMeanOne_{2})$ are the (linear) functions relating the return shocks to the IID shocks.

3Computational Results¶

Figure 1 compares the computed optimal portfolio share for a numerical solution using the built-in numerical optimizer and maximization functions (the lowest, black, locus), the Campbell-Viceira solution (the highest, red locus) and an equiprobable approximation using 20 approximation points (green, middle) as well as the solution using the equiprobable approximation at an evenly-spaced grid of points (blue dots).

Careful examination indicates that the numerical approximation is quite close to the full numerical solution, while the CV approximation diverges substantially from the numerical answer. The tradeoff is that the equiprobable solution is about 2000 times slower than the CV approximation, while the direct solution is more than 100 times slower than the equiprobable solution. Depending on the requirements of the problem being examined, these differences in efficiency can make a tremendous difference in the feasibility of a research project.^[2]