Bootstrap Wielomianowy regresji w R

Question

Próbuję bootstrap proste wielomianu regresji w R, a ja otrzymuję błąd:

Error in is.data.frame(data) : object 'd' not found

Co jest naprawdę dziwne, że używam tego samego kodu (dostosowane do tego szczególny problem) jak w tutorialu dla pakietu startowego at Quick-R, a ten sam kod zadziałał, gdy używam innej funkcji (jak lm()). Na pewno robię coś głupiego, ale nie widzę co. Proszę, jeśli ktokolwiek może pomóc, byłbym bardzo wdzięczny.

To jest przykład:

require(foreign) 
require(nnet) 
require(boot) 

# an example for multinomial logistic regression 
ml = read.dta('http://www.ats.ucla.edu/stat/data/hsbdemo.dta') 
ml = ml[,c(5,7,3)] 

bs <- function(formula, data, indices) { 
    d = data[indices,] # allows boot to select sample 
    fit = multinom(formula, data=d) 
    s = summary(fit) 
    return(list(fit$coefficients, fit$standard.errors)) 
} 

# 5 replications 
results = list() 
results <- boot(
    data=ml, statistic=bs, R=5, parallel='multicore', 
    formula=prog~write 
) 

Answer 1

Błąd dzieje się w summary() części, również obiekt zwrócony przez multinom() nie posiada coefficients i standard.errors. Wygląda na to, że summary.multinom() z kolei oblicza hessian z twoich danych, d, który z jakiegoś powodu (prawdopodobnie problem z zasięgiem) nie może zostać znaleziony. Szybkie ustalenie jest dodanie Hess = TRUE:

bs <- function(formula, data, indices) { 
    d = data[indices,] # allows boot to select sample 
    fit = multinom(formula, data=d, Hess = TRUE) 
    s = summary(fit) 
    return(cbind(s$coefficients, s$standard.errors)) 
} 

# 5 replications 
results = list() 
results <- boot(
    data=ml, statistic=bs, R=5, parallel='multicore', 
    formula=prog~write 
) 

Answer 2

Wielomianowy regresji logistycznej zwraca macierz współczynników z wykorzystaniem funkcji coef(). Różni się to od modelu lm lub glm, który zwraca wektor współczynników.

library(foreign)  # read.dta() 
library(nnet)  # multinom() 
require(boot)  # boot() 

# an example for multinomial logistic regression 
ml = read.dta('http://www.ats.ucla.edu/stat/data/hsbdemo.dta') 
ml = ml[,c(5,7,3)] 

names(ml) 

bs <- function(formula, data, indices) { 
    d = data[indices,] # allows boot to select sample 
    fit = multinom(formula, data=d, maxit=1000, trace=FALSE) 
    #s = summary(fit) 
    #return(list(fit$coefficients, fit$standard.errors)) 

    estimates <- coef(fit) 
    return(t(estimates)) 
} 

# enable parallel 

library(parallel) 
cl <- makeCluster(2) 
clusterExport(cl, "multinom") 

# 10000 replications 
set.seed(1984) 

results <- boot(
    data=ml, statistic=bs, R=10000, parallel = "snow", ncpus=2, cl=cl, 
    formula=prog~write 
) 

# label the estimates 

subModelNames <- colnames(results$t0) 
varNames <- rownames(results$t0) 

results$t0 

estNames <- apply(expand.grid(varNames,subModelNames),1,function(x) paste(x,collapse="_")) 

estNames 

colnames(results$t) <- estNames 

# summary of results 

library(car) 

summary(results) 

confint(results, level=0.95, type="norm") 
confint(results, level=0.95, type="perc") 
confint(results, level=0.95, type="bca") 

# plot the results 

hist(results, legend="separate")