ARIMA estimation by hand - Cross Validated most recent 30 from stats.stackexchange.com 2019-10-21T16:04:31Z https://stats.stackexchange.com/feeds/question/77663 https://creativecommons.org/licenses/by-sa/4.0/rdf https://stats.stackexchange.com/q/77663 15 ARIMA estimation by hand forecaster https://stats.stackexchange.com/users/29137 2013-11-25T20:12:38Z 2017-11-28T04:14:40Z <p>I'm trying to understand how the parameters are estimated in ARIMA modeling/Box Jenkins (BJ). Unfortunately none of the books that I have encountered describes the estimation procedure such as Log-Likelihood estimation procedure in detail. I found the <a href="http://spia.uga.edu/faculty_pages/monogan/teaching/ts/Barima.pdf" rel="nofollow noreferrer">website/teaching material</a> that was very helpful. Following is the equation from the source referenced above. $$LL(\theta)=-\frac{n}{2}\log(2\pi) - \frac{n}{2}\log(\sigma^2) - \sum\limits_{t=1}^n\frac{e_t^2}{2\sigma^2}$$</p> <p>I want to learn the ARIMA/BJ estimation by doing it myself. So I used $R$ to write a code for estimating ARMA by hand. Below is what I did in $R$,</p> <ul> <li>I simulated ARMA (1,1)</li> <li>Wrote the above equation as a function</li> <li>Used the simulated data and the optim function to estimate AR and MA parameters.</li> <li>I also ran the ARIMA in the stats package and compared the ARMA parameters from what I did by hand. <strong><em>Below is the comparison:</em></strong></li> </ul> <p><img src="https://i.stack.imgur.com/UUigG.png" alt="Comparison"></p> <p>**Below are my questions:</p> <ul> <li><strong>Why is there a slight difference between the estimated and calculated variables ?</strong></li> <li><strong>Does ARIMA function in R backcasts or does the estimation procedure differently than what is outlined below in my code?</strong></li> <li><strong>I have assigned e1 or error at observation 1 as 0, is this correct ?</strong></li> <li><strong>Also is there a way to estimate confidence bounds of forecasts using the hessian of the optimization ?</strong></li> </ul> <p>Thanks so much for your help as always.</p> <p>Below is the code:</p> <pre><code>## Load Packages library(stats) library(forecast) set.seed(456) ## Simulate Arima y &lt;- arima.sim(n = 250, list(ar = 0.3, ma = 0.7), mean = 5) plot(y) ## Optimize Log-Likelihood for ARIMA n = length(y) ## Count the number of observations e = rep(1, n) ## Initialize e logl &lt;- function(mx){ g &lt;- numeric mx &lt;- matrix(mx, ncol = 4) mu &lt;- mx[,1] ## Constant Term sigma &lt;- mx[,2] rho &lt;- mx[,3] ## AR coeff theta &lt;- mx[,4] ## MA coeff e = 0 ## Since e1 = 0 for (t in (2 : n)){ e[t] = y[t] - mu - rho*y[t-1] - theta*e[t-1] } ## Maximize Log-Likelihood Function g1 &lt;- (-((n)/2)*log(2*pi) - ((n)/2)*log(sigma^2+0.000000001) - (1/2)*(1/(sigma^2+0.000000001))*e%*%e) ##note: multiplying Log-Likelihood by "-1" in order to maximize in the optimization ## This is done becuase Optim function in R can only minimize, "X"ing by -1 we can maximize ## also "+"ing by 0.000000001 sigma^2 to avoid divisible by 0 g &lt;- -1 * g1 return(g) } ## Optimize Log-Likelihood arimopt &lt;- optim(par=c(10,0.6,0.3,0.5), fn=logl, gr = NULL, method = c("L-BFGS-B"),control = list(), hessian = T) arimopt ############# Output Results############### ar1_calculated = arimopt$par ma1_calculated = arimopt$par sigmasq_calculated = (arimopt$par)^2 logl_calculated = arimopt$val ar1_calculated ma1_calculated sigmasq_calculated logl_calculated ############# Estimate Using Arima############### est &lt;- arima(y,order=c(1,0,1)) est </code></pre> https://stats.stackexchange.com/questions/77663/-/78026#78026 3 Answer by mpiktas for ARIMA estimation by hand mpiktas https://stats.stackexchange.com/users/2116 2013-11-29T07:31:06Z 2013-11-29T07:31:06Z <p>Did you read the the help page of <code>arima</code> function? Here is the relevant excerpt:</p> <blockquote> <p>The exact likelihood is computed via a state-space representation of the ARIMA process, and the innovations and their variance found by a Kalman filter. The initialization of the differenced ARMA process uses stationarity and is based on Gardner et al. (1980). For a differenced process the non-stationary components are given a diffuse prior (controlled by kappa). Observations which are still controlled by the diffuse prior (determined by having a Kalman gain of at least 1e4) are excluded from the likelihood calculations. (This gives comparable results to arima0 in the absence of missing values, when the observations excluded are precisely those dropped by the differencing.)</p> </blockquote> <p>Also relevant is a parameter <code>method=c("CSS-ML", "ML", "CSS")</code>:</p> <blockquote> <p>Fitting method: maximum likelihood or minimize conditional sum-of-squares. The default (unless there are missing values) is to use conditional-sum-of-squares to find starting values, then maximum likelihood.</p> </blockquote> <p>Your results do not differ that much from the ones produced by <code>arima</code> function, so you definitely got everything right. </p> <p>Remember that if you want to compare results of two optimisation procedures, you need to make sure, that the starting values are the same, and the same optimisation method is used, otherwise the results might differ.</p> https://stats.stackexchange.com/questions/77663/-/78858#78858 6 Answer by Cagdas Ozgenc for ARIMA estimation by hand Cagdas Ozgenc https://stats.stackexchange.com/users/20980 2013-12-07T14:39:46Z 2013-12-07T14:39:46Z <p>There is the concept of exact likelihood. It requires the knowledge of initial parameters such as the fist value of the MA error (one of your questions). Implementations usually differ regarding how they treat the initial values. What I usually do is (which is not mentioned in many books) is to also maximize ML w.r.t. the initial values as well.</p> <p>Please take a look at the following from Tsay, it is not covering all cases, but was quite helpful for me: </p> <p><a href="http://faculty.chicagobooth.edu/ruey.tsay/teaching/uts/lec8-08.pdf" rel="noreferrer">http://faculty.chicagobooth.edu/ruey.tsay/teaching/uts/lec8-08.pdf</a> </p>