CRAN Package Check Results for Package robustbase

Last updated on 2020-03-27 00:49:35 CET.

Flavor Version Tinstall Tcheck Ttotal Status Flags
r-devel-linux-x86_64-debian-clang 0.93-6 26.38 300.10 326.48 OK
r-devel-linux-x86_64-debian-gcc 0.93-6 20.37 222.68 243.05 OK
r-devel-linux-x86_64-fedora-clang 0.93-6 377.09 OK
r-devel-linux-x86_64-fedora-gcc 0.93-6 319.10 OK
r-devel-windows-ix86+x86_64 0.93-5 49.00 512.00 561.00 ERROR
r-devel-windows-ix86+x86_64-gcc8 0.93-5 53.00 566.00 619.00 OK
r-patched-linux-x86_64 0.93-6 20.22 273.56 293.78 OK
r-patched-solaris-x86 0.93-6 509.80 OK
r-release-linux-x86_64 0.93-6 18.99 269.82 288.81 OK
r-release-windows-ix86+x86_64 0.93-6 57.00 542.00 599.00 OK
r-release-osx-x86_64 0.93-6 OK
r-oldrel-windows-ix86+x86_64 0.93-6 50.00 513.00 563.00 OK
r-oldrel-osx-x86_64 0.93-6 OK

Check Details

