Fix predict with confidence intervals

[nalimilan]

Use keyword arguments, which are passed through by the DataFrameRegressionModel method.

While we're at it, I wonder whether we should change the API a bit. interval=:confint and interval=:predint could be changed to interval=:confidence and interval=:prediction, which are the names used by R (since "int" is redundant with interval).

Cc: @mkborregaard

[codecov-io]

Codecov Report

Merging #253 into master will decrease coverage by 0.27%.
The diff coverage is 28.57%.

@@            Coverage Diff            @@
##           master    #253      +/-   ##
=========================================
- Coverage   51.97%   51.7%   -0.28%     
=========================================
  Files           6       6              
  Lines         583     588       +5     
=========================================
+ Hits          303     304       +1     
- Misses        280     284       +4

Impacted Files	Coverage Δ
src/lm.jl	49.46% <28.57%> (-1.68%)	⬇️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 4966a05...c2a5bd2. Read the comment docs.

[mkborregaard]

Wonderful. I also agree wholeheartedly with the name change.

[mkborregaard]

Does this mean this will work with DataFrames now? There was a PR on DataFrames I never got merged...

[nalimilan]

Does this mean this will work with DataFrames now? There was a PR on DataFrames I never got merged...

Do you have a link? It must have worked at some point, since it's in the examples, but I wonder since when it's been broken.

[mkborregaard]

JuliaData/DataFrames.jl#1160

I don't think it's ever worked from DataFrames...

[mkborregaard]

I'll rebase and write the test if that's still doable, but I'm guessing the interface has changed a lot?

[nalimilan]

Thanks. Actually I now realize I hadn't tested the example based on DataFrames correctly. The test in runtests.jl only covers matrices. So this PR shouldn't be merged without fixing DataFrameRegressionModel first.

The code now lives in StatsModels.jl. If we make interval-related arguments keyword arguments rather than positional, I think we can avoid defining a separate function as in JuliaData/DataFrames.jl#1160. Instead, you can just check whether yp is a vector or a matrix, and slightly adapt the code to suit each possibility.

[mkborregaard]

Here? https://github.com/JuliaStats/StatsModels.jl/blob/33f7b9cceac7e6624bf9e4802df3c4c8b2f5b83d/src/statsmodel.jl#L100

[mkborregaard]

Or here in GLM?

[nalimilan]

Yes, in StatsModels. GLM should be ready with this PR.

[mkborregaard]

By the way @nalimilan I just rememembered that the way I implemented it like it is, and not with a keyword as in this PR, is because the keyword version breaks type stability of predict. With confidence intervals it returns a three-column Matrix, without it's a Vector. Making the interval an optional positional argument allows the compiler to figure out the return type from the method signature.

[nalimilan]

AFAIK type-stability isn't a problem anymore with keyword arguments:

julia> f(; x=nothing) = x === nothing ? 1 : 1.0
f (generic function with 1 method)

julia> @code_warntype f()
Body::Int64
[...]

julia> @code_warntype f(x=nothing)
Body::Int64
[...]

julia> @code_warntype f(x=1)
Body::Float64
[...]

[mkborregaard]

sweet 🎉

[nalimilan]

I've pushed a commit to change the names. Is this good to go? IIRC a separate change is needed in StatsModels to support DataFrames, but that's already broken currently.

[mkborregaard]

Yes, I'll have a look a StatsModels some time during the week

[mkborregaard]

@nalimilan I can't see that the name has been changed? It still says :confint but I agree that :confidence makes more sense.

Here's the StatsModels PR:
JuliaStats/StatsModels.jl#77

I have a question - we currently return a three-row Matrix as [lower prediction upper]. Would it make more sense / be more useful to return a Tuple, e.g. (prediction, (lower, upper)) or maybe a NamedTuple? Allows things like pred, int = predict(mymodel, newdata, interval = :confidence); plot(pred, ribbon = int) etc.?

[mkborregaard]

Oh, I see we already discussed the tuple thing once: #171 (comment)
No reason to flog a dead horse :-)

[nalimilan]

Woops, I had just forgotten to push the new commits.

[mkborregaard]

I've added prediction intervals in a PR against this branch