dx: Diagnostics for binomial regression
In LogisticDx: Diagnostic Tests for Models with a Binomial Response
Description
Usage
Arguments
Value
Note
See Also
Examples
Description
Diagnostics for binomial regression
Returns diagnostic measures for a binary regression
model by covariate pattern
Usage
1
2
3
4
Arguments
x
A regression model with class glm and
x$family$family == "binomial".
...
Additional arguments
which can be passed to:
?stats::model.matrix
e.g. contrasts.arg which can be used for factor coding.
byCov
Return values by covariate pattern,
rather than by individual observation.
Value
A data.table, with rows sorted by dBhat.
If byCov==TRUE, there is one row per covariate pattern
with at least one observation.
The initial columns give the predictor variables 1 ... p.
Subsequent columns are labelled as follows:
y
The actual number of observations with
y=1 in the model data.
P
Probability of this covariate pattern.
This is given by the inverse of the link function,
x$family$linkinv. See:
?stats::family
n
Number of observations with these covariates.
If byCov=FALSE then this will be =1 for all observations.
yhat
The predicted number of observations having
a response of y=1, according to the model.
This is:
yhat[i] = n[i] * P[i]
h
Leverage, the diagonal of the hat matrix used to
generate the model:
H = V^0.5 X (X'VX)^-1 X'V^0.5
Here ^-1 is the inverse and
' is the transpose of a matrix.
X is the matrix of predictors, given by stats::model.matrix.
V is an N * N sparse matrix. All elements are
=0 except for the diagonal, which is:
v[i][i] = n[i]P[i] * (1 - P[i])
Leverage H is also the estimated covariance matrix of
Bhat.
Leverage is measure of the influence of this
covariate pattern on the model and is approximately
h[i] = x[i] - mean(x), 0.1 < P[i] < 0.9
That is, leverage is approximately equal to the distance of
the covariate pattern i from the mean mean(x).
For values of p which are large (>0.9) or
small (<0.1) this relationship no longer holds.
Pr
The Pearson residual, a measure of influence. This is:
Pr[i] = (y[i] - ybar) / SD[y]
where ybar and SD[y] refer
to the mean and standard deviation of a binomial distribution.
SD^2 = Var, is the variance.
E(y=1) = ybar = yhat = nP and SE[y] = (nP(1-P))^0.5
Thus:
Pr[i] = (y[i] - n[i]P[i]) / (n[i]P[i](1 - P[i])^0.5)
dr
The deviance residual, a measure of influence:
dr[i] = sign(y[i] - yhat[i]) * d[i]^0.5
d[i] is the contribution of observation i
to the model deviance.
The sign above is:
-
y[i] > yhat --> sign(i) = 1
-
y[i] = yhat --> sign(i) = 0
-
y[i] < yhat --> sign(i) = -1
In logistic regression this is:
y[i] = 1 --> dr[i] = (2 * log (1 + e^f(x) - f(x)))^0.5
y[i] = 0 --> dr[i] = (2 * log (1 + e^f(x)))
where f(x) is the linear function of the predictors
1 ... p:
f(x) = B[0] + B[1][i] * x[1][i] + ... + B[p][i] * x[p][i]
this is also:
dr[i] = sign(y[i] - yhat[i])
[2 * (y[i] * log(y[i] / n[i] * p[i])) +
(n[i] - y[i]) * log((n[i] - y[i]) / (n[i] * (1 - p[i])))]^0.5
To avoid the problem of division by zero:
y[i] = 0 --> dr[i] = (2 * n[i] * | log(1 - P[i]) |)^0.5
Similarly to avoid log(Inf):
y[i] = n[i] --> dr[i] = (2 * n[i] * | log(P[i]) |)^0.5
The above equations are used when calculating dr[i] by covariate group.
sPr
The standardized Pearson residual.
The residual is standardized by the leverage h[i]:
sPr[i] = Pr[i] / (1 - h[i])^0.5
sdr
The standardized deviance residual.
The residual is standardized by the leverage, as above:
sdr[i] = dr[i] / (1 - h[i])^0.5
dChisq
The change in the Pearson chi-square statistic
with observation i removed.
Given by:
dChi^2 = sPr[i]^2 = Pr[i]^2 / (1 - h[i])
where sPr[i] is the standardized Pearson residual,
Pr[i] is the Pearson residual and
h[i] is the leverage.
dChisq[i] should be <4
if the observation has little influence on the model.
dDev
The change in the deviance statistic
D = SUM dr[i]
with observation i excluded.
It is scaled by the leverage h[i] as above:
dDev[i] = sdr[i]^2 = dr[i]^2 / (1 - h[i])
dBhat
The change in Bhat
with observation i excluded.
This is scaled by the leverage as above:
dBhat = h[i] * sPr[i]^2 / (1 - h[i])
where sPR[i] is the standardized Pearson residual.
dBhat[i] should be <1
if the observation has little influence on the model coefficients.
Note
By default, values for the statistics are calculated by
covariate pattern.
Different values may be obtained if
calculated for each individual
obervation (e.g. rows in a data.frame).
Generally, the values calculated by
covariate pattern are preferred,
particularly where the number of observations in a group is >5.
In this case Pearsons chi-squared and the deviance statistic
should follow a chi-squared distribution with i - p degrees of freedom.
See Also
Examples
1
2
3
4
5
6
7
8
9
10
11
12
## H&L 2nd ed. Table 5.8. Page 182.
## Pattern nos. 31, 477, 468
data(uis)
uis <- within(uis, {
NDRGFP1 <- 10 / (NDRGTX + 1)
NDRGFP2 <- NDRGFP1 * log((NDRGTX + 1) / 10)
})
(d1 <- dx(g1 <- glm(DFREE ~ AGE + NDRGFP1 + NDRGFP2 + IVHX +
RACE + TREAT + SITE +
AGE:NDRGFP1 + RACE:SITE,
family=binomial, data=uis)))
d1[519:521, ]
Example output
(Intercept) AGE NDRGFP1 NDRGFP2 IVHXprevious IVHXrecent RACEother
1: 1 27 0.7142857 0.24033731 0 1 0
2: 1 26 0.3125000 0.36348463 0 1 0
3: 1 33 5.0000000 -8.04718956 0 0 0
4: 1 28 0.9090909 0.08664562 0 1 0
5: 1 33 0.3225806 0.36496842 0 1 0
---
517: 1 46 0.6250000 0.29375227 0 1 1
518: 1 40 10.0000000 -23.02585093 1 0 0
519: 1 40 10.0000000 -23.02585093 0 1 1
520: 1 41 10.0000000 -23.02585093 0 1 1
521: 1 24 0.4761905 0.35330350 1 0 0
TREATlong SITEB AGE:NDRGFP1 RACEother:SITEB y P n yhat
1: 0 0 19.28571 0 0 0.03236867 1 0.03236867
2: 1 0 8.12500 0 0 0.02857757 1 0.02857757
3: 0 1 165.00000 0 1 0.46433729 2 0.92867458
4: 0 0 25.45455 0 0 0.04242218 1 0.04242218
5: 0 0 10.64516 0 0 0.04048063 1 0.04048063
---
517: 1 1 28.75000 1 1 0.22337160 1 0.22337160
518: 1 1 400.00000 0 1 0.22004414 1 0.22004414
519: 0 0 400.00000 0 1 0.16759820 1 0.16759820
520: 0 0 410.00000 0 1 0.16262781 1 0.16262781
521: 0 1 11.42857 0 1 0.03262591 1 0.03262591
Pr dr h sPr sdr dChisq
1: -0.1828974 -0.2565312 0.004923582 -0.1833493 -0.2571650 0.03361696
2: -0.1715176 -0.2408064 0.006861302 -0.1721090 -0.2416368 0.02962152
3: 0.1011269 0.1009980 0.019260759 0.1021151 0.1019849 0.01042750
4: -0.2104793 -0.2944428 0.004987819 -0.2110062 -0.2951799 0.04452363
5: -0.2053983 -0.2874814 0.005968712 -0.2060140 -0.2883432 0.04244177
---
517: 1.8646300 1.7314263 0.053366647 1.9164688 1.7795620 3.67285273
518: 1.8826956 1.7400731 0.052449879 1.9341004 1.7875837 3.74074442
519: 2.2285985 1.8900719 0.043542719 2.2787624 1.9326158 5.19275812
520: 2.2691430 1.9059334 0.047027078 2.3244574 1.9523939 5.40310206
521: 5.4452261 2.6163519 0.009167194 5.4703578 2.6284273 29.92481422
dDev dBhat
1: 0.06613386 0.0001663348
2: 0.05838834 0.0002046463
3: 0.01040093 0.0002047859
4: 0.08713117 0.0002231890
5: 0.08314180 0.0002548438
---
517: 3.16684074 0.2070578164
518: 3.19545552 0.2070619666
519: 3.73500393 0.2364003227
520: 3.81184196 0.2666309753
521: 6.90862989 0.2768646411
(Intercept) AGE NDRGFP1 NDRGFP2 IVHXprevious IVHXrecent RACEother
1: 1 40 10.0000000 -23.0258509 0 1 1
2: 1 41 10.0000000 -23.0258509 0 1 1
3: 1 24 0.4761905 0.3533035 1 0 0
TREATlong SITEB AGE:NDRGFP1 RACEother:SITEB y P n yhat
1: 0 0 400.00000 0 1 0.16759820 1 0.16759820
2: 0 0 410.00000 0 1 0.16262781 1 0.16262781
3: 0 1 11.42857 0 1 0.03262591 1 0.03262591
Pr dr h sPr sdr dChisq dDev dBhat
1: 2.228599 1.890072 0.043542719 2.278762 1.932616 5.192758 3.735004 0.2364003
2: 2.269143 1.905933 0.047027078 2.324457 1.952394 5.403102 3.811842 0.2666310
3: 5.445226 2.616352 0.009167194 5.470358 2.628427 29.924814 6.908630 0.2768646
LogisticDx documentation built on May 2, 2019, 6:15 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Description Usage Arguments Value Note See Also Examples
Description
Diagnostics for binomial regression
Returns diagnostic measures for a binary regression model by covariate pattern
Usage
1 2 3 4 |
Arguments
x |
A regression model with class |
... |
Additional arguments
which can be passed to:
|
byCov |
Return values by covariate pattern, rather than by individual observation. |
Value
A data.table, with rows sorted by dBhat.
If byCov==TRUE, there is one row per covariate pattern
with at least one observation.
The initial columns give the predictor variables 1 ... p.
Subsequent columns are labelled as follows:
y |
The actual number of observations with y=1 in the model data. |
P |
Probability of this covariate pattern.
|
n |
Number of observations with these covariates.
|
yhat |
The predicted number of observations having
a response of y=1, according to the model.
yhat[i] = n[i] * P[i] |
h |
Leverage, the diagonal of the hat matrix used to generate the model: H = V^0.5 X (X'VX)^-1 X'V^0.5 Here ^-1 is the inverse and
' is the transpose of a matrix.
v[i][i] = n[i]P[i] * (1 - P[i]) Leverage H is also the estimated covariance matrix of
Bhat.
h[i] = x[i] - mean(x), 0.1 < P[i] < 0.9 That is, leverage is approximately equal to the distance of
the covariate pattern i from the mean mean(x).
|
Pr |
The Pearson residual, a measure of influence. This is: Pr[i] = (y[i] - ybar) / SD[y] where ybar and SD[y] refer
to the mean and standard deviation of a binomial distribution.
E(y=1) = ybar = yhat = nP and SE[y] = (nP(1-P))^0.5 Thus: Pr[i] = (y[i] - n[i]P[i]) / (n[i]P[i](1 - P[i])^0.5) |
dr |
The deviance residual, a measure of influence: dr[i] = sign(y[i] - yhat[i]) * d[i]^0.5 d[i] is the contribution of observation i
to the model deviance.
In logistic regression this is: y[i] = 1 --> dr[i] = (2 * log (1 + e^f(x) - f(x)))^0.5 y[i] = 0 --> dr[i] = (2 * log (1 + e^f(x))) where f(x) is the linear function of the predictors 1 ... p: f(x) = B[0] + B[1][i] * x[1][i] + ... + B[p][i] * x[p][i] this is also: dr[i] = sign(y[i] - yhat[i]) [2 * (y[i] * log(y[i] / n[i] * p[i])) + (n[i] - y[i]) * log((n[i] - y[i]) / (n[i] * (1 - p[i])))]^0.5 To avoid the problem of division by zero: y[i] = 0 --> dr[i] = (2 * n[i] * | log(1 - P[i]) |)^0.5 Similarly to avoid log(Inf): y[i] = n[i] --> dr[i] = (2 * n[i] * | log(P[i]) |)^0.5 The above equations are used when calculating dr[i] by covariate group. |
sPr |
The standardized Pearson residual.
sPr[i] = Pr[i] / (1 - h[i])^0.5 |
sdr |
The standardized deviance residual.
sdr[i] = dr[i] / (1 - h[i])^0.5 |
dChisq |
The change in the Pearson chi-square statistic with observation i removed. Given by: dChi^2 = sPr[i]^2 = Pr[i]^2 / (1 - h[i]) where sPr[i] is the standardized Pearson residual,
Pr[i] is the Pearson residual and
h[i] is the leverage.
|
dDev |
The change in the deviance statistic
D = SUM dr[i]
with observation i excluded.
dDev[i] = sdr[i]^2 = dr[i]^2 / (1 - h[i]) |
dBhat |
The change in Bhat
with observation i excluded.
dBhat = h[i] * sPr[i]^2 / (1 - h[i]) where sPR[i] is the standardized Pearson residual.
|
Note
By default, values for the statistics are calculated by
covariate pattern.
Different values may be obtained if
calculated for each individual
obervation (e.g. rows in a data.frame).
Generally, the values calculated by
covariate pattern are preferred,
particularly where the number of observations in a group is >5.
In this case Pearsons chi-squared and the deviance statistic
should follow a chi-squared distribution with i - p degrees of freedom.
See Also
Examples
1 2 3 4 5 6 7 8 9 10 11 12 | ## H&L 2nd ed. Table 5.8. Page 182.
## Pattern nos. 31, 477, 468
data(uis)
uis <- within(uis, {
NDRGFP1 <- 10 / (NDRGTX + 1)
NDRGFP2 <- NDRGFP1 * log((NDRGTX + 1) / 10)
})
(d1 <- dx(g1 <- glm(DFREE ~ AGE + NDRGFP1 + NDRGFP2 + IVHX +
RACE + TREAT + SITE +
AGE:NDRGFP1 + RACE:SITE,
family=binomial, data=uis)))
d1[519:521, ]
|
Example output
(Intercept) AGE NDRGFP1 NDRGFP2 IVHXprevious IVHXrecent RACEother
1: 1 27 0.7142857 0.24033731 0 1 0
2: 1 26 0.3125000 0.36348463 0 1 0
3: 1 33 5.0000000 -8.04718956 0 0 0
4: 1 28 0.9090909 0.08664562 0 1 0
5: 1 33 0.3225806 0.36496842 0 1 0
---
517: 1 46 0.6250000 0.29375227 0 1 1
518: 1 40 10.0000000 -23.02585093 1 0 0
519: 1 40 10.0000000 -23.02585093 0 1 1
520: 1 41 10.0000000 -23.02585093 0 1 1
521: 1 24 0.4761905 0.35330350 1 0 0
TREATlong SITEB AGE:NDRGFP1 RACEother:SITEB y P n yhat
1: 0 0 19.28571 0 0 0.03236867 1 0.03236867
2: 1 0 8.12500 0 0 0.02857757 1 0.02857757
3: 0 1 165.00000 0 1 0.46433729 2 0.92867458
4: 0 0 25.45455 0 0 0.04242218 1 0.04242218
5: 0 0 10.64516 0 0 0.04048063 1 0.04048063
---
517: 1 1 28.75000 1 1 0.22337160 1 0.22337160
518: 1 1 400.00000 0 1 0.22004414 1 0.22004414
519: 0 0 400.00000 0 1 0.16759820 1 0.16759820
520: 0 0 410.00000 0 1 0.16262781 1 0.16262781
521: 0 1 11.42857 0 1 0.03262591 1 0.03262591
Pr dr h sPr sdr dChisq
1: -0.1828974 -0.2565312 0.004923582 -0.1833493 -0.2571650 0.03361696
2: -0.1715176 -0.2408064 0.006861302 -0.1721090 -0.2416368 0.02962152
3: 0.1011269 0.1009980 0.019260759 0.1021151 0.1019849 0.01042750
4: -0.2104793 -0.2944428 0.004987819 -0.2110062 -0.2951799 0.04452363
5: -0.2053983 -0.2874814 0.005968712 -0.2060140 -0.2883432 0.04244177
---
517: 1.8646300 1.7314263 0.053366647 1.9164688 1.7795620 3.67285273
518: 1.8826956 1.7400731 0.052449879 1.9341004 1.7875837 3.74074442
519: 2.2285985 1.8900719 0.043542719 2.2787624 1.9326158 5.19275812
520: 2.2691430 1.9059334 0.047027078 2.3244574 1.9523939 5.40310206
521: 5.4452261 2.6163519 0.009167194 5.4703578 2.6284273 29.92481422
dDev dBhat
1: 0.06613386 0.0001663348
2: 0.05838834 0.0002046463
3: 0.01040093 0.0002047859
4: 0.08713117 0.0002231890
5: 0.08314180 0.0002548438
---
517: 3.16684074 0.2070578164
518: 3.19545552 0.2070619666
519: 3.73500393 0.2364003227
520: 3.81184196 0.2666309753
521: 6.90862989 0.2768646411
(Intercept) AGE NDRGFP1 NDRGFP2 IVHXprevious IVHXrecent RACEother
1: 1 40 10.0000000 -23.0258509 0 1 1
2: 1 41 10.0000000 -23.0258509 0 1 1
3: 1 24 0.4761905 0.3533035 1 0 0
TREATlong SITEB AGE:NDRGFP1 RACEother:SITEB y P n yhat
1: 0 0 400.00000 0 1 0.16759820 1 0.16759820
2: 0 0 410.00000 0 1 0.16262781 1 0.16262781
3: 0 1 11.42857 0 1 0.03262591 1 0.03262591
Pr dr h sPr sdr dChisq dDev dBhat
1: 2.228599 1.890072 0.043542719 2.278762 1.932616 5.192758 3.735004 0.2364003
2: 2.269143 1.905933 0.047027078 2.324457 1.952394 5.403102 3.811842 0.2666310
3: 5.445226 2.616352 0.009167194 5.470358 2.628427 29.924814 6.908630 0.2768646
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.