Version: 0.93-5
Check: running tests for arch ‘i386’
Result: ERROR
     Running 'LTS-specials.R' [1s]
     Running 'MCD-specials.R' [1s]
     Comparing 'MCD-specials.Rout' to 'MCD-specials.Rout.save' ... OK
     Running 'MT-tst.R' [9s]
     Running 'NAcoef.R' [3s]
     Comparing 'NAcoef.Rout' to 'NAcoef.Rout.save' ... OK
     Running 'OGK-ex.R' [1s]
     Comparing 'OGK-ex.Rout' to 'OGK-ex.Rout.save' ... OK
     Running 'Rsquared.R' [1s]
     Comparing 'Rsquared.Rout' to 'Rsquared.Rout.save' ... OK
     Running 'binom-ni-small.R' [0s]
     Comparing 'binom-ni-small.Rout' to 'binom-ni-small.Rout.save' ... OK
     Running 'binom-no-x.R' [0s]
     Running 'comedian-tst.R' [2s]
     Running 'exact-fit-categorical.R' [0s]
     Running 'glmrob-1.R' [14s]
     Running 'glmrob-specials.R' [0s]
     Running 'huber-etc.R' [0s]
     Comparing 'huber-etc.Rout' to 'huber-etc.Rout.save' ... OK
     Running 'lmrob-data.R' [19s]
     Running 'lmrob-ex12.R' [5s]
     Running 'lmrob-methods.R' [1s]
     Comparing 'lmrob-methods.Rout' to 'lmrob-methods.Rout.save' ... OK
     Running 'lmrob-outlierStats.R' [2s]
     Running 'lmrob-psifns.R' [5s]
     Comparing 'lmrob-psifns.Rout' to 'lmrob-psifns.Rout.save' ... OK
     Running 'm-s-estimator.R' [3s]
     Running 'mc-etc.R' [1s]
     Running 'mc-strict.R' [7s]
     Running 'nlregrob-tst.R' [26s]
     Running 'nlrob-tst.R' [3s]
     Running 'plot-ex.R' [1s]
     Running 'poisson-ex.R' [2s]
     Running 'psi-rho-etc.R' [2s]
     Comparing 'psi-rho-etc.Rout' to 'psi-rho-etc.Rout.save' ... OK
     Running 'small-sample.R' [7s]
     Comparing 'small-sample.Rout' to 'small-sample.Rout.save' ... OK
     Running 'subsample.R' [8s]
     Running 'tlts.R' [1s]
     Comparing 'tlts.Rout' to 'tlts.Rout.save' ... OK
     Running 'tmcd.R' [9s]
     Running 'weights.R' [1s]
     Comparing 'weights.Rout' to 'weights.Rout.save' ...442,443c442,443
    < Warning message:
    < In lf.cov(object, complete = complete, ...) :
    ---
    > Warning messages:
    > 1: In lf.cov(object, complete = complete, ...) :
    444a445,446
    > 2: In lf.cov(object, complete = complete, ...) :
    > .vcov.avar1: negative diag(<vcov>) fixed up; consider 'cov=".vcov.w."' instead
     Running 'wgt-himed-xtra.R' [3s]
     Running 'wgt-himed.R' [1s]
     Comparing 'wgt-himed.Rout' to 'wgt-himed.Rout.save' ... OK
    Running the tests in 'tests/m-s-estimator.R' failed.
    Complete output:
     > ## Test implementation of M-S estimator
     > require(robustbase)
     Loading required package: robustbase
     > source(system.file("xtraR/m-s_fns.R", package = "robustbase", mustWork=TRUE))
     > source(system.file("xtraR/ex-funs.R", package = "robustbase", mustWork=TRUE))
     > source(system.file("test-tools-1.R", package = "Matrix", mustWork=TRUE))# assert.EQ
     Loading required package: tools
     >
     > ## dataset with factors and continuous variables:
     > data(education)
     > education <- within(education, Region <- factor(Region))
     > ## for testing purposes:
     > education2 <- within(education, Group <- factor(rep(1:3, length.out=length(Region))))
     >
     > ## Test splitFrame (type fii is the only problematic type)
     > testFun <- function(formula, x1.idx) {
     + obj <- lm(formula, education2)
     + mf <- obj$model
     + ret <- splitFrame(mf, type="fii")
     + if (missing(x1.idx)) {
     + print(ret$x1.idx)
     + return(which(unname(ret$x1.idx)))
     + }
     + stopifnot(identical(x1.idx, which(unname(ret$x1.idx))))
     + }
     > testFun(Y ~ 1, integer(0))
     > testFun(Y ~ X1*X2*X3, integer(0))
     > testFun(Y ~ Region + X1 + X2 + X3, 1:4)
     > testFun(Y ~ 0 + Region + X1 + X2 + X3, 1:4)
     > testFun(Y ~ Region*X1 + X2 + X3, c(1:5, 8:10))
     > testFun(Y ~ Region*X1 + X2 + X3 + Region*Group, c(1:5, 8:18))
     > testFun(Y ~ Region*X1 + X2 + X3 + Region*Group*X2, c(1:6, 8:29))
     > testFun(Y ~ Region*X1 + X2 + Region*Group*X2, 1:28)
     > testFun(Y ~ Region*X1 + X2 + Region:Group:X2, 1:21)
     > testFun(Y ~ Region*X1 + X2*X3 + Region:Group:X2, c(1:6, 8:10, 12:23))
     > testFun(Y ~ (X1+X2+X3+Region)^2, c(1:7,10:12,14:19))
     > testFun(Y ~ (X1+X2+X3+Region)^3, c(1:19, 21:29))
     > testFun(Y ~ (X1+X2+X3+Region)^4, 1:32)
     > testFun(Y ~ Region:X1:X2 + X1*X2, c(1:1, 4:7))
     >
     >
     > control <- lmrob.control()
     > f.lm <- lm(Y ~ Region + X1 + X2 + X3, education)
     > splt <- splitFrame(f.lm$model)
     > y <- education$Y
     >
     > ## test orthogonalizing
     > x1 <- splt$x1
     > x2 <- splt$x2
     > tmp <- lmrob.lar(x1, y, control)
     > y.tilde <- tmp$resid
     > t1 <- tmp$coef
     > x2.tilde <- x2
     > T2 <- matrix(0, nrow=ncol(x1), ncol=ncol(x2))
     > for (i in 1:ncol(x2)) {
     + tmp <- lmrob.lar(x1, x2[,i], control)
     + x2.tilde[,i] <- tmp$resid
     + T2[,i] <- tmp$coef
     + }
     >
     > set.seed(10)
     > mss1 <- m_s_subsample(x1, x2.tilde, y.tilde, control, orth = FALSE)
     > mss1 <- within(mss1, b1 <- drop(t1 + b1 - T2 %*% b2))
     > set.seed(10)
     > mss2 <- m_s_subsample(x1, x2, y, control, orth = TRUE)
     > stopifnot(all.equal(mss1, mss2))
     >
     > res <- vector("list", 100)
     > set.seed(0)
     > time <- system.time(for (i in seq_along(res)) {
     + tmp <- m_s_subsample(x1, x2.tilde, y.tilde, control, FALSE)
     + res[[i]] <- unlist(within(tmp, b1 <- drop(t1 + b1 - T2 %*% b2)))
     + })
     > cat('Time elapsed in subsampling: ', time,'\n')
     Time elapsed in subsampling: 0.62 0 0.66 NA NA
     > ## show a summary of the results {"FIXME": output is platform dependent}
     > summary(res1 <- do.call(rbind, res))
     b11 b12 b13 b14
     Min. :-316.24 Min. :-33.92 Min. :-35.8704 Min. :16.43
     1st Qu.:-223.80 1st Qu.:-23.55 1st Qu.: -8.9558 1st Qu.:30.17
     Median :-163.19 Median :-21.67 Median : -6.9137 Median :32.52
     Mean :-161.83 Mean :-22.06 Mean : -7.4830 Mean :32.40
     3rd Qu.:-103.36 3rd Qu.:-18.47 3rd Qu.: -4.8541 3rd Qu.:36.16
     Max. : 61.83 Max. :-12.03 Max. : 0.7015 Max. :42.23
     b21 b22 b23 scale
     Min. :-0.03808 Min. :0.02111 Min. :0.2555 Min. :29.79
     1st Qu.:-0.01309 1st Qu.:0.03960 1st Qu.:0.4956 1st Qu.:30.42
     Median : 0.02151 Median :0.04717 Median :0.6381 Median :31.08
     Mean : 0.02569 Mean :0.04635 Mean :0.6246 Mean :30.99
     3rd Qu.: 0.05296 3rd Qu.:0.05211 3rd Qu.:0.7506 3rd Qu.:31.35
     Max. : 0.10074 Max. :0.06938 Max. :0.9172 Max. :32.18
     > ## compare with fast S solution
     > fmS <- lmrob(Y ~ Region + X1 + X2 + X3, education, init="S")
     > coef(fmS)
     (Intercept) Region2 Region3 Region4 X1
     -135.72592554 -20.64576283 -9.84881727 24.58013066 0.03405591
     X2 X3
     0.04327562 0.57895741
     > fmS$scale
     [1] 26.40386
     >
     > ### Comparing m-s_descent implementations() {our C and R} : -------------------
     >
     > ctrl <- control
     > #ctrl$trace.lev <- 5
     > ctrl$k.max <- 1
     > mC <- m_s_descent (x1, x2, y, ctrl, mss2$b1, mss2$b2, mss2$scale+10)
     > mR <- m_s_descent_Ronly(x1, x2, y, ctrl, mss2$b1, mss2$b2, mss2$scale+10)
     > nm <- c("b1","b2", "scale", "res")
     > stopifnot(all.equal(mC[nm], mR[nm], check.attributes = FALSE, tolerance = 4e-14))
     > # seen 5.567e-15 in OpenBLAS ^^^^^
     >
     > ## control$k.m_s <- 100
     > res3 <- vector("list", 100)
     > time <- system.time(for (i in seq_along(res3)) {
     + ri <- res[[i]]
     + res3[[i]] <- unlist(m_s_descent(x1, x2, y, control,
     + ri[1:4], ri[5:7], ri[8]))
     + })
     > cat('Time elapsed in descent proc: ', time,'\n')
     Time elapsed in descent proc: 0.12 0 0.13 NA NA
     >
     > ## show a summary of the results {"FIXME": output is platform dependent}
     > res4 <- do.call(rbind, res3)
     > summary(res4[,1:8])
     b11 b12 b13 b14
     Min. :-316.3 Min. :-30.56 Min. :-36.501 Min. :16.43
     1st Qu.:-222.9 1st Qu.:-23.02 1st Qu.: -8.795 1st Qu.:27.55
     Median :-158.4 Median :-20.45 Median : -7.790 Median :29.92
     Mean :-157.4 Mean :-20.18 Mean : -8.850 Mean :30.44
     3rd Qu.:-102.3 3rd Qu.:-15.99 3rd Qu.: -6.848 3rd Qu.:32.59
     Max. : 101.7 Max. :-12.24 Max. : -4.032 Max. :42.23
     b21 b22 b23 scale
     Min. :-0.02141 Min. :0.01459 Min. :0.2034 Min. :29.79
     1st Qu.: 0.02142 1st Qu.:0.03836 1st Qu.:0.4956 1st Qu.:30.35
     Median : 0.04687 Median :0.04154 Median :0.6352 Median :30.57
     Mean : 0.04320 Mean :0.04315 Mean :0.6228 Mean :30.72
     3rd Qu.: 0.06124 3rd Qu.:0.04776 3rd Qu.:0.7355 3rd Qu.:31.05
     Max. : 0.09102 Max. :0.06367 Max. :0.9172 Max. :31.84
     >
     > stopifnot(all.equal( # 'test', not only plot:
     + res1[, "scale"], res4[,"scale"], tol = 0.03),
     + res1[, "scale"] >= res4[,"scale"] - 1e-7 ) # 1e-7 just in case
     > plot(res1[, "scale"], res4[,"scale"])
     > abline(0,1, col=adjustcolor("gray", 0.5))
     >
     > ## Test lmrob.M.S
     > x <- model.matrix(fmS)
     > control$trace.lev <- 3
     > ## --------- --
     > set.seed(1003)
     > fMS <- lmrob.M.S(x, y, control, fmS$model)
     lmrob_M_S(n = 50, nRes = 500, (p1,p2)=(4,3), (orth,subs,desc)=(1,1,1))
     Starting subsampling procedure.. [setup Ok]
     Sample[ 0]: new candidate with sc = 43.227 in 46 iter
     Sample[ 1]: new candidate with sc = 33.125 in 49 iter
     Sample[ 18]: new candidate with sc = 31.856 in 42 iter
     Sample[136]: new candidate with sc = 31.741 in 35 iter
     Sample[338]: new candidate with sc = 31.346 in 39 iter
     Finished M-S subsampling with scale = 31.34552
     b1: -5.445072 7.825257 0.178215 7.656800
     b2: 0.048442 0.037633 0.591759
     Starting descent procedure...
     Scale: 31.34552
     Refinement step 1: no improvement, scale: 31.528
     Refinement step 2: no improvement, scale: 31.520
     Refinement step 3: no improvement, scale: 31.525
     Refinement step 4: no improvement, scale: 31.456
     Refinement step 5: no improvement, scale: 31.373
     Refinement step 6: better fit, scale: 31.343
     Refinement step 7: better fit, scale: 31.337
     Refinement step 8: no improvement, scale: 31.345
     Refinement step 9: no improvement, scale: 31.359
     Refinement step 10: no improvement, scale: 31.377
     Refinement step 11: no improvement, scale: 31.397
     Refinement step 12: no improvement, scale: 31.420
     Refinement step 13: no improvement, scale: 31.443
     Refinement step 14: no improvement, scale: 31.468
     Refinement step 15: no improvement, scale: 31.493
     Refinement step 16: no improvement, scale: 31.520
     Refinement step 17: no improvement, scale: 31.547
     Refinement step 18: no improvement, scale: 31.575
     Refinement step 19: no improvement, scale: 31.605
     Refinement step 20: no improvement, scale: 31.635
     Refinement step 21: no improvement, scale: 31.666
     Refinement step 22: no improvement, scale: 31.698
     Refinement step 23: no improvement, scale: 31.732
     Refinement step 24: no improvement, scale: 31.766
     Refinement step 25: no improvement, scale: 31.802
     Refinement step 26: no improvement, scale: 31.839
     Refinement step 27: no improvement, scale: 31.877
     Descent procedure: not converged (best scale: 31.337, last step: 31.877)
     The procedure stopped after 28 steps because there was no improvement in the last 20 steps.
     To increase this number, use the control parameter 'k.m_s'.
     b1: -111.930390 -14.841307 -9.456777 27.431135
     b2: 0.067235 0.036026 0.536591
     > resid <- drop(y - x %*% fMS$coef)
     > assert.EQ(resid, fMS$resid, check.attributes=FALSE, tol = 1e-12)
     >
     > ## Test direct call to lmrob
     > ## 1. trace_lev output:
     > set.seed(17)
     > fMS <- lmrob(Y ~ Region + X1 + X2 + X3, education, init = "M-S", trace.lev=2)
     lmrob_M_S(n = 50, nRes = 500, (p1,p2)=(4,3), (orth,subs,desc)=(1,1,1))
     Starting subsampling procedure.. [setup Ok]
     Sample[ 0]: new candidate with sc = 720.84 in 42 iter
     Sample[ 1]: new candidate with sc = 32.106 in 48 iter
     Sample[ 69]: new candidate with sc = 31.856 in 40 iter
     Sample[128]: new candidate with sc = 31.122 in 47 iter
     Finished M-S subsampling with scale = 31.12236
     Starting descent procedure...
     Refinement step 1: better fit, scale: 30.963
     Refinement step 4: better fit, scale: 30.888
     Descent procedure: not converged (best scale: 30.888, last step: 31.134)
     The procedure stopped after 25 steps because there was no improvement in the last 20 steps.
     To increase this number, use the control parameter 'k.m_s'.
     init converged (remaining method = "M") -> coef=
     [1] -236.33813174 -22.50174637 -7.69929975 25.25388208 0.05220654
     [6] 0.04769691 0.78984802
     lmrob_MM(): rwls():
     rwls() used 19 it.; last ||b0 - b1||_1 = 1.76872e-05, L(b1) = 0.146343135734; convergence
     lmrob..MM..fit(*, obj) --> updating .. qr(x * rweights) -> rank= 7
     step "M" -> new coef=
     (Intercept) Region2 Region3 Region4 X1
     -150.21716170 -12.68711627 -10.64944153 21.92534635 0.04153037
     X2 X3
     0.04337025 0.61139401
     >
     > set.seed(13)
     > fiMS <- lmrob(Y ~ Region + X1 + X2 + X3, education, init = "M-S")
     > out2 <- capture.output(summary(fiMS))
     > writeLines(out2)
    
     Call:
     lmrob(formula = Y ~ Region + X1 + X2 + X3, data = education, init = "M-S")
     \--> method = "M-SM"
     Residuals:
     Min 1Q Median 3Q Max
     -62.729 -15.529 -1.572 23.392 174.750
    
     Coefficients:
     Estimate Std. Error t value Pr(>|t|)
     (Intercept) -150.07630 143.09300 -1.049 0.30013
     Region2 -12.76767 16.63758 -0.767 0.44704
     Region3 -10.63954 15.92865 -0.668 0.50774
     Region4 21.95445 16.96484 1.294 0.20253
     X1 0.04146 0.05040 0.823 0.41525
     X2 0.04337 0.01373 3.159 0.00289 **
     X3 0.61106 0.35153 1.738 0.08932 .
     ---
     Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
    
     Robust residual standard error: 30.82
     Multiple R-squared: 0.5695, Adjusted R-squared: 0.5095
     Convergence in 19 IRWLS iterations
    
     Robustness weights:
     observation 50 is an outlier with |weight| = 0 ( < 0.002);
     7 weights are ~= 1. The remaining 42 ones are summarized as
     Min. 1st Qu. Median Mean 3rd Qu. Max.
     0.2884 0.8904 0.9508 0.8890 0.9867 0.9985
     Algorithmic parameters:
     tuning.chi bb tuning.psi rel.tol
     1.548e+00 5.000e-01 4.685e+00 1.000e-07
     scale.tol solve.tol eps.outlier eps.x
     1.000e-10 1.000e-07 2.000e-03 1.071e-08
     warn.limit.reject warn.limit.meanrw
     5.000e-01 5.000e-01
     nResample max.it k.max maxit.scale k.m_s
     500 50 200 200 20
     trace.lev mts compute.rd fast.s.large.n
     0 1000 0 2000
     psi subsampling cov
     "bisquare" "nonsingular" ".vcov.w"
     split.type compute.outlier.stats
     "f" "SM"
     seed : int(0)
     >
     > set.seed(13)
     > fiM.S <- lmrob(Y ~ Region + X1 + X2 + X3, education, init=lmrob.M.S)
     > out3 <- capture.output(summary(fiM.S))
     >
     > ## must be the same {apart from the "init=" in the call}:
     > i <- 3
     > stopifnot(identical(out2[-i], out3[-i]))
     > ## the difference:
     > c(rbind(out2[i], out3[i]))
     [1] "lmrob(formula = Y ~ Region + X1 + X2 + X3, data = education, init = \"M-S\")"
     [2] "lmrob(formula = Y ~ Region + X1 + X2 + X3, data = education, init = lmrob.M.S)"
     >
     >
     > ### "Skipping design matrix equilibration" warning can arise for reasonable designs -----
     > set.seed(1)
     > x2 <- matrix(rnorm(2*30), 30, 2)
     > data <- data.frame(y = rnorm(30), group = rep(letters[1:3], each=10), x2)
     >
     > obj <- lmrob(y ~ ., data, init="M-S", trace.lev=1)
     lmrob_S(n = 30, nRes = 500): fast_s() [non-large n]:
     Subsampling to find candidate betas:
     Now refine() to convergence for 2 very best ones:
     Best[0]: convergence (87 iter.): -> improved scale to 0.921076968958332
     Best[1]: convergence (41 iter.)
     lmrob.S(): scale = 0.921077; coeff.=
     [1] 0.055853342 -0.395949198 0.068785230 -0.002081381 0.239242595
     init converged (remaining method = "M") -> coef=
     (Intercept) groupb groupc X1 X2
     0.055853342 -0.395949198 0.068785230 -0.002081381 0.239242595
     lmrob_MM(): rwls():
     rwls() used 12 it.; last ||b0 - b1||_1 = 6.71467e-08, L(b1) = 0.0919204728697; convergence
     lmrob..MM..fit(*, obj) --> updating .. qr(x * rweights) -> rank= 5, outlierStats()
     step "M" -> new coef=
     (Intercept) groupb groupc X1 X2
     0.41184656 -0.74163256 -0.36361231 0.02399208 0.51678384
     Warning message:
     In lmrob.M.S(x, y, control, mf = mf) :
     No categorical variables found in model. Reverting to S-estimator.
     >
     > ## illustration: the zero row is introduced during the orthogonalization of x2 wrt x1
     > ## l1 regression always produces p zero residuals
     > ## by chance, the zero residuals of multiple columns happen to be on the same row
     > sf <- splitFrame(obj$model)
     > x1 <- sf$x1
     > x2 <- sf$x2
     > control <- obj$control
     >
     > ## orthogonalize
     > x2.tilde <- x2
     >
     > for(i in 1:ncol(x2)) {
     + tmp <- lmrob.lar(x1, x2[,i], control)
     + x2.tilde[,i] <- tmp$resid
     + }
     Error in lmrob.lar(x1, x2[, i], control) : p > 0 is not TRUE
     Calls: lmrob.lar -> stopifnot
     Execution halted
