Документация Engee

Справка по API

Типы, определенные в пакете

LinearModel

Сочетание LmResp и LinPred

Члены

  • rr: объект LmResp

  • pp: объект LinPred

DensePredChol{T}

Тип LinPred с плотным разложением Холецкого для X’X

Члены

  • X: модельная матрица размером n × p, где n ≥ p. Должна быть матрицей полного столбцового ранга.

  • beta0: вектор базовых коэффициентов длиной p

  • delbeta: приращение до вектора коэффициентов, также длиной p

  • scratchbeta: вспомогательный вектор длиной p, используется в методе linpred!

  • chol: объект Cholesky, созданный в результате операции X’X, возможно, с использованием весов строк.

  • scratchm1: вспомогательная матрица Matrix{T} того же размера, что и X

  • scratchm2: вспомогательная матрица Matrix{T} того же размера, что и X’X

DensePredQR

Тип LinPred с плотным QR-разложением без выбора ведущего элемента для X

Члены

  • X: модельная матрица размером n × p, где n ≥ p. Должна быть матрицей полного столбцового ранга.

  • beta0: вектор базовых коэффициентов длиной p

  • delbeta: приращение до вектора коэффициентов, также длиной p

  • scratchbeta: вспомогательный вектор длиной p, используется в методе linpred!

  • qr: объект QRCompactWY, созданный в результате операции X, возможно, с использованием весов строк.

LmResp

Инкапсулирует отклик для линейной модели

Члены

  • mu: текущее значение вектора среднего отклика или подобранного значения

  • offset: необязательное смещение, добавляемое к линейному предиктору для образования mu

  • wts: необязательный вектор весов априорных частот для наблюдений

  • y: вектор наблюдаемого отклика

Любой из членов offset или wts либо оба они могут иметь длину 0

GlmResp

Вектор отклика и различные производные векторы в обобщенной линейной модели.

LinPred

Абстрактный тип, представляющий линейный предиктор

ModResp

Абстрактный тип, представляющий вектор отклика модели

Конструкторы моделей

Наиболее общий подход к подгонке модели — с помощью функции fit, как в

julia> using Random

julia> fit(LinearModel, hcat(ones(10), 1:10), randn(MersenneTwister(12321), 10))
LinearModel{GLM.LmResp{Vector{Float64}}, GLM.DensePredChol{Float64, LinearAlgebra.CholeskyPivoted{Float64, Matrix{Float64}, Vector{Int64}}}}:

Coefficients:
────────────────────────────────────────────────────────────────
        Coef.  Std. Error      t  Pr(>|t|)  Lower 95%  Upper 95%
────────────────────────────────────────────────────────────────
x1   0.717436    0.775175   0.93    0.3818  -1.07012    2.50499
x2  -0.152062    0.124931  -1.22    0.2582  -0.440153   0.136029
────────────────────────────────────────────────────────────────

Эту модель можно также подогнать как

julia> using Random

julia> lm(hcat(ones(10), 1:10), randn(MersenneTwister(12321), 10))
LinearModel{GLM.LmResp{Vector{Float64}}, GLM.DensePredChol{Float64, LinearAlgebra.CholeskyPivoted{Float64, Matrix{Float64}, Vector{Int64}}}}:

Coefficients:
────────────────────────────────────────────────────────────────
        Coef.  Std. Error      t  Pr(>|t|)  Lower 95%  Upper 95%
────────────────────────────────────────────────────────────────
x1   0.717436    0.775175   0.93    0.3818  -1.07012    2.50499
x2  -0.152062    0.124931  -1.22    0.2582  -0.440153   0.136029
────────────────────────────────────────────────────────────────
lm(formula, data, allowrankdeficient=false;
   [wts::AbstractVector], dropcollinear::Bool=true)
lm(X::AbstractMatrix, y::AbstractVector;
   wts::AbstractVector=similar(y, 0), dropcollinear::Bool=true)

Подгоняет линейную модель к данным. Псевдоним для fit(LinearModel, X, y; wts=wts, dropcollinear=dropcollinear).

