# Simulating Stochastic Matrix Differential Equation with Arbitrary Autocorrelation Function

I want to numerically simulate a matrix differential equation that includes a stochastic (vector) Gaussian noise $$\mathbf{n}$$, where the different vector components are independent, and each component has an arbitrary autocorrelation function $$R_i(\tau)$$, so:

$$\mathbb{E} [n_i(t) n_j(t+\tau )] \equiv R_i(\tau) \ \delta_{ij}$$

(Let's say the process is zero mean). As a minimal example, consider the following matrix DE:

$$\dot{X}(t) = \Big({F_0}(t)+{n_1}(t) A + n_2(t) B \Big) X(t)$$ $$X(0) = \mathrm{id}_{2\times 2}$$

with

$$F_0(t) = \begin{bmatrix} \sin(t) & t \\ 0 & e^{-t} \end{bmatrix}$$

$$A = \begin{bmatrix} 1 & 0 \\ 0 & -1 \end{bmatrix}$$

$$B = \begin{bmatrix} 0 & -1 \\ 1 & 0 \end{bmatrix}$$

And let's say that I want the noise autocorrelation function to be: $$R_j(\tau) = e^{-j \tau^2} \qquad j \in \{1,2 \}$$

I looked at the documentation and I couldn't really find anything I could use to generate a continuous-time process with an arbitrary autocorrelation (I thought about using TransformedProcess for instance, but I couldn't figure out how exactly it can be done). Most other stochastic DEs solved in other questions are also done with more standard processes such as the Wiener process (see this and this answer). What would be an efficient way of simulating such systems where the noise has some non-standard spectral density?