Generate *p*-dimensional random vectors according to Stahel's
Barrow Wheel Distribution.

1 2 3 4 5 6 |

`n` |
integer, specifying the sample size. |

`p` |
integer, specifying the dimension (aka number of variables). |

`frac` |
numeric, the proportion of outliers.
The default, |

`sig1` |
thickness of the “wheel”, ( |

`sig2` |
thickness of the “axis” (compared to 1). |

`rGood` |
function; the generator for “good” observations. |

`rOut` |
function, generating the outlier observations. |

`U1` |
p-vector to which |

`scaleAfter` |
logical indicating if the matrix is re-scaled |

`scaleBefore` |
logical indicating if the matrix is re-scaled before
rotation (via |

`spherize` |
logical indicating if the matrix is to be
“spherized”, i.e., rotated and scaled to have empirical
covariance |

`fullResult` |
logical indicating if in addition to the |

....

By default (when `fullResult`

is `FALSE`

), an
*n x p* matrix of *n* sample vectors of the
*p* dimensional barrow wheel distribution, with an attribute,
`n1`

specifying the exact number of “good” observations,
*n1 ~= (1-f)*n*, *f = *`frac`

.

If `fullResult`

is `TRUE`

, a list with components

`X` |
the |

`X0` |
......... |

`A` |
the |

`n1` |
the number of “good” observations, see above. |

`n2` |
the number of “outlying” observations, |

Werner Stahel and Martin Maechler

http://stat.ethz.ch/people/maechler/robustness

Stahel, W.~A. and Mächler, M. (2009).
Comment on “invariant co-ordinate selection”,
*Journal of the Royal Statistical Society B* **71**, 584–586.

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 31 32 33 34 35 36 37 38 39 | ```
set.seed(17)
rX8 <- rbwheel(1000,8, fullResult = TRUE, scaleAfter=FALSE)
with(rX8, stopifnot(all.equal(X, X0 %*% A, tol = 1e-15),
all.equal(X0, X %*% t(A), tol = 1e-15)))
##--> here, don't need to keep X0 (nor A, since that is Qrot(p))
## for n = 100, you don't see "it", but may guess .. :
n <- 100
pairs(r <- rbwheel(n,6))
n1 <- attr(r,"n1") ; pairs(r, col=1+((1:n) > n1))
## for n = 500, you *do* see it :
n <- 500
pairs(r <- rbwheel(n,6))
## show explicitly
n1 <- attr(r,"n1") ; pairs(r, col=1+((1:n) > n1))
## but increasing sig2 does help:
pairs(r <- rbwheel(n,6, sig2 = .2))
## show explicitly
n1 <- attr(r,"n1") ; pairs(r, col=1+((1:n) > n1))
set.seed(12)
pairs(X <- rbwheel(n, 7, spherize=TRUE))
colSums(X) # already centered
if(require("ICS")) {
# ICS: Compare M-estimate [Max.Lik. of t_{df = 2}] with high-breakdown :
stopifnot(require("MASS"))
X.paM <- ics(X, S1 = cov, S2 = function(.) cov.trob(., nu=2)$cov, stdKurt = FALSE)
X.paM.<- ics(X, S1 = cov, S2 = function(.) tM(., df=2)$V, stdKurt = FALSE)
X.paR <- ics(X, S1 = cov, S2 = function(.) covMcd(.)$cov, stdKurt = FALSE)
plot(X.paM) # not at all clear
plot(X.paM.)# ditto
plot(X.paR)# very clear
}
## Similar such experiments ---> demo(rbwheel_d) and demo(rbwheel_ics)
## -------------- -----------------
``` |

Questions? Problems? Suggestions? Tweet to @rdrrHQ or email at ian@mutexlabs.com.

All documentation is copyright its authors; we didn't write any of that.