Flavor: r-devel-windows-ix86+x86_64

Version: 0.93-5
Check: running tests for arch ‘x64’
Result: ERROR
     Running 'LTS-specials.R' [1s]
     Running 'MCD-specials.R' [1s]
     Comparing 'MCD-specials.Rout' to 'MCD-specials.Rout.save' ... OK
     Running 'MT-tst.R' [9s]
     Running 'NAcoef.R' [3s]
     Comparing 'NAcoef.Rout' to 'NAcoef.Rout.save' ... OK
     Running 'OGK-ex.R' [1s]
     Comparing 'OGK-ex.Rout' to 'OGK-ex.Rout.save' ... OK
     Running 'Rsquared.R' [1s]
     Comparing 'Rsquared.Rout' to 'Rsquared.Rout.save' ... OK
     Running 'binom-ni-small.R' [1s]
     Comparing 'binom-ni-small.Rout' to 'binom-ni-small.Rout.save' ... OK
     Running 'binom-no-x.R' [0s]
     Running 'comedian-tst.R' [2s]
     Running 'exact-fit-categorical.R' [1s]
     Running 'glmrob-1.R' [14s]
     Running 'glmrob-specials.R' [0s]
     Running 'huber-etc.R' [0s]
     Comparing 'huber-etc.Rout' to 'huber-etc.Rout.save' ... OK
     Running 'lmrob-data.R' [18s]
     Running 'lmrob-ex12.R' [4s]
     Running 'lmrob-methods.R' [1s]
     Comparing 'lmrob-methods.Rout' to 'lmrob-methods.Rout.save' ... OK
     Running 'lmrob-outlierStats.R' [2s]
     Running 'lmrob-psifns.R' [5s]
     Comparing 'lmrob-psifns.Rout' to 'lmrob-psifns.Rout.save' ... OK
     Running 'm-s-estimator.R' [3s]
     Running 'mc-etc.R' [1s]
     Running 'mc-strict.R' [10s]
     Running 'nlregrob-tst.R' [28s]
     Running 'nlrob-tst.R' [3s]
     Running 'plot-ex.R' [1s]
     Running 'poisson-ex.R' [3s]
     Running 'psi-rho-etc.R' [1s]
     Comparing 'psi-rho-etc.Rout' to 'psi-rho-etc.Rout.save' ... OK
     Running 'small-sample.R' [9s]
     Comparing 'small-sample.Rout' to 'small-sample.Rout.save' ... OK
     Running 'subsample.R' [8s]
     Running 'tlts.R' [1s]
     Comparing 'tlts.Rout' to 'tlts.Rout.save' ... OK
     Running 'tmcd.R' [9s]
     Running 'weights.R' [1s]
     Comparing 'weights.Rout' to 'weights.Rout.save' ... OK
     Running 'wgt-himed-xtra.R' [4s]
     Running 'wgt-himed.R' [0s]
     Comparing 'wgt-himed.Rout' to 'wgt-himed.Rout.save' ... OK
    Running the tests in 'tests/m-s-estimator.R' failed.
    Complete output:
     > ## Test implementation of M-S estimator
     > require(robustbase)
     Loading required package: robustbase
     > source(system.file("xtraR/m-s_fns.R", package = "robustbase", mustWork=TRUE))
     > source(system.file("xtraR/ex-funs.R", package = "robustbase", mustWork=TRUE))
     > source(system.file("test-tools-1.R", package = "Matrix", mustWork=TRUE))# assert.EQ
     Loading required package: tools
     >
     > ## dataset with factors and continuous variables:
     > data(education)
     > education <- within(education, Region <- factor(Region))
     > ## for testing purposes:
     > education2 <- within(education, Group <- factor(rep(1:3, length.out=length(Region))))
     >
     > ## Test splitFrame (type fii is the only problematic type)
     > testFun <- function(formula, x1.idx) {
     + obj <- lm(formula, education2)
     + mf <- obj$model
     + ret <- splitFrame(mf, type="fii")
     + if (missing(x1.idx)) {
     + print(ret$x1.idx)
     + return(which(unname(ret$x1.idx)))
     + }
     + stopifnot(identical(x1.idx, which(unname(ret$x1.idx))))
     + }
     > testFun(Y ~ 1, integer(0))
     > testFun(Y ~ X1*X2*X3, integer(0))
     > testFun(Y ~ Region + X1 + X2 + X3, 1:4)
     > testFun(Y ~ 0 + Region + X1 + X2 + X3, 1:4)
     > testFun(Y ~ Region*X1 + X2 + X3, c(1:5, 8:10))
     > testFun(Y ~ Region*X1 + X2 + X3 + Region*Group, c(1:5, 8:18))
     > testFun(Y ~ Region*X1 + X2 + X3 + Region*Group*X2, c(1:6, 8:29))
     > testFun(Y ~ Region*X1 + X2 + Region*Group*X2, 1:28)
     > testFun(Y ~ Region*X1 + X2 + Region:Group:X2, 1:21)
     > testFun(Y ~ Region*X1 + X2*X3 + Region:Group:X2, c(1:6, 8:10, 12:23))
     > testFun(Y ~ (X1+X2+X3+Region)^2, c(1:7,10:12,14:19))
     > testFun(Y ~ (X1+X2+X3+Region)^3, c(1:19, 21:29))
     > testFun(Y ~ (X1+X2+X3+Region)^4, 1:32)
     > testFun(Y ~ Region:X1:X2 + X1*X2, c(1:1, 4:7))
     >
     >
     > control <- lmrob.control()
     > f.lm <- lm(Y ~ Region + X1 + X2 + X3, education)
     > splt <- splitFrame(f.lm$model)
     > y <- education$Y
     >
     > ## test orthogonalizing
     > x1 <- splt$x1
     > x2 <- splt$x2
     > tmp <- lmrob.lar(x1, y, control)
     > y.tilde <- tmp$resid
     > t1 <- tmp$coef
     > x2.tilde <- x2
     > T2 <- matrix(0, nrow=ncol(x1), ncol=ncol(x2))
     > for (i in 1:ncol(x2)) {
     + tmp <- lmrob.lar(x1, x2[,i], control)
     + x2.tilde[,i] <- tmp$resid
     + T2[,i] <- tmp$coef
     + }
     >
     > set.seed(10)
     > mss1 <- m_s_subsample(x1, x2.tilde, y.tilde, control, orth = FALSE)
     > mss1 <- within(mss1, b1 <- drop(t1 + b1 - T2 %*% b2))
     > set.seed(10)
     > mss2 <- m_s_subsample(x1, x2, y, control, orth = TRUE)
     > stopifnot(all.equal(mss1, mss2))
     >
     > res <- vector("list", 100)
     > set.seed(0)
     > time <- system.time(for (i in seq_along(res)) {
     + tmp <- m_s_subsample(x1, x2.tilde, y.tilde, control, FALSE)
     + res[[i]] <- unlist(within(tmp, b1 <- drop(t1 + b1 - T2 %*% b2)))
     + })
     > cat('Time elapsed in subsampling: ', time,'\n')
     Time elapsed in subsampling: 0.45 0 0.46 NA NA
     > ## show a summary of the results {"FIXME": output is platform dependent}
     > summary(res1 <- do.call(rbind, res))
     b11 b12 b13 b14
     Min. :-316.24 Min. :-33.92 Min. :-35.8704 Min. :16.43
     1st Qu.:-223.37 1st Qu.:-23.66 1st Qu.: -8.9603 1st Qu.:29.92
     Median :-163.19 Median :-21.37 Median : -7.2929 Median :32.20
     Mean :-161.11 Mean :-22.02 Mean : -8.0727 Mean :32.15
     3rd Qu.:-103.36 3rd Qu.:-18.42 3rd Qu.: -5.9055 3rd Qu.:35.72
     Max. : 61.83 Max. :-12.03 Max. : 0.7015 Max. :42.23
     b21 b22 b23 scale
     Min. :-0.03808 Min. :0.02111 Min. :0.2555 Min. :29.79
     1st Qu.:-0.00397 1st Qu.:0.03927 1st Qu.:0.4956 1st Qu.:30.43
     Median : 0.02160 Median :0.04716 Median :0.6381 Median :30.91
     Mean : 0.02734 Mean :0.04618 Mean :0.6219 Mean :30.95
     3rd Qu.: 0.05561 3rd Qu.:0.05224 3rd Qu.:0.7473 3rd Qu.:31.35
     Max. : 0.10427 Max. :0.06938 Max. :0.9172 Max. :32.18
     > ## compare with fast S solution
     > fmS <- lmrob(Y ~ Region + X1 + X2 + X3, education, init="S")
     > coef(fmS)
     (Intercept) Region2 Region3 Region4 X1
     -135.72592554 -20.64576283 -9.84881727 24.58013066 0.03405591
     X2 X3
     0.04327562 0.57895741
     > fmS$scale
     [1] 26.40386
     >
     > ### Comparing m-s_descent implementations() {our C and R} : -------------------
     >
     > ctrl <- control
     > #ctrl$trace.lev <- 5
     > ctrl$k.max <- 1
     > mC <- m_s_descent (x1, x2, y, ctrl, mss2$b1, mss2$b2, mss2$scale+10)
     > mR <- m_s_descent_Ronly(x1, x2, y, ctrl, mss2$b1, mss2$b2, mss2$scale+10)
     > nm <- c("b1","b2", "scale", "res")
     > stopifnot(all.equal(mC[nm], mR[nm], check.attributes = FALSE, tolerance = 4e-14))
     > # seen 5.567e-15 in OpenBLAS ^^^^^
     >
     > ## control$k.m_s <- 100
     > res3 <- vector("list", 100)
     > time <- system.time(for (i in seq_along(res3)) {
     + ri <- res[[i]]
     + res3[[i]] <- unlist(m_s_descent(x1, x2, y, control,
     + ri[1:4], ri[5:7], ri[8]))
     + })
     > cat('Time elapsed in descent proc: ', time,'\n')
     Time elapsed in descent proc: 0.11 0 0.11 NA NA
     >
     > ## show a summary of the results {"FIXME": output is platform dependent}
     > res4 <- do.call(rbind, res3)
     > summary(res4[,1:8])
     b11 b12 b13 b14
     Min. :-316.3 Min. :-30.56 Min. :-36.501 Min. :16.43
     1st Qu.:-222.4 1st Qu.:-23.09 1st Qu.: -8.960 1st Qu.:27.96
     Median :-160.7 Median :-21.00 Median : -7.868 Median :30.46
     Mean :-158.2 Mean :-20.60 Mean : -9.046 Mean :30.84
     3rd Qu.:-102.7 3rd Qu.:-17.40 3rd Qu.: -6.842 3rd Qu.:32.75
     Max. : 101.7 Max. :-12.24 Max. : -4.032 Max. :42.23
     b21 b22 b23 scale
     Min. :-0.02141 Min. :0.01459 Min. :0.2034 Min. :29.79
     1st Qu.: 0.02048 1st Qu.:0.03873 1st Qu.:0.5007 1st Qu.:30.37
     Median : 0.04169 Median :0.04271 Median :0.6381 Median :30.57
     Mean : 0.03911 Mean :0.04359 Mean :0.6270 Mean :30.70
     3rd Qu.: 0.06102 3rd Qu.:0.04798 3rd Qu.:0.7460 3rd Qu.:30.96
     Max. : 0.09102 Max. :0.06367 Max. :0.9172 Max. :31.84
     >
     > stopifnot(all.equal( # 'test', not only plot:
     + res1[, "scale"], res4[,"scale"], tol = 0.03),
     + res1[, "scale"] >= res4[,"scale"] - 1e-7 ) # 1e-7 just in case
     > plot(res1[, "scale"], res4[,"scale"])
     > abline(0,1, col=adjustcolor("gray", 0.5))
     >
     > ## Test lmrob.M.S
     > x <- model.matrix(fmS)
     > control$trace.lev <- 3
     > ## --------- --
     > set.seed(1003)
     > fMS <- lmrob.M.S(x, y, control, fmS$model)
     lmrob_M_S(n = 50, nRes = 500, (p1,p2)=(4,3), (orth,subs,desc)=(1,1,1))
     Starting subsampling procedure.. [setup Ok]
     Sample[ 0]: new candidate with sc = 43.227 in 46 iter
     Sample[ 1]: new candidate with sc = 33.125 in 49 iter
     Sample[ 18]: new candidate with sc = 31.856 in 42 iter
     Sample[136]: new candidate with sc = 31.741 in 35 iter
     Sample[338]: new candidate with sc = 31.512 in 37 iter
     Finished M-S subsampling with scale = 31.51188
     b1: -5.445072 7.825257 1.756708 7.656800
     b2: 0.048442 0.037633 0.591759
     Starting descent procedure...
     Scale: 31.51188
     Refinement step 1: better fit, scale: 31.486
     Refinement step 2: no improvement, scale: 31.496
     Refinement step 3: no improvement, scale: 31.498
     Refinement step 4: better fit, scale: 31.388
     Refinement step 5: better fit, scale: 31.343
     Refinement step 6: better fit, scale: 31.329
     Refinement step 7: no improvement, scale: 31.331
     Refinement step 8: no improvement, scale: 31.342
     Refinement step 9: no improvement, scale: 31.358
     Refinement step 10: no improvement, scale: 31.377
     Refinement step 11: no improvement, scale: 31.398
     Refinement step 12: no improvement, scale: 31.421
     Refinement step 13: no improvement, scale: 31.445
     Refinement step 14: no improvement, scale: 31.469
     Refinement step 15: no improvement, scale: 31.495
     Refinement step 16: no improvement, scale: 31.521
     Refinement step 17: no improvement, scale: 31.549
     Refinement step 18: no improvement, scale: 31.577
     Refinement step 19: no improvement, scale: 31.607
     Refinement step 20: no improvement, scale: 31.637
     Refinement step 21: no improvement, scale: 31.668
     Refinement step 22: no improvement, scale: 31.700
     Refinement step 23: no improvement, scale: 31.734
     Refinement step 24: no improvement, scale: 31.768
     Refinement step 25: no improvement, scale: 31.804
     Refinement step 26: no improvement, scale: 31.841
     Descent procedure: not converged (best scale: 31.329, last step: 31.841)
     The procedure stopped after 27 steps because there was no improvement in the last 20 steps.
     To increase this number, use the control parameter 'k.m_s'.
     b1: -113.528805 -14.858187 -9.426368 27.485299
     b2: 0.066846 0.036134 0.540701
     > resid <- drop(y - x %*% fMS$coef)
     > assert.EQ(resid, fMS$resid, check.attributes=FALSE, tol = 1e-12)
     >
     > ## Test direct call to lmrob
     > ## 1. trace_lev output:
     > set.seed(17)
     > fMS <- lmrob(Y ~ Region + X1 + X2 + X3, education, init = "M-S", trace.lev=2)
     lmrob_M_S(n = 50, nRes = 500, (p1,p2)=(4,3), (orth,subs,desc)=(1,1,1))
     Starting subsampling procedure.. [setup Ok]
     Sample[ 0]: new candidate with sc = 720.84 in 42 iter
     Sample[ 1]: new candidate with sc = 32.106 in 48 iter
     Sample[ 69]: new candidate with sc = 31.856 in 40 iter
     Sample[128]: new candidate with sc = 31.122 in 47 iter
     Finished M-S subsampling with scale = 31.12236
     Starting descent procedure...
     Refinement step 1: better fit, scale: 30.963
     Refinement step 4: better fit, scale: 30.888
     Descent procedure: not converged (best scale: 30.888, last step: 31.134)
     The procedure stopped after 25 steps because there was no improvement in the last 20 steps.
     To increase this number, use the control parameter 'k.m_s'.
     init converged (remaining method = "M") -> coef=
     [1] -236.33813174 -22.50174637 -7.69929975 25.25388208 0.05220654
     [6] 0.04769691 0.78984802
     lmrob_MM(): rwls():
     rwls() used 19 it.; last ||b0 - b1||_1 = 1.76872e-05, L(b1) = 0.146343135734; convergence
     lmrob..MM..fit(*, obj) --> updating .. qr(x * rweights) -> rank= 7
     step "M" -> new coef=
     (Intercept) Region2 Region3 Region4 X1
     -150.21716170 -12.68711627 -10.64944153 21.92534635 0.04153037
     X2 X3
     0.04337025 0.61139401
     >
     > set.seed(13)
     > fiMS <- lmrob(Y ~ Region + X1 + X2 + X3, education, init = "M-S")
     > out2 <- capture.output(summary(fiMS))
     > writeLines(out2)
    
     Call:
     lmrob(formula = Y ~ Region + X1 + X2 + X3, data = education, init = "M-S")
     \--> method = "M-SM"
     Residuals:
     Min 1Q Median 3Q Max
     -62.729 -15.529 -1.572 23.392 174.750
    
     Coefficients:
     Estimate Std. Error t value Pr(>|t|)
     (Intercept) -150.07630 143.09300 -1.049 0.30013
     Region2 -12.76767 16.63758 -0.767 0.44704
     Region3 -10.63954 15.92865 -0.668 0.50774
     Region4 21.95445 16.96484 1.294 0.20253
     X1 0.04146 0.05040 0.823 0.41525
     X2 0.04337 0.01373 3.159 0.00289 **
     X3 0.61106 0.35153 1.738 0.08932 .
     ---
     Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
    
     Robust residual standard error: 30.82
     Multiple R-squared: 0.5695, Adjusted R-squared: 0.5095
     Convergence in 19 IRWLS iterations
    
     Robustness weights:
     observation 50 is an outlier with |weight| = 0 ( < 0.002);
     7 weights are ~= 1. The remaining 42 ones are summarized as
     Min. 1st Qu. Median Mean 3rd Qu. Max.
     0.2884 0.8904 0.9508 0.8890 0.9867 0.9985
     Algorithmic parameters:
     tuning.chi bb tuning.psi rel.tol
     1.548e+00 5.000e-01 4.685e+00 1.000e-07
     scale.tol solve.tol eps.outlier eps.x
     1.000e-10 1.000e-07 2.000e-03 1.071e-08
     warn.limit.reject warn.limit.meanrw
     5.000e-01 5.000e-01
     nResample max.it k.max maxit.scale k.m_s
     500 50 200 200 20
     trace.lev mts compute.rd fast.s.large.n
     0 1000 0 2000
     psi subsampling cov
     "bisquare" "nonsingular" ".vcov.w"
     split.type compute.outlier.stats
     "f" "SM"
     seed : int(0)
     >
     > set.seed(13)
     > fiM.S <- lmrob(Y ~ Region + X1 + X2 + X3, education, init=lmrob.M.S)
     > out3 <- capture.output(summary(fiM.S))
     >
     > ## must be the same {apart from the "init=" in the call}:
     > i <- 3
     > stopifnot(identical(out2[-i], out3[-i]))
     > ## the difference:
     > c(rbind(out2[i], out3[i]))
     [1] "lmrob(formula = Y ~ Region + X1 + X2 + X3, data = education, init = \"M-S\")"
     [2] "lmrob(formula = Y ~ Region + X1 + X2 + X3, data = education, init = lmrob.M.S)"
     >
     >
     > ### "Skipping design matrix equilibration" warning can arise for reasonable designs -----
     > set.seed(1)
     > x2 <- matrix(rnorm(2*30), 30, 2)
     > data <- data.frame(y = rnorm(30), group = rep(letters[1:3], each=10), x2)
     >
     > obj <- lmrob(y ~ ., data, init="M-S", trace.lev=1)
     lmrob_S(n = 30, nRes = 500): fast_s() [non-large n]:
     Subsampling to find candidate betas:
     Now refine() to convergence for 2 very best ones:
     Best[0]: convergence (87 iter.): -> improved scale to 0.921076968958331
     Best[1]: convergence (41 iter.)
     lmrob.S(): scale = 0.921077; coeff.=
     [1] 0.055853342 -0.395949198 0.068785230 -0.002081381 0.239242595
     init converged (remaining method = "M") -> coef=
     (Intercept) groupb groupc X1 X2
     0.055853342 -0.395949198 0.068785230 -0.002081381 0.239242595
     lmrob_MM(): rwls():
     rwls() used 12 it.; last ||b0 - b1||_1 = 6.71467e-08, L(b1) = 0.0919204728697; convergence
     lmrob..MM..fit(*, obj) --> updating .. qr(x * rweights) -> rank= 5, outlierStats()
     step "M" -> new coef=
     (Intercept) groupb groupc X1 X2
     0.41184656 -0.74163256 -0.36361231 0.02399208 0.51678384
     Warning message:
     In lmrob.M.S(x, y, control, mf = mf) :
     No categorical variables found in model. Reverting to S-estimator.
     >
     > ## illustration: the zero row is introduced during the orthogonalization of x2 wrt x1
     > ## l1 regression always produces p zero residuals
     > ## by chance, the zero residuals of multiple columns happen to be on the same row
     > sf <- splitFrame(obj$model)
     > x1 <- sf$x1
     > x2 <- sf$x2
     > control <- obj$control
     >
     > ## orthogonalize
     > x2.tilde <- x2
     >
     > for(i in 1:ncol(x2)) {
     + tmp <- lmrob.lar(x1, x2[,i], control)
     + x2.tilde[,i] <- tmp$resid
     + }
     Error in lmrob.lar(x1, x2[, i], control) : p > 0 is not TRUE
     Calls: lmrob.lar -> stopifnot
     Execution halted
Flavor: r-devel-windows-ix86+x86_64