Dear R-experts,
I am trying to get the bootstrapped confidence intervals of Mean squared error
(MSE) for a (G)AM model. I get an error message.
Here below the reproducible R code. Many thanks for your response.
####################
install.packages("ISLR")
library(ISLR)
install.packages("mgcv")
library(mgcv)
install.packages("boot")
library(boot)
#MSE calculation
n=dim(Wage)[1]
p=0.667
GAM1<-gam(wage ~education+s(age,bs="ps")+year,data=Wage)
sam=sample(1?:n,floor(p*n),replace=FALSE)
?
Training =Wage [sam,]
Testing = Wage [-sam,]
?
ypred=predict(GAM1,newdata=Testing)
y=Testing$wage
MSE = mean((y-ypred)^2)
# Bootstrap 95% CI for MSE
# function to obtain MSE
MSE <- function(formula, data, indices) {
??d <- data[indices,] # allows boot to select sample
??fit <- gam(formula, data=d)
??return(MSE(fit))
} # bootstrapping with 1000 replications
results <- boot(data=Wage, statistic=MSE,
?? R=1000, formula=gam(wage ~education+s(age,bs="ps")+year,data=Wage)
)
# get 95% confidence intervals
boot.ci(results, type="bca")
##########################
Hello,
There were several errors with your code. The following works but with
the other CI types.
library(ISLR)
library(mgcv)
library(boot)
# function to obtain MSE
MSE <- function(data, indices, formula) {
d <- data[indices, ] # allows boot to select sample
fit <- gam(formula, data = d)
ypred <- predict(fit)
mean((d[["wage"]] - ypred)^2)
}
data(Wage)
# Make the results reproducible
set.seed(1234)
# bootstrapping with 1000 replications
results <- boot(data = Wage, statistic = MSE,
R = 1000, formula = wage ~ education + s(age, bs =
"ps") + year))
# get 95% confidence intervals
# type = "bca" is throwing an error
ci.type <- c("norm","basic", "stud",
"perc")
boot.ci(results, type = ci.type)
Hope this helps,
Rui Barradas
?s 20:36 de 22/11/2018, varin sacha via R-help escreveu:> Dear R-experts,
>
> I am trying to get the bootstrapped confidence intervals of Mean squared
error (MSE) for a (G)AM model. I get an error message.
> Here below the reproducible R code. Many thanks for your response.
>
> ####################
>
>
> install.packages("ISLR")
>
> library(ISLR)
>
> install.packages("mgcv")
>
> library(mgcv)
>
> install.packages("boot")
>
> library(boot)
>
> #MSE calculation
>
> n=dim(Wage)[1]
>
> p=0.667
>
> GAM1<-gam(wage ~education+s(age,bs="ps")+year,data=Wage)
>
> sam=sample(1?:n,floor(p*n),replace=FALSE)
>
>
>
> Training =Wage [sam,]
>
> Testing = Wage [-sam,]
>
>
>
> ypred=predict(GAM1,newdata=Testing)
>
> y=Testing$wage
>
> MSE = mean((y-ypred)^2)
>
>
> # Bootstrap 95% CI for MSE
>
> # function to obtain MSE
> MSE <- function(formula, data, indices) {
> ??d <- data[indices,] # allows boot to select sample
> ??fit <- gam(formula, data=d)
> ??return(MSE(fit))
>
> } # bootstrapping with 1000 replications
> results <- boot(data=Wage, statistic=MSE,
> ?? R=1000, formula=gam(wage
~education+s(age,bs="ps")+year,data=Wage)
>
> )
> # get 95% confidence intervals
> boot.ci(results, type="bca")
>
> ##########################
>
> ______________________________________________
> 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.
>
Hello, Sorry, there's a close parenthesis too many in the call to boot, the very last one. Delete it and it runs with no errors. Rui Barradas ?s 21:55 de 22/11/2018, Rui Barradas escreveu:> Hello, > > There were several errors with your code. The following works but with > the other CI types. > > > library(ISLR) > library(mgcv) > library(boot) > > # function to obtain MSE > MSE <- function(data, indices, formula) { > ? d <- data[indices, ] # allows boot to select sample > ? fit <- gam(formula, data = d) > ? ypred <- predict(fit) > ? mean((d[["wage"]] - ypred)^2) > } > > data(Wage) > > # Make the results reproducible > set.seed(1234) > > # bootstrapping with 1000 replications > results <- boot(data = Wage, statistic = MSE, > ??????????????? R = 1000, formula = wage ~ education + s(age, bs = > "ps") + year)) > > # get 95% confidence intervals > # type = "bca" is throwing an error > ci.type <- c("norm","basic", "stud", "perc") > boot.ci(results, type = ci.type) > > > > Hope this helps, > > Rui Barradas > > ?s 20:36 de 22/11/2018, varin sacha via R-help escreveu: >> Dear R-experts, >> >> I am trying to get the bootstrapped confidence intervals of Mean >> squared error (MSE) for a (G)AM model. I get an error message. >> Here below the reproducible R code. Many thanks for your response. >> >> #################### >> >> >> install.packages("ISLR") >> >> library(ISLR) >> >> install.packages("mgcv") >> >> library(mgcv) >> >> install.packages("boot") >> >> library(boot) >> >> #MSE calculation >> >> n=dim(Wage)[1] >> >> p=0.667 >> >> GAM1<-gam(wage ~education+s(age,bs="ps")+year,data=Wage) >> >> sam=sample(1?:n,floor(p*n),replace=FALSE) >> >> >> Training =Wage [sam,] >> >> Testing = Wage [-sam,] >> >> >> ypred=predict(GAM1,newdata=Testing) >> >> y=Testing$wage >> >> MSE = mean((y-ypred)^2) >> >> >> # Bootstrap 95% CI for MSE >> >> # function to obtain MSE >> MSE <- function(formula, data, indices) { >> ???d <- data[indices,] # allows boot to select sample >> ???fit <- gam(formula, data=d) >> ???return(MSE(fit)) >> >> } # bootstrapping with 1000 replications >> results <- boot(data=Wage, statistic=MSE, >> ??? R=1000, formula=gam(wage ~education+s(age,bs="ps")+year,data=Wage) >> >> ) >> # get 95% confidence intervals >> boot.ci(results, type="bca") >> >> ########################## >> >> ______________________________________________ >> 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.