*
* @CusumTests( options ) resids start end
* Performs CUSUM and CUSUMQ tests of the input series, which should be a
* series of recursive residuals.
*
* Parameters:
*   resids     Series of recursive residuals
*   start end  Range to use. By default, the defined range of <<resids>>
*
* Variables Defined:
*   %CDSTAT    value of the CUSUM statistic (maximum absolute value)
*   %SIGNIF    (approximate) significance level of the CUSUM
*   %%CUSUMQ   value of CUSUMQ statistic (maximum absolute gap)
*   %%CUSUMQCV .05 critical value for CUSUMQ
*
* References:
*
* For the general technique:
*
*   Brown, Durbin & Evans(1975), "Techniques for Testing the Constancy of
*    Regression Relationships over Time", JRSS-B, vol 37, 149-192.
*
* For the significance level on the CUSUM test:
*
*   Zeileis(2004), "Alternative Boundaries for CUSUM tests", Statistical
*    Papers, vol 45, no 1, 123-131.
*
* For the critical boundaries on the CUSUMQ test:
*
*   Edgerton & Wells(1994), "Critical Values for the Cusumsq Statistic in
*    Medium and Large Sized Samples," Oxford Bulletin of Economics and
*    Statistics, vol. 56, no 3, 355-365.
*
* Revision Schedule:
*   09/2009 Written by Tom Doan, Estima
*   07/2010 Added HEADER and FOOTER options
*
procedure CusumTests resids start end
type 		series  	resids
type 		integer 	start end
*
option 	string	header
option   string	footer
*
local 	integer 	startl endl
local 	series  	cusum  cusumsq  cusumstd residssq  target lower upper
local 	real    	tsigma n
*
declare real %%cusumq %%cusumqcv
*
if .not.%defined(resids) {
   disp "Syntax: @cusumtests RESIDS start end"
   return
}
inquire(series=resids) startl<<start endl<<end
*
* Compute the sum of squared deviations from the mean for resids
*
stats(noprint) resids startl endl
compute tsigma=%variance*(%nobs-1)
*
* Compute the cumulated sum of resids scaled to make what should be a
* Brownian motion on [0,1] under the null.
*
acc resids   startl endl cusum
set cusum    startl endl = cusum/sqrt(tsigma)
*
* Compute the cumulated sum of squared residuals pinned to anchors at 0,1
*
set cusumsq  startl endl = resids^2
acc cusumsq  startl endl
set cusumsq  startl endl = cusumsq/cusumsq(endl)
*
* target is a straight line from 0 to 1 (the normalized entry number).
*
set target   startl endl = (t-startl+1.0)/(endl-startl+1.0)
*
* This is the cusum divided by the test boundary, used for computing the
* significance level.
*
set cusumstd startl endl = cusum/(1+2.0*target)
*

if .not.%defined(header).and..not.%defined(footer)
	spgraph(vfields=2,footer="CUSUM Tests")
else
	spgraph(vfields=2,footer=footer,header=header)
*
* Compute the cusum test statistic and its approximate critical value
*
sstats(max) startl endl abs(cusumstd)>>%cdstat
compute %signif=2*(1-%cdf(3*%cdstat)+exp(-4*%cdstat^2)*(%cdf(%cdstat)+%cdf(5*%cdstat)-1)-$
  exp(-16*%cdstat^2)*(1-%cdf(%cdstat)))
if %cdstat<.3
   compute %signif=1-.1465*%cdstat
*
* Graph with the upper and lower asymptotic .05 boundaries. Put the
* significance level in the top left corner.
*
set upper startl endl = .948*(1+2.0*target)
set lower startl endl = -upper
*
graph(hlabel="CUSUM test",grid=(t==%maxent)) 3
# cusum startl endl
# upper startl endl 2
# lower startl endl 2
grtext(position=upleft) "Sig Level "+%strval(%signif,"#.###")
*
* Compute the .05 critical value
*
compute n =0.5*%nobs-1
compute %%cusumqcv=1.3580150*n^-0.5-0.6701218*n^-1.0-0.8858694*n^(-1.5)
*
* Compute the cusumq statistic
*
sstats(max) startl endl abs(cusumsq-target)>>%%cusumq
*
* Graph with the target and the boundaries
*
set upper startl endl = %min(target+%%cusumqcv,1.0)
set lower startl endl = %max(target-%%cusumqcv,0.0)
graph(hlabel="CUSUMSQ test",grid=(t==%maxent)) 4
# cusumsq startl endl
# upper   startl endl 2
# lower   startl endl 2
# target  startl endl 3
spgraph(done)
end

