# simByTransition

Simulate `CIR`

sample paths with transition
density

## Description

`[`

simulates `Paths`

,`Times`

] = simByTransition(`MDL`

,`NPeriods`

)`NTrials`

sample paths of `NVars`

independent state variables driven by the Cox-Ingersoll-Ross (CIR) process sources
of risk over `NPeriods`

consecutive observation periods.
`simByTransition`

approximates a continuous-time CIR model
using an approximation of the transition density function.

`[`

specifies options using one or more name-value pair arguments in addition to the
input arguments in the previous syntax.`Paths`

,`Times`

] = simByTransition(___,`Name,Value`

)

You can perform quasi-Monte Carlo simulations using the name-value arguments for
`MonteCarloMethod`

, `QuasiSequence`

, and
`BrownianMotionMethod`

. For more information, see Quasi-Monte Carlo Simulation.

## Examples

### Simulate Future Term Structures Using a CIR Model

Using the short rate, simulate the rate dynamics and term structures in the future using a CIR model. The CIR model is expressed as

$$dr(t)=\alpha (b-r(t))dt+\sigma \sqrt{r(t)}dW(t)$$

The exponential affine form of the bond price is

$$B(t,T)={e}^{-A(t,T)r(t)+C(t,T)}$$

where

$$A(t,T)=\frac{2({e}^{\gamma (T-t)}-1)}{(\gamma +\alpha )({e}^{\gamma (T-t)}-1)+2\gamma}$$

$$B(t,T)=\frac{2\alpha b}{{\sigma}^{2}}\mathrm{log}\left(\frac{2\gamma {e}^{(\alpha +\gamma )(T-t)/2}}{(\gamma +\alpha )({e}^{\gamma (T-t)}-1)+2\gamma}\right)$$

and

$$\gamma =\sqrt{{\alpha}^{2}+2{\sigma}^{2}}$$

Define the parameters for the `cir`

object.

alpha = .1; b = .05; sigma = .05; r0 = .04;

Define the function for bond prices.

gamma = sqrt(alpha^2 + 2*sigma^2); A_func = @(t, T) ... 2*(exp(gamma*(T-t))-1)/((alpha+gamma)*(exp(gamma*(T-t))-1)+2*gamma); C_func = @(t, T) ... (2*alpha*b/sigma^2)*log(2*gamma*exp((alpha+gamma)*(T-t)/2)/((alpha+gamma)*(exp(gamma*(T-t))-1)+2*gamma)); P_func = @(t,T,r_t) exp(-A_func(t,T)*r_t+C_func(t,T));

Create a `cir`

object.

`obj = cir(alpha,b,sigma,'StartState',r0)`

obj = Class CIR: Cox-Ingersoll-Ross ---------------------------------------- Dimensions: State = 1, Brownian = 1 ---------------------------------------- StartTime: 0 StartState: 0.04 Correlation: 1 Drift: drift rate function F(t,X(t)) Diffusion: diffusion rate function G(t,X(t)) Simulation: simulation method/function simByEuler Sigma: 0.05 Level: 0.05 Speed: 0.1

Define the simulation parameters.

nTrials = 100; nPeriods = 5; % Simulate future short over the next five years nSteps = 12; % Set intermediate steps to improve the accuracy

Simulate the short rates. The returning path is a (NPeriods + 1)-by-`NVars`

-by-`NTrials `

three-dimensional time-series array. For this example, the size of the output is `6`

-by-`1`

-by-`100`

.

rng('default'); % Reproduce the same result rPaths = simByTransition(obj,nPeriods,'nTrials',nTrials,'nSteps',nSteps); size(rPaths)

`ans = `*1×3*
6 1 100

rPathsExp = mean(rPaths,3);

Determine the term structure over the next 30 years.

