Skip to content
Prev 312329 / 398506 Next

Problem with Stationary Bootstrap

Hello,

You should post to the list, not just send your questions to me. If you 
Cc the list the odds of getting more and better answers are bigger.

As for your question,

1. You forgot a variable X1, like this your code doesn't run. In what 
follows I assume it's also a runif(40, 0, 20).
2. Some simplifications in your code could be made. I've left the 
original instructions, commented out, and the corrections below them.
3. When you say that the bias and the se's are very small what do you 
mean exactly? What kind of values are you expecting and what are the 
results you're getting?



fun <- function(A) {
     model <- lm(Y ~ X1 + X2, data = A)
     e <- resid(model)
     mylag <- function(e, d = 1) {
         n <- length(e)
         c(rep(NA, d), e)[1:n]
     }
     n <- length(e)
     e1 <- mylag(e)
     modele <- lm(e ~ e1 - 1)
     rho <- coef(modele)
     reps <- 20
     for (i in 1:reps) {
         n <- length(e)
         x1star <- c(X1[1] * (1 - rho), X1[2:n] - rho * X1[1:(n - 1)])
         x2star <- c(X2[1] * (1 - rho), X2[2:n] - rho * X2[1:(n - 1)])
         ystar <- c(Y[1] * (1 - rho), Y[2:n] - rho * Y[1:(n - 1)])
         #cr <- (1 - rho)
         #cr[1:n] <- cr
         cr <- rep(1- rho, n)
         #W <- 1
         #W[1:n] <- W
         W <- rep(1, n)
         W[1] = (1 - rho^2)/((1 - rho)^2)
         #W <- c(W[1], W[2:n])
         Model <- lm(ystar ~ cr + x1star + x2star - 1, weights = W)
         bstar <- coef(Model)
         #B0 <- (bstar[[1]][[1]])
         #B1 <- bstar[[2]][[1]]
         #B2 <- bstar[[3]][[1]]
         #Yhat <- B0 + B1 * X1 + B2 * X2
         #u <- Y - Yhat
         u <- resid(Model)
         myu <- function(u, d = 1) {
             n <- length(u)
             c(rep(NA, d), u)[1:n]
         }
         u1 <- myu(u)
         modelu <- lm(u ~ u1 - 1)
         Rho <- coef(modelu)
         if (abs(Rho) > 1)
             break
         diff <- abs(Rho - rho)
         if (diff < 1e-05)
             break else rho <- Rho
     }
     return(coef(Model))
}

l <- runif(40, 0, 20)
X1 <- X2 <- runif(40, 0, 20)
U <- rnorm(1, 0, 2)
for (k in 2:40) {
     U[k] = 0.9 * U[k - 1] + rnorm(1, 0, 1)
}
Y <- l + 2 * X1 + 5 * X2 + U
A <- data.frame(X1 = X1, X2 = X2, Y = Y)

result <- boot::tsboot(A, statistic = fun, R = 1000, l = nrow(A), sim = 
"geom", orig.t = TRUE)
result
Call:
boot::tsboot(tseries = A, statistic = fun, R = 1000, l = nrow(A),
     sim = "geom", orig.t = TRUE)


Bootstrap Statistics :
      original      bias    std. error
t1* 12.563106  0.12260192  0.34212328
t2*  6.980546 -0.01316851  0.03659207
WARNING: All values of t3* are NA


Are these results ok?

Rui Barradas

Em 29-11-2012 04:34, Hock Ann Lim escreveu: