Set of functions to compute diagnostics and goodness of fit of filtered series:
cross validation (cv()) and cross validate estimate (cve()),
leave-one-out cross validation estimate (loocve),
CP statistic (cp()) and Rice's T statistics (rt()).
Usage
cve(x, coef, ...)
cv(x, coef, ...)
loocve(x, coef, ...)
rt(x, coef, ...)
cp(x, coef, var, ...)Arguments
- x
input time series.
- coef
vector of coefficients or a moving-average (
moving_average()).- ...
other arguments passed to the function
moving_average()to convertcoefto a"moving_average"object.- var
variance used to compute the CP statistic (
cp()).
Details
Let \((\theta_i)_{-p\leq i \leq q}\) be a moving average of length \(p+q+1\) used to filter a time series \((y_i)_{1\leq i \leq n}\). Let denote \(\hat{\mu}_t\) the filtered series computed at time \(t\) as: $$ \hat{\mu}_t = \sum_{i=-p}^q \theta_i y_{t+i}. $$
The cross validation estimate (cve()) is defined as the time series \(Y_t-\hat{\mu}_{-t}\) where
\(\hat{\mu}_{-t}\) is the leave-one-out cross validation estimate (loocve()) defined as the filtered series
computed deleting the observation \(t\) and remaining all the other points.
The cross validation statistics (cv()) is defined as:
$$
CV=\frac{1}{n-(p+q)}
\sum_{t=p+1}^{n-q} \left(y_t - \hat{\mu}_{-t}\right)^2.
$$
In the case of filtering with a moving average, we can show that:
$$
\hat{\mu}_{-t}= \frac{\hat{\mu}_t - \theta_0 y_t}{1-\theta_0}
$$
and
$$
CV=\frac{1}{n-(p+q)}
\sum_{t=p+1}^{n-q} \left(\frac{y_t - \hat{\mu}_{t}}{1-\theta_0}\right)^2.
$$
In the case of filtering with a moving average,
the CP estimate of risk (introduced by Mallows (1973); cp()) can be defined as:
$$
CP=\frac{1}{\sigma^2}
\sum_{t=p+1}^{n-q} \left(y_t - \hat{\mu}_{t}\right)^2
-(n-(p+q))(1-2\theta_0).
$$
The CP method requires an estimate of \(\sigma^2\) (var parameter).
The usual use of CP is to compare several different fits (for example different bandwidths):
one should use the same estimate of \(\hat{\sigma}^2\) for all fits (using for example var_estimator()).
The recommendation of Cleveland and Devlin (1988) is to compute \(\hat{\sigma}^2\)
from a fit at the smallest bandwidth under consideration,
at which one should be willing to assume that bias is negligible.
The Rice's T statistic (rt()) is defined as:
$$
\frac{1}{n-(p+q)}
\sum_{t=p+1}^{n-q}
\frac{
\left(y_t - \hat{\mu}_{t}\right)^2
}{
1-2\theta_0
}
$$
References
Loader, Clive. 1999. Local regression and likelihood. New York: Springer-Verlag.
Mallows, C. L. (1973). Some comments on Cp. Technometrics 15, 661– 675.
Cleveland, W. S. and S. J. Devlin (1988). Locally weighted regression: An approach to regression analysis by local fitting. Journal of the American Statistical Association 83, 596–610.
Examples
x <- retailsa$DrinkingPlaces
h13 <- lp_filter(6)@sfilter
cv(x, h13)
#> [1] 407.4261
cve(x, h13)
#> Jan Feb Mar Apr May
#> 1992 NA NA NA NA NA
#> 1993 -5.14618919 -17.36539019 -22.66017086 19.48627848 6.38676264
#> 1994 -14.94755369 1.31816849 16.08166194 2.99065561 15.39861048
#> 1995 -1.53508398 -4.48791513 -18.81871480 11.35536743 18.76040896
#> 1996 -11.33393568 5.45864719 -8.62268474 -6.19187107 0.95078248
#> 1997 18.77818119 -18.61322499 -17.17605235 -4.12577780 -4.27674801
#> 1998 10.59258011 28.16293545 8.31659995 -22.05390503 -17.51712183
#> 1999 -48.17323628 17.53856794 38.38939220 18.54077982 -22.70706791
#> 2000 -7.02610985 8.02570873 36.22036597 -8.13817046 -11.42449982
#> 2001 27.33400965 5.80412776 6.99879327 4.78163034 -30.22752227
#> 2002 -11.34762754 -0.20243116 -10.07480443 -25.98719900 6.73661480
#> 2003 11.46806052 -29.37016568 0.62875491 8.51260766 29.26958198
#> 2004 27.23957619 30.30472892 -21.95477831 -16.60558278 14.86214050
#> 2005 12.27255423 9.29866299 -21.40464112 19.84056059 3.56401766
#> 2006 18.54990122 3.33008511 -15.86831521 0.12489661 20.36885048
#> 2007 19.90360577 -27.82185673 -2.48409391 -5.48784718 -0.27915873
#> 2008 -23.53117762 -31.69747196 23.02579057 -0.03672121 -34.76937258
#> 2009 -3.39345458 -29.38214322 11.69201316 3.12072701 -5.46458948
#> 2010 -25.71132511 58.38731382 9.11448335 -20.02745162 -9.70022829
#> Jun Jul Aug Sep Oct
#> 1992 NA -7.45163005 1.83150696 -10.63985277 -1.81096445
#> 1993 -6.91124669 17.38961154 7.21934587 -14.90866741 5.25405342
#> 1994 7.48706942 -26.49465414 -0.61391181 -0.55328564 -3.43809917
#> 1995 -8.77479238 -11.81560024 -18.60041024 29.35430793 -3.10710680
#> 1996 6.99704008 3.42028972 -5.28833799 -1.13947586 -12.35225098
#> 1997 -3.63767869 13.72028750 18.44547779 -10.50813630 19.01695738
#> 1998 5.30638008 21.56486306 -3.08180164 -2.88325027 12.30211219
#> 1999 8.98878086 -11.22572665 -19.97397123 2.67316574 16.61960698
#> 2000 -9.58312781 8.03088417 10.89665600 7.67915962 -22.27612789
#> 2001 -8.71909653 -11.83960284 17.45698711 11.97476138 -8.24814025
#> 2002 14.21159421 40.41760670 -16.44857081 -24.57738713 -21.08680019
#> 2003 -22.09374844 -33.43274719 38.81279522 -8.79340528 9.19923890
#> 2004 7.57149298 -22.90992623 -8.19251285 19.21134252 17.72941598
#> 2005 -13.84431449 -19.53383934 4.49824921 33.49423202 0.22197059
#> 2006 -9.69489499 1.51751027 8.60055368 -21.70330845 11.25841667
#> 2007 -29.52231179 102.55610067 -65.42302394 3.23086048 -3.29366590
#> 2008 -16.49353525 79.73225912 19.58813744 -24.25576408 -29.71545803
#> 2009 43.56059230 7.56011540 -41.59033165 -8.30191350 3.06300827
#> 2010 -0.11265854 NA NA NA NA
#> Nov Dec
#> 1992 18.09364541 15.76603944
#> 1993 8.47467410 -21.28816011
#> 1994 7.03649340 9.10556156
#> 1995 6.18802494 7.38605250
#> 1996 1.88234334 23.30457282
#> 1997 -16.73061099 -31.60062082
#> 1998 -8.66668807 -1.98694389
#> 1999 -7.20028076 -12.11272420
#> 2000 -22.20734469 2.22191852
#> 2001 1.92347751 28.38857900
#> 2002 8.69885196 20.90120091
#> 2003 2.12312265 -41.34684151
#> 2004 -3.47820081 -26.28682206
#> 2005 -11.34082468 -20.55652496
#> 2006 -6.64974968 9.68338896
#> 2007 19.85589339 16.47456313
#> 2008 -14.06707239 31.84197281
#> 2009 4.08447074 -7.70501444
#> 2010 NA NA
loocve(x, h13)
#> Jan Feb Mar Apr May Jun Jul
#> 1992 NA NA NA NA NA NA 1001.4188
#> 1993 1015.1594 1020.3532 1026.8713 1024.3485 1035.0342 1042.5584 1034.8609
#> 1994 988.7961 990.2744 992.1258 997.4044 991.6407 986.1583 985.6242
#> 1995 988.5031 996.0819 1007.7043 1010.0722 1019.6715 1039.2706 1054.0698
#> 1996 1085.0587 1074.7321 1073.9555 1072.4712 1072.1337 1073.6848 1078.4236
#> 1997 1094.6668 1099.3398 1097.6220 1098.7211 1109.4822 1124.5497 1134.4864
#> 1998 1147.8621 1147.0683 1155.6695 1166.6834 1169.6807 1168.2666 1169.3165
#> 1999 1208.4481 1206.2939 1211.1326 1218.6411 1224.1920 1208.3887 1206.5105
#> 2000 1284.3230 1295.6725 1297.4652 1310.2377 1306.8928 1297.0525 1280.5238
#> 2001 1259.6144 1267.0683 1265.9439 1265.6477 1277.8026 1284.3451 1305.4978
#> 2002 1382.4145 1373.5619 1372.4694 1375.9839 1369.0936 1366.9831 1357.6295
#> 2003 1395.0378 1418.9726 1423.0347 1430.7793 1434.8607 1456.6359 1467.2549
#> 2004 1470.8433 1468.0161 1476.9507 1470.2130 1458.2565 1458.6102 1467.6244
#> 2005 1464.2372 1465.9738 1476.3027 1471.7152 1483.7125 1499.3042 1514.4133
#> 2006 1571.5614 1585.6833 1602.3406 1609.4686 1612.8793 1625.8015 1626.8222
#> 2007 1637.9071 1652.7135 1653.8260 1665.2150 1675.7851 1691.9900 1664.8237
#> 2008 1627.1504 1625.0269 1617.3124 1640.4628 1672.1715 1688.3681 1676.8539
#> 2009 1666.6162 1681.0205 1679.3326 1686.2126 1687.6141 1669.3136 1667.3008
#> 2010 1695.6635 1684.4013 1698.9724 1705.3744 1701.9965 1702.1428 NA
#> Aug Sep Oct Nov Dec
#> 1992 998.4183 1005.0462 1007.6766 1006.3912 1008.9190
#> 1993 1028.4304 1019.9247 1000.5621 988.3539 989.4471
#> 1994 973.4110 972.0057 974.8608 976.1581 980.4700
#> 1995 1069.4671 1070.1204 1086.7905 1088.5547 1086.4344
#> 1996 1085.8888 1091.6483 1100.5698 1101.1735 1096.5632
#> 1997 1141.5656 1150.6398 1143.0091 1151.2018 1155.7109
#> 1998 1178.8301 1178.7093 1173.1592 1179.4854 1185.0623
#> 1999 1207.6861 1209.0833 1219.3172 1243.5352 1265.6620
#> 2000 1267.9811 1259.9709 1262.6220 1262.8058 1260.8995
#> 2001 1323.9274 1350.2080 1374.7824 1382.9176 1377.8403
#> 2002 1367.1471 1366.8675 1367.7855 1367.5832 1377.1793
#> 2003 1455.3307 1471.4263 1472.0933 1477.2598 1488.5324
#> 2004 1467.5947 1463.7731 1464.8130 1469.3895 1473.9520
#> 2005 1523.2259 1529.8537 1549.5539 1561.8692 1572.2058
#> 2006 1628.4739 1638.6933 1633.7511 1640.8886 1639.5967
#> 2007 1703.2466 1677.3435 1664.2483 1642.2551 1628.5558
#> 2008 1692.6857 1696.2875 1686.3758 1672.0949 1655.7922
#> 2009 1666.5111 1650.3396 1648.7069 1658.4430 1676.4276
#> 2010 NA NA NA NA NA
rt(x, h13)
#> [1] 452.5878
cp(x, h13, var = var_estimator(x, h13))
#> [1] 44.02421