code: def test_ar1_forecast_simulation(self):
      am = arch_model(self.ar1, mean="AR", vol="GARCH", lags=[1])
      res = am.fit(disp="off")
  
      with preserved_state(self.rng):
          forecast = res.forecast(
              horizon=5,
              start=0,
              method="simulation",
          )
          forecast_reindex = res.forecast(
              horizon=5, start=10, method="simulation", reindex=True
          )
      assert forecast.simulations.index.shape[0] == self.ar1.shape[0]
      assert (
          forecast.simulations.index.shape[0]
          == forecast.simulations.values.shape[0]
      )
  
      with preserved_state(self.rng):
          forecast_reindex = res.forecast(
              horizon=5, start=10, method="simulation", reindex=True
          )
      assert forecast_reindex.mean.shape[0] == self.ar1.shape[0]
      assert (
          forecast_reindex.simulations.index.shape[0]
          == self.ar1.shape[0]
      )
  
      y = np.asarray(self.ar1)
      index = self.ar1.index
      t = y.shape[0]
      params = np.array(res.params)
      resids = np.asarray(y[1:] - params[0] - params[1] * y[:-1])
      vol = am.volatility
      params = np.array(res.params)
      backcast = vol.backcast(resids)
      var_bounds = vol.variance_bounds(resids)
      rng = am.distribution.simulate([])
      vfcast = vol.forecast(
          params[2:],
          resids,
          backcast,
          var_bounds,
          start=0,
          method="simulation",
          rng=rng,
          horizon=5,
      )
      const, ar = params[0], params[1]
      means = np.zeros((t, 5))
      means[:, 0] = const + ar * y
      for i in range(1, 5):
          means[:, i] = const + ar * means[:, i - 1]
      means = pd.DataFrame(
          means, index=index, columns=[f"h.{j}" for j in range(1, 6)]
      )
      assert_frame_equal(means, forecast.mean)
      var = np.concatenate([[[np.nan] * 5], vfcast.forecasts])
      rv = pd.DataFrame(
          var, index=index, columns=[f"h.{j}" for j in range(1, 6)]
      )
      assert_frame_equal(rv, forecast.residual_variance)
  
      lrv = rv.copy()
      for i in range(5):
          weights = (ar ** np.arange(i + 1)) ** 2
          weights = weights[:, None]
          lrv.iloc[:, i : i + 1] = rv.values[:, : i + 1].dot(
              weights[::-1]
          )
      assert_frame_equal(lrv, forecast.variance)
  
imports: import datetime
  from itertools import product
  import numpy
  from numpy.random import RandomState
  from numpy.testing import assert_allclose
  import pandas
  from pandas.testing import assert_frame_equal
  import pytest
  from arch.data import sp500
  from arch.tests.univariate.test_variance_forecasting import preserved_state
  from arch.univariate import APARCH, ARX, EGARCH, FIGARCH, GARCH, HARCH, HARX, ConstantMean, ConstantVariance, EWMAVariance, MIDASHyperbolic, RiskMetrics2006, ZeroMean, arch_model
  from arch.univariate.mean import _ar_forecast, _ar_to_impulse
label_desc: Use image node_id 13 for calling the TestForecasting obj's underlying member method code with example usage: obj.test_ar1_forecast_simulation() without return types