Assignment 3 Problem 3

library(astsa)

# Fitting AR(2) model to cmort dataset
data(cmort)
fit_ar2 = arima(cmort, order = c(2, 0, 0))

# Extracting Coefficients
ar2_coefs = fit_ar2$coef
print("AR(2) Model Coefficients: ")

## [1] "AR(2) Model Coefficients: "

print(ar2_coefs)

##        ar1        ar2  intercept 
##  0.4300634  0.4424108 88.8538358

# Forecast over an eight-week horizon
horizon = 8
forecast_ar2 = predict(fit_ar2, n.ahead = horizon)

# Extracting forecasts and standard errors
forecasts = forecast_ar2$pred
standard_errs = forecast_ar2$se

# Computing 95% Prediction intervals
upperbnd = forecasts + 1.96*standard_errs
lowerbnd = forecasts - 1.96*standard_errs

# Displaying forecasts and intervals
forecast_df = data.frame(
  Week = 1:horizon,
  Forecast = forecasts,
  Lower = lowerbnd,
  Upper = upperbnd
)
 print("Forecasts and 95% Prediction Intervals : ")

## [1] "Forecasts and 95% Prediction Intervals : "

 print(forecast_df)

##   Week Forecast    Lower     Upper
## 1    1 87.66207 76.51057  98.81358
## 2    2 86.85311 74.71407  98.99214
## 3    3 87.46615 73.45539 101.47690
## 4    4 87.37190 72.45134 102.29246
## 5    5 87.60258 71.76813 103.43703
## 6    6 87.66009 71.18495 104.13524
## 7    7 87.78688 70.75928 104.81448
## 8    8 87.86685 70.40654 105.32716

 # Estimate AR(2) model parameters using Yule-Walker Equations
 fit_yw = ar.yw(cmort, order = 2)
 yw_coefs = fit_yw$ar
 
print("Yule-Walker AR(2) Coefficients : ")

## [1] "Yule-Walker AR(2) Coefficients : "

print(yw_coefs)

## [1] 0.4339481 0.4375768

# Compare Yule-Walker estimates with AR(2) MLE estimates
comparison_df = data.frame(
  "MLE Estimates" = ar2_coefs[1:2],
  "Yule-Walker Estimates" = yw_coefs
)

print("Comparison of AR(2) Coefficients : ")

## [1] "Comparison of AR(2) Coefficients : "

print(comparison_df)

##     MLE.Estimates Yule.Walker.Estimates
## ar1     0.4300634             0.4339481
## ar2     0.4424108             0.4375768