JuliaStats / GLM.jl

Generalized linear models in Julia
Other
595 stars 114 forks source link

inexact error when computing aic on proportional dependent variable and Binomial() #402

Open jerlich opened 3 years ago

jerlich commented 3 years ago

In Julia 1.5.1 or 1.6. Here py is a vector of probabilities and y are "samples" drawn from py. py isn't really Binomial() but glm still fits the model. You can get the deviance of this model, modpy, but the aic, bic, loglikelihood all give inexact error

using GLM, StatsBase

begin
    x = randn(100).*3
    py = (1 ./ (1 .+ exp.(-0.5 .* x .+ 2.4 .+ rand(100))))
    y = rand(100) .< py
    mody = glm(@formula(y ~ x), Dict(:x=>x, :y=>y), Binomial())
    modpy = glm(@formula(y ~ x), Dict(:x=>x, :y=>py), Binomial())
    mody, modpy
        aic(mody) # OK! 
        aic(modpy) # error
        deviance.([modpy, mody]) # OK!
        loglikelihood(modpy) # Error
end
nalimilan commented 3 years ago

I'm not sure we can do anything about this, as Binomial(1.5, 1) throws an error. Contrary to us, R gives an AIC when estimating a binomial model with non-integer counts, but it prints a warning. So I'm not sure their AIC can be trusted. Do you have a reason to think that non-integer counts would make any statistical sense?

andreasnoack commented 1 week ago

Looks like the error is thrown from this package

julia> aic(modpy) # error
ERROR: InexactError: Float64(0.006358248591326209)
Stacktrace:
 [1] _safe_int(x::Float64)
   @ GLM ~/.julia/dev/GLM/src/glmtools.jl:510
 [2] loglik_obs
   @ ~/.julia/dev/GLM/src/glmtools.jl:528 [inlined]
 [3] loglikelihood(m::GeneralizedLinearModel{GLM.GlmResp{…}, GLM.DensePredChol{…}})
   @ GLM ~/.julia/dev/GLM/src/glmfit.jl:313
 [4] aic(model::GeneralizedLinearModel{GLM.GlmResp{…}, GLM.DensePredChol{…}})
   @ StatsAPI ~/.julia/packages/StatsAPI/zqIEd/src/statisticalmodel.jl:189
 [5] top-level scope
   @ REPL[270]:1