Nothing
### R code from vignette source 'meboot.Rnw'
###################################################
### code chunk number 1: meboot.Rnw:49-50
###################################################
options(prompt = "R> ", digits = 4, continue = "+ ")
###################################################
### code chunk number 2: meboot.Rnw:391-401
###################################################
library("meboot")
library("car")
library("lmtest")
library("dynlm")
data("USconsum")
USconsum <- log(USconsum)
lmcf <- dynlm(consum ~ L(consum, 1) + L(dispinc, 1), data = USconsum)
coeftest(lmcf)
set.seed(135)
durbinWatsonTest(model = lmcf, max.lag = 4)
###################################################
### code chunk number 3: meboot.Rnw:461-466
###################################################
theta <- function(y, x) {
reg <- dynlm(y ~ L(y, 1) + L(x, 1))
thet <- coef(reg)[3]
return(thet)
}
###################################################
### code chunk number 4: meboot.Rnw:488-493
###################################################
theta.cobbd <- function(y, x1, x2) {
reg <- lm(y ~ x1 + x2)
thet <- coef(reg)[2] + coef(reg)[3]
return(thet)
}
###################################################
### code chunk number 5: meboot.Rnw:523-539
###################################################
bstar.consu <- function(y, x, theta,
level = 0.95, bigJ = 999, seed1 = 135) {
set.seed(seed1)
semy <- meboot(x = y, reps = bigJ)$ensemble
semx <- meboot(x = x, reps = bigJ)$ensemble
n <- NROW(y)
m <- length(theta(y, x))
if(m!=1) stop("too many parameters in theta")
bb <- matrix(NA, bigJ)
for(j in 1:bigJ) {
yy <- semy[,j]
xx <- semx[,j]
bb[j] <- theta(yy, xx)
}
return(bb)
}
###################################################
### code chunk number 6: meboot.Rnw:558-576
###################################################
bstar.cobbd <- function(y, x1, x2, theta = theta.cobbd,
level = 0.95, bigJ = 999, seed1 = 135) {
set.seed(seed1)
semy <- meboot(x = y, reps = bigJ)$ensemble
semx1 <- meboot(x = x1, reps = bigJ)$ensemble
semx2 <- meboot(x = x2, reps = bigJ)$ensemble
n <- NROW(y)
m <- length(theta.cobbd(y, x1, x2))
if(m!=1) stop("too many parameters in theta")
bb <- matrix(NA, bigJ)
for(j in 1:bigJ) {
yy <- semy[,j]
xx1 <- semx1[,j]
xx2 <- semx2[,j]
bb[j] <- theta.cobbd(yy, xx1, xx2)
}
return(bb)
}
###################################################
### code chunk number 7: meboot.Rnw:593-610
###################################################
y <- USconsum[,2]
x <- USconsum[,1]
reg <- dynlm(y ~ L(y, 1) + L(x, 1))
su <- summary(reg)
se <- su$coefficients[3,2]
t0 <- theta(y, x)
b3s <- bstar.consu(y, x, theta)
simple.percentile <- quantile(b3s, c(0.025, 0.975), type = 8)
asymmetric.around.0 <- null.ci(b3s)
out <- list(t = b3s, t0 = t0, var.t0 = se^2, R = 999)
class(out) <- "boot"
library("boot")
boot.percentile <- boot.ci(out, type = "perc")$percent[4:5]
boot.norm <- boot.ci(out, type = "norm")$normal[2:3]
boot.basic <- boot.ci(out, type = "basic")$basic[4:5]
rbind(simple.percentile, asymmetric.around.0, boot.percentile,
boot.norm, boot.basic)
###################################################
### code chunk number 8: hdrcde
###################################################
library("hdrcde")
hdr.den(b3s, main = expression(Highest ~ density ~ region ~
of ~ beta [3] ~ estimates : ~ Hall ~ model))
###################################################
### code chunk number 9: meboot.Rnw:752-760
###################################################
library("plm")
data("ullwan")
attach(ullwan)
LMV <- log(MktVal)
ullwan[,3] <- LMV
names(ullwan)[3] <- "LMV"
head(ullwan)
tail(ullwan)
###################################################
### code chunk number 10: meboot.Rnw:775-776
###################################################
summary(lm(Price ~ LMV + Tb3))
###################################################
### code chunk number 11: meboot.Rnw:818-823
###################################################
jboot <- 999
set.seed(567)
LMV.ens <- meboot(x = ullwan, reps = jboot, colsubj = 1, coldata = 3)
Price.ens <- meboot(x = ullwan, reps = jboot, colsubj = 1, coldata = 4)
Tb3.ens <- meboot(x = ullwan, reps = jboot, colsubj = 1, coldata = 6)
###################################################
### code chunk number 12: meboot.Rnw:830-844
###################################################
slopeTb3 <- slopeLMV <- rep(0, jboot)
for(j in 1:jboot) {
frm <- data.frame(Subj = ullwan[,1], Time = ullwan[,2],
Price = Price.ens[,j], Tb3 = Tb3.ens[,j], LMV = LMV.ens[,j])
frm <- pdata.frame(frm, 7)
gip <- coef(plm(Price ~ LMV + Tb3, model = "pooling", data = frm))
slopeTb3[j] <- gip[3]
slopeLMV[j] <- gip[2]
}
Percentile.Tb3 <- quantile(slopeTb3, c(0.025, 0.975), type = 8)
Refined.Tb3 <- null.ci(slopeTb3)
Percentile.LMV <- quantile(slopeLMV, c(0.025, 0.975), type = 8)
Refined.LMV <- null.ci(slopeLMV)
rbind(Percentile.Tb3, Refined.Tb3, Percentile.LMV, Refined.LMV)
###################################################
### code chunk number 13: meboot.Rnw:864-876
###################################################
thetp <- plm(Price ~ LMV + Tb3, model = "pooling", data = ullwan)
varTb3 <- thetp$vcov[3,3]
plm1 <- coef(thetp)
t0Tb3 <- plm1[3]
t0LMV <- plm1[2]
out2 <- list(t = as.matrix(slopeTb3), t0 = t0Tb3,
var.t0 = varTb3, R = 999)
class(out2) <- "boot"
boot.percentile <- boot.ci(out2, type = "perc")$percent[4:5]
boot.norm <- boot.ci(out2, type = "norm")$normal[2:3]
boot.basic <- boot.ci(out2, type = "basic")$basic[4:5]
rbind(Percentile.Tb3, boot.percentile, boot.norm, boot.basic)
###################################################
### code chunk number 14: meboot.Rnw:882-885
###################################################
znp <- pvcm(Price ~ LMV + Tb3, data = ullwan, model = "within")
zplm <- plm(Price ~ LMV + Tb3, data = ullwan)
pooltest(zplm, znp)
###################################################
### code chunk number 15: meboot.Rnw:901-903
###################################################
gir <- plm(Price ~ LMV + Tb3, data = ullwan, model = "random")
coef(gir)
###################################################
### code chunk number 16: meboot.Rnw:910-924
###################################################
slopeTb3 <- slopeLMV <- rep(0, jboot)
for(j in 1:jboot) {
frm <- data.frame(Subj = ullwan[,1], Tim = ullwan[,2],
Price = Price.ens[,j], Tb3 = Tb3.ens[,j], LMV = LMV.ens[,j])
frm <- pdata.frame(frm, 7)
gip <- coef(plm(Price ~ LMV + Tb3, model = "random", data = frm))
slopeTb3[j] <- gip[3]
slopeLMV[j] <- gip[2]
}
Percentile.Tb3 <- quantile(slopeTb3, c(0.025, 0.975), type = 8)
Refined.Tb3 <- null.ci(slopeTb3)
Percentile.LMV <- quantile(slopeLMV, c(0.025, 0.975), type = 8)
Refined.LMV <- null.ci(slopeLMV)
rbind(Percentile.Tb3, Refined.Tb3, Percentile.LMV, Refined.LMV)
###################################################
### code chunk number 17: meboot.Rnw:933-943
###################################################
thetr <- plm(Price ~ LMV + Tb3, model = "random", data = ullwan)
varTb3 <- thetr$vcov[3,3]
plm1 <- coef(thetr)
t0Tb3 <- plm1[3]
out3 <- list(t = as.matrix(slopeTb3), t0 = t0Tb3, var.t0 = varTb3, R = 999)
class(out3) <- "boot"
boot.percentile <- boot.ci(out3, type = "perc")$percent[4:5]
boot.norm <- boot.ci(out3, type = "norm")$normal[2:3]
boot.basic <- boot.ci(out3, type = "basic")$basic[4:5]
rbind(boot.percentile, boot.norm, boot.basic)
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.