Geek asks and answers

What function / package for reliable linear regression works with glmulti (i.e. behaves like glm)?

Background: multi-model output with glmulti

glmulti is an R function / package for automatic model selection for common linear models, which builds all possible common linear models taking into account the dependent variable and a set of predictors, is suitable for them using the classical glm and then allows for outputting several models (for example, using weighting coefficients models derived from AICc, BIC). glmulti also works in theory with any other function that returns coefficients, model log-likelihood and the number of free parameters (and possibly other information?) in the same format as glm.

My goal: multi-model output with reliable errors

I would like to use glmulti with robust error modeling of a quantitative dependent variable to protect against the outliers effect.

For example, I could assume that errors in a linear model are distributed as t distribution instead of the usual distribution. With its parameter kurtosis, the distribution of t can have heavy tails and, therefore, is more resistant to emissions (compared to the normal distribution).

However, I do not intend to use the t distribution approach. I am pleased with any approach that returns log verisimilitude and therefore works with the multimodal approach in glmulti. But this means that, unfortunately, I cannot use well-known linear models in R (for example, lmRob from robust or lmrob from robustbase ), since they do not work within the framework of the log-likelihood and, therefore, cannot work with glmulti.

Problem: I cannot find a reliable regression function that works with glmulti

R I, -, - heavyLm ( heavy); t. , heavyLm glmulti ( , ), S3 loglik (, , ).

:

library(glmulti)
library(heavy)

stackloss

head(stackloss)

:

summary(glm(stack.loss ~ ., data = stackloss))

glmulti glm

stackloss.glmulti <- glmulti(stack.loss ~ ., data = stackloss, level=1, crit=bic)
print(stackloss.glmulti)
plot(stackloss.glmulti)

t ( df = 4)

summary(heavyLm(stack.loss ~ ., data = stackloss))

glmulti, heavyLm

stackloss.heavyLm.glmulti <- glmulti(stack.loss ~ .,
 data = stackloss, level=1, crit=bic, fitfunction=heavyLm)

:

    Initialization...
    Error in UseMethod("logLik") :
    no applicable method for 'logLik' applied to an object of class "heavyLm".

,

logLik.heavyLm <- function(x){x$logLik}

glmulti -, :

    Initialization...
    Error in .jcall(molly, "V", "supplyErrorDF", 
    as.integer(attr(logLik(fitfunc(as.formula(paste(y,  : 
    method supplyErrorDF with signature ([I)V not found

: / glmulti (.. glm)?

, , heavyLm glmulti, , , -

  • , () - () glm (, , glmulti out -of ).
  • heavyLm, glmulti.

!

+5
1

, heavyLm. , , , , heavyLm ( Error in .jcall(molly, "V", "supplyErrorDF"…).

, glmulti , , , logLik.heavyLm; . logLik. , , , , (AIC, BIC,...). nobs.heavyLm .

:

nobs.heavyLm <- function(mdl) mdl$dims[1] # the sample size (number of data points)

logLik.heavyLm <- function(mdl) {
    res <- mdl$logLik
    attr(res, "nobs") <- nobs.heavyLm(mdl) # this is not really needed for 'glmulti', but is included to adhere to the format of 'logLik'
    attr(res, "df") <- length(mdl$coefficients) + 1 + 1 # I am also considering the scale parameter that is estimated; see mdl$family
    class(res) <- "logLik"
    res
}

, , :

Initialization...
TASK: Exhaustive screening of candidate set.
Fitting...
Completed.

> print(stackloss.glmulti)
glmulti.analysis
Method: h / Fitting: glm / IC used: bic
Level: 1 / Marginality: FALSE
From 8 models:
Best IC: 117.892471265874
Best model:
[1] "stack.loss ~ 1 + Air.Flow + Water.Temp"
Evidence weight: 0.709174196998897
Worst IC: 162.083142797858
2 models within 2 IC units.
1 models to reach 95% of evidence weight.

2 2 BIC.

: , . p + 1, p - , (+ 1) - "" ( ). logLik.heavyLm , " ", heavyLm , , , p + 1 + 1, , . , , , heavyLm ( ., 1980). - , ( ) , "" . , .

+1

More articles:


All Articles