Last updated on 2024-05-27 12:51:45 CEST.
| Flavor | Version | Tinstall | Tcheck | Ttotal | Status | Flags |
|---|---|---|---|---|---|---|
| r-devel-linux-x86_64-debian-clang | 2.7-8 | 56.36 | 63.83 | 120.19 | ERROR | |
| r-devel-linux-x86_64-debian-gcc | 2.7-8 | 22.99 | 46.09 | 69.08 | OK | |
| r-devel-linux-x86_64-fedora-clang | 2.7-8 | 149.78 | ERROR | |||
| r-devel-linux-x86_64-fedora-gcc | 2.7-8 | 115.85 | OK | |||
| r-devel-windows-x86_64 | 2.7-8 | 25.00 | 92.00 | 117.00 | OK | |
| r-patched-linux-x86_64 | 2.7-8 | 30.44 | 58.82 | 89.26 | OK | |
| r-release-linux-x86_64 | 2.7-8 | 27.51 | 59.47 | 86.98 | OK | |
| r-release-macos-arm64 | 2.7-8 | 37.00 | OK | |||
| r-release-windows-x86_64 | 2.7-8 | 25.00 | 91.00 | 116.00 | OK | |
| r-oldrel-macos-arm64 | 2.7-8 | 40.00 | OK | |||
| r-oldrel-macos-x86_64 | 2.7-8 | 93.00 | OK | |||
| r-oldrel-windows-x86_64 | 2.7-8 | 27.00 | 102.00 | 129.00 | OK |
Version: 2.7-8
Check: examples
Result: ERROR
Running examples in ‘robeth-Ex.R’ failed
The error most likely occurred in:
> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: robeth-package
> ### Title: Interface for the FORTRAN programs developed at the ETH-Zuerich
> ### and then at IUMSP-Lausanne
> ### Aliases: robeth-package robeth
> ### Keywords: package robust
>
> ### ** Examples
>
> library(robeth)
>
> #
> # ------------- Examples of Chapter 1: Location problems ------------------------------
> #
> y <- c(6.0,7.0,5.0,10.5,8.5,3.5,6.1,4.0,4.6,4.5,5.9,6.5)
> n <- 12
> dfvals()
NULL
> #-----------------------------------------------------------------------
> # M-estimate (tm) of location and confidence interval (tl,tu)
> #
> dfrpar(as.matrix(y),"huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> libeth()
$bta
[1] 0.3550823
attr(,"Csingle")
[1] TRUE
> s <- lilars(y); t0 <- s$theta; s0 <- s$sigma
> s <- lyhalg(y=y,theta=t0,sigmai=s0)
> tm <- s$theta; vartm <- s$var
> s <- quant(0.975)
> tl <- tm-s$x*sqrt(vartm)
> tu <- tm+s$x*sqrt(vartm)
> #-----------------------------------------------------------------------
> # Hodges and Lehmann estimate (th) and confidence interval (zl,zu)
> #
> m <- n*(n+1)/2 # n even
> k1 <- m/2; k2 <- k1+1
> z1 <- lyhdle(y=y,k=k1); z2 <- lyhdle(y=y,k=k2)
> th <- (z1$hdle+z2$hdle)/2.
> ku <- liindh(0.95,n); kl <- liindh(0.05,n)
> zu <- lyhdle(y=y,k=ku$k); zl <- lyhdle(y=y,k=kl$k)
> #.......................................................................
> {
+ cat(" tm, tl, tu : ",round(c(tm,tl,tu),3),"\n")
+ cat(" th, zl, zu : ",round(c(th,zl$hdle,zu$hdle),3),"\n")
+ }
tm, tl, tu : 5.809 4.748 6.87
th, zl, zu : 5.85 5 7
> # tm, tl, tu : 5.809 4.748 6.87
> # th, zl, zu : 5.85 5 7
> #=======================================================================
> #
> # Two sample problem
> #
> y <- c(17.9,13.3,10.6,7.6,5.7,5.6,5.4,3.3,3.1,0.9)
> n <- 10
> x <- c(7.7,5.0,1.7,0.0,-3.0,-3.1,-10.5)
> m <- 7
> #-----------------------------------------------------------------------
> # Estimate (dm) and confidence interval [dl,du] based on Mann-Whitney
> #
> k1 <- m*n/2; k2 <- k1+1
> s1 <- lymnwt(x=x,y=y,k=k1); s2 <- lymnwt(x=x,y=y,k=k2)
> dm <- (s1$tmnwt+s2$tmnwt)/2.0
> sl <- liindw(0.05,m,n); kl <- sl$k
> s <- lymnwt(x=x,y=y,k=kl); dl <- s$tmnwt
> s <- lymnwt(x=x,y=y,k=m*n-kl+1); du <- s$tmnwt
> #-----------------------------------------------------------------------
> # Tau-test . P-value (P)
> #
> z <- c(x,y)
> dfrpar(as.matrix(z),"huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> libeth()
$bta
[1] 0.3550823
attr(,"Csingle")
[1] TRUE
> s <- lytau2(z=z,m=m,n=n)
> P <- s$p
> #
> # estimate (tm) and confidence interval (tl,tu)
> #
> tm <- s$deltal
> c22<- s$cov[3]
> s <- quant(0.975)
> tl <- tm-s$x*sqrt(c22)
> tu <- tm+s$x*sqrt(c22)
> #.......................................................................
> {
+ cat("dm, dl, du:",round(c(dm,dl,du),3),"\n")
+ cat("P, tm, tl, tu:",round(c(P,tm,tl,tu),3),"\n")
+ }
dm, dl, du: 6.35 2.9 12.9
P, tm, tl, tu: 0.014 6.918 1.562 12.273
> # dm, dl, du: 6.35 2.9 12.9
> # P, tm, tl, tu: 0.014 6.918 1.562 12.273
>
> #
> # Examples of Chapter 2: M-estimates of coefficients and scale in linear regression
> #
> # Read data; declare the vector wgt; load defaults
> #
> z <- c(-1, -2, 0, 35, 1, 0, -3, 20,
+ -1, -2, 0, 30, 1, 0, -3, 39,
+ -1, -2, 0, 24, 1, 0, -3, 16,
+ -1, -2, 0, 37, 1, 0, -3, 27,
+ -1, -2, 0, 28, 1, 0, -3, -12,
+ -1, -2, 0, 73, 1, 0, -3, 2,
+ -1, -2, 0, 31, 1, 0, -3, 31,
+ -1, -2, 0, 21, 1, 0, -1, 26,
+ -1, -2, 0, -5, 1, 0, -1, 60,
+ -1, 0, 0, 62, 1, 0, -1, 48,
+ -1, 0, 0, 67, 1, 0, -1, -8,
+ -1, 0, 0, 95, 1, 0, -1, 46,
+ -1, 0, 0, 62, 1, 0, -1, 77,
+ -1, 0, 0, 54, 1, 0, 1, 57,
+ -1, 0, 0, 56, 1, 0, 1, 89,
+ -1, 0, 0, 48, 1, 0, 1, 103,
+ -1, 0, 0, 70, 1, 0, 1, 129,
+ -1, 0, 0, 94, 1, 0, 1, 139,
+ -1, 0, 0, 42, 1, 0, 1, 128,
+ -1, 2, 0, 116, 1, 0, 1, 89,
+ -1, 2, 0, 105, 1, 0, 1, 86,
+ -1, 2, 0, 91, 1, 0, 3, 140,
+ -1, 2, 0, 94, 1, 0, 3, 133,
+ -1, 2, 0, 130, 1, 0, 3, 142,
+ -1, 2, 0, 79, 1, 0, 3, 118,
+ -1, 2, 0, 120, 1, 0, 3, 137,
+ -1, 2, 0, 124, 1, 0, 3, 84,
+ -1, 2, 0, -8, 1, 0, 3, 101)
> xx <- matrix(z,ncol=4, byrow=TRUE)
> dimnames(xx) <- list(NULL,c("z2","xS","xT","y"))
> z2 <- xx[,"z2"]; xS <- xx[,"xS"]; xT <- xx[,"xT"]
> x <- cbind(1, z2, xS+xT, xS-xT, xS^2+xT^2, xS^2-xT^2, xT^3)
> y <- xx[,"y"]
> wgt <- vector("numeric",length(y))
> n <- 56; np <- 7
> dfvals()
NULL
> # Set parameters for Huber estimate
> dfrpar(x, "huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> # Compute the constants beta, bet0, epsi2 and epsip
> ribeth(wgt)
$d
[1] 1.345
attr(,"Csingle")
[1] TRUE
$bta
[1] 0.3550823
attr(,"Csingle")
[1] TRUE
> ribet0(wgt)
$bt0
[1] 0.6741892
attr(,"Csingle")
[1] TRUE
> s <- liepsh()
> epsi2 <- s$epsi2; epsip <- s$epsip
> #
> # Least squares solution (theta0,sigma0)
> #
> z <- riclls(x, y)
> theta0<- z$theta; sigma0 <- z$sigma
> # Preliminary estimate of the covariance matrix of the coefficients
> cv <- kiascv(z$xt, fu=epsi2/epsip^2, fb=0.)
> cov <- cv$cov
> #-----------------------------------------------------------------------
> # Solution (theta1,sigma1) by means of RYHALG.
> #
> zr <- ryhalg(x,y,theta0,wgt,cov,sigmai=sigma0,ic=0)
> theta1<- zr$theta[1:np]; sigma1 <- zr$sigmaf; nit1 <- zr$nit
> #-----------------------------------------------------------------------
> # Solution (theta2,sigma2) by means of RYWALG (recompute cov)
> #
> cv <- ktaskv(x, f=epsi2/epsip^2)
> zr <- rywalg(x, y, theta0, wgt, cv$cov, sigmai=sigma0)
> theta2 <- zr$theta[1:np]; sigma2 <- zr$sigmaf; nit2 <- zr$nit
> #-----------------------------------------------------------------------
> # Solution (theta3,sigma3) by means of RYSALG with ISIGMA=2.
> #
> zr <- rysalg(x,y, theta0, wgt, cv$cov, sigma0, isigma=2)
> theta3 <- zr$theta[1:np]; sigma3 <- zr$sigmaf; nit3 <- zr$nit
> #-----------------------------------------------------------------------
> # Solution (theta4,sigma4) by means of RYNALG with ICNV=2 and ISIGMA=0.
> #
> # Invert cov
> covm1 <- cv$cov
> zc <- mchl(covm1,np)
> zc <- minv(zc$a, np)
> zc <- mtt1(zc$r,np)
> covm1 <- zc$b
> zr <- rynalg(x,y, theta0, wgt, covm1, sigmai=sigma3,
+ iopt=1, isigma=0, icnv=2)
> theta4 <- zr$theta[1:np]; sigma4 <- zr$sigmaf; nit4 <- zr$nit
> #.......................................................................
> {
+ cat("theta0 : ",round(theta0[1:np],3),"\n")
+ cat("sigma0 : ",round(sigma0,3),"\n")
+ cat("theta1 : ",round(theta1,3),"\n")
+ cat("sigma1, nit1 : ",round(sigma1,3),nit1,"\n")
+ cat("theta2 : ",round(theta2,3),"\n")
+ cat("sigma2, nit2 : ",round(sigma2,3),nit2,"\n")
+ cat("theta3 : ",round(theta3,3),"\n")
+ cat("sigma3, nit3 : ",round(sigma3,3),nit3,"\n")
+ cat("theta4 : ",round(theta4,3),"\n")
+ cat("sigma4, nit4 : ",round(sigma4,3),nit4,"\n")
+ }
theta0 : 68.634 3.634 24.081 -8.053 -0.446 -0.179 -1.634
sigma0 : 26.635
theta1 : 70.006 5.006 24.742 -6.246 -0.079 0.434 -1.487
sigma1, nit1 : 23.564 7
theta2 : 70.006 5.006 24.742 -6.245 -0.079 0.434 -1.487
sigma2, nit2 : 23.563 7
theta3 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
sigma3, nit3 : 22.249 3
theta4 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
sigma4, nit4 : 22.249 3
> # theta0 : 68.634 3.634 24.081 -8.053 -0.446 -0.179 -1.634
> # sigma0 : 26.635
> # theta1 : 70.006 5.006 24.742 -6.246 -0.079 0.434 -1.487
> # sigma1, nit1 : 23.564 7
> # theta2 : 70.006 5.006 24.742 -6.245 -0.079 0.434 -1.487
> # sigma2, nit2 : 23.563 7
> # theta3 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
> # sigma3, nit3 : 22.249 3
> # theta4 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
> # sigma4, nit4 : 22.249 3
>
>
> #
> # ---- Examples of Chapter 3: Weights for bounded influence regression ------
> #
>
> #=======================================================================
> rbmost <- function(x,y,cc,usext=userfd) {
+ n <- nrow(x); np <- ncol(x); dfcomn(xk=np)
+ .dFvPut(1,"itw")
+ z <- wimedv(x)
+ z <- wyfalg(x, z$a, y, exu=usext); nitw <- z$nit
+ wgt <- 1/z$dist; wgt[wgt>1.e6] <- 1.e6
+ z <- comval()
+ bto <- z$bt0; ipso <- z$ipsi; co <- z$c
+ z <- ribet0(wgt, itype=2, isqw=0)
+ xt <- x*wgt; yt <- y * wgt
+ z <- rilars(xt, yt)
+ theta0 <- z$theta; sigma0 <- z$sigma
+ rs <- z$rs/wgt; r1 <- rs/sigma0
+ dfcomn(ipsi=1,c=cc)
+ z <- liepsh(cc)
+ den <- z$epsip
+ g <- Psp(r1)/den # (see Psi in Chpt. 14)
+ dfcomn(ipsi=ipso, c=co, bet0=bto)
+ list(theta=theta0, sigma=sigma0, rs=rs, g=g, nitw=nitw)
+ }
> #-----------------------------------------------------------------------
> # Mallows-standard estimate (with wyfalg and rywalg)
> #
> Mal.Std <- function(x, y, b2=-1, cc=-1, isigma=2) {
+ n <- length(y); np <- ncol(x)
+ dfrpar(x, "Mal-Std", b2, cc); .dFv <- .dFvGet()
+ if (isigma==1) {dfcomn(d=.dFv$ccc); .dFvPut(1,"isg")}
+ # Weights
+ z <- wimedv(x)
+ z <- wyfalg(x, z$a, y); nitw <- z$nit
+ wgt <- Www(z$dist) # See Www in Chpt. 14
+ # Initial cov. matrix of coefficient estimates
+ z <- kiedch(wgt)
+ cov <- ktaskw(x, z$d, z$e, f=1/n)
+ # Initial theta and sigma
+ z <- rbmost(x,y,1.5,userfd)
+ theta0 <- z$theta; sigma0 <- z$sigma; nitw0 <- z$nitw
+ # Final theta and sigma
+ if (isigma==1) ribeth(wgt) else ribet0(wgt)
+ z <- rywalg(x, y,theta0,wgt,cov$cov, sigmai=sigma0)
+ theta1 <- z$theta[1:np]; sigma1 <- z$sigmaf; nit1 <- z$nit
+ # Covariance matrix of coefficient estimates
+ z <- kfedcc(wgt, z$rs, sigma=sigma1)
+ cov <- ktaskw(x, z$d, z$e, f=(sigma1^2)/n)
+ sd1 <- NULL; for (i in 1:np) { j <- i*(i+1)/2
+ sd1 <- c(sd1,cov$cov[j]) }
+ sd1 <- sqrt(sd1)
+ #.......................................................................
+ {
+ cat("rbmost: theta0 : ",round(theta0[1:np],3),"\n")
+ cat("rbmost: sigma0, nitw : ",round(sigma0,3),nitw0,"\n")
+ cat("Mallows-Std: theta1 : ",round(theta1,3),"\n")
+ cat("Mallows-Std: sd1 : ",round(sd1,3),"\n")
+ cat("Mallows-Std: sigma1, nitw, nit1 : ",round(sigma1,3),nitw,nit1,"\n")
+ }
+
+ #.......................................................................
+ list(theta0=theta0[1:np], sigma0=sigma0, nitw=nitw,
+ theta1=theta1, sigma1=sigma1, nit1=nit1, sd1=sd1)}
> #-----------------------------------------------------------------------
> # Krasker-Welsch estimate (with wynalg and rynalg)
> #
> Kra.Wel <- function(x, y, ckw=-1, isigma=2) {
+ n <- length(y); np <- ncol(x)
+ dfrpar(x, "Kra-Wel", ckw); .dFv <- .dFvGet()
+ if (isigma==1) {dfcomn(d=.dFv$ccc); .dFvPut(1,"isg")}
+ # Weights
+ z <- wimedv(x)
+ z <- wynalg(x, z$a); nitw <- z$nit
+ wgt <- Www(z$dist) # See Www in Chpt. 14
+ # Initial cov. matrix of coefficient estimates
+ z <- kiedch(wgt)
+ cov <- ktaskw(x, z$d, z$e, f=1/n)
+ # Initial theta and sigma
+ z <- rbmost(x, y, cc=1.5)
+ theta0 <- z$theta; sigma0 <- z$sigma; nitw0 <- z$nitw
+ # Final theta and sigma
+ if (isigma==1) ribeth(wgt) else ribet0(wgt)
+ z <- rynalg(x, y,theta0,wgt,cov$cov, sigmai=sigma0)
+ theta2 <- z$theta[1:np]; sigma2 <- z$sigma; nit2 <- z$nit
+ #
+ # Covariance matrix of coefficient estimates
+ #
+ z <- kfedcc(wgt, z$rs, sigma=sigma2)
+ cov <- ktaskw(x, z$d, z$e, f=(sigma2^2)/n)
+ sd2 <- NULL; for (i in 1:np) { j <- i*(i+1)/2
+ sd2 <- c(sd2,cov$cov[j]) }
+ sd2 <- sqrt(sd2)
+ #.......................................................................
+ {
+ cat("rbmost: theta0 : ",round(theta0[1:np],3),"\n")
+ cat("rbmost: sigma0, nitw : ",round(sigma0,3),nitw0,"\n")
+ cat("Krasker-Welsch: theta2 : ",round(theta2,3),"\n")
+ cat("Krasker-Welsch: sd2 : ",round(sd2,3),"\n")
+ cat("Krasker-Welsch: sigma2, nitw, nit2 : ",round(sigma2,3),nitw,nit2,"\n")
+ }
+ #.......................................................................
+ list(theta0=theta0[1:np], sigma0=sigma0, nitw=nitw,
+ theta2=theta2, sigma2=sigma2, nit2=nit2, sd2=sd2)}
> #-----------------------------------------------------------------------
> # Read data; load defaults
> #
> z <- c( 8.2, 4, 23.005, 1, 7.6, 5, 23.873, 1,
+ 4.6, 0, 26.417, 1, 4.3, 1, 24.868, 1,
+ 5.9, 2, 29.895, 1, 5.0, 3, 24.200, 1,
+ 6.5, 4, 23.215, 1, 8.3, 5, 21.862, 1,
+ 10.1, 0, 22.274, 1, 13.2, 1, 23.830, 1,
+ 12.6, 2, 25.144, 1, 10.4, 3, 22.430, 1,
+ 10.8, 4, 21.785, 1, 13.1, 5, 22.380, 1,
+ 13.3, 0, 23.927, 1, 10.4, 1, 33.443, 1,
+ 10.5, 2, 24.859, 1, 7.7, 3, 22.686, 1,
+ 10.0, 0, 21.789, 1, 12.0, 1, 22.041, 1,
+ 12.1, 4, 21.033, 1, 13.6, 5, 21.005, 1,
+ 15.0, 0, 25.865, 1, 13.5, 1, 26.290, 1,
+ 11.5, 2, 22.932, 1, 12.0, 3, 21.313, 1,
+ 13.0, 4, 20.769, 1, 14.1, 5, 21.393, 1)
> x <- matrix(z, ncol=4, byrow=TRUE)
> y <- c( 7.6, 7.7, 4.3, 5.9, 5.0, 6.5, 8.3, 8.2, 13.2, 12.6,
+ 10.4, 10.8, 13.1, 12.3, 10.4, 10.5, 7.7, 9.5, 12.0, 12.6,
+ 13.6, 14.1, 13.5, 11.5, 12.0, 13.0, 14.1, 15.1)
> dfvals()
NULL
> dfcomn(xk=4)
$ipsi
[1] -9
$iucv
[1] -1
$iwww
[1] -1
> cat("Results\n")
Results
> z1 <- Mal.Std(x, y)
rbmost: theta0 : 0.674 -0.171 -0.678 20.043
rbmost: sigma0, nitw : 0.846 20
Mallows-Std: theta1 : 0.721 -0.174 -0.654 18.868
Mallows-Std: sd1 : 0.052 0.145 0.166 4.285
Mallows-Std: sigma1, nitw, nit1 : 0.774 38 7
> z2 <- Kra.Wel(x, y)
rbmost: theta0 : 0.674 -0.171 -0.678 20.043
rbmost: sigma0, nitw : 0.846 20
Krasker-Welsch: theta2 : 0.68 -0.177 -0.679 19.974
Krasker-Welsch: sd2 : 0.047 0.067 0.141 3.664
Krasker-Welsch: sigma2, nitw, nit2 : 0.732 21 4
>
>
> #
> # ---- Examples of Chapter 4: Covariance matrix of the coefficient estimates ------
> #
>
> #
> # Read x[1:4] and then set x[,4] <- 1
> #
> z <- c(80, 27, 89, 1, 80, 27, 88, 1, 75, 25, 90, 1,
+ 62, 24, 87, 1, 62, 22, 87, 1, 62, 23, 87, 1,
+ 62, 24, 93, 1, 62, 24, 93, 1, 58, 23, 87, 1,
+ 58, 18, 80, 1, 58, 18, 89, 1, 58, 17, 88, 1,
+ 58, 18, 82, 1, 58, 19, 93, 1, 50, 18, 89, 1,
+ 50, 18, 86, 1, 50, 19, 72, 1, 50, 19, 79, 1,
+ 50, 20, 80, 1, 56, 20, 82, 1, 70, 20, 91, 1)
> x <- matrix(z, ncol=4, byrow=TRUE)
> n <- 21; np <- 4; ncov <- np*(np+1)/2
> dfvals()
NULL
> # Cov. matrix of Huber-type estimates
> dfrpar(x, "huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> s <- liepsh()
> epsi2 <- s$epsi2; epsip <- s$epsip
> z <- rimtrf(x)
> xt <- z$x; sg <- z$sg; ip <- z$ip
> zc <- kiascv(xt, fu=epsi2/epsip^2, fb=0.)
> covi <- zc$cov # Can be used in ryhalg with ic=0
> zc <- kfascv(xt, covi, f=1, sg=sg, ip=ip)
> covf <- zc$cov
> #.......................................................................
> str <- rep(" ", ncov); str[cumsum(1:np)] <- "\n"
> {
+ cat("covf:\n")
+ cat(round(covf,6),sep=str)
+ }
covf:
0.00182
-0.003653 0.013553
-0.000715 1e-06 0.002444
0.028778 -0.065222 -0.167742 14.16076
>
>
> #
> # ---- Examples of Chapter 5: Asymptotic relative efficiency ------
> #
> #-----------------------------------------------------------------------
> # Huber
> #
> dfcomn(ipsi=1,c=1.345,d=1.345)
$ipsi
[1] 1
$iucv
[1] -1
$iwww
[1] -1
> .dFvPut(1,"ite")
NULL
> z <- airef0(mu=3, ialfa=1, sigmx=1)
> #.......................................................................
> {
+ cat(" airef0 : Huber\n reff, alfa, beta, nit: ")
+ cat(round(c(z$reff,z$alfa,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airef0 : Huber
reff, alfa, beta, nit: 0.95, 0, 0, 0
> #-----------------------------------------------------------------------
> # Schweppe: Krasker-Welsch
> #
> dfcomn(ipsi=1,c=3.76,iucv=3,ckw=3.76,iwww=1)
$ipsi
[1] 1
$iucv
[1] 3
$iwww
[1] 1
> .dFvPut(3,"ite")
NULL
> z <- airef0(mu=3, ialfa=1, sigmx=1)
> #.......................................................................
> {
+ cat(" airef0 : Krasker-Welsch\n reff, alfa, beta, nit: ")
+ cat(round(c(z$reff,z$alfa,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airef0 : Krasker-Welsch
reff, alfa, beta, nit: 0.95, 1.091, 1.17, 6
> #-----------------------------------------------------------------------
> # Mallows-Standard
> #
> dfcomn(ipsi=1,c=1.5,iucv=1,a2=0,b2=6.16,iww=3)
$ipsi
[1] 1
$iucv
[1] 1
$iwww
[1] 3
> .dFvPut(2,"ite")
NULL
> z <- airef0(mu=3, ialfa=1, sigmx=1)
> #.......................................................................
> {
+ cat(" airef0 : Mallows-Std \n reff, alfa, beta, nit: ")
+ cat(round(c(z$reff,z$alfa,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airef0 : Mallows-Std
reff, alfa, beta, nit: 0.95, 1.031, 1.09, 5
> #=======================================================================
> z <- c(1, 0, 0,
+ 1, 0, 0,
+ 1, 0, 0,
+ 1, 0, 0,
+ 0, 1, 0,
+ 0, 1, 0,
+ 0, 1, 0,
+ 0, 1, 0,
+ 0, 0, 1,
+ 0, 0, 1,
+ 0, 0, 1,
+ 0, 0, 1)
> tt <- matrix(z,ncol=3,byrow=TRUE)
> n <- nrow(tt); mu <- 2
> nu <- ncol(tt)
>
> #-----------------------------------------------------------------------
> # Huber
> #
> dfrpar(tt,"Huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> z <- airefq(tt, mu=mu, sigmx=1)
> #.......................................................................
> {
+ cat(" airefq : Huber\n reff, beta, nit: ")
+ cat(round(c(z$reff,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airefq : Huber
reff, beta, nit: 0.95, 0, 0
> #-----------------------------------------------------------------------
> # Krasker-Welsch
> #
> dfrpar(tt,"kra-wel",upar=3.755)
$itypw
[1] 1
$itype
[1] 3
$isigma
[1] 2
> z <- airefq(tt, mu=mu, sigmx=1,init=1)
Floating point exception
Flavor: r-devel-linux-x86_64-debian-clang
Version: 2.7-8
Check: examples
Result: ERROR
Running examples in ‘robeth-Ex.R’ failed
The error most likely occurred in:
> ### Name: robeth-package
> ### Title: Interface for the FORTRAN programs developed at the ETH-Zuerich
> ### and then at IUMSP-Lausanne
> ### Aliases: robeth-package robeth
> ### Keywords: package robust
>
> ### ** Examples
>
> library(robeth)
>
> #
> # ------------- Examples of Chapter 1: Location problems ------------------------------
> #
> y <- c(6.0,7.0,5.0,10.5,8.5,3.5,6.1,4.0,4.6,4.5,5.9,6.5)
> n <- 12
> dfvals()
NULL
> #-----------------------------------------------------------------------
> # M-estimate (tm) of location and confidence interval (tl,tu)
> #
> dfrpar(as.matrix(y),"huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> libeth()
$bta
[1] 0.3550823
attr(,"Csingle")
[1] TRUE
> s <- lilars(y); t0 <- s$theta; s0 <- s$sigma
> s <- lyhalg(y=y,theta=t0,sigmai=s0)
> tm <- s$theta; vartm <- s$var
> s <- quant(0.975)
> tl <- tm-s$x*sqrt(vartm)
> tu <- tm+s$x*sqrt(vartm)
> #-----------------------------------------------------------------------
> # Hodges and Lehmann estimate (th) and confidence interval (zl,zu)
> #
> m <- n*(n+1)/2 # n even
> k1 <- m/2; k2 <- k1+1
> z1 <- lyhdle(y=y,k=k1); z2 <- lyhdle(y=y,k=k2)
> th <- (z1$hdle+z2$hdle)/2.
> ku <- liindh(0.95,n); kl <- liindh(0.05,n)
> zu <- lyhdle(y=y,k=ku$k); zl <- lyhdle(y=y,k=kl$k)
> #.......................................................................
> {
+ cat(" tm, tl, tu : ",round(c(tm,tl,tu),3),"\n")
+ cat(" th, zl, zu : ",round(c(th,zl$hdle,zu$hdle),3),"\n")
+ }
tm, tl, tu : 5.809 4.748 6.87
th, zl, zu : 5.85 5 7
> # tm, tl, tu : 5.809 4.748 6.87
> # th, zl, zu : 5.85 5 7
> #=======================================================================
> #
> # Two sample problem
> #
> y <- c(17.9,13.3,10.6,7.6,5.7,5.6,5.4,3.3,3.1,0.9)
> n <- 10
> x <- c(7.7,5.0,1.7,0.0,-3.0,-3.1,-10.5)
> m <- 7
> #-----------------------------------------------------------------------
> # Estimate (dm) and confidence interval [dl,du] based on Mann-Whitney
> #
> k1 <- m*n/2; k2 <- k1+1
> s1 <- lymnwt(x=x,y=y,k=k1); s2 <- lymnwt(x=x,y=y,k=k2)
> dm <- (s1$tmnwt+s2$tmnwt)/2.0
> sl <- liindw(0.05,m,n); kl <- sl$k
> s <- lymnwt(x=x,y=y,k=kl); dl <- s$tmnwt
> s <- lymnwt(x=x,y=y,k=m*n-kl+1); du <- s$tmnwt
> #-----------------------------------------------------------------------
> # Tau-test . P-value (P)
> #
> z <- c(x,y)
> dfrpar(as.matrix(z),"huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> libeth()
$bta
[1] 0.3550823
attr(,"Csingle")
[1] TRUE
> s <- lytau2(z=z,m=m,n=n)
> P <- s$p
> #
> # estimate (tm) and confidence interval (tl,tu)
> #
> tm <- s$deltal
> c22<- s$cov[3]
> s <- quant(0.975)
> tl <- tm-s$x*sqrt(c22)
> tu <- tm+s$x*sqrt(c22)
> #.......................................................................
> {
+ cat("dm, dl, du:",round(c(dm,dl,du),3),"\n")
+ cat("P, tm, tl, tu:",round(c(P,tm,tl,tu),3),"\n")
+ }
dm, dl, du: 6.35 2.9 12.9
P, tm, tl, tu: 0.014 6.918 1.562 12.273
> # dm, dl, du: 6.35 2.9 12.9
> # P, tm, tl, tu: 0.014 6.918 1.562 12.273
>
> #
> # Examples of Chapter 2: M-estimates of coefficients and scale in linear regression
> #
> # Read data; declare the vector wgt; load defaults
> #
> z <- c(-1, -2, 0, 35, 1, 0, -3, 20,
+ -1, -2, 0, 30, 1, 0, -3, 39,
+ -1, -2, 0, 24, 1, 0, -3, 16,
+ -1, -2, 0, 37, 1, 0, -3, 27,
+ -1, -2, 0, 28, 1, 0, -3, -12,
+ -1, -2, 0, 73, 1, 0, -3, 2,
+ -1, -2, 0, 31, 1, 0, -3, 31,
+ -1, -2, 0, 21, 1, 0, -1, 26,
+ -1, -2, 0, -5, 1, 0, -1, 60,
+ -1, 0, 0, 62, 1, 0, -1, 48,
+ -1, 0, 0, 67, 1, 0, -1, -8,
+ -1, 0, 0, 95, 1, 0, -1, 46,
+ -1, 0, 0, 62, 1, 0, -1, 77,
+ -1, 0, 0, 54, 1, 0, 1, 57,
+ -1, 0, 0, 56, 1, 0, 1, 89,
+ -1, 0, 0, 48, 1, 0, 1, 103,
+ -1, 0, 0, 70, 1, 0, 1, 129,
+ -1, 0, 0, 94, 1, 0, 1, 139,
+ -1, 0, 0, 42, 1, 0, 1, 128,
+ -1, 2, 0, 116, 1, 0, 1, 89,
+ -1, 2, 0, 105, 1, 0, 1, 86,
+ -1, 2, 0, 91, 1, 0, 3, 140,
+ -1, 2, 0, 94, 1, 0, 3, 133,
+ -1, 2, 0, 130, 1, 0, 3, 142,
+ -1, 2, 0, 79, 1, 0, 3, 118,
+ -1, 2, 0, 120, 1, 0, 3, 137,
+ -1, 2, 0, 124, 1, 0, 3, 84,
+ -1, 2, 0, -8, 1, 0, 3, 101)
> xx <- matrix(z,ncol=4, byrow=TRUE)
> dimnames(xx) <- list(NULL,c("z2","xS","xT","y"))
> z2 <- xx[,"z2"]; xS <- xx[,"xS"]; xT <- xx[,"xT"]
> x <- cbind(1, z2, xS+xT, xS-xT, xS^2+xT^2, xS^2-xT^2, xT^3)
> y <- xx[,"y"]
> wgt <- vector("numeric",length(y))
> n <- 56; np <- 7
> dfvals()
NULL
> # Set parameters for Huber estimate
> dfrpar(x, "huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> # Compute the constants beta, bet0, epsi2 and epsip
> ribeth(wgt)
$d
[1] 1.345
attr(,"Csingle")
[1] TRUE
$bta
[1] 0.3550823
attr(,"Csingle")
[1] TRUE
> ribet0(wgt)
$bt0
[1] 0.6741892
attr(,"Csingle")
[1] TRUE
> s <- liepsh()
> epsi2 <- s$epsi2; epsip <- s$epsip
> #
> # Least squares solution (theta0,sigma0)
> #
> z <- riclls(x, y)
> theta0<- z$theta; sigma0 <- z$sigma
> # Preliminary estimate of the covariance matrix of the coefficients
> cv <- kiascv(z$xt, fu=epsi2/epsip^2, fb=0.)
> cov <- cv$cov
> #-----------------------------------------------------------------------
> # Solution (theta1,sigma1) by means of RYHALG.
> #
> zr <- ryhalg(x,y,theta0,wgt,cov,sigmai=sigma0,ic=0)
> theta1<- zr$theta[1:np]; sigma1 <- zr$sigmaf; nit1 <- zr$nit
> #-----------------------------------------------------------------------
> # Solution (theta2,sigma2) by means of RYWALG (recompute cov)
> #
> cv <- ktaskv(x, f=epsi2/epsip^2)
> zr <- rywalg(x, y, theta0, wgt, cv$cov, sigmai=sigma0)
> theta2 <- zr$theta[1:np]; sigma2 <- zr$sigmaf; nit2 <- zr$nit
> #-----------------------------------------------------------------------
> # Solution (theta3,sigma3) by means of RYSALG with ISIGMA=2.
> #
> zr <- rysalg(x,y, theta0, wgt, cv$cov, sigma0, isigma=2)
> theta3 <- zr$theta[1:np]; sigma3 <- zr$sigmaf; nit3 <- zr$nit
> #-----------------------------------------------------------------------
> # Solution (theta4,sigma4) by means of RYNALG with ICNV=2 and ISIGMA=0.
> #
> # Invert cov
> covm1 <- cv$cov
> zc <- mchl(covm1,np)
> zc <- minv(zc$a, np)
> zc <- mtt1(zc$r,np)
> covm1 <- zc$b
> zr <- rynalg(x,y, theta0, wgt, covm1, sigmai=sigma3,
+ iopt=1, isigma=0, icnv=2)
> theta4 <- zr$theta[1:np]; sigma4 <- zr$sigmaf; nit4 <- zr$nit
> #.......................................................................
> {
+ cat("theta0 : ",round(theta0[1:np],3),"\n")
+ cat("sigma0 : ",round(sigma0,3),"\n")
+ cat("theta1 : ",round(theta1,3),"\n")
+ cat("sigma1, nit1 : ",round(sigma1,3),nit1,"\n")
+ cat("theta2 : ",round(theta2,3),"\n")
+ cat("sigma2, nit2 : ",round(sigma2,3),nit2,"\n")
+ cat("theta3 : ",round(theta3,3),"\n")
+ cat("sigma3, nit3 : ",round(sigma3,3),nit3,"\n")
+ cat("theta4 : ",round(theta4,3),"\n")
+ cat("sigma4, nit4 : ",round(sigma4,3),nit4,"\n")
+ }
theta0 : 68.634 3.634 24.081 -8.053 -0.446 -0.179 -1.634
sigma0 : 26.635
theta1 : 70.006 5.006 24.742 -6.246 -0.079 0.434 -1.487
sigma1, nit1 : 23.564 7
theta2 : 70.006 5.006 24.742 -6.245 -0.079 0.434 -1.487
sigma2, nit2 : 23.563 7
theta3 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
sigma3, nit3 : 22.249 3
theta4 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
sigma4, nit4 : 22.249 3
> # theta0 : 68.634 3.634 24.081 -8.053 -0.446 -0.179 -1.634
> # sigma0 : 26.635
> # theta1 : 70.006 5.006 24.742 -6.246 -0.079 0.434 -1.487
> # sigma1, nit1 : 23.564 7
> # theta2 : 70.006 5.006 24.742 -6.245 -0.079 0.434 -1.487
> # sigma2, nit2 : 23.563 7
> # theta3 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
> # sigma3, nit3 : 22.249 3
> # theta4 : 69.993 5.002 24.766 -6.214 -0.055 0.44 -1.48
> # sigma4, nit4 : 22.249 3
>
>
> #
> # ---- Examples of Chapter 3: Weights for bounded influence regression ------
> #
>
> #=======================================================================
> rbmost <- function(x,y,cc,usext=userfd) {
+ n <- nrow(x); np <- ncol(x); dfcomn(xk=np)
+ .dFvPut(1,"itw")
+ z <- wimedv(x)
+ z <- wyfalg(x, z$a, y, exu=usext); nitw <- z$nit
+ wgt <- 1/z$dist; wgt[wgt>1.e6] <- 1.e6
+ z <- comval()
+ bto <- z$bt0; ipso <- z$ipsi; co <- z$c
+ z <- ribet0(wgt, itype=2, isqw=0)
+ xt <- x*wgt; yt <- y * wgt
+ z <- rilars(xt, yt)
+ theta0 <- z$theta; sigma0 <- z$sigma
+ rs <- z$rs/wgt; r1 <- rs/sigma0
+ dfcomn(ipsi=1,c=cc)
+ z <- liepsh(cc)
+ den <- z$epsip
+ g <- Psp(r1)/den # (see Psi in Chpt. 14)
+ dfcomn(ipsi=ipso, c=co, bet0=bto)
+ list(theta=theta0, sigma=sigma0, rs=rs, g=g, nitw=nitw)
+ }
> #-----------------------------------------------------------------------
> # Mallows-standard estimate (with wyfalg and rywalg)
> #
> Mal.Std <- function(x, y, b2=-1, cc=-1, isigma=2) {
+ n <- length(y); np <- ncol(x)
+ dfrpar(x, "Mal-Std", b2, cc); .dFv <- .dFvGet()
+ if (isigma==1) {dfcomn(d=.dFv$ccc); .dFvPut(1,"isg")}
+ # Weights
+ z <- wimedv(x)
+ z <- wyfalg(x, z$a, y); nitw <- z$nit
+ wgt <- Www(z$dist) # See Www in Chpt. 14
+ # Initial cov. matrix of coefficient estimates
+ z <- kiedch(wgt)
+ cov <- ktaskw(x, z$d, z$e, f=1/n)
+ # Initial theta and sigma
+ z <- rbmost(x,y,1.5,userfd)
+ theta0 <- z$theta; sigma0 <- z$sigma; nitw0 <- z$nitw
+ # Final theta and sigma
+ if (isigma==1) ribeth(wgt) else ribet0(wgt)
+ z <- rywalg(x, y,theta0,wgt,cov$cov, sigmai=sigma0)
+ theta1 <- z$theta[1:np]; sigma1 <- z$sigmaf; nit1 <- z$nit
+ # Covariance matrix of coefficient estimates
+ z <- kfedcc(wgt, z$rs, sigma=sigma1)
+ cov <- ktaskw(x, z$d, z$e, f=(sigma1^2)/n)
+ sd1 <- NULL; for (i in 1:np) { j <- i*(i+1)/2
+ sd1 <- c(sd1,cov$cov[j]) }
+ sd1 <- sqrt(sd1)
+ #.......................................................................
+ {
+ cat("rbmost: theta0 : ",round(theta0[1:np],3),"\n")
+ cat("rbmost: sigma0, nitw : ",round(sigma0,3),nitw0,"\n")
+ cat("Mallows-Std: theta1 : ",round(theta1,3),"\n")
+ cat("Mallows-Std: sd1 : ",round(sd1,3),"\n")
+ cat("Mallows-Std: sigma1, nitw, nit1 : ",round(sigma1,3),nitw,nit1,"\n")
+ }
+
+ #.......................................................................
+ list(theta0=theta0[1:np], sigma0=sigma0, nitw=nitw,
+ theta1=theta1, sigma1=sigma1, nit1=nit1, sd1=sd1)}
> #-----------------------------------------------------------------------
> # Krasker-Welsch estimate (with wynalg and rynalg)
> #
> Kra.Wel <- function(x, y, ckw=-1, isigma=2) {
+ n <- length(y); np <- ncol(x)
+ dfrpar(x, "Kra-Wel", ckw); .dFv <- .dFvGet()
+ if (isigma==1) {dfcomn(d=.dFv$ccc); .dFvPut(1,"isg")}
+ # Weights
+ z <- wimedv(x)
+ z <- wynalg(x, z$a); nitw <- z$nit
+ wgt <- Www(z$dist) # See Www in Chpt. 14
+ # Initial cov. matrix of coefficient estimates
+ z <- kiedch(wgt)
+ cov <- ktaskw(x, z$d, z$e, f=1/n)
+ # Initial theta and sigma
+ z <- rbmost(x, y, cc=1.5)
+ theta0 <- z$theta; sigma0 <- z$sigma; nitw0 <- z$nitw
+ # Final theta and sigma
+ if (isigma==1) ribeth(wgt) else ribet0(wgt)
+ z <- rynalg(x, y,theta0,wgt,cov$cov, sigmai=sigma0)
+ theta2 <- z$theta[1:np]; sigma2 <- z$sigma; nit2 <- z$nit
+ #
+ # Covariance matrix of coefficient estimates
+ #
+ z <- kfedcc(wgt, z$rs, sigma=sigma2)
+ cov <- ktaskw(x, z$d, z$e, f=(sigma2^2)/n)
+ sd2 <- NULL; for (i in 1:np) { j <- i*(i+1)/2
+ sd2 <- c(sd2,cov$cov[j]) }
+ sd2 <- sqrt(sd2)
+ #.......................................................................
+ {
+ cat("rbmost: theta0 : ",round(theta0[1:np],3),"\n")
+ cat("rbmost: sigma0, nitw : ",round(sigma0,3),nitw0,"\n")
+ cat("Krasker-Welsch: theta2 : ",round(theta2,3),"\n")
+ cat("Krasker-Welsch: sd2 : ",round(sd2,3),"\n")
+ cat("Krasker-Welsch: sigma2, nitw, nit2 : ",round(sigma2,3),nitw,nit2,"\n")
+ }
+ #.......................................................................
+ list(theta0=theta0[1:np], sigma0=sigma0, nitw=nitw,
+ theta2=theta2, sigma2=sigma2, nit2=nit2, sd2=sd2)}
> #-----------------------------------------------------------------------
> # Read data; load defaults
> #
> z <- c( 8.2, 4, 23.005, 1, 7.6, 5, 23.873, 1,
+ 4.6, 0, 26.417, 1, 4.3, 1, 24.868, 1,
+ 5.9, 2, 29.895, 1, 5.0, 3, 24.200, 1,
+ 6.5, 4, 23.215, 1, 8.3, 5, 21.862, 1,
+ 10.1, 0, 22.274, 1, 13.2, 1, 23.830, 1,
+ 12.6, 2, 25.144, 1, 10.4, 3, 22.430, 1,
+ 10.8, 4, 21.785, 1, 13.1, 5, 22.380, 1,
+ 13.3, 0, 23.927, 1, 10.4, 1, 33.443, 1,
+ 10.5, 2, 24.859, 1, 7.7, 3, 22.686, 1,
+ 10.0, 0, 21.789, 1, 12.0, 1, 22.041, 1,
+ 12.1, 4, 21.033, 1, 13.6, 5, 21.005, 1,
+ 15.0, 0, 25.865, 1, 13.5, 1, 26.290, 1,
+ 11.5, 2, 22.932, 1, 12.0, 3, 21.313, 1,
+ 13.0, 4, 20.769, 1, 14.1, 5, 21.393, 1)
> x <- matrix(z, ncol=4, byrow=TRUE)
> y <- c( 7.6, 7.7, 4.3, 5.9, 5.0, 6.5, 8.3, 8.2, 13.2, 12.6,
+ 10.4, 10.8, 13.1, 12.3, 10.4, 10.5, 7.7, 9.5, 12.0, 12.6,
+ 13.6, 14.1, 13.5, 11.5, 12.0, 13.0, 14.1, 15.1)
> dfvals()
NULL
> dfcomn(xk=4)
$ipsi
[1] -9
$iucv
[1] -1
$iwww
[1] -1
> cat("Results\n")
Results
> z1 <- Mal.Std(x, y)
rbmost: theta0 : 0.674 -0.171 -0.678 20.043
rbmost: sigma0, nitw : 0.846 20
Mallows-Std: theta1 : 0.721 -0.174 -0.654 18.868
Mallows-Std: sd1 : 0.052 0.145 0.166 4.285
Mallows-Std: sigma1, nitw, nit1 : 0.774 38 7
> z2 <- Kra.Wel(x, y)
rbmost: theta0 : 0.674 -0.171 -0.678 20.043
rbmost: sigma0, nitw : 0.846 20
Krasker-Welsch: theta2 : 0.68 -0.177 -0.679 19.974
Krasker-Welsch: sd2 : 0.047 0.067 0.141 3.664
Krasker-Welsch: sigma2, nitw, nit2 : 0.732 21 4
>
>
> #
> # ---- Examples of Chapter 4: Covariance matrix of the coefficient estimates ------
> #
>
> #
> # Read x[1:4] and then set x[,4] <- 1
> #
> z <- c(80, 27, 89, 1, 80, 27, 88, 1, 75, 25, 90, 1,
+ 62, 24, 87, 1, 62, 22, 87, 1, 62, 23, 87, 1,
+ 62, 24, 93, 1, 62, 24, 93, 1, 58, 23, 87, 1,
+ 58, 18, 80, 1, 58, 18, 89, 1, 58, 17, 88, 1,
+ 58, 18, 82, 1, 58, 19, 93, 1, 50, 18, 89, 1,
+ 50, 18, 86, 1, 50, 19, 72, 1, 50, 19, 79, 1,
+ 50, 20, 80, 1, 56, 20, 82, 1, 70, 20, 91, 1)
> x <- matrix(z, ncol=4, byrow=TRUE)
> n <- 21; np <- 4; ncov <- np*(np+1)/2
> dfvals()
NULL
> # Cov. matrix of Huber-type estimates
> dfrpar(x, "huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> s <- liepsh()
> epsi2 <- s$epsi2; epsip <- s$epsip
> z <- rimtrf(x)
> xt <- z$x; sg <- z$sg; ip <- z$ip
> zc <- kiascv(xt, fu=epsi2/epsip^2, fb=0.)
> covi <- zc$cov # Can be used in ryhalg with ic=0
> zc <- kfascv(xt, covi, f=1, sg=sg, ip=ip)
> covf <- zc$cov
> #.......................................................................
> str <- rep(" ", ncov); str[cumsum(1:np)] <- "\n"
> {
+ cat("covf:\n")
+ cat(round(covf,6),sep=str)
+ }
covf:
0.00182
-0.003653 0.013553
-0.000715 1e-06 0.002444
0.028778 -0.065222 -0.167742 14.16076
>
>
> #
> # ---- Examples of Chapter 5: Asymptotic relative efficiency ------
> #
> #-----------------------------------------------------------------------
> # Huber
> #
> dfcomn(ipsi=1,c=1.345,d=1.345)
$ipsi
[1] 1
$iucv
[1] -1
$iwww
[1] -1
> .dFvPut(1,"ite")
NULL
> z <- airef0(mu=3, ialfa=1, sigmx=1)
> #.......................................................................
> {
+ cat(" airef0 : Huber\n reff, alfa, beta, nit: ")
+ cat(round(c(z$reff,z$alfa,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airef0 : Huber
reff, alfa, beta, nit: 0.95, 0, 0, 0
> #-----------------------------------------------------------------------
> # Schweppe: Krasker-Welsch
> #
> dfcomn(ipsi=1,c=3.76,iucv=3,ckw=3.76,iwww=1)
$ipsi
[1] 1
$iucv
[1] 3
$iwww
[1] 1
> .dFvPut(3,"ite")
NULL
> z <- airef0(mu=3, ialfa=1, sigmx=1)
> #.......................................................................
> {
+ cat(" airef0 : Krasker-Welsch\n reff, alfa, beta, nit: ")
+ cat(round(c(z$reff,z$alfa,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airef0 : Krasker-Welsch
reff, alfa, beta, nit: 0.95, 1.091, 1.17, 6
> #-----------------------------------------------------------------------
> # Mallows-Standard
> #
> dfcomn(ipsi=1,c=1.5,iucv=1,a2=0,b2=6.16,iww=3)
$ipsi
[1] 1
$iucv
[1] 1
$iwww
[1] 3
> .dFvPut(2,"ite")
NULL
> z <- airef0(mu=3, ialfa=1, sigmx=1)
> #.......................................................................
> {
+ cat(" airef0 : Mallows-Std \n reff, alfa, beta, nit: ")
+ cat(round(c(z$reff,z$alfa,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airef0 : Mallows-Std
reff, alfa, beta, nit: 0.95, 1.031, 1.09, 5
> #=======================================================================
> z <- c(1, 0, 0,
+ 1, 0, 0,
+ 1, 0, 0,
+ 1, 0, 0,
+ 0, 1, 0,
+ 0, 1, 0,
+ 0, 1, 0,
+ 0, 1, 0,
+ 0, 0, 1,
+ 0, 0, 1,
+ 0, 0, 1,
+ 0, 0, 1)
> tt <- matrix(z,ncol=3,byrow=TRUE)
> n <- nrow(tt); mu <- 2
> nu <- ncol(tt)
>
> #-----------------------------------------------------------------------
> # Huber
> #
> dfrpar(tt,"Huber")
$itypw
[1] 0
$itype
[1] 1
$isigma
[1] 1
> z <- airefq(tt, mu=mu, sigmx=1)
> #.......................................................................
> {
+ cat(" airefq : Huber\n reff, beta, nit: ")
+ cat(round(c(z$reff,z$beta,z$nit),3),sep=c(", ",", ",", ","\n"))
+ }
airefq : Huber
reff, beta, nit: 0.95, 0, 0
> #-----------------------------------------------------------------------
> # Krasker-Welsch
> #
> dfrpar(tt,"kra-wel",upar=3.755)
$itypw
[1] 1
$itype
[1] 3
$isigma
[1] 2
> z <- airefq(tt, mu=mu, sigmx=1,init=1)
Flavor: r-devel-linux-x86_64-fedora-clang