Selecting Bootstrap Statistics in the boot package

classic Classic list List threaded Threaded
3 messages Options
Reply | Threaded
Open this post in threaded view
|

Selecting Bootstrap Statistics in the boot package

Jorge I Velez
Dear all,

Let's say I have the following:

# Loading the boot package
# install.packages(boot)
library(boot)

# Generating data
set.seed(123)
x <- rnorm(100)

# Bootstrap for the sample mean
bmean <- boot(x, function(x,d) mean(x[d]), R=1000)
bmean
#
#ORDINARY NONPARAMETRIC BOOTSTRAP
#
#
#Call:
#boot(data = x, statistic = function(x, d) mean(x[d]), R = 1000)
#
#
#Bootstrap Statistics :
#     original      bias    std. error
#t1* 0.0904059 0.004641537  0.09239923


and I would like to get just this:


Bootstrap Statistics :
     original      bias    std. error
t1* 0.0904059 0.004641537  0.09239923


How can I do that?  I'm running R 2.8.1 Patched on XP SP2. Here is my
sessionInfo():

> sessionInfo()
R version 2.8.1 Patched (2009-03-09 r48093)
i386-pc-mingw32

locale:
LC_COLLATE=English_United States.1252;LC_CTYPE=English_United States.1252;
LC_MONETARY=English_United States.1252;LC_NUMERIC=C;LC_TIME=English_United
States.1252

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base

other attached packages:
[1] RODBC_1.2-5 boot_1.2-35


As always, thank you so much for all your help.

Regards,

Jorge

Jorge Ivan Velez
jorgeivanvelez [at] gmail [dot] com

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Selecting Bootstrap Statistics in the boot package

David Winsemius
Greetings Jorge;

There certainly did not seem to be something readily apparent with  
str(bmean), so the next logical place to look would be for a print  
method. If you look at print.boot with:

getAywhere(print.boot)

  ... you see that the first portion sets up an internal structure  
called "op" and the last of the function prints this with optional  
headers. If you strip out the sections with headers you get:

lim.boot <- function (x, digits = getOption("digits"), index =  
1L:ncol(boot.out$t),
     ...)
{
     boot.out <- x
     sim <- boot.out$sim
     cl <- boot.out$call
     t <- matrix(boot.out$t[, index], nrow = nrow(boot.out$t))
     allNA <- apply(t, 2L, function(t) all(is.na(t)))
     ind1 <- index[allNA]
     index <- index[!allNA]
     t <- matrix(t[, !allNA], nrow = nrow(t))
     rn <- paste("t", index, "*", sep = "")
     if (length(index) == 0L)
         op <- NULL
     else if (is.null(t0 <- boot.out$t0)) {
         if (is.null(boot.out$call$weights))
             op <- cbind(apply(t, 2L, mean, na.rm = TRUE), sqrt(apply(t,
                 2L, function(t.st) var(t.st[!is.na(t.st)]))))
         else {
             op <- NULL
             for (i in index) op <- rbind(op, imp.moments(boot.out,
                 index = i)$rat)
             op[, 2L] <- sqrt(op[, 2])
         }
         dimnames(op) <- list(rn, c("mean", "std. error"))
     }
     else {
         t0 <- boot.out$t0[index]
         if (is.null(boot.out$call$weights)) {
             op <- cbind(t0, apply(t, 2L, mean, na.rm = TRUE) -
                 t0, sqrt(apply(t, 2L, function(t.st) var(t.st[!
is.na(t.st)]))))
             dimnames(op) <- list(rn, c("original", " bias  ",
                 " std. error"))
         }
         else {
             op <- NULL
             for (i in index) op <- rbind(op, imp.moments(boot.out,
                 index = i)$rat)
             op <- cbind(t0, op[, 1L] - t0, sqrt(op[, 2L]), apply(t,
                 2L, mean, na.rm = TRUE))
             dimnames(op) <- list(rn, c("original", " bias  ",
                 " std. error", " mean(t*)"))
         }
     }

     cat("\n\nBootstrap Statistics :\n")
     if (!is.null(op))
         print(op, digits = digits)
     if (length(ind1) > 0L)
         for (j in ind1) cat(paste("WARNING: All values of t",
             j, "* are NA\n", sep = ""))
     invisible(boot.out)
}
#--------------------------------
lim.boot(bmean)


Bootstrap Statistics :
      original      bias    std. error
t1* 0.0904059 0.004641537  0.09239923

Probably too much in that function, but it does what was requested.

--
David Winsemius


On Mar 31, 2009, at 5:49 PM, Jorge Ivan Velez wrote:

>
> Let's say I have the following:
>
> # Loading the boot package
> # install.packages(boot)
> library(boot)
>
> # Generating data
> set.seed(123)
> x <- rnorm(100)
>
> # Bootstrap for the sample mean
> bmean <- boot(x, function(x,d) mean(x[d]), R=1000)
> bmean
> #
> #ORDINARY NONPARAMETRIC BOOTSTRAP
> #
> #
> #Call:
> #boot(data = x, statistic = function(x, d) mean(x[d]), R = 1000)
> #
> #
> #Bootstrap Statistics :
> #     original      bias    std. error
> #t1* 0.0904059 0.004641537  0.09239923
>
>
> and I would like to get just this:
>
>
> Bootstrap Statistics :
>     original      bias    std. error
> t1* 0.0904059 0.004641537  0.09239923
>
>
> How can I do that?  I'm running R 2.8.1 Patched on XP SP2. Here is my
> sessionInfo():
>
>
David Winsemius, MD
Heritage Laboratories
West Hartford, CT

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
Reply | Threaded
Open this post in threaded view
|

