# Adding shocks to arima simulation does not work

I am simulating US oil prices, here is the data and a simple model:

oil <- c(20, 18, 20, 23, 24, 23, 25, 24, 20, 39, 32, 23, 23, 25, 22, 23,
27, 27, 25, 25, 24, 23, 19, 20, 26, 24, 24, 26, 25, 25, 26, 29, 28, 33, 34,
29, 27, 27, 26, 21, 20, 23, 18, 22, 26, 35, 38, 44, 47, 50, 43, 42, 42, 41,
31, 39, 39, 47, 46, 49, 41, 39, 42, 49, 52, 59, 53, 66, 70, 77, 69, 73, 79,
71, 71, 71, 72, 82, 92, 106, 136, 110, 51, 61, 78, 75, 79, 83, 78, 78, 90,
100, 94, 86, 101, 106, 85, 92, 87, 95, 99, 110, 104, 112, 115, 103, 68, 60,
74, 60, 51, 50, 63, 59, 69, 66, 60)

lambda <- BoxCox.lambda(oil)
model <- auto.arima(ts(data = oil, frequency = 4),  stepwise=FALSE,
approximation=FALSE, lambda=lambda)


Now, I want to simulate a path:

1. This is what simulate function uses in its actual code, everything is fine:

innovs <- rnorm(5, 0, sqrt(model$sigma2)) simulate(model, innov = innovs, future=TRUE, lambda=lambda) 66 79 68 64 67  2. Now add a simple shock: innovs <- rnorm(5, 0, sqrt(model$sigma2)) + c(0,0,10,0,0)
simulate(model, innov = innovs, future=TRUE, lambda=lambda)
66   59   NA   NA   NA


It is strange that adding lambda causes this problem. If I run this model without lambda, then there is no problem. Any insights?

EDIT:

Based on the great answer from Rob, I tried adding one standard deviation shock and everything works fine:

innovs <- rnorm(5, 0, sqrt(model$sigma2)) + c(0,0,1,0,0)*sqrt(model$sigma2)
simulate(model, innov = innovs, future=TRUE, lambda=lambda)


The lambda value is -0.47 which is an extremely strong Box-Cox transformation. When you add the innovation, and then back-transform, the resulting values are essentially infinite. If you set lambda=0, you will see the problem:

library(forecast)
library(ggplot2)

oil <- c(20, 18, 20, 23, 24, 23, 25, 24, 20, 39, 32, 23, 23, 25, 22, 23,
27, 27, 25, 25, 24, 23, 19, 20, 26, 24, 24, 26, 25, 25, 26, 29, 28, 33, 34,
29, 27, 27, 26, 21, 20, 23, 18, 22, 26, 35, 38, 44, 47, 50, 43, 42, 42, 41,
31, 39, 39, 47, 46, 49, 41, 39, 42, 49, 52, 59, 53, 66, 70, 77, 69, 73, 79,
71, 71, 71, 72, 82, 92, 106, 136, 110, 51, 61, 78, 75, 79, 83, 78, 78, 90,
100, 94, 86, 101, 106, 85, 92, 87, 95, 99, 110, 104, 112, 115, 103, 68, 60,
74, 60, 51, 50, 63, 59, 69, 66, 60) %>%
ts(frequency=4)

lambda <- 0
model <- auto.arima(oil, stepwise=FALSE, approximation=FALSE, lambda=lambda)
innovs <- rnorm(5, 0, sqrt(model\$sigma2)) + c(0,0,10,0,0)
cbind(oil, simulate(model, innov = innovs) %>%
autoplot()