Lazy RecurrenceMatrix to compute certain rqa metrics

[felixcremer]

Describe the feature you'd like to have
We are using RecurrenceAnalysis.jl for computing the RQA Trend to detect changes in satellite time series.
When we upscaled this to larger areas ( meaning more time series) we found that a lot of time was spent in garbage collection and that a lot of the allocations are coming from the construction of the sparse recurrence matrix.
When we switched the algorithm to compute the tau_recurrence directly from our input time series we got a huge improvement in the number of allocations and thereby also in the timing of the RQA Trend computation.
This is the performance and memory consumption from saving the RecurrenceMatrix and then computing the Trend versus computing the tau_recurrence in one sweep and then computing the trend from there.
These are the compared functions:

"""
rqatrend(xin, xout, thresh)

Compute the RQA trend metric for the non-missing time steps of xin, and save it to xout. 
`thresh` specifies the epsilon threshold of the Recurrence Plot computation
"""
function rqatrend(pix_trend, pix, thresh=2)
    #replace!(pix, -9999 => missing)
    ts = collect(skipmissing(pix))
    #@show length(ts)
    tau_pix = tau_recurrence(ts, thresh)
    pix_trend .= RA._trend(tau_pix)
end

function rqatrend_matrix(pix_trend, pix, thresh=2)
    #replace!(pix, -9999 => missing)
    ts = collect(skipmissing(pix))
    rm = RecurrenceMatrix(ts, thresh)
    pix_trend .= RA.trend(rm)
end

And this is the number of allocations:

If possible, sketch out an implementation strategy

This is our implementation of the tau_recurrence function:

function tau_recurrence(ts::AbstractVector, thresh, metric=Euclidean())
    n = length(ts)
    rr_τ = zeros(n)
    for col in 1:n
        for row in 1:(col-1)
            d = evaluate(metric, ts[col], ts[row])
            #@show row, col, d
            rr_τ[col-row+1] += d <= thresh
        end
    end
    rr_τ[1] = n
    rr_τ ./ (n:-1:1)
    #rr_τ
end

I am envisioning to have a LazyRecurrenceMatrix type which would just keep the input vector and the arguments for the RecurrenceMatrix computation but where the actual matrix would not be immidiately computed but where the tau_recurrence would be computed like above.

I would be happy to provide a PR if you think this would be a useful addition to RecurrenceAnalysis.jl

[Datseris]

Hi,

So if I understand correctly you use the same _trend function and just alter the tau_recurrence.

In any case, what you propose here makes perfect sense to me.

This LazyRecurrenceMatrix type is totally fine and then you would dispatch some specific functions on it like trend. Then you would say in the docs: "for functions that are missing simply convert this explicitly to RecurrenceMatrix by calling RecurrenceMatrix(lazy). "

If I understood correctly however, your implementation for the tau_recurrence only works in the case of a fixed and constant recurrence threshold, right? Generally RecurrenceMatrix has several possibilities for how to specify the recurrence threshold. How would we handle this in the lazy type? Would you say the lazy type only works for a fixed recurrence threshold, or would you check at the function dispatch level (so your tau_recurrence would throw an error saying to convert to RecurrenceMatrix ?)

[pucicu]

I can confirm George's guess. It should only work for fixed threshold. Therefore, I would not call it an "improved" or "refined" algorithm. It is just a different algorithm which works fine in specific use cases (like the other implementation is fine for other use cases).
Nevertheless, make sure that it gives the same results for the same parameters (norm, fixed threshold, etc). Use unit tests for this.

Thanks for your contribution!

[felixcremer]

Yes, our current implementation would only work for a fixed threshold, and I haven't looked into the other threshold computations whether they could also be done lazily. I would throw the error at the initialisation of the LazyRecurrenceMatrix if the indicated threshold could not be computed lazily and hint to using RecurrenceMatrix instead.

We had some inital tests that our version of tau_recurrence gives the same tau_recurrence values as when we computed the RecurrencePlot eagerly and take the tau_recurrence then.

The "improved/refined algorithm" was only for our very specific use case and only in terms of memory and computation costs and I had the plots still laying around...

[Datseris]

that sounds fine to me, feel free to make a PR whenver you have the time!