Statistics
Class FinMath

java.lang.Object
  Statistics.FinMath

public class FinMath

extends java.lang.Object

This class contains static methods to compute functions or solve eqations * useful in basic financial mathematics. * * @author Michael J. Meyer

Constructor Summary

FinMath()

Method Summary

static double	BisectionSolveBSF(double Q, double k, double y) Same as NewtonSolveBSF(double, double, double) but uses continued bisection instead of * Newton's algorithm.
static double	blackScholesFunction(double Q, double k, double Sigma) Computes the function QN(h_+)-kN(h_-).
static double	bsDiscountedCallPrice(double S, double K, double tau, double sigma, double B) Discounted Black-Scholes call price (note that current time may * not be zero).
static double	bsDiscountedPutPrice(double S, double K, double tau, double sigma, double B) Discounted Black-Scholes put price (note that current time may * not be zero).
static void	CholeskyRoot(double[][] C, double[][] L) Let C be symmetric positive definite n by n matrix.
static double	d_minus(double Q, double k, double Sigma) The quantity d_- in Margrabe's formula for the option to * exchange assets.
static double	d_plus(double Q, double k, double Sigma) The quantity d_+ in Margrabe's formula for the option to * exchange assets S_1, S_2 (receive S_1 for kS_2).
static double	dBSF(double Q, double k, double Sigma) Derivative of the Black-Scholes function QN(h_+)-kN(h_-) with * respect to Sigma.
static double	dN(double x) Derivative of N'(x) of N(double).
static cern.colt.matrix.DoubleMatrix1D	linearSystemSolution(cern.colt.matrix.DoubleMatrix2D A, cern.colt.matrix.DoubleMatrix1D y) Returns cern.colt.matrix.impl.DenseDoubleMatrix1D solution vector x * of the linear equation Ax=y.
static void	main(java.lang.String[] args) TEST PROGRAM
static double	N_Inverse(double x) The Inverse of the cumulative normal distribution function N(x).
static double	N_Solve(double y) Solves the equation N(x)=y for x=x(y) using Newton's algorithm.
static double	N_Solve1(double y) Solves the equation N(x)=y for x=x(y) using continued bisection.
static double	N(double x) The cumulative distribution function N(x)=Prob(X<=x) of a standard normal * variable X.
static double	NewtonSolveBSF(double Q, double k, double y) Given b > 0 solves the equation * bLackScholesFunction(Q,k,Sigma)=y for * Sigma> 0 using Newton's algorithm.
double	normalDensity(double[] z) The density of the standard multinormal distribution N(0,I).
static void	PlusEquals(cern.colt.matrix.DoubleMatrix1D w, double k, cern.colt.matrix.DoubleMatrix2D A, boolean transposeA, cern.colt.matrix.DoubleMatrix1D z) Implements the vector operation w+=kAz, or w+=kA'z, * where z,w are vectors and A a matrix.
static double	round(double x, int n) rounds off the double x to n decimals.

Methods inherited from class java.lang.Object

clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait

Constructor Detail

FinMath

public FinMath()

Method Detail

round

public static double round(double x,
                           int n)

rounds off the double x to n decimals.

* * @param x number to be rounded off * @param n number of decimals

normalDensity

public double normalDensity(double[] z)

The density of the standard multinormal distribution N(0,I).

N

public static double N(double x)

The cumulative distribution function N(x)=Prob(X<=x) of a standard normal * variable X.

* *

It is derived from the error function erf(x)=Prob(-sqrt(2)x<=X<=sqrt(2)x) * included in the class VisualNumerics.math.Sfun.

* * @param x Any real number.

d_plus

public static double d_plus(double Q,
                            double k,
                            double Sigma)

The quantity d_+ in Margrabe's formula for the option to * exchange assets S_1, S_2 (receive S_1 for kS_2). In this context:

* * d_+=log(Q/k)/Sigma+Sigma/2 with Q=S_1(t)/S_2(t) * *

and Sigma is the aggregate volatility of log(Q) from current * time t to the horizon T, that is, the quadratic variation

* * Sigma = ‹ log(Q) ›_[t,T]. * *

If Q has constant annual volatility sigma this becomes

* * Sigma=sigma*sqrt(tau) * *

where tau=T-t is time to expiry. In the case of the European call on * S_1 with strike k the asset S_2 * is the zero coupon bond maturing at call expiration.

* * @param Q The quotient S_1(t)/S_2(t) * @param k Ratio of exchange between assets (receive S_1 for kS_2). * @param Sigma Quadratic variation ‹ log(Q) ›_[t,T].

d_minus

public static double d_minus(double Q,
                             double k,
                             double Sigma)

The quantity d_- in Margrabe's formula for the option to * exchange assets.

*

See d_plus(double, double, double) and replace "+Sigma/2" with "-Sigma/2".

blackScholesFunction

public static double blackScholesFunction(double Q,
                                          double k,
                                          double Sigma)

Computes the function QN(h_+)-kN(h_-).

* *

Useful in application of the Black-Scholes formula or Margrabe's formula * to price calls or options to exchange assets. For the option to * exchange assets S_1,S_2 we have Q=S_1(t)/S_2(t) and

* * Price(t,S_1,S_2)=QN(h_+)-kN(h_-) * *

See d_plus(double, double, double) for information on the parameter Sigma. * In the special case of a call on S we have S_1=S, S_2 is the zero coupon * bond maturing at call expiry, Q the forward price of S and the forward price * of the call can be written as

* *

Call.ForwardPrice(t,S)=S_2*(QN(h_+)-kN(h_-))

