# Series simulator to make pretty pictures from dynardl
library(dynamac)

set.seed(1)
x.error <- rnorm(500, 0, 0.5)
x.full <- x.error
y.error <- rnorm(500, 0, 2)
y.full <- y.error
y.diff.full <- rep(NA, length(y.full))
phi <- -0.8 # coef on ldv
beta.lx <- 1 # coef on l.x
beta.diff.x <- -2 # coef on x


for(i in 2:length(x.full)) {
	x.full[i] <- x.full[i - 1] + x.error[i]
	y.diff.full[i] <- phi*y.full[i - 1] + beta.diff.x*(x.full[i] - x.full[i - 1]) + beta.lx*x.full[i - 1] + y.error[i] # rnorm(1, 0, 1) + const*rnorm(1, 1, 0.25) +
	y.full[i] <- y.full[i - 1] + y.diff.full[i]
}

# Trim the burn in
x <- x.full[101:500]
y.diff <- y.diff.full[101:500]
y <- y.full[101:500]

set.seed(1)
model.ec <- dynardl(y ~ x, lags = list("x" = c(1)), diffs = c("x"), 
	simulate = TRUE, fullsims = TRUE, 
	shockvar = "x",
	sims = 10000, range = 20,
	ec = TRUE)
	
summary(model.ec)

# This is NEW Figure 1.
dynardl.all.plots(model.ec, bw = TRUE, tol = 0.05)