В первом методе formula должно быть объектом Formula из StatsModels.jl, а data — таблицей (из определения Tables.jl, например фреймом данных). Во втором методе X должно быть матрицей, содержащей значения независимых переменных в столбцах (включая, если применимо, свободный коэффициент), а y — вектором, содержащим значения зависимых переменных.

Именованный аргумент wts может быть объектом Vector, в котором указываются веса частот для наблюдений. Вес эквивалентен повторению каждого наблюдения количество раз, равное весу. Обратите внимание, что при такой интерпретации получаются те же точечные оценки, но иные среднеквадратичные погрешности по сравнению с аналитическими весами (весами обратной дисперсии) и весами вероятности (выборочными весами), которые применяются по умолчанию в некоторых программах.

dropcollinear определяет, принимает ли lm модульную матрицу неполного ранга. При значении true (по умолчанию) используется только первый столбец из каждого набора линейно зависимых столбцов. Коэффициент для избыточных линейно зависимых столбцов равен 0.0, а все связанные статистические характеристики задаются равными NaN.

glm(formula, data,
    distr::UnivariateDistribution, link::Link = canonicallink(distr); <keyword arguments>)
glm(X::AbstractMatrix, y::AbstractVector,
    distr::UnivariateDistribution, link::Link = canonicallink(distr); <keyword arguments>)

Подгоняет обобщенную линейную модель к данным. Псевдоним для fit(GeneralizedLinearModel, ...).

В первом методе formula должно быть объектом Formula из StatsModels.jl, а data — таблицей (из определения Tables.jl, например фреймом данных). Во втором методе X должно быть матрицей, содержащей значения независимых переменных в столбцах (включая, если применимо, свободный коэффициент), а y — вектором, содержащим значения зависимых переменных. В обоих случаях в distr должно указываться распределение, а в link может указываться связывающая функция (если она не задана, используется каноническая связь для distr; в описании типа Link приведен список встроенных связей).

Именованные аргументы

  • dropcollinear::Bool=true: определяет, принимает ли lm модульную матрицу неполного ранга. При значении true (по умолчанию) коэффициент для избыточных линейно зависимых столбцов равен 0.0, а все связанные статистические характеристики задаются равными NaN. Обычно в наборе линейно зависимых столбцов избыточными считаются последние (хотя сказать с точностью, какие столбцы будут выбраны как избыточные, невозможно).

  • dofit::Bool=true: определяет, будет ли подогнана модель

  • wts::Vector=similar(y,0): веса априорных частот для наблюдений. Вес эквивалентен повторению каждого наблюдения количество раз, равное весу. Обратите внимание, что при такой интерпретации получаются те же точечные оценки, но иные среднеквадратичные погрешности по сравнению с аналитическими весами (весами обратной дисперсии) и весами вероятности (выборочными весами), которые применяются по умолчанию в некоторых программах. Нулевая длина означает, что взвешивание не производится (по умолчанию).

  • offset::Vector=similar(y,0): смещение, добавляемое к для образования eta. Может быть нулевой длины

  • verbose::Bool=false: для каждой итерации выводится информация о сходимости

  • maxiter::Integer=30: максимальное допустимое количество итераций, за которое должна достигаться сходимость

  • atol::Real=1e-6: сходимость достигается, когда относительное изменение отклонения меньше max(rtol*dev, atol).

  • rtol::Real=1e-6: сходимость достигается, когда относительное изменение отклонения меньше max(rtol*dev, atol).

  • minstepfac::Real=0.001: минимальная доля линейного шага. Значение должно быть в диапазоне от 0 до 1.

  • start::AbstractVector=nothing: начальные значения для беты. Длина должна совпадать с количеством столбцов в модельной матрице.

negbin(formula, data, [link::Link];
       <keyword arguments>)
negbin(X::AbstractMatrix, y::AbstractVector, [link::Link];
       <keyword arguments>)

Подгоняет отрицательную биномиальную обобщенную линейную модель к данным, одновременно оценивая параметр формы θ. Дополнительные аргументы и именованные аргументы передаются в glm.