maturity = 30; T = 1:maturity; futuresTimes = 1:nPeriods+1; % Preallocate simTermStruc simTermStructure = zeros(nPeriods+1,30); for i = futuresTimes for t = T bondPrice = P_func(i,i+t,rPathsExp(i)); simTermStructure(i,t) = -log(bondPrice)/t; end end plot(simTermStructure') legend('Current','1-year','2-year','3-year','4-year','5-year') title('Projected Term Structure for Next 5 Years') ylabel('Long Rate Maturity R(t,T)') xlabel('Time')

### Quasi-Monte Carlo Simulation Using a CIR Model

The Cox-Ingersoll-Ross (CIR) short rate class derives directly from SDE with mean-reverting drift (`SDEMRD`

): $d{X}_{t}=S(t)[L(t)-{X}_{t}]dt+D(t,{X}_{t}^{\frac{1}{2}})V(t)dW$

where $$D$$ is a diagonal matrix whose elements are the square root of the corresponding element of the state vector.

Create a `cir`

object to represent the model: $d{X}_{t}=0.2(0.1-{X}_{t})dt+0.05{X}_{t}^{\frac{1}{2}}dW$.

`cir_obj = cir(0.2, 0.1, 0.05) % (Speed, Level, Sigma)`

cir_obj = Class CIR: Cox-Ingersoll-Ross ---------------------------------------- Dimensions: State = 1, Brownian = 1 ---------------------------------------- StartTime: 0 StartState: 1 Correlation: 1 Drift: drift rate function F(t,X(t)) Diffusion: diffusion rate function G(t,X(t)) Simulation: simulation method/function simByEuler Sigma: 0.05 Level: 0.1 Speed: 0.2

Define the quasi-Monte Carlo simulation using the optional name-value arguments for `'MonteCarloMethod'`

,`'QuasiSequence'`

, and `'BrownianMotionMethod'`

.

[paths,time] = simByTransition(cir_obj,10,'ntrials',4096,'montecarlomethod','quasi','quasisequence','sobol','BrownianMotionMethod','principal-components');

## Input Arguments

`MDL`

— Stochastic differential equation model

object

Stochastic differential equation model, specified as a
`cir`

object. For more information on creating a
`CIR`

object, see `cir`

.

**Data Types: **`object`

`NPeriods`

— Number of simulation periods

positive scalar integer

Number of simulation periods, specified as a positive scalar integer. The
value of `NPeriods`

determines the number of rows of the
simulated output series.

**Data Types: **`double`

### Name-Value Arguments

Specify optional pairs of arguments as
`Name1=Value1,...,NameN=ValueN`

, where `Name`

is
the argument name and `Value`

is the corresponding value.
Name-value arguments must appear after other arguments, but the order of the
pairs does not matter.

*
Before R2021a, use commas to separate each name and value, and enclose*
`Name`

*in quotes.*

**Example: **```
[Paths,Times] =
simByTransition(CIR,NPeriods,'DeltaTime',dt)
```

`NTrials`

— Simulated trials (sample paths)

`1`

(single path of correlated state
variables) (default) | positive integer

Simulated trials (sample paths) of `NPeriods`

observations each, specified as the comma-separated pair consisting of
`'NTrials'`

and a positive scalar integer.

**Data Types: **`double`

`DeltaTime`

— Positive time increments between observations

`1`

(default) | scalar | column vector

Positive time increments between observations, specified as the
comma-separated pair consisting of `'DeltaTime'`

and a
scalar or a `NPeriods`

-by-`1`

column
vector.

`DeltaTime`

represents the familiar
*dt* found in stochastic differential equations,
and determines the times at which the simulated paths of the output
state variables are reported.

**Data Types: **`double`

`NSteps`

— Number of intermediate time steps

`1`

(indicating no intermediate
evaluation) (default) | positive integer

Number of intermediate time steps within each time increment
*dt* (defined as `DeltaTime`

),
specified as the comma-separated pair consisting of
`'NSteps'`

and a positive scalar integer.

The `simByTransition`

function partitions each time
increment *dt* into `NSteps`

subintervals of length *dt*/`NSteps`

,
and refines the simulation by evaluating the simulated state vector at
`NSteps − 1`

intermediate points. Although
`simByTransition`

does not report the output state
vector at these intermediate points, the refinement improves accuracy by
enabling the simulation to more closely approximate the underlying
continuous-time process.

**Data Types: **`double`

`StorePaths`

— Flag for storage and return method

`True`

(default) | logical with values `True`

or
`False`

Flag for storage and return method that indicates how the output array
`Paths`

is stored and returned, specified as the
comma-separated pair consisting of `'StorePaths'`

and a
scalar logical flag with a value of `True`

or
`False`

.

If

`StorePaths`

is`True`

(the default value) or is unspecified, then`simByTransition`

returns`Paths`

as a three-dimensional time series array.If

`StorePaths`

is`False`

(logical`0`

), then`simByTransition`

returns the`Paths`

output array as an empty matrix.

**Data Types: **`logical`

`MonteCarloMethod`

— Monte Carlo method to simulate stochastic processes

`"standard"`

(default) | string with values `"standard"`

,
`"quasi"`

, or
`"randomized-quasi"`

| character vector with values `'standard'`

,
`'quasi'`

, or
`'randomized-quasi'`

Monte Carlo method to simulate stochastic processes, specified as the
comma-separated pair consisting of `'MonteCarloMethod'`

and a string or character vector with one of the following values:

`"standard"`

— Monte Carlo using pseudo random numbers that has a convergence rate of*O*(*N*^{-½}).`"quasi"`

— Quasi-Monte Carlo rate of convergence is faster than standard Monte Carlo with errors approaching size of*O*(*N*^{-1}).`"randomized-quasi"`

— Quasi-random sequences, also called low-discrepancy sequences, are deterministic uniformly distributed sequences which are specifically designed to place sample points as uniformly as possible.

**Data Types: **`string`

| `char`

`QuasiSequence`

— Low discrepancy sequence to drive stochastic processes

`"sobol"`

(default) | string with value `"sobol"`

| character vector with value `'sobol'`

Low discrepancy sequence to drive the stochastic processes, specified
as the comma-separated pair consisting of
`'QuasiSequence'`

and a string or character vector
with the following value:

`"sobol"`

— Quasi-random low-discrepancy sequences that use a base of two to form successively finer uniform partitions of the unit interval and then reorder the coordinates in each dimension

**Note**

If `MonteCarloMethod`

option is not specified
or specified as`"standard"`

,
`QuasiSequence`

is ignored.

**Data Types: **`string`

| `char`

`BrownianMotionMethod`

— Brownian motion construction method

`"standard"`

(default) | string with value `"brownian-bridge"`

or
`"principal-components"`

| character vector with value `'brownian-bridge'`

or
`'principal-components'`

Brownian motion construction method, specified as the comma-separated
pair consisting of `'BrownianMotionMethod'`

and a
string or character vector with one of the following values:

`"standard"`

— The Brownian motion path is found by taking the cumulative sum of the Gaussian variates.`"brownian-bridge"`

— The last step of the Brownian motion path is calculated first, followed by any order between steps until all steps have been determined.`"principal-components"`

— The Brownian motion path is calculated by minimizing the approximation error.

The starting point for a Monte Carlo simulation is the construction of a Brownian motion sample path (or Wiener path). Such paths are built from a set of independent Gaussian variates, using either standard discretization, Brownian-bridge construction, or principal components construction.

Both standard discretization and Brownian-bridge construction share
the same variance and, therefore, the same resulting convergence when
used with the `MonteCarloMethod`

using pseudo random
numbers. However, the performance differs between the two when the
`MonteCarloMethod`

option
`"quasi"`

is introduced, with faster convergence
seen for `"brownian-bridge"`

construction option and
the fastest convergence when using the
`"principal-components"`

construction
option.

**Data Types: **`string`

| `char`

`Processes`

— Sequence of end-of-period processes or state vector adjustments

`simByEuler`

makes no adjustments and
performs no processing (default) | function | cell array of functions

Sequence of end-of-period processes or state vector adjustments,
specified as the comma-separated pair consisting of
`'Processes'`

and a function or cell array of
functions of the form

$${X}_{t}=P(t,{X}_{t})$$

.

`simByTransition`

applies processing functions at the
end of each observation period. The processing functions accept the
current observation time *t* and the current state
vector *X*_{t},
and return a state vector that may adjust the input state.

If you specify more than one processing function,
`simByTransition`

invokes the functions in the
order in which they appear in the cell array.

**Data Types: **`cell`

| `function`

## Output Arguments

`Paths`

— Simulated paths of correlated state variables

array

Simulated paths of correlated state variables, returned as an
```
(NPeriods +
1)
```

-by-`NVars`

-by-`NTrials`

three-dimensional time series array.

For a given trial, each row of `Paths`

is the transpose
of the state vector
*X*_{t} at time
*t*. When the input flag
`StorePaths`

= `False`

,
`simByTransition`

returns `Paths`

as
an empty matrix.

`Times`

— Observation times associated with simulated paths

column vector

Observation times associated with the simulated paths, returned as an
`(NPeriods + 1)`

-by-`1`

column vector.
Each element of `Times`

is associated with the
corresponding row of `Paths`

.

## More About

### Transition Density Simulation

The SDE has no solution such that
*r*(*t*) =
*f*(*r*(0),⋯).

In other words, the equation is not explicitly solvable. However, the transition density for the process is known.

The exact simulation for the distribution of
*r*(*t*_1
),⋯,*r*(*t*_*n*) is that of
the process at times
*t*_1,⋯,*t*_*n* for the same
value of *r*(0). The transition density for this process is known
and is expressed as

$$\begin{array}{l}r(t)=\frac{{\sigma}^{2}(1-{e}^{-\alpha (t-u)}}{4\alpha}{x}_{d}^{2}\left(\frac{4\alpha {e}^{-\alpha (t-u)}}{{\sigma}^{2}(1-{e}^{-\alpha (t-u)})}r(u)\right),t>u\\ \text{where}\\ d\equiv \frac{4b\alpha}{{\sigma}^{2}}\end{array}$$

## Algorithms

Use the `simByTransition`

function to simulate any vector-valued CIR
process of the form

$$d{X}_{t}=S(t)[L(t)-{X}_{t}]dt+D(t,{X}_{t}^{\frac{1}{2}})V(t)d{W}_{t}$$

where

*X*is an_{t}`NVars`

-by-`1`

state vector of process variables.*S*is an`NVars`

-by-`NVars`

matrix of mean reversion speeds (the rate of mean reversion).*L*is an`NVars`

-by-`1`

vector of mean reversion levels (long-run mean or level).*D*is an`NVars`

-by-`NVars`

diagonal matrix, where each element along the main diagonal is the square root of the corresponding element of the state vector.*V*is an`NVars`

-by-`NBrowns`

instantaneous volatility rate matrix.*dW*is an_{t}`NBrowns`

-by-`1`

Brownian motion vector.

## References

[1] Glasserman, P. *Monte Carlo Methods in Financial
Engineering*, New York: Springer-Verlag, 2004.

## Version History

**Introduced in R2018b**

### R2022b: Perform Brownian bridge and principal components construction

Perform Brownian bridge and principal components construction using the name-value
argument `BrownianMotionMethod`

.

### R2022a: Perform Quasi-Monte Carlo simulation

Perform Quasi-Monte Carlo simulation using the name-value arguments
`MonteCarloMethod`

and
`QuasiSequence`

.

## Open Example

You have a modified version of this example. Do you want to open this example with your edits?

## MATLAB Command

You clicked a link that corresponds to this MATLAB command:

Run the command by entering it in the MATLAB Command Window. Web browsers do not support MATLAB commands.

Select a Web Site

Choose a web site to get translated content where available and see local events and offers. Based on your location, we recommend that you select: .

You can also select a web site from the following list:

## How to Get Best Site Performance

Select the China site (in Chinese or English) for best site performance. Other MathWorks country sites are not optimized for visits from your location.

### Americas

- América Latina (Español)
- Canada (English)
- United States (English)

### Europe

- Belgium (English)
- Denmark (English)
- Deutschland (Deutsch)
- España (Español)
- Finland (English)
- France (Français)
- Ireland (English)
- Italia (Italiano)
- Luxembourg (English)

- Netherlands (English)
- Norway (English)
- Österreich (Deutsch)
- Portugal (English)
- Sweden (English)
- Switzerland
- United Kingdom (English)