Generating significant amount of FBM instances

[alexander-myltsev]

Thanks for the library. I needed to generate multiple FBM instances at the same time. I used your original implementation, but it appeared too slow since there is too much preliminary work.

My implementation than generates [m x n] matrix:

    def _hosking(self, gn):
        """Generate a fGn realization using Hosking's method.

        Method of generation is Hosking's method (exact method) from his paper:
        Hosking, J. R. (1984). Modeling persistence in hydrological time series
        using fractional differencing. Water resources research, 20(12),
        1898-1908.
        """
        fgn = np.zeros((self.m, self.n))
        phi = np.zeros(self.n)
        psi = np.zeros(self.n)
        # Monte carlo consideration
        if self._cov is None or self._changed:
            self._cov = np.array([self._autocovariance(i) for i in range(self.n)])
            self._changed = False

        # First increment from stationary distribution
        fgn[:, 0] = gn[:, 0]
        v = 1
        phi[0] = 0

        # Generate fgn realization with n increments of size 1
        for i in range(1, self.n):
            phi[i - 1] = self._cov[i]
            for j in range(i - 1):
                psi[j] = phi[j]
                phi[i - 1] -= psi[j] * self._cov[i - j - 1]
            phi[i - 1] /= v
            for j in range(i - 1):
                phi[j] = psi[j] - phi[i - 1] * psi[i - j - 2]
            v *= 1 - phi[i - 1] * phi[i - 1]
            for j in range(i):
                fgn[:, i] += phi[j] * fgn[:, i - j - 1]
            fgn[:, i] += np.sqrt(v) * gn[:, i]

        return fgn

Feel free to include it to your codebase if you wish.

[crflynn]

Hosking's method is O(N^2). You will get must faster realizations using the Davies-Harte method which is O(NlogN). Have you tried using it instead?

The code looks good to me. I maintain a more general package https://github.com/crflynn/stochastic which also support fBm. I don't plan on maintaining the fBm package since it's superceded by the stochastic package.

[alexander-myltsev]

Hi @crflynn , Thanks for your reply.

[alexander-myltsev]

@crflynn , I have a task as follows. I have an FBM trajectory of, say, 100 points each. Then I need to take, say, 80 values and generate from scratch the rest of 20 values. In other words, I need to generate a new FBM trajectory that matches at first 80 values from the first trajectory. It's obvious for me how to achieve it with Hosking method: I store first 80 values, then add 20 more at the last nested for loop.

On the other hand, Davies-Harte is much faster. I don't understand how to do it with Fourie transform. Do you think it's possible?

[crflynn]

AFAIK this would only be possible using Hosking's method.