R help - May 2018 - Bootstrap and average median squared error

I forgot, you should also set.seed() before calling boot() to make the 
results reproducible.

Rui Barradas

On 5/22/2018 10:00 AM, Rui Barradas wrote:
> Hello,
> 
> If you want to bootstrap a statistic, I suggest you use base package boot.
> You would need the data in a data.frame, see how you could do it.
> 
> 
> library(boot)
> 
> bootMedianSE <- function(data, indices){
>  ??? d <- data[indices, ]
>  ??? fit <- rq(crp ~ bmi + glucose, tau = 0.5, data = d)
>  ??? ypred <- predict(fit)
>  ??? y <- d$crp
>  ??? median(y - ypred)^2
> }
> 
> dat <- data.frame(crp, bmi, glucose)
> nboot <- 100
> 
> medse <- boot(dat, bootMedianSE, R = nboot)
> 
> medse$t0
> mean(medse$t)??? # This is the value you want
> 
> 
> Hope this helps,
> 
> Rui Barradas
> 
> 
> 
> On 5/22/2018 12:19 AM, varin sacha via R-help wrote:
>> Dear R-experts,
>>
>> I am trying to bootstrap (and average) the median squared error 
>> evaluation metric for a robust regression. I can't get it. What is 
>> going wrong ?
>>
>> Here is the reproducible example.
>>
>> #############################
>> install.packages( "quantreg" )
>> library(quantreg)
>>
>> crp <-c(12,14,13,24,25,34,45,56,25,34,47,44,35,24,53,44,55,46,36,67)
>> bmi <-c(34,32,12,76,54,34,21,18,92,32,11,13,45,46,56,57,67,87,12,13)
>> glucose
<-c(23,54,11,12,13,21,32,12,45,54,65,87,21,23,12,12,23,23,43,54)
>>
>> # Create a list to store the results
>> lst<-list()
>>
>> # Numbers of bootstrap samples
>> nboot=100
>> bootstrap.MedAESQ =rep(NA,nboot)
>>
>> for(i in 1?:nboot)
>> {
>>
>> fit <- rq( crp ~ bmi+glucose, tau = 0.5)
>> ypred=predict(fit)
>> y=new$crp
>> bootstrap.MedAESQ [i]=median(y-ypred)^2
>> lst[i]<-bootstrap.MedAESQ
>>
>> }
>> mean(unlist(lst))
>> ###################################
>>
>>
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> 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.
>>
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.

On 5/22/2018 2:32 AM, Rui Barradas wrote:
> bootMedianSE <- function(data, indices){
>  ???? d <- data[indices, ]
>  ???? fit <- rq(crp ~ bmi + glucose, tau = 0.5, data = d)
>  ???? ypred <- predict(fit)
>  ???? y <- d$crp
>  ???? median(y - ypred)^2
> }
since the OP is looking for the "median squared error", shouldn't
the
final line of the function be

         median((y - ypred)^2)


Dan

-- 
Daniel Nordlund
Port Townsend, WA  USA

Hello,

Right!
I copied from the OP's question without thinking about it.
Corrected would be


bootMedianSE <- function(data, indices){
     d <- data[indices, ]
     fit <- rq(crp ~ bmi + glucose, tau = 0.5, data = d)
     ypred <- predict(fit)
     y <- d$crp
     median((y - ypred)^2)
}

Sorry,

rui Barradas

On 5/22/2018 11:32 AM, Daniel Nordlund wrote:
> On 5/22/2018 2:32 AM, Rui Barradas wrote:
>> bootMedianSE <- function(data, indices){
>> ????? d <- data[indices, ]
>> ????? fit <- rq(crp ~ bmi + glucose, tau = 0.5, data = d)
>> ????? ypred <- predict(fit)
>> ????? y <- d$crp
>> ????? median(y - ypred)^2
>> }
> 
> since the OP is looking for the "median squared error",
shouldn't the
> final line of the function be
> 
>  ??????? median((y - ypred)^2)
> 
> 
> Dan
>