[Robast-commits] r237 - pkg/RobRex/inst/scripts

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Sat Nov 29 18:35:46 CET 2008 

Author: stamats
Date: 2008-11-29 18:35:46 +0100 (Sat, 29 Nov 2008)
New Revision: 237

Added:
   pkg/RobRex/inst/scripts/ExamplesEstimation.R
   pkg/RobRex/inst/scripts/ExamplesIC.R
Removed:
   pkg/RobRex/inst/scripts/example.R
Log:
added new script, changed name of script

Added: pkg/RobRex/inst/scripts/ExamplesEstimation.R
===================================================================
--- pkg/RobRex/inst/scripts/ExamplesEstimation.R	                        (rev 0)
+++ pkg/RobRex/inst/scripts/ExamplesEstimation.R	2008-11-29 17:35:46 UTC (rev 237)
@@ -0,0 +1,91 @@
+###############################################################################
+## Optimally Robust Estimation
+###############################################################################
+
+###############################################################################
+## Example 1
+###############################################################################
+require(MASS)
+require(robustbase)
+library(ISwR)
+data(thuesen)
+attach(thuesen)
+
+## LS-estimator
+fit.LS <- lm(short.velocity ~ blood.glucose)
+## M-estimator
+fit.M <- rlm(short.velocity ~ blood.glucose)
+## MM-estimator
+fit.MM <- lmrob(short.velocity ~ blood.glucose)
+## LTS-estimator
+fit.LTS <- ltsReg(short.velocity ~ blood.glucose)
+
+## AL-estimators
+## regressor distribution: design measure
+require(RobRex)
+K <- DiscreteMVDistribution(cbind(1,blood.glucose))
+## ALc-estimator: conditional neighborhood; i.e., only y is contaminated
+## about 10 sec. on my system
+system.time(IC1 <- rgsOptIC.ALc(r = 0.5, K = K, theta = fit.M$coeff, scale = fit.M$s))
+ALc1 <- oneStepEstimator(cbind(1,as.matrix(thuesen)), IC1, c(fit.M$coeff, fit.M$s))
+## AL-estimator: unconditional neighborhood; i.e., x and y are contaminated
+## about 85 sec. on my system
+system.time(IC2 <- rgsOptIC.AL(r = 0.5, K = K, theta = fit.MM$coeff, scale = fit.MM$s))
+AL1 <- oneStepEstimator(cbind(1,as.matrix(thuesen)), IC2, c(fit.MM$coeff, fit.MM$s))
+
+## Plot
+require(RColorBrewer)
+myCol <- brewer.pal(6, "Set1")
+plot(short.velocity ~ blood.glucose, ylab = "fasting blood glucose [mmol/l]", 
+     xlab = "mean circumferential shortening velocity [%/s]", 
+     main = "Ventricular shortening velocity", pch = 20)
+abline(fit.LS, lwd = 2, col = myCol[1])
+abline(fit.M, lwd = 2, col = myCol[2])
+abline(fit.MM, lwd = 2, col = myCol[3])
+abline(fit.LTS, lwd = 2, col = myCol[4])
+lines(c(1, c(blood.glucose,25)), ALc1[1] + ALc1[2]*c(1,c(blood.glucose,25)), 
+      col = myCol[5], lwd = 2)
+lines(c(1, c(blood.glucose,25)), AL1[1] + AL1[2]*c(1,c(blood.glucose,25)), 
+      col = myCol[1], lwd = 2)
+legend("topleft", legend = c("LS", "M", "MM", "LTS", "ALc", "AL"), 
+       fill = myCol, ncol = 2)
+detach(thuesen)
+
+###############################################################################
+## Example 2
+###############################################################################
+data(phones)
+attach(phones)
+
+## LS estimator
+fit2.LS <- lm(calls ~ year)
+## M estimator
+fit2.M <- rlm(calls ~ year, maxit = 50)
+## MM estimator
+fit2.MM <- lmrob(calls ~ year)
+## LTS estimator
+fit2.LTS <- ltsReg(calls ~ year)
+
+## AL estimators
+## regressor distribution: design measure
+K <- DiscreteMVDistribution(cbind(1,year))
+## ALc estimator: only y contaminated
+system.time(IC1 <- rgsOptIC.ALc(r = 0.5, K = K, theta = fit2.M$coeff, scale = fit2.M$s))
+## takes about 9 sec. on my system
+ALc2 <- oneStepEstimator(cbind(1,year,calls), IC1, c(fit2.M$coeff, fit2.M$s))
+## AL estimator: x and y contaminated
+system.time(IC2 <- rgsOptIC.AL(r = 0.5, K = K, theta = fit2.MM$coeff, scale = fit2.MM$s))
+## takes about 80 sec. on my system
+AL2 <- oneStepEstimator(cbind(1,year,calls), IC2, c(fit2.MM$coeff, fit2.MM$s))
+
+## Plot
+plot(calls ~ year, ylab = "phone calls [Mio.]", xlab = "year", 
+     main = "Belgium Phone Calls 1950-1973", pch = 20)
+abline(fit2.LS, lwd = 2, col = "black")
+abline(fit2.M, lwd = 2, col = 1)
+abline(fit2.MM, lwd = 2, col = 3)
+abline(fit2.LTS, lwd = 2, col = 2)
+lines(c(1, c(year,75)), ALc2[1] + ALc2[2]*c(1,c(year,75)), col = 4, lwd = 2)
+lines(c(1, c(year,75)), AL2[1] + AL2[2]*c(1,c(year,75)), col = 5, lwd = 2)
+legend("topleft", legend = c("LS", "M", "MM", "LTS", "ALc", "AL"), 
+       fill = c("black", 1, 3, 2, 4, 5), ncol = 2)


