General Project Objectives:
Testing is needed for methods for computing multi-dimensional integrals
(here denoted by integral(f,[a,b])
for the integral over function f of
d-dimensional vector argument
x=[x(i)]_{d×1}
over a HyperRectangle given by opposite corner vectors
[a,b]) based on better efficiency and accuracy.
An appropriate sample d-dimensional integral is given as
the expectation or mean,
I_d = E[f(x)|a<x<b]=integral(f(x)*p(x),[a,b]),
where, in MATLAB notation,
f(x) = max(sum(exp(0.06325*s.*x))-K,0)
denotes roughly the payoff for a basket option with
i=1:d stocks,
p(x) = (2*pi)^(-d/2)*exp(-0.5*sum((x-c).*(x-c)./(s.^2)))/prod(s),
(Changes: matched ")" and replaced "sum(s.*s)" by "prod(s)", 01dec04)
denotes the d-dimensional normal density with
vector mean c
and scalar standard deviation s,
where pi is area of a circle of unit radius.
Integration Methods:
The objective is to compute the the integral, using either MATLAB or
C/F programming languages, by several integration methods:
- Multi-dimensional Midpoint Rule (MPR) using for this example
M=16
midpoints within each evenly spaced
subinterval of vector length h=(b-a)/M along each dimension
so that there are N_{mpr}=M^{d} total approximation points.
- Pseudo Monte Carlo (PMC) estimates with variance based on both
- normal deviates with the method of acceptance/rejection and
- uniform deviates
using a total of N_{pmc} points
in the d-dimension, i.e,
so that N_{pmc}=N_{mpr}.
Caution: for the first part you will need to scale the integral so that
the density p(x)
is transformed to the standard zero mean and unit variance used by
MATLAB's "randn"
See class MATLAB examples
Also, noting that since the intervals
are finite, the method of acceptance/rejection must be used to exclude normal
pseudo-random deviates outside the interval of integration when computing your
sample average estimate for f(x), but be sure to use the
full point count and be sure your scaling is correct.}
- Quasi-Monte Carlo (QMC) estimates based on both uniform quasi-random
- halton sequences
- Sobol' sequences
in multi-dimensions using a the closest number of points
N(qmc)
in dimension d so the order of the error is comparable to
that used in MPR. See
that used in MPR (for the revision, then make
the N the same as for MPR).
using same total of number of points N points in MPR.
in the d-dimension.
{ Code can be converted to MATLAB or perhaps you can find a MATLAB
version on line, which are strongly recommended.}
Validation and Performance Testing:
- Random Number Testing: Test any random number generator you use
in the state/seed that you use it in (e.g., state 3 of randn is known to
generate poor values). Check the means and standard
deviations and do a statistical test, either chisquares or
KS or spectral (plotting samples x(k+1) versus x(k)) tests.
For PMC, compute the 95% confidence intervals using the sample
variance corresponding to the sample mean of f(x)
- Time each method using a timer function with sufficient resolution.
- Compare accuracy of the results for each integral against a more
refined Sobol' calculation with 4 times the number of points used.
- Test accuracy also by using p(x)f(x) = 1 for
uniform RNG methods or f(x) = 1 for
normal RNG type methods.
{MATAB-C++ Conversion Hints: DP=DoublePrecision=double (default for
MATLAB); Vec_*=Vector; ULNG=UnsignedLong; Mat_*= Matrix;
regn.size=SizeDetector (dimensions are optional
in MATLAB as long a properly referenced; size is a built-in function);
SQR=SquaringFunction; MAX=max; MIN=min;}