mad() is broken -- should promote integers to floats before calling mad!()

[dcarrera]

The mad function should promote integers to floats before calling mad!

julia> mad([1,2,3,4])
ERROR: InexactError()
 in mad! at /home/daniel/.julia/v0.4/StatsBase/src/scalarstats.jl:184
 in mad! at /home/daniel/.julia/v0.4/StatsBase/src/scalarstats.jl:183
 in mad at /home/daniel/.julia/v0.4/StatsBase/src/scalarstats.jl:179

[ararslan]

Would it make sense to define it like this?

function proposed_mad{T<:Real}(v::AbstractArray{T}, center::Real=median!(v); constant::Real=1.4826)
    S = typeof(constant * (one(T) - center))
    n = ndims(v)
    mad!(Array{S,n}(size(v)...), center; constant=constant)
end

Then I get

julia> proposed_mad([1,2,3,4])
3.7064999999999997

[ararslan]

cc @andreasnoack

[andreasnoack]

I think we should remove the optional center and constant from the mad but not the mad! function. I don't think anybody really uses them and they complicate matters a bit here. Your version modifies v through the center computation. If you remove the ! then v is copied, sorted and thrown away in the median computation just before we'd have to copy and sort it again in the body of the mad function. Even though the promotion rule you use is used (by me) all over the linear algebra code in base it's flawed so it should either be something like

S = promote_type(T,typeof(constant * (one(T) - center)))

or

promote_type(Float32,T)

Then size(v)... should also just be v, right? If you make a PR for this then it might also be a good idea to remove the method for ranges. It should be covered by the new one for AbstractArray.

[ararslan]

+1 for removing center and constant (especially constant; I'm still not clear on what it's actually for).

Sure, I'll make a PR for this and will use promote_type as you suggest.

[johnmyleswhite]

Just to be clear, constant is what make this function produce a consistent estimator of the standard deviation for normally distributed data. So if you intend to remove it, you need to either decide that (a) your estimator is inconsistent or (b) you don't provide any mechanism for producing the number some people believe MAD should generate.

[ararslan]

I think Andreas just meant remove constant from the mad API but leave it in the unexported mad!; it wouldn't be removed entirely.

[andreasnoack]

Yes. I just meant to remove it as a keyword argument for mad. It would still be applied to the results of median(abs(v - median(v))) so users that just call mad(x) wouldn't see a difference. I don't think anybody sets it and if they do, it's a just a multiplicative factor anyway to it would be simpler to just multiply than using a keyword.

That reminds me that we should avoid the promotion to Float64 that the multiplication with 1.4826 causes for e.g. Float32 input now that we try to be generic.

[ararslan]

Not sure how to get around multiplication with a hard-coded Float64 literal since we can't compute normal quantiles in StatsBase.

[andreasnoack]

1/(-sqrt(2*one(T))*erfcinv(3*one(T)/2)).

[ararslan]

Just when I think you can't get any more awesome 😎