# stability: Calculates stability values for results of 'lmInfl', 'lmMult'... In reverseR: Linear Regression Stability to Significance Reversal

## Description

This function calculates stability values for LOO (lmInfl), LMO (lmMult) and response value shifting/addition (lmThresh).

## Usage

 1 stability(x, pval = FALSE, ...) 

## Arguments

 x a result of either lmInfl, lmMult or lmThresh. pval logical. If TRUE, for lmThresh, objects an exact p-value is calculated for a future response to reverse significance. ... other parameters, not yet implemented.

## Details

For results of lmInfl:
A [0, 1]-bounded stability measure S = 1-\frac{n}{N}, with n = number of influencers (significance reversers) and N = total number of response values.

For results of lmMult:
For each 1...max, the percentage of all resamples that did *NOT* result in significance reversal.

For results of lmThresh:
A [0, 1]-bounded stability measure S = 1-\frac{n}{N}, with n = number of response values where one of the ends of the significance region is within the prediction interval and N = total number of response values.
If pval = TRUE, the exact p-value is calculated in the following manner:

1) Mean square error (MSE) and prediction standard error (se) are calculated from the linear model:

\mathrm{MSE} = ∑_{i=1}^n \frac{(y_i - \hat{y}_i)^2}{n-2} \quad\quad \mathrm{se}_i = √{\mathrm{MSE} \cdot ≤ft(1 + \frac{1}{n} + \frac{(x_i - \bar{x}_i)^2}{∑_{i=1}^n (x_i - \bar{x}_i)^2}\right)}

2) Upper and lower prediction intervals boundaries are calculated for each \hat{y}_i:

\hat{y}_i \pm Q_t(α/2, n-2) \cdot \rm{se}_i

The prediction interval around \hat{y}_i is a scaled/shifted t-distribution with density function

P_{tss}(y, n-2) = \frac{1}{\rm{se}_i} \cdot P_t≤ft(\frac{y - \hat{y}_i}{\rm{se}_i}, n-2\right)

, where P_t is the density function of the central, unit-variance t-distribution.
3) The probability of either shifting the response value (if lmThresh(..., newobs = FALSE)) or including a future response value y_{2i} (if lmThresh(..., newobs = TRUE)) to reverse the significance of the linear model is calculated as the integral between the end of the significance region (eosr) and the upper/lower α/2, 1-α/2 prediction interval:

## Value

The stability value.

## Author(s)

Andrej-Nikolai Spiess

## Examples

  1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 ## See examples in 'lmInfl' and 'lmThresh'. ## The implemented strategy of calculating the ## probability of significance reversal, as explained above ## and compared to 'stabPlot'. set.seed(125) a <- 1:20 b <- 5 + 0.08 * a + rnorm(length(a), 0, 1) LM1 <- lm(b ~ a) res1 <- lmThresh(LM1, newobs = TRUE) st1 <- stability(res1, pval = TRUE) ## Let's check that the prediction interval encompasses 95%: dt.scaled <- function(x, df, mu, s) 1/s * dt((x - mu)/s, df) integrate(dt.scaled, lower = st1$stats[1, "lower"], st1$stats[1, "upper"], df = 18, mu = st1$stats[1, "fitted"], s = st1$stats[1, "se"]) ## => 0.95 with absolute error < 8.4e-09 ## This is the interval between "end of significance region" and upper ## prediction boundary: integrate(dt.scaled, lower = st1$stats[1, "eosr.2"], st1$stats[1, "upper"], df = 18, mu = st1$stats[1, "fitted"], s = st1$stats[1, "se"]) ## => 0.09264124 with absolute error < 1e-15 ## We can recheck this value by P(B) - P(A): pt.scaled <- function(x, df, mu, s) pt((x - mu)/s, df) pA <- pt.scaled(x = st1$stats[1, "eosr.2"], df = 18, mu = st1$stats[1, "fitted"], s = st1\$stats[1, "se"]) 0.975 - pA ## => 0.09264124 as above 

reverseR documentation built on May 2, 2019, 10:59 a.m.