File: psi_functions.R

package info (click to toggle)
robustbase 0.99-7-1
  • links: PTS, VCS
  • area: main
  • in suites: sid
  • size: 4,596 kB
  • sloc: fortran: 3,245; ansic: 3,243; sh: 15; makefile: 2
file content (122 lines) | stat: -rw-r--r-- 4,285 bytes parent folder | download | duplicates (3)
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
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
### R code from vignette source 'psi_functions.Rnw'

###################################################
### code chunk number 1: init
###################################################
# set margins for plots
options(SweaveHooks=list(fig=function() par(mar=c(3,3,1.4,0.7),
                         mgp=c(1.5, 0.5, 0))))
## x axis for plots:
x. <- seq(-5, 10, length.out = 1501)
require(robustbase)


###################################################
### code chunk number 2: source-p-psiFun
###################################################
source(system.file("xtraR/plot-psiFun.R", package = "robustbase", mustWork=TRUE))


###################################################
### code chunk number 3: Huber
###################################################
getOption("SweaveHooks")[["fig"]]()
plot(huberPsi, x., ylim=c(-1.4, 5), leg.loc="topright", main=FALSE)


###################################################
### code chunk number 4: lmrob-psi
###################################################
names(.Mpsi.tuning.defaults)


###################################################
### code chunk number 5: tuning-defaults
###################################################
print(c(k.M = .Mpsi.tuning.default("bisquare"),
        k.S = .Mchi.tuning.default("bisquare")), digits = 10)


###################################################
### code chunk number 6: note-p-psiFun (eval = FALSE)
###################################################
getOption("SweaveHooks")[["fig"]]()
## source(system.file("xtraR/plot-psiFun.R", package = "robustbase", mustWork=TRUE))


###################################################
### code chunk number 7: bisquare
###################################################
getOption("SweaveHooks")[["fig"]]()
p.psiFun(x., "biweight", par = 4.685)


###################################################
### code chunk number 8: Hampel
###################################################
getOption("SweaveHooks")[["fig"]]()
## see also hampelPsi
p.psiFun(x., "Hampel", par = ## Default, but rounded:
                             round(c(1.5, 3.5, 8) * 0.9016085, 1))


###################################################
### code chunk number 9: GGW-const
###################################################
cT <- rbind(cc1 = .psi.ggw.findc(ms = -0.5, b = 1.5, eff = 0.95        ),
            cc2 = .psi.ggw.findc(ms = -0.5, b = 1.5,          bp = 0.50)); cT


###################################################
### code chunk number 10: rhoInf-ggw
###################################################
ipsi.ggw <- .psi2ipsi("GGW") # = 5
ccc <- c(0, cT[1, 2:4], 1)
integrate(.Mpsi, 0, Inf, ccc=ccc, ipsi=ipsi.ggw)$value # = rho(Inf)


###################################################
### code chunk number 11: GGW
###################################################
getOption("SweaveHooks")[["fig"]]()
p.psiFun(x., "GGW", par = c(-.5, 1, .95, NA))


###################################################
### code chunk number 12: lqq-const
###################################################
cT <- rbind(cc1 = .psi.lqq.findc(ms= -0.5, b.c = 1.5, eff=0.95, bp=NA ),
            cc2 = .psi.lqq.findc(ms= -0.5, b.c = 1.5, eff=NA , bp=0.50))
colnames(cT) <- c("b", "c", "s"); cT


###################################################
### code chunk number 13: LQQ
###################################################
getOption("SweaveHooks")[["fig"]]()
p.psiFun(x., "LQQ", par = c(-.5,1.5,.95,NA))


###################################################
### code chunk number 14: optimal
###################################################
getOption("SweaveHooks")[["fig"]]()
p.psiFun(x., "optimal", par = 1.06, leg.loc="bottomright")


###################################################
### code chunk number 15: Welsh-GGW
###################################################
ccc <- c(0, a = 2.11^2, b = 2, c = 0, 1)
(ccc[5] <- integrate(.Mpsi, 0, Inf, ccc=ccc, ipsi = 5)$value) # = rho(Inf)
stopifnot(all.equal(Mpsi(x., ccc,  "GGW"),   ## psi[ GGW ](x; a=k^2, b=2, c=0) ==
                    Mpsi(x., 2.11, "Welsh")))## psi[Welsh](x; k)


###################################################
### code chunk number 16: Welsh
###################################################
getOption("SweaveHooks")[["fig"]]()
p.psiFun(x., "Welsh", par = 2.11)