Re: Selecting Bootstrap Statistics in the boot package

Jorge I Velez
Many, many thanks to Phil Spector and David Winsemius for their helpful
advice regarding this.

Phil Spector suggested me to use

library(boot)
xx = getAnywhere('print.boot')
myprint = xx$objs[[1]]
myprint = edit(myprint,editor='notepad')
vi(myprint)

to edit the function deleting the original function from line 45 ( if
(cl[[1]] == "boot") {  )  all the way down to line 154
(dput(cl, control = NULL) ) and then use

myprint(bemean)

to get what I asked.

David Winsemius' function does also the trick and gets me what I asked for
in my post this afternoon. Thanks David.


Best regards,

Jorge



On Tue, Mar 31, 2009 at 10:29 PM, David Winsemius <[hidden email]>wrote:

> Greetings Jorge;
>
> There certainly did not seem to be something readily apparent with
> str(bmean), so the next logical place to look would be for a print method.
> If you look at print.boot with:
>
> getAywhere(print.boot)
>
>  ... you see that the first portion sets up an internal structure called
> "op" and the last of the function prints this with optional headers. If you
> strip out the sections with headers you get:
>
> lim.boot <- function (x, digits = getOption("digits"), index =
> 1L:ncol(boot.out$t),
>    ...)
> {
>    boot.out <- x
>    sim <- boot.out$sim
>    cl <- boot.out$call
>    t <- matrix(boot.out$t[, index], nrow = nrow(boot.out$t))
>    allNA <- apply(t, 2L, function(t) all(is.na(t)))
>    ind1 <- index[allNA]
>    index <- index[!allNA]
>    t <- matrix(t[, !allNA], nrow = nrow(t))
>    rn <- paste("t", index, "*", sep = "")
>    if (length(index) == 0L)
>        op <- NULL
>    else if (is.null(t0 <- boot.out$t0)) {
>        if (is.null(boot.out$call$weights))
>            op <- cbind(apply(t, 2L, mean, na.rm = TRUE), sqrt(apply(t,
>                2L, function(t.st) var(t.st[!is.na(t.st)]))))
>        else {
>            op <- NULL
>            for (i in index) op <- rbind(op, imp.moments(boot.out,
>                index = i)$rat)
>            op[, 2L] <- sqrt(op[, 2])
>        }
>        dimnames(op) <- list(rn, c("mean", "std. error"))
>    }
>    else {
>        t0 <- boot.out$t0[index]
>        if (is.null(boot.out$call$weights)) {
>            op <- cbind(t0, apply(t, 2L, mean, na.rm = TRUE) -
>                t0, sqrt(apply(t, 2L, function(t.st) var(t.st[!is.na(t.st
> )]))))
>            dimnames(op) <- list(rn, c("original", " bias  ",
>                " std. error"))
>        }
>        else {
>            op <- NULL
>            for (i in index) op <- rbind(op, imp.moments(boot.out,
>                index = i)$rat)
>            op <- cbind(t0, op[, 1L] - t0, sqrt(op[, 2L]), apply(t,
>                2L, mean, na.rm = TRUE))
>            dimnames(op) <- list(rn, c("original", " bias  ",
>                " std. error", " mean(t*)"))
>        }
>    }
>
>    cat("\n\nBootstrap Statistics :\n")
>    if (!is.null(op))
>        print(op, digits = digits)
>    if (length(ind1) > 0L)
>        for (j in ind1) cat(paste("WARNING: All values of t",
>            j, "* are NA\n", sep = ""))
>    invisible(boot.out)
> }
> #--------------------------------
> lim.boot(bmean)
>
>
> Bootstrap Statistics :
>     original      bias    std. error
> t1* 0.0904059 0.004641537  0.09239923
>
> Probably too much in that function, but it does what was requested.
>
> --
> David Winsemius
>
>
>
> On Mar 31, 2009, at 5:49 PM, Jorge Ivan Velez wrote:
>
>>
>> Let's say I have the following:
>>
>> # Loading the boot package
>> # install.packages(boot)
>> library(boot)
>>
>> # Generating data
>> set.seed(123)
>> x <- rnorm(100)
>>
>> # Bootstrap for the sample mean
>> bmean <- boot(x, function(x,d) mean(x[d]), R=1000)
>> bmean
>> #
>> #ORDINARY NONPARAMETRIC BOOTSTRAP
>> #
>> #
>> #Call:
>> #boot(data = x, statistic = function(x, d) mean(x[d]), R = 1000)
>> #
>> #
>> #Bootstrap Statistics :
>> #     original      bias    std. error
>> #t1* 0.0904059 0.004641537  0.09239923
>>
>>
>> and I would like to get just this:
>>
>>
>> Bootstrap Statistics :
>>    original      bias    std. error
>> t1* 0.0904059 0.004641537  0.09239923
>>
>>
>> How can I do that?  I'm running R 2.8.1 Patched on XP SP2. Here is my
>> sessionInfo():
>>
>>
>>  David Winsemius, MD
> Heritage Laboratories
> West Hartford, CT
>
>

        [[alternative HTML version deleted]]

______________________________________________
[hidden email] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.