christophsax / tempdisagg

Methods for Temporal Disaggregation and Interpolation of Time Series
http://cran.r-project.org/web/packages/tempdisagg
37 stars 5 forks source link

benchmarking Santos-Silva-Cardoso #25

Open christophsax opened 8 years ago

christophsax commented 8 years ago

The code below replicates some basic Santos-Silva-Cardoso estimation produced by a GAUSS procedure by Tommaso Di Fonzo. Contains all data and can be run in R as it is.


# install.packages("devtools")  # if not installed yet

# install dev version of tempdisagg
devtools::install_github("christophsax/tempdisagg")

library(tempdisagg)

# benchmark input data (from SECO, provided by by Tommaso Di Fonzo)

X <- structure(c(562, 1284, 1336.4, 1185.4, 617.6, 1318.8, 1403.7, 
1260.4, 790, 1432.2, 1467.2, 1257.4, 874.6, 1612.9, 1540, 1412.6, 
912.1, 1537, 1492.4, 1242.9, 759.2, 1396.1, 1330.5, 1207.6, 752.2, 
1229.1, 1229.3, 1023.2, 697.5, 1127.2, 1138, 1037.7, 826.8, 1216.1, 
1190.4, 1109.4, 781.4, 1154.6, 1081, 988.4, 670.4, 1023.9, 982.7, 
858.6, 633.5, 986, 951.4, 879.9, 682.9, 965.3, 983.8, 814.96, 
616.7, 1001.2, 1011.7, 873.63, 742.3, 1017.42, 999.82, 956.37, 
799.01, 1084.22, 1060.1, 984.43, -12.83, -12.83, 28.32, 45.74, 
12.14, 45.74, 33.26, 41.58, 81.26, 35.22, 72.82, 36.41, 58.21, 
72.71, 74.84, 89.81, 60.68, 4.68, 38.49, -34.21, -83.16, -74.84, 
-72.82, -80.91, -106.92, -89.81, -83.83, -124.74, -98.37, -123.27, 
-89.81, -46.78, 4.54, 9.43, 26.35, 42.81, -37.67, -41.94, -51.29, 
-51.57, -56.01, -66.76, -49.89, -50.62, -32.04, -31, -29.7, -19.98, 
-23.83, 15.47, 21.9, 11.9, 11, 22.1, 21.8, 35.5, 29.6, 31.7, 
38.5, 44.5, 27.9, 20, 6.2, 11.7), .Dim = c(64L, 2L))

y0 <- structure(c(31590.3, 32753.9, 35317.5, 38234, 38375, 36224, 35084, 
34595, 37099.4, 36100.5, 32984.4, 32061.3, 32880.7, 32021.3, 
32009, 31082.4), .Dim = c(16L, 1L))

# benchmark output data (produced by the GAUSS package DYNCHOW, by Tommaso Di Fonzo)

out_maxlog <- c(7434.25625248, 7804.49642733, 8122.40980904, 8229.13751115, 
7824.98234147, 8065.48951243, 8350.64918663, 8512.77895947, 8346.94837991, 
8706.24874703, 9059.11792316, 9205.18494989, 9079.02486863, 9514.94871531, 
9779.40289394, 9860.62352212, 9488.76411553, 9633.05645893, 9710.42888886, 
9542.75053668, 9002.98100167, 9070.69864662, 9098.18616122, 9052.1341905, 
8696.26271061, 8780.08511207, 8854.75806683, 8752.89411049, 8440.25185484, 
8554.69701601, 8725.37787429, 8874.67325486, 8943.56556838, 9253.64837922, 
9437.37003085, 9464.81602155, 9134.4648286, 9122.16131066, 9019.55094686, 
8824.32291388, 8372.33631192, 8290.81845968, 8221.60597858, 8099.63924982, 
7871.3226164, 7980.05268821, 8073.38139176, 8136.54330362, 8085.73379074, 
8236.69188099, 8332.20231728, 8226.07201099, 7925.4108585, 7989.67462186, 
8064.1173843, 8042.09713534, 7945.98517736, 8038.22567347, 8052.47635024, 
7972.31279893, 7708.37300254, 7744.57806824, 7795.04800574, 7834.40092348
)