Property changes on: pkg/RobRex/inst/scripts/ExamplesEstimation.R
___________________________________________________________________
Name: svn:mergeinfo
   + 

Copied: pkg/RobRex/inst/scripts/ExamplesIC.R (from rev 230, pkg/RobRex/inst/scripts/example.R)
===================================================================
--- pkg/RobRex/inst/scripts/ExamplesIC.R	                        (rev 0)
+++ pkg/RobRex/inst/scripts/ExamplesIC.R	2008-11-29 17:35:46 UTC (rev 237)
@@ -0,0 +1,156 @@
+###############################################################################
+## Optimal Influence Curves
+###############################################################################
+
+library(RobRex)
+options("newDevice"=TRUE)
+
+###############################################################################
+## Example 1 (1-dim., discrete Regressor)
+###############################################################################
+K1 <- DiscreteDistribution(supp = 1:5)
+
+# AL-estimator
+system.time(IC.AL1 <- rgsOptIC.AL(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
+distrExOptions(ErelativeTolerance, 1e-10)
+checkIC(IC.AL1)
+distrExOptions(ErelativeTolerance, .Machine$double.eps^0.25)
+Risks(IC.AL1)
+
+# M-estimator
+system.time(IC.M1 <- rgsOptIC.M(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
+checkIC(IC.M1)
+Risks(IC.M1)
+
+# MK-estimator
+system.time(IC.MK1 <- rgsOptIC.MK(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
+checkIC(IC.MK1)
+Risks(IC.MK1)
+
+# ALc-estimator
+system.time(IC.ALc1 <- rgsOptIC.ALc(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
+checkIC(IC.ALc1)
+Risks(IC.ALc1)
+
+# Mc-estimator
+system.time(IC.Mc1 <- rgsOptIC.Mc(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
+checkIC(IC.Mc1)
+Risks(IC.Mc1)
+
+# ALs-estimator
+system.time(IC.ALs1 <- rgsOptIC.ALs(r = 0.1, scale = 2, K = K1, check = TRUE), gcFirst = TRUE)
+checkIC(IC.ALs1)
+Risks(IC.ALs1)
+
+# Ms-estimator
+system.time(IC.Ms1 <- rgsOptIC.Ms(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
+checkIC(IC.Ms1)
+Risks(IC.Ms1)
+
+# BM-estimator
+system.time(IC.BM1 <- rgsOptIC.BM(r = 0.1, K = K1), gcFirst = TRUE)
+checkIC(IC.BM1)
+Risks(IC.BM1)
+
+# asymptotic MSEs
+Risks(IC.AL1)$asMSE
+Risks(IC.M1)$asMSE
+Risks(IC.MK1)$asMSE
+Risks(IC.ALc1)$asMSE
+Risks(IC.Mc1)$asMSE
+Risks(IC.ALs1)$asMSE
+Risks(IC.Ms1)$asMSE
+Risks(IC.BM1)$asMSE
+
+
+###############################################################################
+## Example 2 (1-dim., abs. cont. Regressor)
+###############################################################################
+K2 <- Unif(Min = 0, Max = 1)
+
+# AL-estimator
+IC.AL2 <- rgsOptIC.AL(r = 0.1, K = K2, check = TRUE)
+#checkIC(IC.AL2, eval(CallL2Fam(IC.AL2)))
+Risks(IC.AL2)
+
+# M-estimator
+IC.M2 <- rgsOptIC.M(r = 0.1, K = K2, check = TRUE)
+#checkIC(IC.M2, eval(CallL2Fam(IC.M2)))
+Risks(IC.M2)
+
+# MK-estimator
+IC.MK2 <- rgsOptIC.MK(r = 0.1, K = K2, check = TRUE)
+#checkIC(IC.MK2, eval(CallL2Fam(IC.MK2)))
+Risks(IC.MK2)
+
+# ALc-estimator
+# only implemented for discrete distributions
+
+# Mc-estimator
+# only implemented for 1-dim. discrete distributions
+
+# ALs-estimator
+IC.ALs2 <- rgsOptIC.ALs(r = 0.1, K = K2, check = TRUE)
+#checkIC(IC.ALs2, eval(CallL2Fam(IC.ALs2)))
+Risks(IC.ALs2)
+
+# Ms-estimator
+# only implemented for 1-dim. discrete distributions
+
+# BM-estimator
+# only defined for discrete distributions
+
+
+# asymptotic MSEs
+Risks(IC.AL2)$asMSE
+Risks(IC.M2)$asMSE
+Risks(IC.MK2)$asMSE
+Risks(IC.ALs2)$asMSE
+
+
+###############################################################################
+## Example 3 (2-dim., discrete Regressor)
+###############################################################################
+K3 <- DiscreteMVDistribution(supp = matrix(c(0,1,1,1,0,2,2,0,2,1), ncol=2, byrow = TRUE))
+
+# AL-estimator
+IC.AL3 <- rgsOptIC.AL(r = 0.1, K = K3, check = TRUE)
+checkIC(IC.AL3)
+Risks(IC.AL3)
+
+# M-estimator
+IC.M3 <- rgsOptIC.M(r = 0.1, K = K3, check = TRUE)
+checkIC(IC.M3)
+Risks(IC.M3)
+
+# MK-estimator
+IC.MK3 <- rgsOptIC.MK(r = 0.1, K = K3, check = TRUE)
+checkIC(IC.MK3)
+Risks(IC.MK3)
+
+# ALc-estimator
+IC.ALc3 <- rgsOptIC.ALc(r = 0.1, K = K3, check = TRUE)
+checkIC(IC.ALc3)
+Risks(IC.ALc3)
+
+# Mc-estimator
+# only implemented for 1-dim. discrete distributions
+
+# ALs-estimator
+IC.ALs3 <- rgsOptIC.ALs(r = 0.1, K = K3, check = TRUE)
+checkIC(IC.ALs3)
+Risks(IC.ALs3)
+
+# Ms-estimator
+# only implemented for 1-dim. discrete distributions
+
+# BM-estimator
+# only defined for discrete distributions
+# and only implemented for 1-dim. discrete distributions
+
+# asymptotic MSEs
+Risks(IC.AL3)$asMSE
+Risks(IC.M3)$asMSE
+Risks(IC.MK3)$asMSE
+Risks(IC.ALc3)$asMSE
+Risks(IC.ALs3)$asMSE


Property changes on: pkg/RobRex/inst/scripts/ExamplesIC.R
___________________________________________________________________
Name: svn:mergeinfo
   + 

Deleted: pkg/RobRex/inst/scripts/example.R
===================================================================
--- pkg/RobRex/inst/scripts/example.R	2008-11-29 13:40:11 UTC (rev 236)
+++ pkg/RobRex/inst/scripts/example.R	2008-11-29 17:35:46 UTC (rev 237)
@@ -1,163 +0,0 @@
-library(RobRex)
-options("newDevice"=TRUE)
-
-###############################################################################
-## start of tests
-###############################################################################
-
-###############################################################################
-## Example 1 (1-dim., discrete Regressor)
-###############################################################################
-K1 <- DiscreteDistribution(supp = 1:5)
-
-# AL-estimator
-system.time(IC.AL1 <- rgsOptIC.AL(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
-distrExOptions(ErelativeTolerance, 1e-10)
-checkIC(IC.AL1)
-distrExOptions(ErelativeTolerance, .Machine$double.eps^0.25)
-Risks(IC.AL1)
-
-# M-estimator
-system.time(IC.M1 <- rgsOptIC.M(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
-checkIC(IC.M1)
-Risks(IC.M1)
-
-# MK-estimator
-system.time(IC.MK1 <- rgsOptIC.MK(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
-checkIC(IC.MK1)
-Risks(IC.MK1)
-
-# ALc-estimator
-system.time(IC.ALc1 <- rgsOptIC.ALc(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
-checkIC(IC.ALc1)
-Risks(IC.ALc1)
-
-# Mc-estimator
-system.time(IC.Mc1 <- rgsOptIC.Mc(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
-checkIC(IC.Mc1)
-Risks(IC.Mc1)
-
-# ALs-estimator
-system.time(IC.ALs1 <- rgsOptIC.ALs(r = 0.1, scale = 2, K = K1, check = TRUE), gcFirst = TRUE)
-checkIC(IC.ALs1)
-Risks(IC.ALs1)
-
-# Ms-estimator
-system.time(IC.Ms1 <- rgsOptIC.Ms(r = 0.1, K = K1, check = TRUE), gcFirst = TRUE)
-checkIC(IC.Ms1)
-Risks(IC.Ms1)
-
-# BM-estimator
-system.time(IC.BM1 <- rgsOptIC.BM(r = 0.1, K = K1), gcFirst = TRUE)
-checkIC(IC.BM1)
-Risks(IC.BM1)
-
-# asymptotic MSEs
-Risks(IC.AL1)$asMSE
-Risks(IC.M1)$asMSE
-Risks(IC.MK1)$asMSE
-Risks(IC.ALc1)$asMSE
-Risks(IC.Mc1)$asMSE
-Risks(IC.ALs1)$asMSE
-Risks(IC.Ms1)$asMSE
-Risks(IC.BM1)$asMSE
-
-
-###############################################################################
-## Example 2 (1-dim., abs. cont. Regressor)
-###############################################################################
-K2 <- Unif(Min = 0, Max = 1)
-
-# AL-estimator
-IC.AL2 <- rgsOptIC.AL(r = 0.1, K = K2, check = TRUE)
-#checkIC(IC.AL2, eval(CallL2Fam(IC.AL2)))
-Risks(IC.AL2)
-
-# M-estimator
-IC.M2 <- rgsOptIC.M(r = 0.1, K = K2, check = TRUE)
-#checkIC(IC.M2, eval(CallL2Fam(IC.M2)))
-Risks(IC.M2)
-
-# MK-estimator
-IC.MK2 <- rgsOptIC.MK(r = 0.1, K = K2, check = TRUE)
-#checkIC(IC.MK2, eval(CallL2Fam(IC.MK2)))
-Risks(IC.MK2)
-
-# ALc-estimator
-# only implemented for discrete distributions
-
-# Mc-estimator
-# only implemented for 1-dim. discrete distributions
-
-# ALs-estimator
-IC.ALs2 <- rgsOptIC.ALs(r = 0.1, K = K2, check = TRUE)
-#checkIC(IC.ALs2, eval(CallL2Fam(IC.ALs2)))
-Risks(IC.ALs2)
-
-# Ms-estimator
-# only implemented for 1-dim. discrete distributions
-
-# BM-estimator
-# only defined for discrete distributions
-
-
-# asymptotic MSEs
-Risks(IC.AL2)$asMSE
-Risks(IC.M2)$asMSE
-Risks(IC.MK2)$asMSE
-Risks(IC.ALs2)$asMSE
-
-
-###############################################################################
-## Example 3 (2-dim., discrete Regressor)
-###############################################################################
-K3 <- DiscreteMVDistribution(supp = matrix(c(0,1,1,1,0,2,2,0,2,1), ncol=2, byrow = TRUE))
-
-# AL-estimator
-IC.AL3 <- rgsOptIC.AL(r = 0.1, K = K3, check = TRUE)
-checkIC(IC.AL3)
-Risks(IC.AL3)
-
-# M-estimator
-IC.M3 <- rgsOptIC.M(r = 0.1, K = K3, check = TRUE)
-checkIC(IC.M3)
-Risks(IC.M3)
-
-# MK-estimator
-IC.MK3 <- rgsOptIC.MK(r = 0.1, K = K3, check = TRUE)
-checkIC(IC.MK3)
-Risks(IC.MK3)
-
-# ALc-estimator
-IC.ALc3 <- rgsOptIC.ALc(r = 0.1, K = K3, check = TRUE)
-checkIC(IC.ALc3)
-Risks(IC.ALc3)
-
-# Mc-estimator
-# only implemented for 1-dim. discrete distributions
-
-# ALs-estimator
-IC.ALs3 <- rgsOptIC.ALs(r = 0.1, K = K3, check = TRUE)
-checkIC(IC.ALs3)
-Risks(IC.ALs3)
-
-# Ms-estimator
-# only implemented for 1-dim. discrete distributions
-
-# BM-estimator
-# only defined for discrete distributions
-# and only implemented for 1-dim. discrete distributions
-
-# asymptotic MSEs
-Risks(IC.AL3)$asMSE
-Risks(IC.M3)$asMSE
-Risks(IC.MK3)$asMSE
-Risks(IC.ALc3)$asMSE
-Risks(IC.ALs3)$asMSE
-
-
-###############################################################################
-## end of tests
-###############################################################################
-
-q("no")

More information about the Robast-commits mailing list