Date: Mon, 5 Apr 2004 15:29:25 -0400
Reply-To: malamill@engmail.uwaterloo.ca
Sender: "SAS(r) Discussion" <SAS-L@LISTSERV.UGA.EDU>
From: Mauricio Alamillo <malamillo@MYREALBOX.COM>
Subject: Re: SAS & Tukey 4253H- twice method
In-Reply-To: <200404051528.i35FSe414737@listserv.cc.uga.edu>
Content-Type: text/plain; charset="us-ascii"
Hello,
Thanks for the code, however, I've found for this particular problem that
SPSS turned out to be a quicker solution than SAS. Why? Well, SPSS has a
built in function to perform Tukey 4253H and I could save a lot of
programming, like as you demonstrate, it is a bit complicated for a beginner
like me in the middle of hurries. However I appreciate your code and will
keep it for studying, most likely others will do so. Still I desire to
learn how to program it in SAS and become good at it.
Mauricio
-----Mensaje original-----
De: SAS(r) Discussion [mailto:SAS-L@LISTSERV.UGA.EDU] En nombre de Chang Y.
Chung
Enviado el: Lunes, 05 de Abril de 2004 11:29 a.m.
Para: SAS-L@LISTSERV.UGA.EDU
Asunto: Re: [SAS-L] SAS & Tukey 4253H- twice method
On Fri, 2 Apr 2004 10:36:38 -0500, Mauricio Alamillo
<malamillo@MAILUP.NET> wrote:
>Hello everyone,
>
>I'm new on the list and full of questions. My first and most urgent one
>is regarding how to smooth data using Tukey 4253H twice method. Any
>hint will be appreciated. Thanks,
>
>Cheers,
>
Hi, Mauricio A.,
Found a S implementation of 4253H, twice at
http://lib.stat.cmu.edu/S/smoothers . It seems legit, so I translated it
into a sas data step and got the same result on a test case. Hope this
helps.
Cheers,
Chang
p.s. The sas code can be much improved in terms of speed. For example, the
current implementation copies out[*] array into in[*] before every call to
the running median routine. One way to speed things up here is to use a
macro variable IN and OUT for the array name every where and switch the
names (instead of copying).
<sasl:code>
data xes;
input x @@;
call symput('n', put(_n_, best12.-l)); * I am so lazy!;
cards;
569.0 416.0 422.0 565.0 484.0 520.0 573.0 518.0 501.0 505.0
468.0 382.0 310.0 334.0 359.0 372.0 439.0 446.0 349.0 395.0
461.0 511.0 583.0 590.0 620.0 578.0 534.0 631.0 600.0 438.0
516.0 534.0 467.0 457.0 392.0 467.0 500.0 493.0 410.0 412.0
416.0 403.0 422.0 459.0 467.0 512.0 534.0 552.0 545.0
run;
%put n=&n.;
proc transpose data=xes out=xes_transposed;
run;
data xes_4253Htwiced;
set xes_transposed;
array xes [1:&n.] col:;
array in [1:&n.] _temporary_;
array out [1:&n.] _temporary_;
array once [1:&n.] _temporary_;
n = &n.;
* main -------------------------------------------------------------------;
do i = 1 to n; in[i] = xes[i]; end;
first = xes[1]; last = xes[n]; * for evenup;
link smooth; * once;
do i = 1 to n; once[i] = out[i]; end;
do i = 1 to n; in[i] = xes[i] - once[i]; end;
first = in[1]; last = in[n]; * for evenup;
link smooth; * twice;
do i = 1 to n; out[i] = once[i] + out[i]; end;
do i = 1 to n; * output;
x = xes[i];
sx = out[i];
keep x sx;
output;
end;
stop;
* subs -------------------------------------------------------------------;
smooth: * smooth in[*] with 4253H (once) and output out[*];
link m4;
link flip; link m2;
link flip; link evenup;
link flip; link m5;
link flip; link m3;
link flip; link evr;
link flip; link hann;
return;
flip: * copy out[*] to in[*];
do i = 1 to n;
in[i] = out[i];
end;
return;
m4: * running medians of 4 with decreased span at ends and end-value copy
on;
do i = 1 to n-1;
select ( i );
when ( 1 ) out[ 1 ] = median(in[ 1 ], in[2]);
when (n-1) out[n-1] = median(in[n-1], in[n]);
otherwise out[ i ] = median(in[i-1], in[i], in[i+1], in[i+2]);
end;
end;
out[n] = .;
return;
m2: * running medians of 2 with decreased span at ends and end-value copy
on;
do i = 1 to n-2;
out[i] = median(in[i], in[i+1]);
end;
out[n-1] = .;
return;
evenup: * restore end values to the out, from first and last;
out[1] = first;
do i = 2 to n-1;
out[i] = in[i-1];
end;
out[n] = last;
return;
m5: * running medians of 5 with decreased span at ends and end-value copy
on;
do i = 1 to n;
select ( i );
when ( 1 ) out[ 1 ] = in[1];
when ( 2 ) out[ 2 ] = median(in[ 1 ], in[ 2 ], in[3]);
when (n-1) out[n-1] = median(in[n-2], in[n-1], in[n]);
when ( n ) out[ n ] = in[n];
otherwise out[ i ] = median(in[i-2], in[i-1], in[i], in[i+1], in
[i+2]);
end;
end;
return;
m3: * running medians of 3 with end-value copy on;
do i = 1 to n;
select (i);
when (1) out[1] = median(in[ 1 ], in[1], in[ 2 ]);
when (n) out[n] = median(in[n-1], in[n], in[ n ]);
otherwise out[i] = median(in[i-1], in[i], in[i+1]);
end;
end;
return;
evr: * end value rule applied to ends;
do i = 1 to n;
out[i] = in[i];
end;
out[1] = median(3 * out[ 2 ] - 2 * out[ 3 ], out[1], out[ 2 ]);
out[n] = median(3 * out[n-1] - 2 * out[n-2], out[n], out[n-1]);
return;
hann: * 1/4, 1/2, 1/4, with lost end value copy on;
out[1] = in[1];
do i = 2 to n-1;
out[i] = (in[i-1] + 2 * in[i] + in[i+1]) / 4;
end;
out[n] = in[n];
return;
run;
proc print data=xes_4253Htwiced noobs;
format sx 8.4;
run;
/*
x sx
569 491.3750
416 491.3750
422 491.3750
565 498.8828
484 514.9375
520 524.6602
573 525.0352
518 521.1602
501 512.5742
505 493.1680
468 449.7422
382 391.6133
310 353.4297
334 343.8438
359 355.1602
372 382.7930
439 405.5469
446 411.8633
349 411.5586
395 420.9375
461 456.1250
511 513.8516
583 565.2422
590 589.4688
620 594.7188
578 594.5625
534 591.8125
631 583.8438
600 569.0313
438 546.3438
516 517.2578
534 489.6445
467 471.2383
457 463.4844
392 464.1875
467 468.4688
500 470.6094
493 462.2617
410 438.5703
412 416.1094
416 408.8711
403 412.1836
422 424.8750
459 448.1445
467 478.7578
512 510.0234
534 534.1250
552 547.0000
545 550.9375
*/
</sasl:code>