out_minrss <- c(7659.28713123, 7818.59224169, 8004.39192547, 8108.02870161, 
7917.98168629, 8089.65492176, 8294.98150935, 8451.2818826, 8448.89900606, 
8703.08684756, 8997.17297681, 9168.34116957, 9199.32038264, 9510.2152571, 
9716.22111684, 9808.24324342, 9597.85937787, 9624.60713696, 9649.42026233, 
9503.11322283, 9117.62248461, 9074.57877766, 9041.97368026, 8989.82505747, 
8767.31699261, 8785.30793751, 8811.04601261, 8720.32905726, 8523.69367166, 
8559.30159411, 8677.08780008, 8834.91693415, 8994.00710228, 9245.23930448, 
9406.18802363, 9453.96556961, 9211.55176156, 9124.08900849, 8981.08097156, 
8783.77825839, 8434.713759, 8289.57597596, 8185.84651316, 8074.26375189, 
7930.45527137, 7979.36800795, 8041.67517791, 8109.80154278, 8122.44882395, 
8236.42157191, 8300.5014401, 8221.32816404, 8000.11636657, 7993.05341739, 
8015.07886455, 8013.0513515, 7979.49829595, 8029.3266962, 8030.02503062, 
7970.14997723, 7784.58875038, 7764.48935583, 7761.90718784, 7771.41470595
)

# Replication with tempdisagg

yts <- ts(y0)
Xts <- ts(X, frequency = 4)

m1 <- td(yts ~ Xts, method = "dynamic-maxlog")
m2 <- td(yts ~ Xts, method = "dynamic-minrss")

out_maxlog - predict(m1)
out_minrss - predict(m2)

summary(m1)
summary(m2)

### GAUSS Output

### MAXLOG

#                   Dynchow estimation results

# First order dynamic model with constant
# Maximum likelihood
# Estimation approach: Santos Silva and Cardoso (2001)
# ----------------------------------------------------------

#       ------------------------------------------------------------------------------
#        Variables         Coeff.            s.e.           t-stat            p-value
#       ------------------------------------------------------------------------------
#         CONSTANT     874.50754906     206.76159845       4.22954531       0.00058447 
#            REG01       0.75714420       0.19292737       3.92450385       0.00100931 
#            REG02       0.30333384       0.57079863       0.53142005       0.30241276 
#         TR. REM.    7620.84018240     549.81454541      13.86074677       0.00000000 
#       ------------------------------------------------------------------------------
#              PHI       0.80543946       0.07735514      10.41222989       0.00000006 
#       ------------------------------------------------------------------------------

# Elapsed time:       0.03100000  seconds.
# Execution stopped in line 184 of C:\Users\christoph\Desktop\dynchow\seco.g

### MINRSS

#                   Dynchow estimation results

# First order dynamic model with constant
# Estimated GLS (min SSR)
# Estimation approach: Santos Silva and Cardoso (2001)
# ----------------------------------------------------------

#       ------------------------------------------------------------------------------
#        Variables         Coeff.            s.e.           t-stat            p-value
#       ------------------------------------------------------------------------------
#         CONSTANT     290.98749192     188.66916988       1.54231607       0.07447261 
#            REG01       0.42791734       0.17597359       2.43171340       0.01581667 
#            REG02       0.76124301       0.52196813       1.45840899       0.08520006 
#         TR. REM.    7809.47851900     588.35595036      13.27339090       0.00000001 
#       ------------------------------------------------------------------------------
#              PHI       0.91396331       0.00000000       0.00000000       0.00000000 
#       ------------------------------------------------------------------------------

# Elapsed time:       0.01600000  seconds.
# Execution stopped in line 184 of C:\Users\christoph\Desktop\dynchow\seco.g
christophsax commented 7 years ago

This has been benchmarked by @petersteiner with the Matlab routines and everything seems fine. Any details?

Still open:

christophsax commented 5 years ago

Method is still marked as 'experimental'