Function: generator - find generator of a continuous-time Markov chain with a given time-1 transition matrix

Calling sequence:

generator(P,digits);

Parameters:

P - the transition matrix: a square matrix with nonnegative entries and row-sums 1

digits - (optional) positive integer

Description:

This function attempts to construct a real matrix Q such that , each off-diagonal element of Q is nonnegative, and all its row-sums are 0.

If the optional second argument digits is included, it determines the initial setting of Digits to use in the calculations, and the precision with which the results are reported. Regardless of the value of digits , the calculations will always be done with at least evalhf(Digits) digits, and the setting may be increased if needed.

The algorithm used should find a generator, if one exists, as long as has all distinct eigenvalues. If the eigenvalues are not all distinct, it may or may not find a generator. A warning is given if the eigenvalues are not all distinct (or if two are extremely close to each other). If (despite having non-distinct eigenvalues) there is a complete set of eigenvectors of , generator should find any generators that can be expressed as polynomial functions of .

Due to roundoff error, some off-diagonal entries of a computed generator may be slightly negative even though the exact generator's entries would be nonnegative. This is particularly likely to happen when the generator has some entries of 0. Therefore a matrix can be considered as a possible generator as long as all its off-diagonal entries are all greater than .

If the environment variable EnvAllGenerators is , generator will return an expression sequence of all possible generators. Otherwise, generator will return at most one generator, the first one found with all off-diagonal entries nonnegative or, if none of those are found, the one whose minimum off-diagonal entry is greatest, as long as that value is greater than . In either case, a return of indicates that no generator was found.

If infolevel[generator] is set to a positive value, generator will provide more information on the solution process.

This code is copyright (c) 2000 by Robert B. Israel. Permission is hereby granted to anyone to copy, use and distribute it, provided that it is not altered in any way and that appropriate credit is given to the author.

Reference: Finding Generators for Markov Chains via Empirical Transition Matrices by Robert B. Israel, Jeffrey S. Rosenthal and Jason Z. Wei, to appear in Mathematical Finance.

Examples:

The function generator must be defined by reading in the file gen4.txt .

The example of Jarrow et al.

> P:= matrix([[0.8910, 0.0963, 0.0078, 0.0019, 0.0030, 0.0000, 0.0000, 0.0000],
[0.0086, 0.9010, 0.0747, 0.0099, 0.0029, 0.0029, 0.0000, 0.0000],
[0.0009, 0.0291, 0.8894, 0.0649, 0.0101, 0.0045, 0.0000, 0.0009],
[0.0006, 0.0043, 0.0656, 0.8427, 0.0644, 0.0160, 0.0018, 0.0045],
[0.0004, 0.0022, 0.0079, 0.0719, 0.7764, 0.1043, 0.0127, 0.0241],
[0.0000, 0.0019, 0.0031, 0.0066, 0.0517, 0.8246, 0.0435, 0.0685],
[0.0000, 0.0000, 0.0116, 0.0116, 0.0203, 0.0754, 0.6493, 0.2319],
[0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 0.0000, 1.0000]]);

> generator(P);

`Warning,   No completely valid generator found`

Note the small negative off-diagonal entries, which are the cause of the warning message. The matrix is a bit easier to see rounded to 4 digits. We also increase infolevel[generator] to see some information about the process.

> infolevel[generator]:= 3;
generator(P,4);

```generator:   Searching for generator with Digits=14
generator/generator:   Eigenvalues    [.98132856250027, .91861155279296, .88452961412060, .85513460895896, .79851236664780, .70457875269079, .63206192880463, 1.0000000000000]
generator/generator:   Approximate generator with qmin =   -.4198317763600e-3
Warning,   No completely valid generator found```

An example from the Israel, Rosenthal and Wei paper where the matrix calculated using the principal branch of has some negative off-diagonal elements. In this case the matrix does not have distinct eigenvalues, and generator does not find a valid generator (even though one exists).

> P:= matrix([
[.284779445, .284035268, .283826586, .147358701],
[.284191780, .284779445, .284035268, .146993507],
[.283487477, .284191780, .284779445, .147541298],
[.284543931, .283487477, .284191780, .147776812]]);

> generator(P);

```generator:   Searching for generator with Digits=14
generator/generator:   Eigenvalues    [.99999999999996, .74635808122890e-6, .10572003209687e-2+.18096310016163e-9*I, .10572003209687e-2-.18096310016163e-9*I]
Warning,   Matrix does not have distinct eigenvalues
generator/generator:   Going to   25   digits
generator/generator:   Eigenvalues    [.9999999999999999999999985, .7463580247075163958855773e-6, .1057200320987646241801978e-2+.1810231689277175342091561e-9*I, .1057200320987646241801978e-2-.1810231689277175342091561e-9*I]
Warning,   Matrix does not have distinct eigenvalues
Warning,   No completely valid generator found```

Another example from Israel, Rosenthal and Wu, having more than one valid generator. Since this example also does not have distinct eigenvalues, generator only finds one valid generator (which is of the paper). We use the environment variable EnvAllGenerators to search for more than one generator.

> EnvAllGenerators:= true;
b:= exp(-4*Pi);
P:= matrix([[(2+3*b)/5, (2-2*b)/5, (1-b)/5],
[(2-2*b)/5, (2+3*b)/5, (1-b)/5],
[(2-2*b)/5, (2-2*b)/5, (1+4*b)/5]]);

> generator(P);

```generator:   Searching for generator with Digits=14
generator/generator:   Eigenvalues    [.348734236e-5, 1.0000000000000, .348734235e-5]
Warning,   Matrix does not have distinct eigenvalues
generator/generator:   Going to   24   digits
generator/generator:   Eigenvalues    [.3487342350000000000e-5, 1.00000000000000000000000, .3487342350000000000e-5]
Warning,   Matrix does not have distinct eigenvalues
generator/generator:   Possible generator with qmin =   2.51327412322792211319233```

Here is an example where the eigenvalues are distinct, and generator finds four valid generators.

> P:= matrix([
[.3333333237, .3333333354, .3333333409],
[.3333333409, .3333333237, .3333333354],
[.3333333354, .3333333409, .3333333237]]);

> generator(P);

```generator:   Searching for generator with Digits=14
generator/generator:   Eigenvalues    [1.0000000000001, -.14450002347608e-7+.47631408489549e-8*I, -.14450002347608e-7-.47631408489549e-8*I]
generator/generator:   Going to   27   digits
generator/generator:   Eigenvalues    [.99999999999999999999999999, -.144499999999999999997367015e-7+.476313972081441255778620119e-8*I, -.144499999999999999997367015e-7-.476313972081441255778620119e-8*I]
generator/generator:   Going to   29   digits
generator/generator:   Eigenvalues    [1.0000000000000000000000000001, -.14450000000000000000004356291e-7+.47631397208144125571998945276e-8*I, -.14450000000000000000004356291e-7-.47631397208144125571998945276e-8*I]
generator/generator:   Possible generator with qmin =   4.3703689664263725349517924462
generator/generator:   Possible generator with qmin =   4.0026973043864609613705104422
generator/generator:   Possible generator with qmin =   .74277023795793683376204721595
generator/generator:   Possible generator with qmin =   .37509857591802526018094667289```

An example where the matrix violates one of the necessary conditions for a generator.

> P:= matrix([[.1,.8,.1],[.1,.1,.8],[.8,.1,.1]]);

> generator(P);

```generator:   Searching for generator with Digits=14
generator/generator:   Determinant exceeds product of diagonal elements```