p: integer, specifying the dimension (aka number of variables).
frac: numeric, the proportion of outliers. The default, 1/p, corresponds to the (asymptotic) breakdown point of M-estimators.
sig1: thickness of the wheel , (=σ
(good[,1])), a non-negative numeric.
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 (1,0,…,0) is rotated.
scaleAfter: logical indicating if the matrix is re-scaled after
rotation (via scale()).. Default TRUE; note that this used to be false by default in the first public version.
scaleBefore: logical indicating if the matrix is re-scaled before rotation (via scale()).
spherize: logical indicating if the matrix is to be spherized , i.e., rotated and scaled to have empirical covariance Ip. This means that the principal components are used (before rotation).
fullResult: logical indicating if in addition to the nxp matrix, some intermediate quantities are returned as well.
Details
....
Returns
By default (when fullResult is FALSE), an nxp 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 nxp matrix of above, X = X0 %*% A, where A <- Qrot(p, u = U1), and X0 is the corresponding matrix before rotation, see below.
Stahel, W.~A. and , M. (2009). Comment on ``invariant co-ordinate selection'', Journal of the Royal Statistical Society B 71 , 584--586. tools:::Rd_expr_doi("10.1111/j.1467-9868.2009.00706.x")
Author(s)
Werner Stahel and Martin Maechler
Examples
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 <-100pairs(r <- rbwheel(n,6))n1 <- attr(r,"n1"); pairs(r, col=1+((1:n)> n1))## for n = 500, you *do* see it :n <-500pairs(r <- rbwheel(n,6))## show explicitlyn1 <- attr(r,"n1"); pairs(r, col=1+((1:n)> n1))## but increasing sig2 does help:pairs(r <- rbwheel(n,6, sig2 =.2))## show explicitlyn1 <- attr(r,"n1"); pairs(r, col=1+((1:n)> n1))set.seed(12)pairs(X <- rbwheel(n,7, spherize=TRUE))colSums(X)# already centeredif(require("ICS")&& require("robustbase")){# 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)## -------------- -----------------