| moments {custom.gauss.quad} | R Documentation |
Moments Computed in Multiple Precision Using the Package Rmpfr
Description
This module computes the r'th moment
\int_{-\infty}^{\infty} x^r f(x) dx,
where f is the weight function (specified by the list which.f), for any nonnegative integer r using nbits bits of precision for its computation, via the R package Rmpfr.
Usage
moments(which.f, r, nbits)
Arguments
which.f |
a list specifying the nonnegative integrable
weight function |
r |
nonnegative integer, specifying that it is the |
nbits |
number of bits in the multiple precision numbers used by the |
Details
Suppose, for example, that we wish to find the Gauss quadrature nodes and weights
for the weight function f that is the probability density function of a random
variable with the same distribution as R/m^{1/2} where R has a
\chi_m distribution (i.e. R^2 has a \chi_m^2 distribution).
In this case, the r'th moment is
\int_{-\infty}^{\infty} x^r f(x) dx
= \left(\frac{2}{m} \right)^{r/2}
\frac{\Gamma((r+m)/2)}{\Gamma(m/2)},
which can be computed to an arbitrary number of bits of precision
nbits using the R package Rmpfr.
In this case, we specify this weight function f by first
assigning the value of m and then using the R command
which.f <- list(name="scaled.chi.pdf", support=c(0, Inf), parameters=m)
The code within the function moments used to compute the
r'th moment, to an arbitrary number of bits of precision
nbits using the package Rmpfr, is listed in the
Examples section.
Value
The r'th moment with number of bits of precision
nbits used in its computation, via the R package Rmpfr
See Also
custom
Examples
# The code for the function moments must include a section
# that computes the r th moment to an arbitrary number of bits
# of precision nbits using the R package Rmpfr for the particular
# weight function f of interest.
# Suppose that the weight function f is the probability density
# function of a random variable with the same probability
# distribution as R divided by the square root of m, where R has a
# chi distribution with m degrees of freedom.
# The code for the function moments includes the following:
#
# if (which.f$name == "scaled.chi.pdf"){
# m <- which.f$parameters
# if (r == 0){
# return(mpfr(1, nbits))
# }
# mp.2 <- mpfr(2, nbits)
# mp.r <- mpfr(r, nbits)
# mp.m <- mpfr(m, nbits)
# term1 <- (mp.r/ mp.2) * log(mp.2 / mp.m)
# term2 <- lgamma((mp.r + mp.m) / mp.2)
# term3 <- lgamma(mp.m / mp.2)
# return(exp(term1 + term2 - term3))
# }