* * @param Q The quotient S_1(t)/S_2(t) * @param k Ratio of exchange between assets (receive S_1 for kS_2). * @param Sigma Quadratic variation ‹ log(Q) ›_[t,T].

dBSF

public static double dBSF(double Q,
                          double k,
                          double Sigma)

Derivative of the Black-Scholes function QN(h_+)-kN(h_-) with * respect to Sigma.

* * @param Q The quotient S_1(t)/S_2(t) * @param k Ratio of exchange between assets (receive S_1 for kS_2). * @param Sigma Quadratic variation ‹log(Q) ›_[t,T].

bsDiscountedCallPrice

public static double bsDiscountedCallPrice(double S,
                                           double K,
                                           double tau,
                                           double sigma,
                                           double B)

Discounted Black-Scholes call price (note that current time may * not be zero). Both the current time t and the time * T to expiry (from time zero) are needed since the formula * will be applied at time other than time zero and S

* * @param S the forward asset price at expiry . * @param K strike price * @param tau time to expiration * @param sigma annual asset volatility * @param B risk free bond B_T at expiry

bsDiscountedPutPrice

public static double bsDiscountedPutPrice(double S,
                                          double K,
                                          double tau,
                                          double sigma,
                                          double B)

Discounted Black-Scholes put price (note that current time may * not be zero). Here we need the interest rate, the current time * t and the time T of expiration explicitly. * The formula will be applied at times other than time zero.

* * @param S the forward asset price at expiry T. * @param K strike price * @param tau time to expiry * @param sigma annual asset volatility * @param B risk free bond B_T at expiry

NewtonSolveBSF

public static double NewtonSolveBSF(double Q,
                                    double k,
                                    double y)
                             throws NoSolutionException

Given b > 0 solves the equation * bLackScholesFunction(Q,k,Sigma)=y for * Sigma> 0 using Newton's algorithm.

* *

A solution Sigma=Sigma(Q,k,y) exists if and only if b < Q. * To be used to compute implied * volatilities from call forward prices y. Note that the solution Sigma * is not the implied annual volatility sigma, indeed these are related as * Sigma=sigma*sqrt(tau), where tau=T-t is time to expiry.

* *

Since Newton's algorithm often fails to converge it is better to rely * on the method using continued bisection below.

Throws:

NoSolutionException

BisectionSolveBSF

public static double BisectionSolveBSF(double Q,
                                       double k,
                                       double y)
                                throws NoSolutionException

Same as NewtonSolveBSF(double, double, double) but uses continued bisection instead of * Newton's algorithm.

* * Search restricted to [0.00001,1000], smaller or larger solutions cut off * to interval endpoints (irrelevant anyway). The blackScholesFunction * f(Q,k,Sigma) increases with Sigma. The algorithm preserves the relation
* f(Q,k,a) <=y<=f(Q,k,b).

Throws:

NoSolutionException

dN

public static double dN(double x)

Derivative of N'(x) of N(double).

N_Solve

public static double N_Solve(double y)

Solves the equation N(x)=y for x=x(y) using Newton's algorithm.

* *

Application: confidence intervals for normal variables. * Newton's algorithm often fails to converge. Better rely on continued * bisection N_Solve1(double) below.

*

N_Solve1

public static double N_Solve1(double y)
                       throws NoSolutionException

Solves the equation N(x)=y for x=x(y) using continued bisection.

* * Search restricted to [-100000,+100000], smaller or larger solutions * throw NoSolutionException (irrelevant anyway). Solution x(y) keeps error in * equation N(x)=y smaller than 0.000001.

Throws:

NoSolutionException

N_Inverse

public static double N_Inverse(double x)

The Inverse of the cumulative normal distribution function N(x). * Algorithm of Peter J. Acklam * http://www.math.uio.no/~jacklam/notes/invnorm.

* *

Adapted from from: Peter Jaeckel, "Monte Carlom Methods in Finance", Wiley * ISBN 0-471-497-41-X

* * @param x must be in the interval (0,1).

CholeskyRoot

public static void CholeskyRoot(double[][] C,
                                double[][] L)

Let C be symmetric positive definite n by n matrix. The routine computes * the lower triangular matrix L satisfying LL'=C. Routine writes to the * memory referenced by L. Does not handle the positive semidefinite case * (all diagonal elements of C must be > 0). Only the lower triangular half * of C is used and only the lower triangular half of L is written.

PlusEquals

public static void PlusEquals(cern.colt.matrix.DoubleMatrix1D w,
                              double k,
                              cern.colt.matrix.DoubleMatrix2D A,
                              boolean transposeA,
                              cern.colt.matrix.DoubleMatrix1D z)

Implements the vector operation w+=kAz, or w+=kA'z, * where z,w are vectors and A a matrix.

* * @param w,k,A,z the players in w+=kAz. * @param wether A should be transposed or not.

linearSystemSolution

public static cern.colt.matrix.DoubleMatrix1D linearSystemSolution(cern.colt.matrix.DoubleMatrix2D A,
                                                                   cern.colt.matrix.DoubleMatrix1D y)

Returns cern.colt.matrix.impl.DenseDoubleMatrix1D solution vector x * of the linear equation Ax=y. Exact solution if A is * square, least squares solution otherwise. * * @param A a cern.colt.matrix.impl.DenseDoubleMatrix2D matrix. * @param y a cern.colt.matrix.impl.DenseDoubleMatrix1D vector.

main

public static void main(java.lang.String[] args)

TEST PROGRAM