Numerical Methods for Stochastic Differential Equations
Lead Research Organisation:
University of Oxford
Department Name: Mathematical Institute
Abstract
I am interested in stochastic analysis, an area of mathematics that describes systems which have randomness. One of the classical tools for modelling such random phenomena is the stochastic differential equation (or SDE) driven by Brownian motion. In 1997, Robert Merton and Myron Scholes famously won the Nobel prize in economics for applying SDEs to the pricing of financial derivatives. In addition to mathematical finance, SDEs have seen widespread use within physics, economics and biology.
Just as for ordinary differential equations (ODEs), it is often necessary to numerically approximate the solution of an SDE. However, due to the highly oscillatory and random nature of Brownian motion, there are several challenges one faces when developing numerical methods for SDEs.
The first challenge is that, in order to produce methods with a high order of convergence, one must first estimate certain iterated integrals of Brownian motion and time. The study of such iterated integrals forms a core part of both stochastic numerics and rough path theory. In fact, a rough path is defined as a continuous process together with a collection of iterated integrals that satisfy both analytic and algebraic conditions. By understanding these Brownian integrals (as well as their probability distributions), one can hope to develop novel and effective numerical approximations for solutions of SDEs.
The second challenge is partly due to the random nature of Brownian motion. Since it has independent and normally distributed increments, there is always a small probability that a Brownian motion will experience large and sudden excursions. Furthermore, many real-world systems are not time-homogenous and can exhibit a range of behaviours depending on the state of the system. Therefore it is often necessary to vary the time steps used to propagate a numerical solution. From a mathematical viewpoint, the added complication is that a good choice of time step should be finalized after information of the Brownian path is generated. This means that there are very few theoretical results on the convergence of SDE methods that use an "adaptive stepsize". So although adaptive stepsize methods are commonly used for ODEs, the vast majority of SDE methods rely on a fixed stepsize.
The third challenge is simply designing numerical methods that are practical (i.e. numerically efficient and/or stable). One approach is the so-called "log-ODE" method, which is a natural technique for approximating an SDE solution with an ODE solution. A particularly appealing feature of the log-ODE method is that gives one the freedom to choose an appropriate ODE solver for a given (stochastic) problem. On the other hand, the high order versions of the method use derivatives of the SDE vector fields, which can be troublesome for high-dimensional or stiff systems.
In the DPhil thesis, I plan to develop mathematics that helps tackle these three challenges. In order to address the first challenge, I hope to identify and explore the relationships between certain iterated integrals of Brownian motion and time. Due to the pathwise nature of SDE simulation, I believe that rough path theory holds the key to understanding the convergence of adaptive stepsize methods. In addition, I am optimistic that one can develop state-of-the-art numerical discretizations for SDEs through creative applications of the log-ODE method.
This project falls within the EPSRC Mathematical aspects of operational research research area
Just as for ordinary differential equations (ODEs), it is often necessary to numerically approximate the solution of an SDE. However, due to the highly oscillatory and random nature of Brownian motion, there are several challenges one faces when developing numerical methods for SDEs.
The first challenge is that, in order to produce methods with a high order of convergence, one must first estimate certain iterated integrals of Brownian motion and time. The study of such iterated integrals forms a core part of both stochastic numerics and rough path theory. In fact, a rough path is defined as a continuous process together with a collection of iterated integrals that satisfy both analytic and algebraic conditions. By understanding these Brownian integrals (as well as their probability distributions), one can hope to develop novel and effective numerical approximations for solutions of SDEs.
The second challenge is partly due to the random nature of Brownian motion. Since it has independent and normally distributed increments, there is always a small probability that a Brownian motion will experience large and sudden excursions. Furthermore, many real-world systems are not time-homogenous and can exhibit a range of behaviours depending on the state of the system. Therefore it is often necessary to vary the time steps used to propagate a numerical solution. From a mathematical viewpoint, the added complication is that a good choice of time step should be finalized after information of the Brownian path is generated. This means that there are very few theoretical results on the convergence of SDE methods that use an "adaptive stepsize". So although adaptive stepsize methods are commonly used for ODEs, the vast majority of SDE methods rely on a fixed stepsize.
The third challenge is simply designing numerical methods that are practical (i.e. numerically efficient and/or stable). One approach is the so-called "log-ODE" method, which is a natural technique for approximating an SDE solution with an ODE solution. A particularly appealing feature of the log-ODE method is that gives one the freedom to choose an appropriate ODE solver for a given (stochastic) problem. On the other hand, the high order versions of the method use derivatives of the SDE vector fields, which can be troublesome for high-dimensional or stiff systems.
In the DPhil thesis, I plan to develop mathematics that helps tackle these three challenges. In order to address the first challenge, I hope to identify and explore the relationships between certain iterated integrals of Brownian motion and time. Due to the pathwise nature of SDE simulation, I believe that rough path theory holds the key to understanding the convergence of adaptive stepsize methods. In addition, I am optimistic that one can develop state-of-the-art numerical discretizations for SDEs through creative applications of the log-ODE method.
This project falls within the EPSRC Mathematical aspects of operational research research area
People |
ORCID iD |
| Terry Lyons (Primary Supervisor) | |
| James Foster (Student) |
Publications
Foster J
(2020)
An optimal polynomial approximation of Brownian motion
in SIAM Journal on Numerical Analysis
Studentship Projects
| Project Reference | Relationship | Related To | Start | End | Student Name |
|---|---|---|---|---|---|
| EP/N509711/1 | 01/10/2016 | 30/09/2021 | |||
| 1789673 | Studentship | EP/N509711/1 | 01/10/2016 | 31/03/2020 | James Foster |
| Description | My research focuses on the numerical analysis of stochastic differential equations (SDEs) driven by Brownian motion. I am particularly interested in the applications of SDEs in mathematical finance and physics, where numerical methods are often required to quantify the behaviour of random systems. Code and presentation slides for my research can be found at http://github.com/james-m-foster. The most well studied numerical methods for SDEs typically use only the discretized increments of the driving Brownian motion. In addition to generating increments, it is also straightforward to generate certain time integrals of Brownian motion. These quantities give extra information about the Brownian path and are known to improve the strong convergence of methods for one-dimensional SDEs. I believe my most important discoveries are primarily theoretical results on Brownian motion and its time integrals. I will briefly outline them as follows: 1. A new polynomial approximation of Brownian motion. I have identified a special class of orthogonal polynomials that can be efficiently used to approximate Brownian motion. The efficiency comes from the fact that the coefficients appearing in the expansion of Brownian motion in this polynomial basis are independent Gaussian random variables. Furthermore, these (random) coefficients can be expressed as time integrals of Brownian motion against certain polynomials. A Matlab demo can be found at: https://www.chebfun.org/examples/stats/RandomPolynomials.html (courtesy of Prof. Nick Trefethen). 2. New Lévy constructions of Brownian motion. The (random) coefficients from the previous result can be easily generated and carry useful information about the Brownian path. This means they can be utilized by high order numerical methods for SDEs. More specifically, they are applicable to methods that use predictable time steps to propagate the numerical solution. A predictable time step is one that is based solely of past behaviour of the Brownian motion. The advantage of having non-predictable time steps is that it helps alleviate the numerical error caused by the driving Brownian motion experiencing large and sudden excursions. In order for time integrals to be used by such variable stepsize methods, I discovered generalizations of Lévy construction of Brownian motion. These generalizations allow one to recursively generate Brownian increments with additional time integrals. Currently, the most significant consequence of this result is that certain high order methods for SDEs can now propagate solutions using non-predictable step sizes. 3. Conditional mean and variance of a third order iterated Brownian time integral. For general one-dimensional SDEs, the largest source of error that high order methods experience comes from approximating a third order iterated integral of Brownian motion and time. I have derived analytic formulae for the mean and variance of this iterated integral conditional on the increment and time integral of the Brownian motion over the same interval. The conditional mean of this third order integral can be used in numerical methods to improve accuracy, whilst the conditional variance can be used to estimate errors within a variable stepsize framework. 4. Conditional moments of multi-dimensional Brownian Lévy area. For general multi-dimensional SDEs, the largest source of error that numerical methods experience comes from approximating the second order iterated integral (or Lévy area) of Brownian motion. Furthermore, there are well established limitations for approximating Lévy area in a "pathwise" sense, e.g. by trapezium rule or Fourier series. This motivates one to further understand the distribution of Lévy area. I have derived analytic formulae for the first four moments of multi-dimensional Lévy area conditional on the increment and time integral of the Brownian motion over the same interval. The hope is that any random variable that matches all of these conditional moments can be used instead of Lévy area within numerical methods. 5. Efficient approximation of two-dimensional Brownian Lévy area. I have developed a method that efficiently generates a random variable conditional on the increment and time integral of the Brownian motion that matches the mean, variance, skewness and kurtosis of two-dimensional Lévy area. Moreover, this approximation exhibits very low numerical error in the 2-Wasserstein metric. I believe that this is a significant breakthrough for the development of numerical methods for SDEs driven by two-dimensional Brownian motion. 6. Convergence of a class of variable step size methods. Using rough path theory, I was able to prove that the new polynomial approximations of Brownian motion can give numerical methods that use non-predictable time steps and still have guaranteed convergence. I have extended this result to a larger class of variable stepsize methods (including Runge-Kutta methods). I believe that this gives a significant improvement to our understanding of numerical approximations for SDEs. 7. High order piecewise linear discretizations of Brownian motion. I have developed piecewise linear discretizations of Brownian motion that result in high order numerical methods for one-dimensional SDEs. Over each interval, the piecewise linear path is constructed to match the increment, time integral and estimated third order integral of the driving Brownian motion. One can then numerical solve the same system but driven instead by this piecewise linear path. A key feature of this approach is that is allows one to use numerical methods for ordinary differential equations (ODEs) to accurately discretize SDEs. 8. A three-stage stochastic Runge-Kutta method for SDEs with additive noise. I have designed an efficient numerical method for SDEs with additive noise where each step is locally close to a conditional expectation of the solution. Therefore, for certain problems, it could achieve a state-of-the-art L2 convergence. One potential disadvantage of this method is that it requires one to estimate the Laplacian of the vector field. This research is ongoing and the method has yet to be numerically tested. As well as the above results, I have investigated a number of specific numerical examples - which are listed below: 9. Inhomogeneous Geometric Brownian Motion (IGBM). In mathematical finance, IGBM is an example of a short rate model that can be both mean-reverting and non-negative. As a result, it is suitable for modelling interest rates, stochastic volatilities and default intensities. By incorporating the Brownian increment, time integral and conditional mean of a third order iterated integral into the so-called "log-ODE" method, I was able to produce a state-of-the-art pathwise discretization of IGBM. 10. Cox-Ingersoll-Ross model (CIR). In mathematical finance, the CIR process has similar applications as IGBM but is much more popular. For example, the famous Heston model uses the CIR process as its stochastic volatility. From a numerical standpoint, the CIR model is a challenging SDE to discretize due to its non-Lipschitz vector field. By combining ideas used to approximate IGBM with a variable stepsize methodology, I was able to produce a state-of-the-art pathwise discretization of the CIR process. 11. The underdamped Langevin equation. The underdamped Langevin equation is an important model within physics that describes a Brownian particle in a potential. More recently, it has appeared in the machine learning literature as it can be used as MCMC method. Using a high order piecewise linear discretization of Brownian motion and a third order symplectic ODE solver, I was able to accurately discretize the Langevin equation without having to compute derivatives of the vector field. Furthermore, due to the specific structure of the equation, this numerical method displays a third order convergence rate. I have recently been able to analyse the convergence of these underdamped Langevin ODEs, which has lead to other (potentially more efficient) ODE methods. |
| Exploitation Route | I have discovered a few general results that can be utilized in numerical methods for solutions of SDEs. Since SDEs are used within mathematical finance, physics, economics and biology to model random systems, I am optimistic that my findings have real-world applications. In any case, methods for the numerical examples I have already investigated can be used directly by practitioners. I am particularly excited about publishing my research on underdamped Langevin MCMC and the convergence of variable step size methods. I can envision several possible directions that my research can be taken forward. Firstly, I believe that new stochastic Runge-Kutta methods can be constructed which correctly use the increment, time integral and conditional mean of a third order integral of the driving Brownian motion. This is perhaps the most significant way in which my research can be improved upon. In addition, I have primarily focused on SDEs driven by a low-dimensional Brownian motion. I hope that aspects my research can be extended to the high-dimensional case, such as my results on Lévy area and piecewise linear discretizations. Similarly, it may also be possible to apply some of the key ideas of my research to numerical methods for stochastic PDEs with space-time white noise. Finally, there should be a generalization of my polynomial approximation theorem for fractional Brownian motion. That said, I am unsure if this would lead to more efficient methodologies for generating fractional Brownian motion. |
| Sectors | Digital/Communication/Information Technologies (including Software),Energy,Financial Services, and Management Consultancy,Other |
| URL | http://github.com/james-m-foster |
| Description | Numerical simulation of the Schramm-Loewner Evolution (SLE) |
| Organisation | New York University Shanghai |
| Country | China |
| Sector | Academic/University |
| PI Contribution | I developed a new "high order" approach for simulating the Schramm-Loewner evolution, which I implemented in Python. My collaborator and I believe this is the first high order numerical method that has been successfully applied to SLE traces. In addition, I helped analyse Taylor approximations of the Loewner equation using rough path theory. |
| Collaborator Contribution | My collaborator (Dr Vlad Margarint) focused on the analysis of Taylor approximations for the Loewner equation. In addition, he has a strong understanding of SLE theory and contributed to the writing of the preprint. |
| Impact | There are two main outputs of the project: 1. We present a new result concerning the convergence of Taylor approximations for the Loewner equation. 2. We present a new "high order" approach for simulating the Schramm-Loewner evolution, which is then implemented in Python. This may be the first high order numerical method that has been successfully applied to SLE traces. The project sits at the intersection of the following three areas of mathematics: 1. Stochastic differential equations (and their numerical analysis) 2. Rough path theory 3. SLE theory |
| Start Year | 2019 |
| Description | 10th Oxford-Berlin Young Researchers Meeting on Applied Stochastic Analysis |
| Form Of Engagement Activity | A talk or presentation |
| Part Of Official Scheme? | No |
| Geographic Reach | International |
| Primary Audience | Postgraduate students |
| Results and Impact | Approximately 40 young researchers (PhD students and postdocs) in applied stochastic analysis attended this meeting and 26 talks were given over the three day period. These meetings take place biannually and help maintain the strong research connections the University of Oxford has with TU and WIAS Berlin. During this event, I presented a high order variable stepsize discretization of the Cox-Ingersoll-Ross model. |
| Year(s) Of Engagement Activity | 2018 |
| URL | http://www.maths.ox.ac.uk/events/conferences/10th-oxford-berlin-young-researchers-meeting-applied-st... |
| Description | 11th Oxford-Berlin Young Researchers Meeting on Applied Stochastic Analysis |
| Form Of Engagement Activity | A talk or presentation |
| Part Of Official Scheme? | No |
| Geographic Reach | International |
| Primary Audience | Postgraduate students |
| Results and Impact | Approximately 40 young researchers (PhD students and postdocs) in applied stochastic analysis attended this meeting and around 25 talks were given over the three day period. These meetings take place biannually and help maintain the strong research connections the University of Oxford has with TU and WIAS Berlin. During this event, I presented a polynomial approximation of Brownian motion and discussed its applications to SDE simulation. |
| Year(s) Of Engagement Activity | 2019 |
| URL | https://www.wias-berlin.de/workshops/YRM2019/ |
| Description | 12th Oxford-Berlin Young Researchers Meeting on Applied Stochastic Analysis |
| Form Of Engagement Activity | A talk or presentation |
| Part Of Official Scheme? | No |
| Geographic Reach | International |
| Primary Audience | Postgraduate students |
| Results and Impact | Approximately 40 young researchers (PhD students and postdocs) in applied stochastic analysis attended this meeting and around 25 talks were given over the three day period. These meetings take place biannually and help maintain the strong research connections the University of Oxford has with TU and WIAS Berlin. During this event, I presented a high order numerical method for the Cox-Ingersoll-Ross process based on ODEs with piecewise linear drivers. |
| Year(s) Of Engagement Activity | 2019 |
| URL | https://www.maths.ox.ac.uk/events/conferences/12th-oxford-berlin-conference |
| Description | Oxford numerical analysis group internal seminar |
| Form Of Engagement Activity | A talk or presentation |
| Part Of Official Scheme? | No |
| Geographic Reach | Local |
| Primary Audience | Postgraduate students |
| Results and Impact | During each academic term, there are seminar talks organized by the Oxford numerical analysis group. In my presentation, I introduced an approximation of Brownian motion based on orthogonal polynomials. |
| Year(s) Of Engagement Activity | 2019 |
| URL | https://www.maths.ox.ac.uk/node/33109 |
| Description | Oxford stochastic analysis group seminar |
| Form Of Engagement Activity | A talk or presentation |
| Part Of Official Scheme? | No |
| Geographic Reach | Local |
| Primary Audience | Postgraduate students |
| Results and Impact | During each academic term, there are seminar talks organized by the Oxford stochastic analysis group. These talks are usually attended by members of the university's stochastic analysis, probability and mathematical finance research groups. In my presentation, I introduced a generalization of Lévy's construction of Brownian motion and discussed its potential impact within numerical methods for SDEs. |
| Year(s) Of Engagement Activity | 2018 |
| URL | https://www.maths.ox.ac.uk/node/27590 |
| Description | Oxford stochastic simulation seminar |
| Form Of Engagement Activity | A talk or presentation |
| Part Of Official Scheme? | No |
| Geographic Reach | Local |
| Primary Audience | Postgraduate students |
| Results and Impact | I have given three presentations at the Oxford stochastic simulation seminar, each to a small audience consisting of graduate students, postdoctoral researchers and faculty at the Oxford mathematics institute. The titles for these talks were 1. High order numerical simulation of the underdamped Langevin diffusion 2. Applications of Brownian polynomials to SDE simulation 3. An optimal polynomial approximation of Brownian motion Each talk sparked questions and discussion afterwards. |
| Year(s) Of Engagement Activity | 2019,2020 |