google-research/torchsde

Correlated Brownian Motions?

mdiamore opened this issue · 2 comments

Hello!

First of all, fantastic work!! This framework is amazing and I've managed to implement some basic SDE's and solve them. My current goal is to go through some of the canonical SDEs used in Finance and implement them, which leads to my question.

Q: I discovered some slides pertaining to the Neural SDE as GANs paper, in which you reference the Heston model, in torchsde how would I link the two Brownian motions together by their correlation factor rho ? I think I managed to do it in a quick and dirty way below

class CorrelatedBrownianIntervals(BrownianInterval):

    def __init__(self, t0, t1, rho, w1, W=None, H=None):
        super().__init__(t0, t1, (batch_size, 1), W, H)
        self.rho = rho  # Correlation Coefficient between the 2 Brownians
        self.w1 = w1  # The first correlated Brownian

    def correlated_brownian(self, ta, tb):
        z1 = super().__call__(ta, tb)  # This generates an uncorrelated Brownian motion Z1
        corr = (self.rho * self.w1_val + torch.sqrt(1 - self.rho ** 2) * z1)  ## Creates the correlated Brownian W2
        return corr

    def __call__(self, ta, tb):
        self.w1_val = self.w1(tb) - self.w1(ta)
        return self.correlated_brownian(ta, tb)

Is this the right approach? Or is there a better / more flexible way to do this?

Thanks again!

Thanks for the positive review!

So I wouldn't try to create two different BrownianIntervals as you seem to be doing (self, w1). It's enough to create a single multidimensional BrownianInterval, and then apply the appropriate correlation. For example:

class CorrelatedBrownianInterval(BrownianInterval):
    def __init__(self, corr, **kwargs):
        super().__init__(**kwargs) 
        self.corr = corr

    def __call__(self, ta, tb):
        dW = super().__call__(ta, tb)
        # trick for computing batch matrix-vector product
        # matrix `corr` has shape (channels, channels)
        # batch of vectors `dW` has shape (batch, channels)
        return dW @ self.corr.t()

rho = ...
corr = torch.tensor([[1, 0], [rho, (1 - rho.pow(2)).sqrt()]])
CorrelatedBrownianInterval(corr=corr, size=(batch_size, 2), ...)

I would however note that the SDE solvers are generally written assuming that you're using Brownian motion -- off the top of my head I'm not sure if some of the higher-order solvers like SRK will necessarily produce the right numerics given an arbitrary noise process. (The simpler ones like Euler will no doubt still be fine.)

The safer/more general option is to put any desired correlation into the diffusion matrix of the SDE you want to solve. That is, instead of having g(Y_t) (corr dW) you have (g(Y_t) corr) dW.

class SDE:
    def __init__(self, corr, ...):
        ...
        self.corr = corr

    def _g(self, t, y):
        ...

    def g(self, t, y):
        return self._g(t, y) @ self.corr

And moreover the efficiency of this can be improved for most solvers by also providing

def g_prod(self, t, y, dW):
    return self._g(t, y) @ (dW @ self.corr.t())

so that one evaluates a matrix-(matrix-vector) product rather than a (matrix-matrix)-vector product.

Does that answer your question?


(All code untested.)

Yes, it does! I figured there was a better way to do it than creating two Brownian instances, and the more general form of adding the correlation directly to the diffusion makes a lot of sense. So I'll go with that approach. Thanks for your help, really great work!