d1_i {qfratio} | R Documentation |
Coefficients in polynomial expansion of generating function—single matrix
Description
These are internal functions to calculate the coefficients in polynomial expansion of generating functions for quadratic forms in multivariate normal variables.
d1_i()
is for standard multivariate normal variables,
\mathbf{x} \sim N_n(\mathbf{0}_n, \mathbf{I}_n)
.
dtil1_i_v()
is for noncentral multivariate normal variables,
\mathbf{x} \sim N_n(\bm{\mu}, \mathbf{I}_n)
.
dtil1_i_m()
is a wrapper for dtil1_i_v()
and takes the argument matrix rather than its eigenvalues.
Usage
d1_i(L, m = 100L, thr_margin = 100)
dtil1_i_v(L, mu = rep.int(0, n), m = 100L, thr_margin = 100)
dtil1_i_m(A, mu = rep.int(0, n), m = 100L, thr_margin = 100)
Arguments
L |
Vector of eigenvalues of the argument matrix |
m |
Integer-alike to specify the order of polynomials |
thr_margin |
Optional argument to adjust the threshold for scaling (see “Details”) |
mu |
Mean vector |
A |
Argument matrix. Assumed to be symmetric in these functions. |
Details
d1_i()
calculates d_k(\mathbf{A})
, and
dtil1_i_v()
and dtil1_i_m()
calculate
\tilde{d}_k(\mathbf{A})
in
Hillier et al. (2009, 2014) and Bao and Kan (2013). The former is
related to the top-order zonal polynomial
C_{[k]}(\mathbf{A})
in the following way:
d_k(\mathbf{A}) = \frac{1}{k!} \left( \frac{1}{2} \right)_k
C_{[k]}(\mathbf{A})
,
where (x)_k = x (x + 1) \dots (x + k - 1)
.
These functions calculate the coefficients based on the super-short recursion algorithm described in Hillier et al. (2014: 3.2, eqs. 28–30).
Scaling
The coefficients described herein (and in d2_ij
and
d3_ijk
) can become very large for higher-order terms,
so there is a practical risk of numerical overflow when applied to
large matrices or matrices with many large eigenvalues
(note that the latter typically arises from those with many small
eigenvalues for the front-end qfrm()
functions). To avoid
numerical overflow, these functions automatically scale
coefficients (and temporary objects used to calculate them) by a large number
(1e10
at present) when any value in the temporary objects exceeds
a threshold, .Machine$double.xmax / thr_margin / n
, where n
is the number of variables. This default value empirically seems to work well
in most conditions, but use a large thr_margin
(e.g., 1e5
)
if you encounter numerical overflow. (The C++ functions use
an equivalent expression,
std::numeric_limits<Scalar>::max() / thr_margin / Scalar(n)
, with
Scalar
being double
or long double
.)
In these R functions, the scaling happens order-wise;
i.e., it influences all the coefficients of the same order in
multidimensional coefficients (in d2_ij
and
d3_ijk
) and the coefficients of the subsequent orders.
These scaling factors are recorded in the attribute "logscale"
of the
return value, which is a vector/matrix/array whose size is identical to the
return value, so that value / exp(attr(value, "logscale"))
equals
the original quantities to be obtained (if there were no overflow).
The qfrm
and qfmrm
functions handle return values of these
functions by first multiplying them with hypergeometric coefficients
(which are typically \ll 1
) and then scaling the products back
to the original scale using the recorded scaling factors. (To be precise,
this typically happens within hgs
functions.) The
C++ functions handle the problem similarly (but by using
separate objects rather than attributes).
However, this procedure does not always mitigate the problem in
multiple series; when there are very large and very small
coefficients in the same order, smaller ones can diminish/underflow to
the numerical 0
after repeated scaling. (The qfrm
and
qfmrm
functions try to detect and warn against
this by examining whether any of the highest-order terms is 0
.)
The present version of this package has implemented two methods to mitigate
this problem, but only through C++ functions. One is to use the
long double
variable type, and the other is to use coefficient-wise
scaling (see qfrm
and qfmrm
).
Value
Vector of length m + 1
, corresponding to
the 0th, 1st, ..., and mth order terms. Hence, the [m + 1]
-th
element should be extracted when the coefficient for the m
th order
term is required.
Has the attribute "logscale"
as described in “Scaling” above.
References
Bao, Y. and Kan, R. (2013) On the moments of ratios of quadratic forms in normal random variables. Journal of Multivariate Analysis, 117, 229–245. doi:10.1016/j.jmva.2013.03.002.
Hillier, G., Kan, R. and Wang, X. (2009) Computationally efficient recursions for top-order invariant polynomials with applications. Econometric Theory, 25, 211–242. doi:10.1017/S0266466608090075.
Hillier, G., Kan, R. and Wang, X. (2014) Generating functions and short recursions, with applications to the moments of quadratic forms in noncentral normal vectors. Econometric Theory, 30, 436–473. doi:10.1017/S0266466613000364.
See Also
qfpm
, qfrm
, and qfmrm
are
major front-end functions that utilize these functions
dtil2_pq
for \tilde{d}
used for moments of a product of quadratic forms
d2_ij
and d3_ijk
for d
, h
,
\tilde{h}
, and \hat{h}
used for moments of ratios
of quadratic forms