1import numpy as np 2 3import pytest 4from numpy.testing import assert_allclose, assert_raises 5 6from statsmodels.tsa.innovations.arma_innovations import arma_innovations 7from statsmodels.tsa.arima.datasets.brockwell_davis_2002 import dowj, lake 8from statsmodels.tsa.arima.estimators.durbin_levinson import durbin_levinson 9 10 11@pytest.mark.low_precision('Test against Example 5.1.1 in Brockwell and Davis' 12 ' (2016)') 13def test_brockwell_davis_example_511(): 14 # Note: this example is primarily tested in 15 # test_yule_walker::test_brockwell_davis_example_511. 16 17 # Difference the series 18 endog = dowj.diff().iloc[1:] 19 20 # Durbin-Levinson 21 dl, _ = durbin_levinson(endog, ar_order=2, demean=True) 22 23 assert_allclose(dl[0].params, np.var(endog)) 24 assert_allclose(dl[1].params, [0.4219, 0.1479], atol=1e-4) 25 assert_allclose(dl[2].params, [0.3739, 0.1138, 0.1460], atol=1e-4) 26 27 28def check_itsmr(lake): 29 # Test against R itsmr::yw; see results/results_yw_dl.R 30 dl, _ = durbin_levinson(lake, 5) 31 32 assert_allclose(dl[0].params, np.var(lake)) 33 assert_allclose(dl[1].ar_params, [0.8319112104]) 34 assert_allclose(dl[2].ar_params, [1.0538248798, -0.2667516276]) 35 desired = [1.0887037577, -0.4045435867, 0.1307541335] 36 assert_allclose(dl[3].ar_params, desired) 37 desired = [1.08425065810, -0.39076602696, 0.09367609911, 0.03405704644] 38 assert_allclose(dl[4].ar_params, desired) 39 desired = [1.08213598501, -0.39658257147, 0.11793957728, -0.03326633983, 40 0.06209208707] 41 assert_allclose(dl[5].ar_params, desired) 42 43 # itsmr::yw returns the innovations algorithm estimate of the variance 44 # we'll just check for p=5 45 u, v = arma_innovations(np.array(lake) - np.mean(lake), 46 ar_params=dl[5].ar_params, sigma2=1) 47 desired_sigma2 = 0.4716322564 48 assert_allclose(np.sum(u**2 / v) / len(u), desired_sigma2) 49 50 51def test_itsmr(): 52 # Note: apparently itsmr automatically demeans (there is no option to 53 # control this) 54 endog = lake.copy() 55 56 check_itsmr(endog) # Pandas series 57 check_itsmr(endog.values) # Numpy array 58 check_itsmr(endog.tolist()) # Python list 59 60 61def test_nonstationary_series(): 62 # Test against R stats::ar.yw; see results/results_yw_dl.R 63 endog = np.arange(1, 12) * 1.0 64 res, _ = durbin_levinson(endog, 2, demean=False) 65 66 desired_ar_params = [0.92318534179, -0.06166314306] 67 assert_allclose(res[2].ar_params, desired_ar_params) 68 69 70@pytest.mark.xfail(reason='Different computation of variances') 71def test_nonstationary_series_variance(): 72 # See `test_nonstationary_series`. This part of the test has been broken 73 # out as an xfail because we compute a different estimate of the variance 74 # from stats::ar.yw, but keeping the test in case we want to also implement 75 # that variance estimate in the future. 76 endog = np.arange(1, 12) * 1.0 77 res, _ = durbin_levinson(endog, 2, demean=False) 78 79 desired_sigma2 = 15.36526603 80 assert_allclose(res[2].sigma2, desired_sigma2) 81 82 83def test_invalid(): 84 endog = np.arange(2) * 1.0 85 assert_raises(ValueError, durbin_levinson, endog, ar_order=2) 86 assert_raises(ValueError, durbin_levinson, endog, ar_order=-1) 87 assert_raises(ValueError, durbin_levinson, endog, ar_order=1.5) 88 89 endog = np.arange(10) * 1.0 90 assert_raises(ValueError, durbin_levinson, endog, ar_order=[1, 3]) 91 92 93def test_misc(): 94 # Test defaults (order = 0, demean=True) 95 endog = lake.copy() 96 res, _ = durbin_levinson(endog) 97 assert_allclose(res[0].params, np.var(endog)) 98 99 # Test that integer input gives the same result as float-coerced input. 100 endog = np.array([1, 2, 5, 3, -2, 1, -3, 5, 2, 3, -1], dtype=int) 101 res_int, _ = durbin_levinson(endog, 2, demean=False) 102 res_float, _ = durbin_levinson(endog * 1.0, 2, demean=False) 103 assert_allclose(res_int[0].params, res_float[0].params) 104 assert_allclose(res_int[1].params, res_float[1].params) 105 assert_allclose(res_int[2].params, res_float[2].params) 106