В первом методе formula должно быть объектом Formula из StatsModels.jl, а data — таблицей (из определения Tables.jl, например фреймом данных). Во втором методе X должно быть матрицей, содержащей значения независимых переменных в столбцах (включая, если применимо, свободный коэффициент), а y — вектором, содержащим значения зависимых переменных. В обоих случаях в link может указываться связывающая функция (если она не задана, используется NegativeBinomial(θ)).

Именованные аргументы

  • initialθ::Real=Inf: начальное значение параметра формы θ. При значении Inf начальное значение оценивается путем подгонки распределения Пуассона.

  • maxiter::Integer=30: см. описание maxiter для glm

  • atol::Real=1.0e-6: см. описание atol для glm

  • rtol::Real=1.0e-6: см. описание rtol для glm

  • verbose::Bool=false: см. описание verbose для glm

fit(LinearModel, formula, data, allowrankdeficient=false;
   [wts::AbstractVector], dropcollinear::Bool=true)
fit(LinearModel, X::AbstractMatrix, y::AbstractVector;
    wts::AbstractVector=similar(y, 0), dropcollinear::Bool=true)

Fit a linear model to data.

In the first method, formula must be a StatsModels.jl Formula object and data a table (in the Tables.jl 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 values of the dependent variable.

The keyword argument wts can be a Vector specifying frequency weights for observations. Such weights are equivalent to repeating each observation a number of times equal to its weight. Do note that this interpretation gives equal point estimates but different standard errors from analytical (a.k.a. inverse variance) weights and from probability (a.k.a. sampling) weights which are the default in some other software.

dropcollinear controls whether or not lm accepts a model matrix which is less-than-full rank. If true (the default), only the first of each set of linearly-dependent columns is used. The coefficient for redundant linearly dependent columns is 0.0 and all associated statistics are set to NaN.

fit(GeneralizedLinearModel, formula, data,
    distr::UnivariateDistribution, link::Link = canonicallink(d); <keyword arguments>)
fit(GeneralizedLinearModel, X::AbstractMatrix, y::AbstractVector,
    distr::UnivariateDistribution, link::Link = canonicallink(d); <keyword arguments>)

Fit a generalized linear model to data.

In the first method, formula must be a StatsModels.jl Formula object and data a table (in the Tables.jl 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 values of the dependent variable. In both cases, distr must specify the distribution, and link may specify the link function (if omitted, it is taken to be the canonical link for distr; see Link for a list of built-in links).

Keyword Arguments

  • dropcollinear::Bool=true: Controls whether or not lm accepts a model matrix which is less-than-full rank. 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).

  • dofit::Bool=true: Determines whether model will be fit

  • wts::Vector=similar(y,0): Prior frequency (a.k.a. case) weights of observations. Such weights are equivalent to repeating each observation a number of times equal to its weight. Do note that this interpretation gives equal point estimates but different standard errors from analytical (a.k.a. inverse variance) weights and from probability (a.k.a. sampling) weights which are the default in some other software. Can be length 0 to indicate no weighting (default).

  • offset::Vector=similar(y,0): offset added to to form eta. Can be of length 0

  • verbose::Bool=false: Display convergence information for each iteration

  • 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 same length as the number of columns in the model matrix.

Методы для моделей

deviance(obj::LinearModel)

For linear models, the deviance is equal to the residual sum of squares (RSS).

dispersion(m::AbstractGLM, sqr::Bool=false)

Возвращает оценочный параметр дисперсии (или масштаба) для распределения модели, который обычно записывается как σ для линейных моделей и ϕ для обобщенных линейных моделей. По определению он равен 1 для семейств распределений Бернулли, биномиальных распределений и распределений Пуассона.

Если sqr имеет значение true, возвращается параметр квадратичной дисперсии.

ftest(mod::LinearModel)

Выполняет F-тест для определения того, подгоняется ли модель mod значительно лучше, чем нулевая модель (то есть такая модель, которая содержит только свободный коэффициент).

julia> dat = DataFrame(Result=[1.1, 1.2, 1, 2.2, 1.9, 2, 0.9, 1, 1, 2.2, 2, 2],
                       Treatment=[1, 1, 1, 2, 2, 2, 1, 1, 1, 2, 2, 2]);


julia> model = lm(@formula(Result ~ 1 + Treatment), dat);


julia> ftest(model.model)
F-test against the null model:
F-statistic: 241.62 on 12 observations and 1 degrees of freedom, p-value: <1e-07
ftest(mod::LinearModel...; atol::Real=0.0)

Для каждой последовательной пары линейных моделей в mod... выполняет F-тест для определения того, подгоняется ли одна модель значительно лучше, чем другая. Модели должны быть подогнаны на основе одних и тех же данных и вложены в прямом или обратном направлении.

Возвращается таблица, содержащая используемые степени свободы (DOF), отличие в DOF от предыдущей модели, сумма квадратов остатков (SSR), отличие в SSR от предыдущей модели, коэффициент детерминации, отличие в коэффициенте детерминации от предыдущей модели, а также F-статистику и p-значение для сравнения двух моделей.

Эту функцию можно использовать для выполнения ANOVA путем проверки относительного соответствия двух моделей данным

Необязательный именованный аргумент atol определяет числовой допуск при проверке вложенности моделей.

Примеры

Предположим, мы хотим сравнить эффект двух или более обработок на некоторый результат. Так как это ANOVA, нулевая гипотеза заключается в том, что Result ~ 1 соответствует данным так же хорошо, как Result ~ 1 + Treatment.

julia> dat = DataFrame(Result=[1.1, 1.2, 1, 2.2, 1.9, 2, 0.9, 1, 1, 2.2, 2, 2],
                       Treatment=[1, 1, 1, 2, 2, 2, 1, 1, 1, 2, 2, 2],
                       Other=categorical([1, 1, 2, 1, 2, 1, 3, 1, 1, 2, 2, 1]));


julia> nullmodel = lm(@formula(Result ~ 1), dat);


julia> model = lm(@formula(Result ~ 1 + Treatment), dat);


julia> bigmodel = lm(@formula(Result ~ 1 + Treatment + Other), dat);


julia> ftest(nullmodel.model, model.model)
F-test: 2 models fitted on 12 observations
─────────────────────────────────────────────────────────────────
     DOF  ΔDOF     SSR     ΔSSR      R²     ΔR²        F*   p(>F)
─────────────────────────────────────────────────────────────────
[1]    2        3.2292           0.0000
[2]    3     1  0.1283  -3.1008  0.9603  0.9603  241.6234  <1e-07
─────────────────────────────────────────────────────────────────

julia> ftest(nullmodel.model, model.model, bigmodel.model)
F-test: 3 models fitted on 12 observations
─────────────────────────────────────────────────────────────────
     DOF  ΔDOF     SSR     ΔSSR      R²     ΔR²        F*   p(>F)
─────────────────────────────────────────────────────────────────
[1]    2        3.2292           0.0000
[2]    3     1  0.1283  -3.1008  0.9603  0.9603  241.6234  <1e-07
[3]    5     2  0.1017  -0.0266  0.9685  0.0082    1.0456  0.3950
─────────────────────────────────────────────────────────────────
installbeta!(p::LinPred, f::Real=1.0)

Устанавливает pbeta0 .+= f * p.delbeta и обнуляет p.delbeta. Возвращает обновленное значение p.beta0.

nobs(obj::LinearModel)
nobs(obj::GLM)

For linear and generalized linear models, returns the number of rows, or, when prior weights are specified, the sum of weights.

nulldeviance(obj::LinearModel)

For linear models, the deviance of the null model is equal to the total sum of squares (TSS).

predict(mm::LinearModel, newx::AbstractMatrix;
        interval::Union{Symbol,Nothing} = nothing, level::Real = 0.95)

If interval is nothing (the default), return a vector with the predicted values for model mm and new data newx. Otherwise, return a vector with the predicted values, as well as vectors with the lower and upper confidence bounds for a given level (0.95 equates alpha = 0.05). Valid values of interval are :confidence delimiting the uncertainty of the predicted relationship, and :prediction delimiting estimated bounds for new data points.

predict(mm::AbstractGLM, newX::AbstractMatrix; offset::FPVector=eltype(newX)[],
        interval::Union{Symbol,Nothing}=nothing, level::Real = 0.95,
        interval_method::Symbol = :transformation)

Return the predicted response of model mm from covariate values newX and, optionally, an offset.

If interval=:confidence, also return upper and lower bounds for a given coverage level. By default (interval_method = :transformation) the intervals are constructed by applying the inverse link to intervals for the linear predictor. If interval_method = :delta, the intervals are constructed by the delta method, i.e., by linearization of the predicted response around the linear predictor. The :delta method intervals are symmetric around the point estimates, but do not respect natural parameter constraints (e.g., the lower bound for a probability could be negative).

Ссылки и применимые методы

Link

Абстрактный тип, подтипы которого ссылаются на связывающие функции.

GLM в настоящее время поддерживает следующие связи: CauchitLink, CloglogLink, IdentityLink, InverseLink, InverseSquareLink, LogitLink, LogLink, NegativeBinomialLink, PowerLink, ProbitLink, SqrtLink.

Подтипы Link должны реализовывать методы GLM.linkfun, GLM.linkinv, GLM.mueta, и GLM.inverselink.

Link01

Абстрактный подтип типа Link, который представляет связи, определенные в интервале (0, 1)

CauchitLink

Подтип Link01, соответствующий стандартному распределению Коши, Distributions.Cauchy.

CloglogLink

Подтип Link01, соответствующий распределению экстремального значения (логарифмическому распределению Вейбулла). Связь представляет собой комплементарное двойное логарифмическое преобразование, log(1 - log(-μ)).

IdentityLink

Каноническая связь Link для распределения Normal, определяемая как η = μ.

InverseLink

Каноническая связь Link для распределения Distributions.Gamma, определяемая как η = inv(μ).

InverseSquareLink

Каноническая связь Link для распределения Distributions.InverseGaussian, определяемая как η = inv(abs2(μ)).

LogitLink

Каноническая связь Link01 для Distributions.Bernoulli и Distributions.Binomial. Обратная связь, linkinv, — это интегральная функция распределения стандартного логистического распределения, Distributions.Logistic.

LogLink

Каноническая связь Link для Distributions.Poisson, определяемая как η = log(μ).

NegativeBinomialLink

Каноническая связь Link для распределения Distributions.NegativeBinomial, определяемая как η = log(μ/(μ+θ)). Для распределения, относящегося к экспоненциальному семейству, параметр формы θ должен быть фиксированным.

PowerLink

Связь Link, определяемая как η = μ^λ при λ ≠ 0 или как η = log(μ) при λ = 0, то есть класс преобразований, использующих степенную или логарифмическую функцию.

Многие другие связи являются особыми случаями PowerLink:

ProbitLink

Связь Link01, linkinv которой — это интегральная функция распределения стандартного нормального распределения, Distributions.Normal().

SqrtLink

Связь Link определяемая как η = √μ

GLM.linkfun(L::Link, μ::Real)

Возвращает η, значение линейного предиктора для связи L при среднем μ.

Примеры

julia> μ = inv(10):inv(5):1
0.1:0.2:0.9

julia> show(linkfun.(LogitLink(), μ))
[-2.197224577336219, -0.8472978603872036, 0.0, 0.8472978603872034, 2.1972245773362196]
GLM.linkinv(L::Link, η::Real)

Возвращает μ (среднее значение) для связи L при значении линейного предиктора η.

Примеры

julia> μ = 0.1:0.2:1
0.1:0.2:0.9

julia> η = logit.(μ);


julia> linkinv.(LogitLink(), η) ≈ μ
true
GLM.mueta(L::Link, η::Real)

Возвращает производную linkinv (dμ/dη) для связи L при значении линейного предиктора η.

Примеры

julia> mueta(LogitLink(), 0.0)
0.25

julia> mueta(CloglogLink(), 0.0) ≈ 0.36787944117144233
true

julia> mueta(LogLink(), 2.0) ≈ 7.38905609893065
true
GLM.inverselink(L::Link, η::Real)

Возвращает кортеж из трех элементов: обратной связи, ее производной и, если применимо, функции дисперсии μ*(1 - μ).

Если диапазон μ отличается от (0, 1), для функции дисперсии возвращается NaN

Примеры

julia> GLM.inverselink(LogitLink(), 0.0)
(0.5, 0.5, 0.25)

julia> μ, oneminusμ, variance = GLM.inverselink(CloglogLink(), 0.0);



julia> μ + oneminusμ ≈ 1
true

julia> μ*(1 - μ) ≈ variance
false

julia> isnan(last(GLM.inverselink(LogLink(), 2.0)))
true
canonicallink(D::Distribution)

Возвращает каноническую связь для распределения D, которое должно относиться к экспоненциальному семейству.

Примеры

julia> canonicallink(Bernoulli())
LogitLink()
GLM.glmvar(D::Distribution, μ::Real)

Возвращает значение функции дисперсии для D при μ

Дисперсия D при μ является произведением параметра дисперсии ϕ, который не зависит от μ, и значения glmvar. Иными словами, glmvar возвращает коэффициент дисперсии, который зависит от μ.

Примеры

julia> μ = 1/6:1/3:1;

julia> glmvar.(Normal(), μ)    # константа для Normal()
3-element Vector{Float64}:
 1.0
 1.0
 1.0

julia> glmvar.(Bernoulli(), μ) ≈ μ .* (1 .- μ)
true

julia> glmvar.(Poisson(), μ) == μ
true

julia> glmvar.(Geometric(), μ) ≈ μ .* (1 .+ μ)
true
GLM.mustart(D::Distribution, y, wt)

Возвращает начальное значение для μ.

Для некоторых распределений можно задать μ = y для инициализации алгоритма IRLS, но для других, в частности для распределений Бернулли, значения y не допускаются в качестве значений μ и должны быть изменены.

Примеры

julia> GLM.mustart(Bernoulli(), 0.0, 1) ≈ 1/4
true

julia> GLM.mustart(Bernoulli(), 1.0, 1) ≈ 3/4
true

julia> GLM.mustart(Binomial(), 0.0, 10) ≈ 1/22
true

julia> GLM.mustart(Normal(), 0.0, 1) ≈ 0
true

julia> GLM.mustart(Geometric(), 4, 1) ≈ 4
true
devresid(D, y, μ::Real)

Возвращает квадрат остатка отклонения μ от y для распределения D

Отклонение GLM можно вычислить как сумму квадратов остатков отклонений. Это основное назначение данных значений. Фактический остаток отклонения, например при построении графиков, — это квадратный корень со знаком из этого значения

sign(y - μ) * sqrt(devresid(D, y, μ))

Примеры

julia> devresid(Normal(), 0, 0.25) ≈ abs2(0.25)
true

julia> devresid(Bernoulli(), 1, 0.75) ≈ -2*log(0.75)
true

julia> devresid(Bernoulli(), 0, 0.25) ≈ -2*log1p(-0.25)
true
GLM.dispersion_parameter(D)

Есть ли у распределения D отдельный параметр дисперсии ϕ?

Возвращает false для распределений Bernoulli, Binomial и Poisson; в противном случае возвращает true.

Примеры

julia> show(GLM.dispersion_parameter(Normal()))
true
julia> show(GLM.dispersion_parameter(Bernoulli()))
false
GLM.loglik_obs(D, y, μ, wt, ϕ)

Возвращает wt * logpdf(D(μ, ϕ), y), где параметры D выводятся из μ и ϕ.

Аргумент wt — это множитель результата, кроме случая с Binomial, в котором wt — это количество испытаний, а μ — доля успешных исходов.

Правдоподобие подогнанной модели — это сумма данных значений для всех наблюдений.

cancancel(r::GlmResp{V,D,L})

Возвращает true, если dμ/dη для связи L является функцией дисперсии для распределения D

Когда L — это каноническая связь для D, производная обратной связи — это величина, кратная функции дисперсии для D. Если они совпадают, числитель и знаменатель в выражении для рабочих весов сокращаются.