# SDEProblem driven by two different NoiseProcess(es)
**URL:** https://discourse.julialang.org/t/sdeproblem-driven-by-two-different-noiseprocess-es/97608
**Category:** Modelling & Simulations
**Tags:** diffeq, sde
**Created:** [April 18, 2023, 9:56am UTC](https://discourse.julialang.org/t/sdeproblem-driven-by-two-different-noiseprocess-es/97608 "2023-04-18T09:56:27Z")
**Posts on this page:** 3
**Page:** 1
### Author:  [@tom-plaa](https://discourse.julialang.org/u/tom-plaa)
#### Post date: [April 18, 2023, 9:56am UTC](https://discourse.julialang.org/t/sdeproblem-driven-by-two-different-noiseprocess-es/97608/1 "2023-04-18T09:56:28Z")
How can we pass two different driving noise processes for an SDEProblem? The goal is to simulate a simple system of SDE’s:
dX(t) = f(X, Y, t)dt + g(X, Y, t)dZ\_1(t) + h(X, Y, t)dZ\_2(t)
dY(t) = \alpha(Y,t)dt + \beta(Y,t)dZ\_2(t)
To exemplify with a concrete example of a process with stochastic volatility:
dX(t) = bdt +Y(t)dZ\_1(t) - \rho dZ\_2(t)
dY(t) = \eta (\mu - Y(t))dt + dZ\_2(t)
Where Z\_2(t) is some subordinator (always positive and increasing).
Let’s assume we are using custom noise processes declared using the interface from `DiffEqNoiseProcess`, so no possibility of exploring the multivariate properties of Wiener processes. I am using my own noise processes here, let’s call them `MyProcess1` for Z\_1 and `MyProcess2` for Z\_2.
Both have completely different properties so I don’t see how to do this without creating a single `NoiseProcess` equivalent to Z(t) = [Z\_1(t), Z\_2(t)], which seems to imply a tremendous amount of juggling to define `dist` and `bridge` functions that act in different ways on each coordinate.
```julia
function f(du,u,p,t)
du[1] = b
du[2] = eta * (mu - u[2])
end
function g(du, u, p ,t)
du[1,1] = u[2]
du[1,2] = -rho
du[2,1] = 0
du[2,2] = 1
end
```
And then we would define our NoiseProcess in the way that I said above (quite troublesome in practice).
Is there any handier way of doing this? Something similar in principle to:
```julia
function g1(du, u, p ,t)
du[1] = u[2]
du[2] = 0
end
function g2(du, u, p ,t)
du[1] = -rho
du[2] = 1
end
SDEProblem(f, g1, g2, u0, tspan; noise=(MyProcess1, MyProcess2))
```
---
### Author:  [@ChrisRackauckas](https://discourse.julialang.org/u/ChrisRackauckas)
#### Post date: [April 19, 2023, 7:35am UTC](https://discourse.julialang.org/t/sdeproblem-driven-by-two-different-noiseprocess-es/97608/2 "2023-04-19T07:35:46Z")
> [@tom-plaa](#):
>
> Is there any handier way of doing this? Something similar in principle to:
There isn’t. That’s a cool idea but right now it cannot be done.
---
### Author:  [@tom-plaa](https://discourse.julialang.org/u/tom-plaa)
#### Post date: [April 19, 2023, 10:37am UTC](https://discourse.julialang.org/t/sdeproblem-driven-by-two-different-noiseprocess-es/97608/3 "2023-04-19T10:37:20Z")
Got you. For anyone curious, here’s a simple MWE on how to implement this the laborious way, based on the system in the first post and using a difference of gamma’s for the first process and a gamma subordinator for the second. I think this a generalization of a variance gamma process with independent gamma-driven OU volatility.
```julia
using DifferentialEquations
using DiffEqNoiseProcess
using Distributions
using Plots
tspan = (0.0, 1.0)
#just to have less boring parameter values
addnoise(x) = x * (1 + 0.1*(rand() - 0.5))
αplus, αminus, αsub = [addnoise(10) for i=1:3]
θplus, θminus, θsub = [addnoise(1) for i=1:3]
function TestDualNoiseProcess!()
function multidist!(dX, X, dt, u, p, t, rng)
dX[1] = rand(rng, Gamma(αplus*dt, θplus)) - rand(rng, Gamma(αminus*dt, θminus))
dX[2] = rand(rng, Gamma(αsub*dt, θsub))
return dX
end
noiseprocess = NoiseProcess{true}(0.0, zeros(2), nothing,
multidist!, nothing;
continuous=false
)
return noiseprocess
end
testnoise = TestDualNoiseProcess!()
testnoise_prob = NoiseProblem(testnoise, tspan)
testnoise_solve = solve(testnoise_prob, alg=EM(), dt=0.001)
# plot(testnoise_solve.t, transpose(reduce(hcat, testnoise_solve.u)))
function f!(du, u, p, t)
du[1] = 0.0 # b
du[2] = 1 - u[2]
end
function g!(du, u, p ,t)
du[1,1] = u[2]
du[1,2] = -0.01
du[2,1] = 0
du[2,2] = 1
end
testsde_prob = SDEProblem(f!, g!, [0.0, 1.0], tspan;
noise_rate_prototype = zeros(2,2),
noise=testnoise
)
testsde_solve = solve(testsde_prob, alg=EM(), dt=0.001)
plot(testsde_solve)
```
Here’s a plot of a realization for this particular system:
To somewhat automatize this, I’m trying to create a function that creates a vector noise from independent marginals of a specific kind, fed as a collection of these marginals. Or making it even more general, by building different NoiseProcesses with independent components by making their dists/bridges act separately on each corresponding coordinate. I’ll try to post a draft here. It would still require the user to define the noise integrand `g(du,u,p,t)` in a consistent manner.