envelope<-function(size,x.res,x.hat)){ ident <-diag(size) # ident is set to be a size x size identity matrix # using the diag function epsilon <- matrix(0,size,100) # epsilon is the matrix to contain one hundred size x 1 vectors of standard # normal variates e <- matrix(0,size,100) # e is the matrix to contain the hundred size x 1 vectors of pseudo- # residuals e1 <- numeric(size) # sets up numeric vector of length size e2 <- numeric(size) for(i in 1:100){ epsilon[,i] <- rnorm(size, 0, 1) # generates 31 standard normal variates e[,i] <- (ident-x.hat)%*%epsilon[,i] e[,i]<- sort(e[,i]) # each column of e now contains 31 sorted pseudo-residuals } for(i in 1:size){ eo <- sort(e[i,]) # eo now contains the sorted elements of the 100 observations # defining the empirical distribution of the ith ranked residual e1[i] <- eo[5] e2[i] <- eo[95] # e1 and e2 define the 95% confidence envelope cbind(e1,e2) }