* * This is program SPF_RTSIC_roll5_08Nov.prg on 11/20/08 * It runs real-time forecasts to see if improvement is possible over SPF * 1. Use SIC on latest available to create new forecasts in quasi real time (QRTSIC) * 2. Use SIC on real time data to create new forecasts (RTSIC) * 3. Compare survey forecast RMSEs to QRTSIC & RTSIC via RMSE * Does this with latest available & calculates RMSEs for subsamples * Adds subsamples to what mod2 version did * Modified 06Nov to do DM test * 5-year rolling window * cal 1947 1 4 alloc 0 2010:4 * compute vtot = 170 ;* change this declare series z declare series z1 pi1 pi1a declare series newfor declare series pinf declare vector[series] x(vtot) dxq(vtot) dxq4(vtot) declare series dxq4_hat foreareq dxq4q_hat source(noecho) c:\data\aicsic.src source(noecho) c:\Data\msetest2.src * * set some parameters * ffirst = first forecast date * flast = last forecast date * afirst = first actual date * alast = last actual date * compute ffirst = 1971:1 ; * first forecast compute flast = 2008:1 ; * last forecast compute afirst = 1965:3 ; * first actual compute alast = 2007:4 ; * last actual compute aflast = 2006:4 ; * last actual minus one year for lag in forecast * open data c:\data\SPFPIE1_08Jun.xls data(format=xls,org=obs) ffirst flast spfpie1 open data c:\data\actuals_p_08Mar.xls data(format=xls,org=obs) afirst alast actual_l_q4 actual_1q_q4 actual_1y_q4 $ actual_3y_q4 actual_pb_q4 * print(dates) afirst flast spfpie1 actual_l_q4 actual_1q_q4 actual_1y_q4 $ actual_3y_q4 actual_pb_q4 * * set some parameters * send = last observation date in last vintage * vend = date of last vintage * vtot = number of vintages * obstot = number of observations from 1947:2 to send compute send = 2007:4 ;* change this compute vend = send+1 compute obstot = vtot+75 display 'send = ' send %datelabel(send) display 'vend = ' vend %datelabel(vend) display 'vtot = ' vtot display 'obstot = ' obstot open data c:\data\RTDSM_p_08Mar.xls ;* change this data(format=xls,org=obs) 1947:1 send $ P65Q4 $ P66Q1 P66Q2 P66Q3 P66Q4 $ P67Q1 P67Q2 P67Q3 P67Q4 $ P68Q1 P68Q2 P68Q3 P68Q4 $ P69Q1 P69Q2 P69Q3 P69Q4 $ P70Q1 P70Q2 P70Q3 P70Q4 $ P71Q1 P71Q2 P71Q3 P71Q4 $ P72Q1 P72Q2 P72Q3 P72Q4 $ P73Q1 P73Q2 P73Q3 P73Q4 $ P74Q1 P74Q2 P74Q3 P74Q4 $ P75Q1 P75Q2 P75Q3 P75Q4 $ P76Q1 P76Q2 P76Q3 P76Q4 $ P77Q1 P77Q2 P77Q3 P77Q4 $ P78Q1 P78Q2 P78Q3 P78Q4 $ P79Q1 P79Q2 P79Q3 P79Q4 $ P80Q1 P80Q2 P80Q3 P80Q4 $ P81Q1 P81Q2 P81Q3 P81Q4 $ P82Q1 P82Q2 P82Q3 P82Q4 $ P83Q1 P83Q2 P83Q3 P83Q4 $ P84Q1 P84Q2 P84Q3 P84Q4 $ P85Q1 P85Q2 P85Q3 P85Q4 $ P86Q1 P86Q2 P86Q3 P86Q4 $ P87Q1 P87Q2 P87Q3 P87Q4 $ P88Q1 P88Q2 P88Q3 P88Q4 $ P89Q1 P89Q2 P89Q3 P89Q4 $ P90Q1 P90Q2 P90Q3 P90Q4 $ P91Q1 P91Q2 P91Q3 P91Q4 $ P92Q1 P92Q2 P92Q3 P92Q4 $ P93Q1 P93Q2 P93Q3 P93Q4 $ P94Q1 P94Q2 P94Q3 P94Q4 $ P95Q1 P95Q2 P95Q3 P95Q4 $ P96Q1 P96Q2 P96Q3 P96Q4 $ P97Q1 P97Q2 P97Q3 P97Q4 $ P98Q1 P98Q2 P98Q3 P98Q4 $ P99Q1 P99Q2 P99Q3 P99Q4 $ P00Q1 P00Q2 P00Q3 P00Q4 $ P01Q1 P01Q2 P01Q3 P01Q4 $ P02Q1 P02Q2 P02Q3 P02Q4 $ P03Q1 P03Q2 P03Q3 P03Q4 $ P04Q1 P04Q2 P04Q3 P04Q4 $ P05Q1 P05Q2 P05Q3 P05Q4 $ P06Q1 P06Q2 P06Q3 P06Q4 $ P07Q1 P07Q2 P07Q3 P07Q4 $ P08Q1 ;* change all vintage names do v = 1965:4, vend compute [label] vintage1 = %datelabel(v) compute [integer] vintage2 = v - 1965:3 /* if v.eq.1965:4 { display @10 'Vintage Index' @25 'Vintage Date' } display @10 vintage2 @25 vintage1 */ end do v **************************************************** * * * Part 1. Index the Vintage and Create * * Annualized Qtr-Over-Qtr Growth Rates * * * **************************************************** compute [integer] i = 1 dofor z = $ P65Q4 $ P66Q1 P66Q2 P66Q3 P66Q4 $ P67Q1 P67Q2 P67Q3 P67Q4 $ P68Q1 P68Q2 P68Q3 P68Q4 $ P69Q1 P69Q2 P69Q3 P69Q4 $ P70Q1 P70Q2 P70Q3 P70Q4 $ P71Q1 P71Q2 P71Q3 P71Q4 $ P72Q1 P72Q2 P72Q3 P72Q4 $ P73Q1 P73Q2 P73Q3 P73Q4 $ P74Q1 P74Q2 P74Q3 P74Q4 $ P75Q1 P75Q2 P75Q3 P75Q4 $ P76Q1 P76Q2 P76Q3 P76Q4 $ P77Q1 P77Q2 P77Q3 P77Q4 $ P78Q1 P78Q2 P78Q3 P78Q4 $ P79Q1 P79Q2 P79Q3 P79Q4 $ P80Q1 P80Q2 P80Q3 P80Q4 $ P81Q1 P81Q2 P81Q3 P81Q4 $ P82Q1 P82Q2 P82Q3 P82Q4 $ P83Q1 P83Q2 P83Q3 P83Q4 $ P84Q1 P84Q2 P84Q3 P84Q4 $ P85Q1 P85Q2 P85Q3 P85Q4 $ P86Q1 P86Q2 P86Q3 P86Q4 $ P87Q1 P87Q2 P87Q3 P87Q4 $ P88Q1 P88Q2 P88Q3 P88Q4 $ P89Q1 P89Q2 P89Q3 P89Q4 $ P90Q1 P90Q2 P90Q3 P90Q4 $ P91Q1 P91Q2 P91Q3 P91Q4 $ P92Q1 P92Q2 P92Q3 P92Q4 $ P93Q1 P93Q2 P93Q3 P93Q4 $ P94Q1 P94Q2 P94Q3 P94Q4 $ P95Q1 P95Q2 P95Q3 P95Q4 $ P96Q1 P96Q2 P96Q3 P96Q4 $ P97Q1 P97Q2 P97Q3 P97Q4 $ P98Q1 P98Q2 P98Q3 P98Q4 $ P99Q1 P99Q2 P99Q3 P99Q4 $ P00Q1 P00Q2 P00Q3 P00Q4 $ P01Q1 P01Q2 P01Q3 P01Q4 $ P02Q1 P02Q2 P02Q3 P02Q4 $ P03Q1 P03Q2 P03Q3 P03Q4 $ P04Q1 P04Q2 P04Q3 P04Q4 $ P05Q1 P05Q2 P05Q3 P05Q4 $ P06Q1 P06Q2 P06Q3 P06Q4 $ P07Q1 P07Q2 P07Q3 P07Q4 $ P08Q1 ;* change all vintage names set x(i) / = z(t) compute i = i + 1 end dofor do i = 1,vtot set dxq(i) / = (((x(i)(t)/x(i)(t-1))**4)-1.0)*100.0 end do i *************************************************** * * * Part 2. Replace Missing Values for Various * * Vintages With Data From Preceding * * Vintage * * * * Everything in this part is * * variable-specific * * * *************************************************** * Part 2.A. Data Sets Feb92 to Nov92 * Problem: levels of PGDP exist from 59Q1-present only * Solution: replace missing growth rates with those in vintage * Nov91 = x(105) set dxq(106) 1947:2 1959:1 = dxq(105)(t) set dxq(107) 1947:2 1959:1 = dxq(105)(t) set dxq(108) 1947:2 1959:1 = dxq(105)(t) set dxq(109) 1947:2 1959:1 = dxq(105)(t) * Part 2.B. Data Sets Feb96 to Feb97 * Problem: levels of PGDP exist from 59Q3-present only * Solution: replace missing growth rates with those in vintage * Nov95 = x(121) set dxq(122) 1947:2 1959:3 = dxq(121)(t) set dxq(123) 1947:2 1959:3 = dxq(121)(t) set dxq(124) 1947:2 1959:3 = dxq(121)(t) set dxq(125) 1947:2 1959:3 = dxq(121)(t) set dxq(126) 1947:2 1959:3 = dxq(121)(t) * Part 2.C. Data Set Feb96 * Problem: missing 95Q4 observation due to Fed Gov shutdown * Solution: replace missing 95Q4 observation with predicted value * from AR(4) specification in growth rates set x(122) 1995:4 1995:4 = 108.5 set x(122) 1995:3 1995:3 = 107.9 set dxq(122) 1995:4 1995:4 = (((x(122)(t)/x(122)(t-1))**4)-1.0)*100.0 set dxq(122) 1995:3 1995:3 = (((x(122)(t)/x(122)(t-1))**4)-1.0)*100.0 * print 1995:3 1995:4 x(122) dxq(122) * Part 2.D. Data Sets Nov99 to Feb00 * Problem: levels of PGDP exist from 59Q1-present only * Solution: replace missing growth rates with those in vintage * Aug99 = x(136) set dxq(137) 1947:2 1959:1 = dxq(136)(t) set dxq(138) 1947:2 1959:1 = dxq(136)(t) *************************************************** * * * Part 3. Construct Indexed Four-Quarter * * Averages * * * *************************************************** do i = 1,vtot * set dxq4(i) / = ((x(i)(t)/x(i)(t-4))-1.0)*100.0 set dxq4(i) / = (((((dxq(i)(t)/100.0)+1.0)**0.25)* $ (((dxq(i)(t-1)/100.0)+1.0)**0.25)* $ (((dxq(i)(t-2)/100.0)+1.0)**0.25)* $ (((dxq(i)(t-3)/100.0)+1.0)**0.25))-1.0)*100.0 * display 'vintage = ' i %datelabel(i+75) * print(dates) / dxq(i) dxq4(i) end do i * * *************************************************** * * * Part 4. Create SIC based forecasts * * * * * *************************************************** * * * * Here's the new part * clear dxq_hat dxq4_hat compute v = 22 do fbeg = ffirst, flast display 'fbeg = ' fbeg ' v = ' v * We will need to do 1, 2, 3, 4, and 5 steps ahead @aicsic(more, suppress) 1947:2 fbeg-1 8 k_aic k_sic # dxq(v) # constant linreg(define=areq,print) dxq(v) 1947:2 fbeg-1 # constant dxq(v){1 to k_sic} forecast(print) 1 5 fbeg # areq foreareq set dxq4_hat fbeg fbeg = (((((foreareq(t+1)/100.0)+1.0)**0.25)* $ (((foreareq(t+2)/100.0)+1.0)**0.25)* $ (((foreareq(t+3)/100.0)+1.0)**0.25)* $ (((foreareq(t+4)/100.0)+1.0)**0.25))-1.0)*100.0 print fbeg fbeg dxq4_hat clear foreareq compute k_aic = -99, k_sic = -99 * display 'Have computed for v = ' v ' for fbeg = ' fbeg compute v = v + 1 end do fbeg *************************************************** * * * Part 4a. Create SIC based forecasts * * Using latest available vintage * * So quasi-real-time forecasts * * * * * *************************************************** * * * * Here's the new part * clear dxq4q_hat foreareq compute v = 170 ;* change this do fbeg = ffirst, flast display 'fbeg = ' fbeg ' v = ' v * We will need to do 1, 2, 3, 4, and 5 steps ahead @aicsic(more, suppress) 1947:2 fbeg-1 8 k_aic k_sic # dxq(v) # constant linreg(define=areq,print) dxq(v) 1947:2 fbeg-1 # constant dxq(v){1 to k_sic} forecast(print) 1 5 fbeg # areq foreareq set dxq4q_hat fbeg fbeg = (((((foreareq(t+1)/100.0)+1.0)**0.25)* $ (((foreareq(t+2)/100.0)+1.0)**0.25)* $ (((foreareq(t+3)/100.0)+1.0)**0.25)* $ (((foreareq(t+4)/100.0)+1.0)**0.25))-1.0)*100.0 print fbeg fbeg dxq4q_hat clear foreareq compute k_aic = -99, k_sic = -99 * display 'Have computed for v = ' v ' for fbeg = ' fbeg * compute v = v + 1 end do fbeg * *************************************************** * * * Part 5. Calculate RMSEs * * * * * *************************************************** compute [vector[string]] actualbl = $ ||'latest','1 quarter','1 year','3 years','pre-bench'|| compute i = 2 * set z1 = actual_1q_q4 * display ' ' display 'SIC vs SPF vs SICLA' display 'Actuals are ' actualbl(i) display ' ' * compute flast5 = flast - 5 * set pi1 / = z1(t+4) set e = pi1 - spfpie1 set newerror = pi1 - dxq4_hat set quasie = pi1 - dxq4q_hat print(dates) ffirst flast5 pi1 spfpie1 dxq4_hat e newerror print(dates) ffirst flast5 pi1 dxq4q_hat quasie * ***descriptive statistics of forecast errors compute naflast = aflast - 20 do fdate = ffirst, naflast smpl fdate fdate+19 * * * * display ' ' * ***Calculate RMSE * compute tse = %nobs compute aee = %dot(e,e) compute rmsee = sqrt(aee/tse) *compute rmse(2,1) = rmsee display 'This is actuals = ' actualbl(i) display ' ' display 'Root Mean Square Forecast Error for survey =' rmsee * compute tsn = %nobs compute aen = %dot(newerror,newerror) compute rmsen = sqrt(aen/tsn) *compute rmse(2,1) = rmsee display 'Root Mean Square Forecast Error for SIC forecast =' rmsen * set rrmse fdate fdate = rmsen/rmsee * @msetest2(nopictures) spfpie1 dxq4_hat pi1 fdate fdate+19 5 p11 p1 set prmse fdate fdate = p1 * * End of new part end do fdate print(dates) ffirst naflast rrmse prmse *