RealFP types too restrictive for signalcorr.jl methods

[danmackinlay]

the methods in signcorr.jl require arrays of RealFP types;
This seems to me to be too restrictive since it, for example, the Dual types used by ForwardDiff.jl, and prevents differentiating through these functions.

Does anyone know a reason that we have this restriction? If I can get some clues about the reasons for this design choice I can propose a patch.

julia> using StatsBase:  autocov
julia> using ForwardDiff
julia> f2(x::Vector) = autocov(x, [0])[1]
julia> g2 = x -> ForwardDiff.gradient(f2, x);
julia> f2([1.0]), g2([1.0])

MethodError: no method matching fptype(::Type{ForwardDiff.Dual{ForwardDiff.Tag{typeof(f2),Float64},Float64,1}})
Closest candidates are:
  fptype(!Matched::Type{Complex{Float32}}) at /home/dan/.julia/packages/StatsBase/yNoh6/src/common.jl:27
  fptype(!Matched::Type{Complex{Float64}}) at /home/dan/.julia/packages/StatsBase/yNoh6/src/common.jl:28
  fptype(!Matched::Type{T<:Union{Bool, Float32, Int16, Int8, UInt16, UInt8}}) where T<:Union{Bool, Float32, Int16, Int8, UInt16, UInt8} at /home/dan/.julia/packages/StatsBase/yNoh6/src/common.jl:25
  ...

[nalimilan]

This code is quite old, it probably just hasn't been updated. As you can see, the error is due to the lack of an appropriate fptype method). I guess we could add a fallback definition fptype(::Type{T}) where {T} = float(T). Does that fix your problem? If so, feel free to make a pull request.

[danmackinlay]

Hm - possibly. Thanks @nalimilan. I will check that for ForwardDiff. It will for sure help with things like BigFloats, which is another limitation the current code has.

Just so I can understand here, though, what we are doing with these fptype methods is trying to find the smallest float type such that there will be no precision lost when the covariance calculation is conducted in float arithmetic. Is that right?

If so, that looks like the promotion methods built in to julia. Is there a reason not to use those built in promotion methods? signalcorr.jl doesn't do anything very exotic, so the custom promotion logic is a bit weird.

[danmackinlay]

To answer my own question, I suppose the answer here is in the docs:

Some loss is tolerated; for example, promote_type(Int64, Float64) returns Float64 even though strictly, not all Int64 values can be represented exactly as Float64 values.

It sounds like this policy is contentious.

So presumably this is our logic in StatsBase, that for operations in signalcorr.jl it is unacceptable to ever lose integer precision. This policy is understandable, but it leads to another question: 😉
why is the strictly-no-precision-loss policy only observed in signalcorr.jl but everywhere else in StatsBase (and AFAICT Julia core) it is acceptable to lose precision in extreme cases unless the user supplies a wide enough output array? fptypes is not used anywhere else.

[danmackinlay]

No wait! we do allow precision loss, since Int128 is "promoted" to Float64. OK I give in, I think @nalimilan is correct - it's probably just old code. I think we can use promote_type to replace it, unless anyone thinks that is foolish for some other reason?

[nalimilan]

Yes, actually fptype allows for more precision loss than promote_type. But here there's no other type with which to promote, so we should just use float(eltype(X)). It generally returns Float64 instead of a smaller type, but at least that's more standard than using custom rules just for one function.

Anyway, even if you use a small input type, you may want to perform computations in Float64 to avoid precision loss. And the standard way to choose a smaller type is to allocate the result manually and call autocov! on it.

[nalimilan]

Fixed by #567.