* * ROBUSTSTAR.RPF * Test for STAR with outlier adjustments. * * Adjustment polynomial for downweighting data points in STAR linearity * tests, from van Dijk, Franses and Lucas(1999), "Testing for smooth * transition nonlinearity in the presence of additive outliers", JBES, * vol 17, no 2, 217-235. * * Weights are 1 below c1 standard errors, 0 above c2 standard errors and * this polynomial in between. * procedure RobustPoly px type vect *px * option real c1 2.576 option real c2 3.291 * local rect a(6,6) local vect b(6) * do j=1,6 compute a(1,j)=%if(j<=2,0.0,(j-1)*(j-2)*c1^(j-3)) compute a(2,j)=%if(j<=2,0.0,(j-1)*(j-2)*c2^(j-3)) compute a(3,j)=%if(j<=1,0.0,(j-1)*c1^(j-2)) compute a(4,j)=%if(j<=1,0.0,(j-1)*c2^(j-2)) compute a(5,j)=c1^(j-1) compute a(6,j)=c2^(j-1) end do j compute b(1)=0.0,b(2)=0.0,b(3)=1.0,b(4)=0.0,b(5)=c1,b(6)=0.0 compute px=%solve(a,b) end * ****************************************************************** * * This is applied to the lynx data contaminated by outliers * cal(a) 1821 open data lynx.dat data(org=cols) 1821:1 1934:1 lynx set x = log(lynx)/log(10) diff(center) x / xc * * Pick contaminated points * seed 532343 set outliers = %uniform(0.0,1.0)<.1 set xc = xc + %if(outliers,%ran(2.0),0.0) * * Get the weight transition polynomial at the default knots * @RobustPoly px * * Start with equal weights. Do 10 iterations on weighted least squares * on the AR model. * set weight = 1.0 do iters=1,10 linreg(weight=weight,noprint) xc # constant xc{1 to 11} compute sigma=sqrt(%seesq) set weight %regstart() %regend() = \$ stdu=abs(%resids)/sigma,%if(stdu<2.576,1.0,%if(stdu>3.291,0.0,%polyvalue(px,stdu)/stdu)) end do iters * * Apply the STAR test with the chosen weights to the contaminated data. * Also apply the equally weighted STAR test. * @startest(p=11,d=1,weights=weight,title="Outlier-Adjusted STAR Test") xc @startest(p=11,d=1) xc