workfile DTR 7 1/1/1960 7/31/2018 smpl 1/1/1960 7/31/2018 'Create trend and monthly dummies genr time = @trend + 1 genr d1=@month=1 genr d2=@month=2 genr d3=@month=3 genr d4=@month=4 genr d5=@month=5 genr d6=@month=6 genr d7=@month=7 genr d8=@month=8 genr d9=@month=9 genr d10=@month=10 genr d11=@month=11 genr d12=@month=12 for %station atl bos bwi cvg dfw dsm dtw las lga msp ord pdx phl slc tus 'Read data and create AVG and DTR read {%station}.txt max{%station} min{%station} genr AVG_{%station}=(max{%station}+min{%station})/2 genr DTR_{%station}=max{%station}-min{%station} 'Kernel density estimates smpl 1960 2017 freeze(a0_{%station}) AVG_{%station}.distplot kernel a0_{%station}.setelem(1) lwidth(3) freeze(a00_{%station}) DTR_{%station}.distplot kernel a00_{%station}.setelem(1) lwidth(3) graph a000_{%station}.merge a0_{%station} a00_{%station} a000_{%station}.align(2,1.5,0) a000_{%station}.axis(bottom) range(-10,100) a000_{%station}.axis(left) range(0,.065) '***** Sequential regressions ***** 'Trend smpl 1960 2017 equation a155_{%station} a155_{%station}.ls(cov=hac) AVG_{%station} c time freeze(taba155_{%station}) a155_{%station}.stats genr AVGtrend_{%station} = AVG_{%station} - resid group groupxxx_{%station} AVG_{%station} AVGtrend_{%station} freeze(a1_{%station}) groupxxx_{%station}.line a1_{%station}.setelem(1) lcolor(gray) a1_{%station}.setelem(2) lwidth(3) lcolor(blue) a1_{%station}.legend -display genr AVGdet_{%station} = resid equation a156_{%station} a156_{%station}.ls(cov=hac) DTR_{%station} c time freeze(taba156_{%station}) a156_{%station}.stats genr DTRtrend_{%station} = DTR_{%station} - resid group groupyyy_{%station} DTR_{%station} DTRtrend_{%station} freeze(a11_{%station}) groupyyy_{%station}.line a11_{%station}.setelem(1) lcolor(gray) a11_{%station}.setelem(2) lwidth(3) lcolor(blue) a11_{%station}.legend -display genr DTRdet_{%station} = resid graph a111_{%station}.merge a1_{%station} a11_{%station} a111_{%station}.align(2,1.5,0) 'Fixed seasonality smpl 1960 2017 equation a444a_{%station} a444a_{%station}.ls(cov=hac) AVGdet_{%station} d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 d12 freeze(taba444a_{%station}) a444a_{%station}.stats genr AVGdetseas_{%station} = AVGdet_{%station} - resid group groupy2a_{%station} avgdet_{%station} avgdetseas_{%station} freeze(a4444a_{%station}) groupy2a_{%station}.line a4444a_{%station}.setelem(1) lcolor(gray) a4444a_{%station}.setelem(2) lwidth(3) lcolor(blue) a4444a_{%station}.options size(20,3) a4444a_{%station}.legend -display genr AVGdetdes_{%station}=resid equation a444_{%station} a444_{%station}.ls(cov=hac) DTRdet_{%station} d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 d12 freeze(taba444_{%station}) a444_{%station}.stats genr DTRdetseas_{%station} = DTRdet_{%station} - resid group groupy2_{%station} DTRdet_{%station} DTRdetseas_{%station} freeze(a4444_{%station}) groupy2_{%station}.line a4444_{%station}.setelem(1) lcolor(gray) a4444_{%station}.setelem(2) lwidth(3) lcolor(blue) a4444_{%station}.options size(20,3) a4444_{%station}.legend -display genr DTRdetdes_{%station}=resid graph a4444b_{%station}.merge a4444a_{%station} a4444_{%station} '12-month fixed seasonal patterns smpl 2017 2017 freeze(a6666db_{%station}) AVGdetseas_{%station}.line a6666db_{%station}.setelem(1) lwidth(3) freeze(a6666d_{%station}) DTRdetseas_{%station}.line a6666d_{%station}.setelem(1) lwidth(3) graph a6666ddd_{%station}.merge a6666db_{%station} a6666d_{%station} a6666ddd_{%station}.align(2,1.5,0) 'Time-varying seasonality: full-sample fits smpl 1960 2017 equation a666AVG_{%station} a666AVG_{%station}.ls(cov=hac) AVGdet_{%station} d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 d12 d1*time d2*time d3*time d4*time d5*time d6*time d7*time d8*time d9*time d10*time d11*time d12*time freeze(taba666AVG_{%station}) a666AVG_{%station}.stats genr AVGdetseas2_{%station} = AVGdet_{%station} - resid genr AVGdetdes2_{%station} = resid equation a666_{%station} a666_{%station}.ls(cov=hac) DTRdet_{%station} d1 d2 d3 d4 d5 d6 d7 d8 d9 d10 d11 d12 d1*time d2*time d3*time d4*time d5*time d6*time d7*time d8*time d9*time d10*time d11*time d12*time freeze(taba666_{%station}) a666_{%station}.stats genr DTRdetseas2_{%station} = DTRdet_{%station} - resid genr DTRdetdes2_{%station} = resid 'Time-Varying Seasonality: Fits at beginning and end smpl 1/1/2017 12/31/2017 genr earlyshiftAVG_{%station} = AVGdetseas2_{%station}(-20820) group patternsAVG_{%station} earlyshiftAVG_{%station} AVGdetseas2_{%station} freeze(help2AVG_{%station}) patternsAVG_{%station}.line help2AVG_{%station}.setelem(1) lwidth(3) lcolor(blue) help2AVG_{%station}.setelem(2) lwidth(3) lcolor(red) help2AVG_{%station}.legend -display genr earlyshiftDTR_{%station} = DTRdetseas2_{%station}(-20820) group patternsDTR_{%station} earlyshiftDTR_{%station} DTRdetseas2_{%station} freeze(help2DTR_{%station}) patternsDTR_{%station}.line help2DTR_{%station}.setelem(1) lwidth(3) lcolor(blue) help2DTR_{%station}.setelem(2) lwidth(3) lcolor(red) help2DTR_{%station}.legend -display graph help2222_{%station}.merge help2AVG_{%station} help2DTR_{%station} help2222_{%station}.align(2,1.5,0) '***** JOINT regressions ***** smpl 1960 2017 'Conditional Mean Regression, AVG equation avg_np1_{%station} avg_np1_{%station}.ls(cov=hac) avg_{%station} c time d1 d2 d3 d4 d5 d6 d8 d9 d10 d11 d12 d1*time d2*time d3*time d4*time d5*time d6*time d8*time d9*time d10*time d11*time d12*time avg_{%station}(-1) genr e_avg_np1_{%station}=resid 'p(nt) freeze(avg_np1_nt_{%station}) avg_np1_{%station}.wald c(2)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(ns) freeze(avg_np1_ns_{%station}) avg_np1_{%station}.wald c(3)=c(4)=c(5)=c(6)=c(7)=c(8)=c(9)=c(10)=c(11)=c(12)=c(13)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(nts) freeze(avg_np1_nts_{%station}) avg_np1_{%station}.wald c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 freeze(taba777cc_{%station}) avg_np1_{%station}.stats 'Conitional Mean Regression, DTR equation dtr_np1_{%station} dtr_np1_{%station}.ls(cov=hac) dtr_{%station} c time d1 d2 d3 d4 d5 d6 d8 d9 d10 d11 d12 d1*time d2*time d3*time d4*time d5*time d6*time d8*time d9*time d10*time d11*time d12*time dtr_{%station}(-1) genr e_dtr_np1_{%station}=resid 'p(nt) freeze(dtr_np1_nt_{%station}) dtr_np1_{%station}.wald c(2)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(ns) freeze(dtr_np1_ns_{%station}) dtr_np1_{%station}.wald c(3)=c(4)=c(5)=c(6)=c(7)=c(8)=c(9)=c(10)=c(11)=c(12)=c(13)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(nts) freeze(dtr_np1_nts_{%station}) dtr_np1_{%station}.wald c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 freeze(taba777dd_{%station}) dtr_np1_{%station}.stats 'Conditional Variance Regression, AVG genr e2_avg_{%station} = (e_avg_np1_{%station})^2 equation avg_np2_{%station} avg_np2_{%station}.ls(cov=hac) e2_avg_{%station} c time d1 d2 d3 d4 d5 d6 d8 d9 d10 d11 d12 d1*time d2*time d3*time d4*time d5*time d6*time d8*time d9*time d10*time d11*time d12*time e2_avg_{%station}(-1) 'p(nt) freeze(avg_np2_nt_{%station}) avg_np2_{%station}.wald c(2)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(ns) freeze(avg_np2_ns_{%station}) avg_np2_{%station}.wald c(3)=c(4)=c(5)=c(6)=c(7)=c(8)=c(9)=c(10)=c(11)=c(12)=c(13)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(nts) freeze(avg_np2_nts_{%station}) avg_np2_{%station}.wald c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 freeze(taba777cccc_{%station}) avg_np2_{%station}.stats 'Conditional Variance Regression, DTR genr e2_dtr_{%station} = (e_dtr_np1_{%station})^2 equation dtr_np2_{%station} dtr_np2_{%station}.ls(cov=hac) e2_dtr_{%station} c time d1 d2 d3 d4 d5 d6 d8 d9 d10 d11 d12 d1*time d2*time d3*time d4*time d5*time d6*time d8*time d9*time d10*time d11*time d12*time e2_dtr_{%station}(-1) 'p(nt) freeze(dtr_np2_nt_{%station}) dtr_np2_{%station}.wald c(2)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(ns) freeze(dtr_np2_ns_{%station}) dtr_np2_{%station}.wald c(3)=c(4)=c(5)=c(6)=c(7)=c(8)=c(9)=c(10)=c(11)=c(12)=c(13)=c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 'p(nts) freeze(dtr_np2_nts_{%station}) dtr_np2_{%station}.wald c(14)=c(15)=c(16)=c(17)=c(18)=c(19)=c(20)=c(21)=c(22)=c(23)=c(24)=0 freeze(taba777dddd_{%station}) dtr_np2_{%station}.stats 'Standardized residual skewness and kurtosis avg_np2_{%station}.fit sigma2_avg_{%station} genr sigma_avg_{%station} = sigma2_avg_{%station}^.5 genr standardized_avg_{%station} = e_avg_np1_{%station} / sigma_avg_{%station} dtr_np2_{%station}.fit sigma2_dtr_{%station} genr sigma_dtr_{%station} = sigma2_dtr_{%station}^.5 genr standardized_dtr_{%station} = e_dtr_np1_{%station} / sigma_dtr_{%station} next stop