I have been reading many papers where confidence intervals for the impulse responses of autoregressive process have been boot-strapped. The question is, is the usual way sufficient? Take as an example the following often used panel AR(2) model, with other independent variables also included:
$$Y_{it}=\mu_i +\delta_t+ \alpha Y_{i,t-1}+\beta Y_{i,t-2}+\sum_{k=0}^2 (\gamma_k D_{i,t-k}+\theta_kD^2_{i,t-k}+\delta_kD^3_{i,t-k})+\epsilon_{it}$$
Where the variables $D$ are some dummies, and also as can be seen there are fixed time and individual effects. So there are a few ways the equation could potentially be boot strapped. The first way, the usual bootstrap for errors, is to take the predicted values $\hat Y$ and a random sampling of the errors $\epsilon^*$, obtaining $Y^*=\hat Y+\epsilon^*$.Then estimate the equation using $Y^*$ and calculate the responses with respect to the variables of interest. Repeat...
Okay, so the bootstrapped estimates are given from:
$$Y^*_{it}=\mu_i +\delta_t+ \alpha Y^*_{i,t-1}+\beta Y_{i,t-2}^*+\large D+\epsilon^*_{it}$$
Where the big $D$ is the sum term of dummies. My concern is that the independent variables were changed… but this won't matter?
Further, there is a technical challenge as data is dropped when the equation is estimated. I presume calculating the fit for the missing values by dropping the AR terms won't bias the result too much and is the way to go? The time effects are set to zero or do I time demean the observations?
Finally, the issue of computation… Is there R package that can accomplish the goal? Particularly, I am having trouble telling R how to extract the predictions for the first values of each individual without having to code a function that extracts all the predicted values from scratch. It seems to be an incredible pain to extract all the fitted values properly. I am using the plm package for now.