Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
8 changes: 4 additions & 4 deletions src/GLM.jl
Original file line number Diff line number Diff line change
Expand Up @@ -93,19 +93,19 @@ abstract type DensePred <: LinPred end # linear predictor with dense
abstract type LinPredModel <: RegressionModel end # model based on a linear predictor

const COMMON_FIT_KWARGS_DOCS = """
- `dropcollinear::Bool`: Controls whether or not a model matrix
- `dropcollinear::Bool=false`: Controls whether or not a model matrix
less-than-full rank is accepted.
If `true` (the default) the coefficient for redundant linearly dependent columns is
`0.0` and all associated statistics are set to `NaN`.
Typically from a set of linearly-dependent columns the last ones are identified as redundant
(however, the exact selection of columns identified as redundant is not guaranteed).
- `method::Symbol`: Controls which decomposition method to use.
- `method::Symbol=:qr`: Controls which decomposition method to use.
If `method=:qr` (the default), then the `QR` decomposition method will be used.
If `method=:cholesky`, then the `Cholesky` decomposition method will be used.
The Cholesky decomposition is faster and more computationally efficient than
QR, but is less numerically stable and thus may fail or produce less accurate
estimates for some models.
- `wts::AbstractWeights`: Weights of observations.
- `wts::AbstractWeights=uweights(0)`: Weights of observations.
The weights can be of type `AnalyticWeights`, `FrequencyWeights`,
`ProbabilityWeights`, or `UnitWeights`. `AnalyticWeights` describe a non-random
relative importance (usually between 0 and 1) for each observation. These weights may
Expand All @@ -114,7 +114,7 @@ const COMMON_FIT_KWARGS_DOCS = """
`ProbabilityWeights` represent the inverse of the sampling probability for each observation,
providing a correction mechanism for under- or over-sampling certain population groups. `UnitWeights`
(default) describes the case in which all weights are equal to 1 (so no weighting takes place)
- `contrasts::AbstractDict{Symbol}`: a `Dict` mapping term names
- `contrasts::AbstractDict{Symbol}=Dict{Symbol,Any}()`: a `Dict` mapping term names
(as `Symbol`s) to term types (e.g., `ContinuousTerm`) or contrasts
(e.g., `HelmertCoding()`, `SeqDiffCoding(; levels=["a", "b", "c"])`,
etc.). If contrasts are not provided for a variable, the appropriate
Expand Down
87 changes: 68 additions & 19 deletions src/glmfit.jl
Original file line number Diff line number Diff line change
Expand Up @@ -381,7 +381,7 @@ end
dof(obj::GeneralizedLinearModel) = linpred_rank(obj) + dispersion_parameter(obj.rr.d)

function _fit!(m::AbstractGLM, maxiter::Integer, minstepfac::Real,
atol::Real, rtol::Real, start)
atol::Real, rtol::Real, start::Union{AbstractVector,Nothing})
# Return early if model has the fit flag set
m.fit && return m

Expand Down Expand Up @@ -463,7 +463,7 @@ function StatsBase.fit!(m::AbstractGLM;
minstepfac::Real=0.001,
atol::Real=1e-6,
rtol::Real=1e-6,
start=nothing,
start::Union{AbstractVector,Nothing}=nothing,
kwargs...)
if haskey(kwargs, :verbose)
Base.depwarn("""`verbose` argument is deprecated, use `ENV["JULIA_DEBUG"]=GLM` instead.""",
Expand All @@ -483,7 +483,7 @@ end

const FIT_GLM_DOC = """
In the first method, `formula` must be a
[StatsModels.jl `Formula` object](https://juliastats.org/StatsModels.jl/stable/formula/)
[StatsModels.jl `FormulaTerm` object](https://juliastats.org/StatsModels.jl/stable/formula/)
and `data` a table (in the [Tables.jl](https://tables.juliadata.org/stable/) definition, e.g., a data frame).
In the second method, `X` must be a matrix holding values of the independent variable(s)
in columns (including if appropriate the intercept), and `y` must be a vector holding
Expand All @@ -494,20 +494,20 @@ const FIT_GLM_DOC = """

# Keyword Arguments
$COMMON_FIT_KWARGS_DOCS
- `offset::Vector=similar(y,0)`: offset added to `Xβ` to form `eta`. Can be of
length 0
- `offset::Union{AbstractVector{<:Real},Nothing}=nothing`: offset added to `Xβ` to form `eta`.
Can be of length 0.
- `maxiter::Integer=30`: Maximum number of iterations allowed to achieve convergence
- `atol::Real=1e-6`: Convergence is achieved when the relative change in
deviance is less than `max(rtol*dev, atol)`.
- `rtol::Real=1e-6`: Convergence is achieved when the relative change in
deviance is less than `max(rtol*dev, atol)`.
- `minstepfac::Real=0.001`: Minimum line step fraction. Must be between 0 and 1.
- `start::AbstractVector=nothing`: Starting values for beta. Should have the
- `start::Union{AbstractVector,Nothing}=nothing`: Starting values for beta. Should have the
same length as the number of columns in the model matrix.
"""

"""
fit(GeneralizedLinearModel, formula, data,
fit(GeneralizedLinearModel, formula::FormulaTerm, data,
distr::UnivariateDistribution, link::Link = canonicallink(d); <keyword arguments>)
fit(GeneralizedLinearModel, X::AbstractMatrix, y::AbstractVector,
distr::UnivariateDistribution, link::Link = canonicallink(d); <keyword arguments>)
Expand All @@ -524,15 +524,22 @@ function fit(::Type{M},
dropcollinear::Bool=true,
method::Symbol=:qr,
wts::Union{AbstractWeights,AbstractVector{<:Real}}=uweights(length(y)),
offset::AbstractVector{<:Real}=similar(y, 0),
fitargs...) where {M<:AbstractGLM}
offset::Union{AbstractVector{<:Real},Nothing}=nothing,
maxiter::Integer=30,
atol::Real=1e-6,
rtol::Real=1e-6,
minstepfac::Real=0.001,
start::Union{AbstractVector,Nothing}=nothing) where {M<:AbstractGLM}
# Check that X and y have the same number of observations
if size(X, 1) != size(y, 1)
throw(DimensionMismatch("number of rows in X and y must match"))
end

_wts = convert_weights(wts, length(y))
rr = GlmResp(y, d, l, offset, _wts)

off = offset === nothing ? similar(y, 0) : offset

rr = GlmResp(y, d, l, off, _wts)

if method === :cholesky
res = M(rr, cholpred(X, dropcollinear, _wts), nothing, false)
Expand All @@ -542,15 +549,26 @@ function fit(::Type{M},
throw(ArgumentError("The only supported values for keyword argument `method` are `:cholesky` and `:qr`."))
end

return fit!(res; fitargs...)
return fit!(res; maxiter, atol, rtol, minstepfac, start)
end

function fit(::Type{M},
X::AbstractMatrix,
y::AbstractVector,
d::UnivariateDistribution,
l::Link=canonicallink(d); kwargs...) where {M<:AbstractGLM}
return fit(M, float(X), float(y), d, l; kwargs...)
l::Link=canonicallink(d);
wts::AbstractVector{<:Real}=uweights(length(y)),
offset::Union{AbstractVector{<:Real},Nothing}=nothing,
dropcollinear::Bool=true,
method::Symbol=:qr,
maxiter::Integer=30,
atol::Real=1e-6,
rtol::Real=1e-6,
minstepfac::Real=0.001,
start::Union{AbstractVector,Nothing}=nothing) where {M<:AbstractGLM}
return fit(M, float(X), float(y), d, l;
offset, wts, dropcollinear, method,
maxiter, atol, rtol, minstepfac, start)
end

function fit(::Type{M},
Expand All @@ -559,13 +577,16 @@ function fit(::Type{M},
d::UnivariateDistribution,
l::Link=canonicallink(d);
offset::Union{AbstractVector,Nothing}=nothing,
wts::Union{AbstractVector,Nothing}=nothing,
wts::AbstractVector{<:Real}=uweights(0),
dropcollinear::Bool=true,
method::Symbol=:qr,
contrasts::AbstractDict{Symbol}=Dict{Symbol,Any}(),
fitargs...) where {M<:AbstractGLM}
maxiter::Integer=30,
atol::Real=1e-6,
rtol::Real=1e-6,
minstepfac::Real=0.001,
start::Union{AbstractVector,Nothing}=nothing) where {M<:AbstractGLM}
f, (y, X) = modelframe(f, data, contrasts, M)
wts = wts === nothing ? uweights(0) : wts
_wts = convert_weights(wts, length(y))
# Check that X and y have the same number of observations
if size(X, 1) != size(y, 1)
Expand All @@ -584,11 +605,11 @@ function fit(::Type{M},
throw(ArgumentError("The only supported values for keyword argument `method` are `:cholesky` and `:qr`."))
end

return fit!(res; fitargs...)
return fit!(res; maxiter, atol, rtol, minstepfac, start)
end

"""
glm(formula, data,
glm(formula::FormulaTerm, data,
distr::UnivariateDistribution, link::Link = canonicallink(distr); <keyword arguments>)
glm(X::AbstractMatrix, y::AbstractVector,
distr::UnivariateDistribution, link::Link = canonicallink(distr); <keyword arguments>)
Expand All @@ -597,7 +618,35 @@ Fit a generalized linear model to data. Alias for `fit(GeneralizedLinearModel, .

$FIT_GLM_DOC
"""
glm(X, y, args...; kwargs...) = fit(GeneralizedLinearModel, X, y, args...; kwargs...)
glm(X::AbstractMatrix, y::AbstractVector,
d::UnivariateDistribution, l::Link=canonicallink(d);
offset::Union{AbstractVector{<:Real},Nothing}=nothing,
wts::AbstractVector{<:Real}=uweights(length(y)),
dropcollinear::Bool=true,
method::Symbol=:qr,
maxiter::Integer=30,
atol::Real=1e-6,
rtol::Real=1e-6,
minstepfac::Real=0.001,
start::Union{AbstractVector,Nothing}=nothing) =
fit(GeneralizedLinearModel, X, y, d, l;
offset, wts, dropcollinear, method,
maxiter, atol, rtol, minstepfac, start)

glm(formula::FormulaTerm, data, d::UnivariateDistribution, l::Link=canonicallink(d);
offset::Union{AbstractVector{<:Real},Nothing}=nothing,
wts::AbstractVector{<:Real}=uweights(0),
dropcollinear::Bool=true,
method::Symbol=:qr,
contrasts::AbstractDict{Symbol}=Dict{Symbol,Any}(),
maxiter::Integer=30,
atol::Real=1e-6,
rtol::Real=1e-6,
minstepfac::Real=0.001,
start::Union{AbstractVector,Nothing}=nothing) =
fit(GeneralizedLinearModel, formula, data, d, l;
offset, wts, dropcollinear, method, contrasts,
maxiter, atol, rtol, minstepfac, start)

GLM.Link(r::GlmResp) = r.link
GLM.Link(m::GeneralizedLinearModel) = Link(m.rr)
Expand Down
32 changes: 17 additions & 15 deletions src/lm.jl
Original file line number Diff line number Diff line change
Expand Up @@ -148,19 +148,13 @@ const FIT_LM_DOC = """
"""

"""
fit(LinearModel, formula::FormulaTerm, data;
wts::AbstractWeights=uweights(0),
dropcollinear::Bool=true, method::Symbol=:qr,
contrasts::AbstractDict{Symbol}=Dict{Symbol,Any}())
fit(LinearModel, X::AbstractMatrix, y::AbstractVector;
wts::AbstractWeights=uweights(length(y)),
dropcollinear::Bool=true, method::Symbol=:qr)
fit(LinearModel, formula::FormulaTerm, data; <keyword arguments>)
fit(LinearModel, X::AbstractMatrix, y::AbstractVector; <keyword arguments>)

Fit a linear model to data.

$FIT_LM_DOC
"""

function fit(::Type{LinearModel}, X::AbstractMatrix{<:Real}, y::AbstractVector{<:Real};
wts::Union{AbstractWeights,AbstractVector{<:Real}}=uweights(length(y)),
dropcollinear::Bool=true, method::Symbol=:qr)
Expand Down Expand Up @@ -191,18 +185,26 @@ function fit(::Type{LinearModel}, f::FormulaTerm, data;
end

"""
lm(formula, data;
[wts::AbstractVector], dropcollinear::Bool=true, method::Symbol=:qr,
contrasts::AbstractDict{Symbol}=Dict{Symbol,Any}())
lm(X::AbstractMatrix, y::AbstractVector;
wts::AbstractVector=similar(y, 0), dropcollinear::Bool=true, method::Symbol=:cholesky)
lm(formula::FormulaTerm, data; <keyword arguments>)
lm(X::AbstractMatrix, y::AbstractVector; <keyword arguments>)

Fit a linear model to data.
An alias for `fit(LinearModel, X, y; wts=wts, dropcollinear=dropcollinear, method=method)`
An alias for `fit(LinearModel, ...)`.

$FIT_LM_DOC
"""
lm(X, y; kwargs...) = fit(LinearModel, X, y; kwargs...)
lm(X::AbstractMatrix, y::AbstractVector;
wts::Union{AbstractWeights,AbstractVector{<:Real}}=uweights(length(y)),
dropcollinear::Bool=true,
method::Symbol=:qr) =
fit(LinearModel, X, y; wts, dropcollinear, method)

lm(f::FormulaTerm, data;
wts::Union{AbstractWeights,AbstractVector{<:Real}}=uweights(0),
dropcollinear::Bool=true,
method::Symbol=:qr,
contrasts::AbstractDict{Symbol}=Dict{Symbol,Any}()) =
fit(LinearModel, f, data; wts, dropcollinear, method, contrasts)

dof(x::LinearModel) = linpred_rank(x.pp) + 1

Expand Down
Loading