Re: [R] boot statictic fn for dual estimation of 2 stats?

[Roger D. Peng]

Occasionaly, the polr returns a zeta element of length equal to 2 rather 
than 3 (which I guess is what you're expecting).  You're bootfn function 
should always check to see that it's returning an object of the same 
length every time.

-roger

Olivia Lau wrote:
Hi, 

I am trying to use boot() to refit an ordinal logit (polr in MASS) model.  
(A very basic bootstrap which samples from the data frame without 
replacement and updates the model.)

I need to extract two statistics per run (the coefficients and zeta) and I 
tried concatenating them into a single vector after fitting, but I get the 
following error: 

Error in "[<-"(*tmp*, r, , value = statistic(data, i[r, ], ...)) : 
	number of items to replace is not a multiple of replacement length

This error goes away if I just return the coefficients, but I need zeta in 
order to calculate predicted probabilities for each outcome category (1:4 
in this case).   Alternatively, if boot() could return a vector of 
predicted probabilities for each categorical outcome, that would 
work as well, but I don't have a clue as to how to start programming this. 

Thanks, Olivia Lau

------------

Here's some sample data and code:  

cost <- c(4, 3, 2, 2, 2, 2, 2, 2, 1, 2, 1, 2, 2, 4, 3, 1, 2, 2, 1, 3, 2, 
1, 4, 1, 1, 2, 2, 2, 1, 3, 1, 1, 1, 1, 2, 1, 1, 1, 1, 2, 2, 2, 2, 2, 1, 1, 
2, 1, 1, 1, 2, 2, 1, 1, 1, 2, 2, 2, 2, 3, 1, 3, 3, 1, 2, 1, 3, 2, 2, 2, 
3, 1, 2, 1, 2, 2, 1, 2)

mil <- c(1, 0, 0, 1, 0, 0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
1, 0, 0, 0, 0, 1)

coop <-  c(4, 2, 1, 1, 1, 1, 2, 1, 2, 2, 1, 1, 3, 3, 3, 1, 4, 3, 1, 3, 4, 
1, 1, 1, 1, 3, 1, 1, 1, 4, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 2, 2, 1, 2, 1, 1, 
1, 2, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 3, 3, 1, 3, 2, 3, 2, 3, 1, 1, 3, 2, 
2, 3, 2, 1, 4, 1, 3)

sanction <- data.frame(cost=cost, mil=mil, coop=coop)

est <- polr(as.factor(cost) ~ mil + coop, data = sanction)

bootfn <- function(data, i, object) {  
  d <- data[i,]
  fit <- update(object, data = d)
  c(fit$coef, fit$zeta)
}

res <- boot(sanction, bootfn, R = 10, object = est)

______________________________________________
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
______________________________________________
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help