Title: | Robust Score Testing in GLMs, by Sign-Flip Contributions |
---|---|
Description: | Provides robust tests for testing in GLMs, by sign-flipping score contributions. The tests are robust against overdispersion, heteroscedasticity and, in some cases, ignored nuisance variables. See Hemerik, Goeman and Finos (2020) <doi:10.1111/rssb.12369>. |
Authors: | Livio Finos [cre, aut] , Jelle J. Goeman [ctb], Jesse Hemerik [ctb], Riccardo De Santis [ctb] |
Maintainer: | Livio Finos <[email protected]> |
License: | GPL-2 |
Version: | 1.3.2 |
Built: | 2024-12-12 06:48:38 UTC |
Source: | CRAN |
This is the anova
method for flipscores
object. Remark: it performs type III deviance decomposition as in car::Anova
.
## S3 method for class 'flipscores' anova(object, model1 = NULL, score_type = NULL, n_flips = 5000, id = NULL, ...)
## S3 method for class 'flipscores' anova(object, model1 = NULL, score_type = NULL, n_flips = 5000, id = NULL, ...)
object |
(the object) |
model1 |
a |
score_type |
The type of score that is computed. It can be "orthogonalized", "effective" or "basic".
Default is "orthogonalized". "effective" and "orthogonalized" take into account the nuisance estimation. The default is |
n_flips |
The number of random flips of the score contributions.
When |
id |
a |
... |
other parameters allowed in |
set.seed(1) dt=data.frame(X=scale(rnorm(50)), Z=factor(rep(LETTERS[1:3],length.out=50))) dt$Y=rpois(n=nrow(dt),lambda=exp(dt$X*(dt$Z=="C"))) mod0=flipscores(Y~Z+X,data=dt,family="poisson") summary(mod0) anova(mod0) mod1=flipscores(Y~Z*X,data=dt,family="poisson") summary(mod1) anova(mod0,model1 = mod1)
set.seed(1) dt=data.frame(X=scale(rnorm(50)), Z=factor(rep(LETTERS[1:3],length.out=50))) dt$Y=rpois(n=nrow(dt),lambda=exp(dt$X*(dt$Z=="C"))) mod0=flipscores(Y~Z+X,data=dt,family="poisson") summary(mod0) anova(mod0) mod1=flipscores(Y~Z*X,data=dt,family="poisson") summary(mod1) anova(mod0,model1 = mod1)
Same usage as anova.glm
.
The parameter id
is used too,
if present in model0
(with priority) or in model1
.
compute_scores(model0, model1, score_type = "standardized", ...)
compute_scores(model0, model1, score_type = "standardized", ...)
model0 |
a |
model1 |
a |
score_type |
The type of score that is computed. It is "orthogonalized", "effective" or "basic". "effective" and "orthogonalized" take into account the nuisance estimation. |
... |
other arguments. |
Jesse Hemerik, Riccardo De Santis, Vittorio Giatti, Jelle Goeman and Livio Finos
set.seed(1) Z=rnorm(20) X=Z+rnorm(20) Y=rpois(n=20,lambda=exp(Z+X)) mod0=glm(Y~Z,family="poisson") X=data.frame(X=X) scr0=compute_scores(model0 = mod0, model1 = X) head(scr0)
set.seed(1) Z=rnorm(20) X=Z+rnorm(20) Y=rpois(n=20,lambda=exp(Z+X)) mod0=glm(Y~Z,family="poisson") X=data.frame(X=X) scr0=compute_scores(model0 = mod0, model1 = X) head(scr0)
Provides robust tests for testing in GLMs, by sign-flipping score contributions. The tests are often robust against overdispersion, heteroscedasticity and, in some cases, ignored nuisance variables.
flipscores(formula, family, data, score_type = "standardized", n_flips = 5000, alternative = "two.sided", id = NULL, seed = NULL, to_be_tested = NULL, flips = NULL, precompute_flips=TRUE, ...)
flipscores(formula, family, data, score_type = "standardized", n_flips = 5000, alternative = "two.sided", id = NULL, seed = NULL, to_be_tested = NULL, flips = NULL, precompute_flips=TRUE, ...)
formula |
see |
family |
see |
data |
see |
score_type |
The type of score that is computed. It can be "standardized" "orthogonalized", "effective" or "basic". Both "orthogonalized" and "effective" take into account the nuisance estimation and they provide the same test statistic. In case of small samples "effective score" might have a slight anti-conservative behaviour. "standardized effective score" gives a solution for this issue. "orthogonalized" has a similar intent, note however that in case of a big model matrix, it may be slow. |
n_flips |
The number of random flips of the score contributions. Overwritten with the |
alternative |
It can be "greater", "less" or "two.sided" (default) |
id |
a |
seed |
|
to_be_tested |
vector of indices or names of coefficients of the glm model to be tested (it is faster than computing every scores and p-values of course). |
flips |
matrix fo +1 or -1, the matrix has |
precompute_flips |
|
... |
see |
flipscores
borrows the same parameters from function glm
(and glm.nb
). See these helps for more details about parameters such as formula
,
data
, family
. Note: in order to use Negative Binomial family, family
reference must have quotes (i.e. family="negbinom"
).
Furthermore, flipscores
object contains two extra elements: scores
– i.e. a matrix of n score contributions, one column for each tested coefficient – and Tspace
– i.e. a matrix of size n_flips
times ncol(scores)
. The fist row of Tspace
contains column-wise the test statistics generated by randomly flipping the score contributions, each column refers to the same column of scores
, the vector of observed test statistics (i.e. no flips) is in the first row of Tspace
.
an object of class flipscores
.
See also its methods (summary.flipscores
, anova.flipscores
, print.flipscores
).
Livio Finos, Riccardo De Santis, Jesse Hemerik and Jelle Goeman
"Robust testing in generalized linear models by sign-flipping score contributions" by J.Hemerik, J.Goeman and L.Finos.
anova.flipscores
, summary.flipscores
, flip
set.seed(1) dt=data.frame(X=rnorm(20), Z=factor(rep(LETTERS[1:3],length.out=20))) dt$Y=rpois(n=20,lambda=exp(dt$Z=="C")) mod=flipscores(Y~Z+X,data=dt,family="poisson",n_flips=1000) summary(mod) # Equivalent to: model=glm(Y~Z+X,data=dt,family="poisson") mod2=flipscores(model) summary(mod2)
set.seed(1) dt=data.frame(X=rnorm(20), Z=factor(rep(LETTERS[1:3],length.out=20))) dt$Y=rpois(n=20,lambda=exp(dt$Z=="C")) mod=flipscores(Y~Z+X,data=dt,family="poisson",n_flips=1000) summary(mod) # Equivalent to: model=glm(Y~Z+X,data=dt,family="poisson") mod2=flipscores(model) summary(mod2)
Methods for flipscores
objects.
The following are methods to extract and manipulate relevant information from
a flipscores
object.
## S3 method for class 'flipscores' print(x, ...) ## S3 method for class 'flipscores' summary(object, ...)
## S3 method for class 'flipscores' print(x, ...) ## S3 method for class 'flipscores' summary(object, ...)
x |
a flipscores object |
... |
additional arguments to be passed |
object |
a flipscores object |
n_flips
xn_obs
matrix of random +1 and -1.
The first row is made by ones (i.e. the observed test statistic is computed)It creates a n_flips
xn_obs
matrix of random +1 and -1.
The first row is made by ones (i.e. the observed test statistic is computed)
make_flips(n_obs, n_flips)
make_flips(n_obs, n_flips)
n_obs |
number of observations |
n_flips |
number of flips |
# example code make_flips(n_obs=10,n_flips=20)
# example code make_flips(n_obs=10,n_flips=20)