|
SAS软件包是一个大型的统计分析系统,其核心是多个用于实现统计分析的实用过程。统计分析离不开操作的数据对象, SAS采用的办法是建立SAS数据集。而实现上述的一切都需要编制SAS引导程序,简称SAS程序。
SAS程序由一系列符合SAS语言语法规则的语句组成, 正如用任何一门计算机语言编制的程序一样。SAS语言不仅提供了一般程序设计语言拥有的语句(如循环控制,条件判断,赋值,输入输出等),而且,其丰富的概率函数、分位数函数、样本统计量函数以及随机数函数更是其他语言无法比拟的。SAS程序的基本组成是: 若干SAS语句组成数据步(DATA步),若干SAS语句组成过程步(PROC步),若干DATA步和若干PROC步组成一个完整的SAS程序,其中,DATA步通常产生SAS数据集, 而PROC步则对SAS数据集内的数据进行处理并输出结果或产生新数据集。
主要介绍如何产生数据集,而且由于数据集大部分由数据步完成,我们把精力也主要集中于数据步上。
第1节 数据步流程
[例1.5.1] 某小学10名9岁男学生6个项目的智力测验得分资料列于下表中。
表1.5.1 某小学10名9岁男学生6项智力测验得分结果
────────────────────────────────────
被测试者 常识 算术 理解 填图 积木 译码
编号 ID X1 X2 X3 X4 X5 X6
────────────────────────────────────
1 14 13 28 14 22 39
2 10 14 15 14 34 35
3 11 12 19 13 24 39
4 7 7 7 9 20 23
5 13 12 24 12 26 38
6 19 14 22 16 23 37
7 20 16 26 21 38 69
8 9 10 14 9 31 46
9 9 8 15 13 14 46
10 9 9 12 10 23 46
────────────────────────────────────
1. 变量和观测
这是SAS数据集的2个基本概念。可以这样看待它们的关系,SAS对各变量的操作都是在各观测内进行的。从每一个观测对象身上观测到n个变量的具体取值,在SAS中, 把这n个数值写在一行上,称为1个观测。如表1.5.1,共有10个观测,每个观测包含7个变量(含编号)的取值。
2. 建立SAS数据集方法(先看下面的第1个程序)
[SAS程序] [EXAMPLE.SAS] [说明] (1)DATA语句指明我们要建立一个名字为A的
DATA a; 数据集,当该段程序执行完以后,我们可以在DOS目录C:\SAS
INPUT id x1-x6; \SASWORK下看到一个文件 A.SSD,其中A为我们在SAS程序中
CARDS; 指定的数据集名,SSD为SAS系统自动提供的扩展名, 表明该
1 14 13 28 14 22 39 文件为SAS数据集(SAS DATA SET), SAS数据集是SAS系统的
2 10 14 15 14 34 35 内部文件,只能被SAS系统正常引用。
3 11 12 19 13 24 39 (2)INPUT语句指明了将要建立的SAS数据集包含的变量
4 7 7 7 9 20 23 ,其中,ID代编号,X1—X6分别代表常识、算术、理解、填
5 13 12 24 12 26 38 图、积木、和译码这6个方面的得分。 值得注意的是: 象
6 19 14 22 16 23 37 X1-X6这样的简写的列表能够被SAS系统接受, 为SAS程序的
7 20 16 26 21 38 69 书写提供了极大的方便,程序也显得十分简洁。
8 9 10 14 9 31 46 (3)CARDS语句表明下面是用于建立SAS数据集的数据流
9 9 8 15 13 14 46 ,数据流以另起一行的一个分号为结束。参见第1篇第4章第
10 9 9 12 10 23 46 3节。
; (4)RUN语句标志着一个DATA步或PROC步的结束, SAS系
RUN; 统一遇到它,即开始建立SAS数据集或调用过程分析数据。
(5)在程序运行的过程中,我们可以在LOG窗口看到这样的信息:
NOTE:The data set work.a has 10 observations and 7 variables.
它告诉我们,数据集A包含10个观测和7个变量。
(6)每个语句以分号为结束。
下面我们再来看看数据步具体是怎样工作的(数据步流程):
每个DATA步开头的DATA语句总是标志这一步的开始,余下的语句被SAS转换为计算机的机器语句,而且每次通过DATA步执行它们──通常对输入数据中的每个观测执行一次。
在DATA步中所有被提到的变量在SAS内部都被变成当前向量的一部分, 也称为PROGRAM DATA VECTOR(程序数据向量)。即来自输入数据的变量同程序语句产生的变量一起,都在程序数据向量里,而且可以通过程序语句对它们进行操作。如: INPUT X; Y=X+1;
先通过INPUT语句从数据流中获得变量X的值,接着通过赋值语句产生一个新变量Y,它的当前值是X的当前值加1,X和Y都在程序数据向量中。如果不特意指定,DATA步产生的数据集将包含程序数据向量中的所有变量。
用INPUT语句读入或用编程语句产生的变量的初始值,在这个DATA步每执行一次之前置为缺失,除非用RETAIN语句指明,在后面的应用中我们将会看到RETAIN语句将会起什么作用。用SET,MERGE或UPDATE语句读入的变量的初始值自动置为保留。
当SAS执行DATA步的最后一个语句时(或让返回到这一步的开头处理新观测的RETURN语句时),它按常规把程序数据向量的当前值输出到正被创建的SAS数据集中。然后SAS返回到DATA语句之后的第一个语句, 在程序数据向量中初始值不被保留的变量置为缺失。接着开始执行语句来建立下一个观测。如果要在DATA步的一次执行中输出多个观测, 那矩须使用OUTPUT语句。
当SAS已经读完来自输入文件的所有数据时,它继续下一个DATA步或PROC步。
SAS通过DATA步的这些语句处理读入的每个观测,虽然每个语句仅出现一次,但SAS对每个观测都执行这些语句。不包含INPUT,SET,MERGE或UPDATE语句的DATA步只执行一次,否则DATA步中的语句被重复执行直到输入的原始数据执行完了,或直到STOP或ABORT语句被执行。程序语句还可能含有让DATA步中的某些语句执行很多次的语句,如DO循环等。
SAS变量_N_是SAS为每个DATA步自动设置的,它在某时刻的取值表示这个DATA步已经执行的次数。在DATA步中可以引用这个变量,但它不会加到正创建的数据集中。
以上这些原则就是我们编制程序的基础。
[EXAMPLE.SAS]显然是一个最基本的DATA步程序,这里DATA步一共执行了10次, _N_ 的值从1变到10,但这并非最佳的方法。容易看出,在10个观测中,变量ID的值正好从1到10,每增加一个观测, ID值增1。于是我们至少还可以写出3个SAS程序使得CARDS语句后的数据流中不含
[EXAMPLE.SAS]中数据的第1列。
程序二: 程序三: 程序四:
DATA a; DATA a; DATA a;
DO id=1 TO 10; id+1; id=_N_;
INPUT x1-x6; INPUT x1-x6; INPUT x1-x6;
OUTPUT; CARDS; CARDS;
END; 数据流 数据流
CARDS; ; ;
数据流 RUN; RUN;
;
RUN;
[程序二说明] ①DO语句是SAS程序中使用频繁的循环控制语句, SAS规定DO语句必须和END语句配合使用,DO和END之间的语句构成1个DO组。本例中, DO组内的语句被执行10次。ID称为下标变量,如果ID值每次递曰是1,而是别的值,如2,则DO语句写成:
DO ID=1 TO 10 BY 2; DO组内语句执行5次。
②本程序中DATA步只执行了一次,_N_的值始终保持为1, 在循环内部使用了OUTPUT语句,故共输出10个观测,程序中创建的变量一起被输出。
[程序三说明] 累加语句ID+1产生这样的效果,DATA步执行第1次时,先对ID赋初值0,然后对ID值加1;DATA步执行第2次时,上次执行结束时保留下来的ID值加1,……,如此往复,也就是说,DATA步第1次执行时ID=1,第2次ID=2,第3次ID=3,……,这正是我们想要达到的目的。
[程序四说明] 显然,程序三和四中DATA步都执行了10次, 程序四中还充分利用了_N_的值。以上4个程序我们都可以加上如下的两条语句: PROC PRINT DATA=A; RUN;
它们调用SAS的打印过程,输出数据集A,都得到下面的结果。
OBS ID X1 X2 X3 X4 X5 X6 OBS ID X1 X2 X3 X4 X5 X6
1 1 14 13 28 14 22 39 6 6 19 14 22 16 23 37
2 2 10 14 15 14 34 35 7 7 20 16 26 21 38 69
3 3 11 12 19 13 24 39 8 8 9 10 14 9 31 46
4 4 7 7 7 9 20 23 9 9 9 8 15 13 14 46
5 5 13 12 24 12 26 38 10 10 9 9 12 10 23 46
(输出结果的第1部分) (输出结果的第2部分)
[说明] 通过这个简单的例子我们可以看到:正确理解数据步流程, 合理运用SAS语句将为我们带来不少方便,使编出来的程序简洁明了。
第2节 创建数据集的途径
SAS 数据集中数据的来源大致有以下几种方式。
1. 将数据行直接写在CARDS语句后
上一节中我们采用的就是这种办法。
2. 从外部文件(文本文件)读取数据
沿用[例1.5.1],如果事先将数据输入计算机并做成文本文件, 取名为example.txt,数据格式和上节中程序二、三、四的相同,那么建立数据集的程序将更加简洁,即:
DATA a; INFILE 'example.txt'; id+1; INPUT x1-x6; RUN;
[说明] INFILE语句将当前目录下的数据文件EXAMPLE.TXT调入。 INFILE语句有很多选择项,其中很有用的两项是RECL=n定义逻辑记录长度,n的值从1至32767; RECFM=N把整个文件作为一个大记录处理。例如,如果example.txt中的数据都写在同一行上, 那么上面程序中的第2句可改写成:infile 'example.txt' lrecl=200; 或infile 'example.txt' recfm=n; 同时第4句改写成input x1-x6 @@;(从一行读入多个观测)。
3. 利用已经创建的数据集产生所需的新数据集
前面已经创建了数据集A,包含7个变量,如果还想增加1个变量,它的值是x1-x6之和,可利用已创建的数据集A,重新建立一个数据集,程序如下:
DATA b; SET a; ss=SUM(OF x1-x6); RUN;
[说明] SET语句告诉SAS系统从1个或几个数据集中读取观测值,SUM是1个SAS函数,作用是对变量求和,它的一般写法是SUM(X1,X2),即变量间用逗号隔开, 但当自变量是相继的一些变量时,有如上的简略写法。类似SUM这样的样本统计量函数很多, 如MAX(最大值)、MEAN(均值)、STD(标准差),这种近乎自然英语的写法为程序的书写带来了极大的方便,可以说这也是SAS的一大特色。
4. 其他软件产生的标准格式文件与SAS数据集之间的互相转换
SAS提供了2个过程,即DBF和DIF,分别用于dBASE数据库和SAS数据集、DIF文件和SAS数据集之间的相互转换。此处仅以DBF过程说明如下:
假设我们已经用dBASEⅢ或FOXBASE做好了一个DBF文件example.dbf,它包含[例1.4.1]中的所有信息,并且,存在当前目录下,则下面的程序可将它转换成SAS数据集。
FILENAME exam 'example.dbf'; PROC DBF DB3=exam OUT=a; RUN;
这里,exam称为文件关联名。运行结束后,新产生的SAS数据集A便在内存之中, 只要在某个过程名之后写上DATA=A,就可调用。如: PROC PRINT DATA=A; RUN;
反之,如果事先做好了SAS数据集,下面的程序可将它转换成DBF文件。
PROC DBF DB3=exam DATA=a; RUN;
运行结束后,在SAS子目录下自动生成1个名为exam.dbf的dBASEⅢ文件。
[说明] FILENAME语句用于建立1个文件名标记, 在一次SAS作业中一旦建立了1个标记,它在本次作业中均有效,直到退出本次作业为止。要引用某文件,只需引用其文件名标记。
5. 某些SAS过程自动产生数据集
很多过程都可以产生新的数据集,有利于数据的再分析。如下面的程序就产生1个包含x1-x6的均值mx1-mx6的数据集M。
PROC MEANS DATA=a; VAR x1-x6; OUTPUT OUT=m MEAN=mx1-mx6; RUN;
第3节 建立数据集的技巧
SAS软件包的功能十分强大,不仅仅表现在它用于统计分析的各个过程, 还在于它内容丰富的SAS语言,为灵活地加工处理数据提供了极大的方便。要掌握SAS语言,当然需要花一些时r间,但正确理解以下几点,也许是关键所在。
1. 含有相同变量集合的两个数据集的纵向(即前后)联接
沿用[例1.5.1],在已建立了数据集A的基础上,如果逾了几名学生的智力得分, 做成数
据文件add1.dat,要把它添加到A中,程序如下。
DATA a; RETURN;
IF last=0 THEN DO; END;
SET a END=last; INFILE 'add1.dat';
OUTPUT; INPUT id x1-x6; RUN;
(程序的第1部分) (程序的第2部分)
[说明] ①END=是SET语句的选择项,跟在它后面的LAST是任选的1个SAS变量名, 它的值一开始被置为0,直到SET语句读完输入数据集(1个或多个)的最后1个观测时, 其值为1。变量
LAST不包含在新创建的数据集中。
②SAS从数据集A中读入1个观测并重新输入A后,本来应该顺序执行下面的INFILE和INPUT语句,但它执行完OUTPUT语句后又遇到一个RETURN语句,RETURN语句的作用就是让SAS重新回到DATA步的开始,于是重新判断LAST是否为0,由于A中还有观测,LAST为0,故又读入A中一个观测并输出,遇到RETURN语句又返回,直到读完A中所有观测,LAST为1,才跳到INFILE语句, 读入外部文件中数据。
③SET语句的选择项还有POINT=,NOBS=,前者是输入数据集中想要处理的观测序号, 后者是输入数据集的观测总数,用方法与END=类似,限于篇幅,不再赘述。
2. 含有观测数目相同的2个数据集的横向(即左右)联接
沿用[例1.5.1],如果对10名学生又增加了两项智力测验,其变量名为x7,x8, 将得分做成文件add2.dat,要把它添加到A中,程序如下。
DATA a; SET a; INFILE 'add2.dat'; INPUT x7 x8; RUN;
运行上面2个程序可分别得到不同的数据集,但它们的共同的特点是数据都来自2个地方,其一是已经创建的数据集;其二是外部数据文件。这可以说是第2节中创建数据集的不同方法的混合使用。理解了数据步流程,这样的程序也就不难理解了,而且可以根据需要任意组合。
3. 将1个大数据集拆成几个小数据集
有时需将1个大数据集按特定的条拣成几个小数据集。沿用[例1.5.1], 假定要由数据集A拆成2个数据集,分别包含编号为奇数亨号为偶数的学生的结果,可用如下的程序实现。
DATA b; IF k=INT(k) THEN OUTPUT a1;
SET a; ELSE OUTPUT a2;
k=id/2; RUN;
(程序的第1部分) (程序的第2部分)
[说明] INT( )为取整函数,当id为偶数时,k=id/2的结果与对它取整的结果一致, 满足此条件的所有观测一一输出到数据集a1中去, 与id为奇数对应的所有观测便输出到数据集a2中去。
4. 运用数组语句和控制语句产生所需的数据集
数组当然是任何一门程序设计语言都缺少不了的一个概念。比较特别的是SAS 的数组分为2类:一类是普通的,每1个数组元素代表1个变量;另一类是临时的, 数组元素不代表任何变量,当然也就不能出现在输出数据集上,优点是计算速度快, 占用内存少。定义临时数组只需加上_TEMPORARY_选择项,如: ARRAY x(10)_TEMPORARY_;便是定义了1个含10个元素的临时数组x。
控制语句包括DO,END,SELECT,IF,GOTO,LINK,RETURN等,最常用的是DO,END和IF语句, 有些语句的作用可被其他语句的组合所代替。下面举例说明。
如果要把在第1节建立的数据集A中的10个观测变成1个观测,而且信息不丢失, 那么新建立的数据集当然需要60个变量(为简单起见,数据集中不再保留变量ID)。实现的程序如下。
DATA b; DO i=1 TO 6;
SET a; y(_N_,i)=x(i);
ARRAY x(6) x1-x6; END;
ARRAY y(10,6); IF _N_=10;
RETAIN y1-y60; KEEP y1-y60; RUN;
(程序的第1部分) (程序的第2部分)
[说明] ①程序中定义了2个数组,前一个数组为一维数组,不但给出了数组名,而且给出了数组元素对应的变量x1-x6,事实上这里可以省略;第2个数组为二维数组,省略了变量名,按规定缺省值为y1-y60。
②RETAIN语句的作用是使变量的值在处理后面观测时得到保留,因为这里数据步执行一次,只对数组Y中6个元素赋值,于是只有到数据步第10次执行时,60个元素才全部赋值完毕; 加上一个子集IF语句,表明当数据步第10次执行时才输出当前的程序数据向量的数据到新数据集b中去;KEEP语句的作用是使新数据集中只保留y1-y60这60个变量。
反之,要从数据集B重新得到A,只需在循环中用OUTPUT语句强制输出观测即可,程序如下:
DATA a; DO id=1 TO 10; OUTPUT;
SET b; DO j=1 TO 6; END;
ARRAY x(6); x(j)=y(id,j); DROP j y1-y60;
ARRAY y(10,6); END; RUN;
(程序的第1部分) (程序的第2部分) (程序的第3部分)
[例1.5.2] 在培养钩端螺旋体时,除已固定若干因素外,拟研究以下4个因素不同水平的效应,求其最佳组合。A:血清种类─兔、胎盘;B:血清浓度─5%、8%; C:基础液─缓冲剂、
蒸馏水、自来水;D:维生素─加、未加。
表1.5.2 2×2×3×2析因试验的钩端螺旋体计数结果
────────────────────────────────────
A与B 钩端螺旋体计数
(血清种类 ──────────────────────────────
与血清浓度) C与D: 缓冲剂(加 未加) 蒸馏水(加 未加) 自来水(加 未加)
────────────────────────────────────
兔血清: 5% 1426 648 684 1763 1182 580
… … … … … …
8% 1260 1144 875 1447 1220 1789
… … … … … …
胎盘血清: 5% 604 830 867 920 1243 1126
… … … … … …
8% 1108 578 1115 933 1283 685
… … … … … …
────────────────────────────────────
把上面表格中的数据做成文本文件ANOVA.DAT存放于当前目录下,数据格式如下。
[ANOVA.DAT] [SAS程序]
1426 648 684 1763 1182 580 DATA anova;
1183 1264 1430 1241 1512 1026 INFILE 'anova.dat';
2000 1398 1165 1381 1450 1026 DO a=1 TO 2;
1621 909 2022 2421 1385 830 DO b=1 TO 2;
1260 1144 875 1447 1220 1789 DO n=1 TO 4;
1599 1877 2250 1883 1095 1215 DO c=1 TO 3;
1410 1671 1871 1896 1700 1434 DO D=1 TO 2;
2416 1845 1962 1926 2372 1651 INPUT x @@;
604 830 867 920 1243 1126 OUTPUT;
1081 853 771 709 1115 1176 END;
487 441 403 848 416 1280 END;
624 1030 370 574 533 1212 END;
1108 578 1115 933 1283 685 END;
886 669 698 1024 1142 546 END;
831 643 791 1092 677 595 RUN;
1159 1002 559 742 534 566
这样做的好处是显而易见的。首先数据格式保持了与原来表格中的一致,其次,可用任何
编辑软件来建立这样的数据文件,还可以在DATA步中以不同的方式来读外部数据文件,用来建立不同的SAS数据集,以适应对同一批数据用不同方法进行分析。
运行右边的程序可得到所需的数据集,调用SAS中的ANOVA或GLM过程可实现方差分析。
[说明] 在SAS的方差分析过程中, 试验的每个因素的不同水平都在用于分析的SAS数据集中以一个变量的不同取值来表示, 这种变量在方差分析中被称为分组变量。例如从表中任取一个数据867,它是血清种类为胎盘血清,血清浓度为5%,基础液为蒸馏水,维生素加上时得到的钩端螺旋体计数,即A=2,B=1,C=2,D=1条件下的观测结果。为对每个数据的分组变量都赋以正确的值,我们在数据步中采用循环嵌套的办法,这样一共形成2×2×4×3×2=96个观测, 其中4为每种组合下的重复试验次数。在循环嵌套中,外层下标变量的值变化慢, 内层下标变量的值变化快。容易看出,DATA步执行完后得到的数据集有如下结构:
OBS A B N C D X ……………………………………
1 1 1 1 1 1 1426 94 2 2 4 2 2 742
2 1 1 1 1 2 648 95 2 2 4 3 1 534
3 1 1 1 2 1 684 96 2 2 4 3 2 566
(数据集的第1部分) (数据集的第2部分)
第4节 宏处理简介
1. 宏变量
宏变量(有时也称符号变量)属于SAS宏语言的范畴,和数据步中的变量的概念是不一样的。除了数据行外,你可以在SAS程序的任何地方定义和使用宏变量。数据步变量是和数据集相联系的,而宏变量是独立于数据集的。数据集变量的值取决于正在处理的观测,而一个宏变量的值总是保郴变,直到被明确改变。定义一个宏变量的最简单的办法是使用宏语句%LET,如%let dsn=Newdata; DSN就是宏变量的名字,Newdata是它的值。宏变量命名遵从一般的SAS命名规则,宏变量的值是一个字符串。
要引用宏变量的值,在宏变量的前面放“&”,如title "Display of Data Set &dsn";
宏处理器用DSN的值去代替&dsn,结果得到TITLE "Display of Data Set Newdata";
值得注意的是这里的标题必须用双引号括起来,而不能用单引号,因为宏处理器只对双引号中引用的宏变量进行这种处理,而把单引号中的所有字符都看作是标题的内容。
在一个SAS程序中,你可以根据需要引用宏变量任意次, 在你改变它之前它的值一直保持不变。下面的程序就引用了DSN2次。
DATA temp; SET &dsn; IF age>=20;
PROC PRINT; TITLE "Subset of Data Set &dsn"; RUN;
每次&dsn一出现,宏处理器就用Newdata代替它,SAS实际执行的语句是DATA TEMP;
SET NEWDATA; IF AGE>=20;
PROC PRINT; TITLE "Subset of Data Set Newdata"; RUN;
你也可以创建一个包含SAS语句的宏变量,如:
%LET plot=%str( 这种情况下,宏变量的值都作为函数%STR的实参, 于是分
PROC PLOT; 号成为宏变量值的一部分,而不会被当做%LET语句的结束。
PLOT income*age; 要改变宏变量的值,只需用%LET语句重定义, 宏变量也可
RUN; 以赋空值,如:
) %LET dsn=Nextdata; %LET plot=;
第1个语句对DSN重新赋值Nextdata,第2个语句对PLOT赋空值。宏变量也可以嵌套引用:
%LET dsn=Olddata; DATA temp;
%LET yvar=Income; SET &dsn;
%LET xvar=Age; IF age>=20;
%LET plot=%str( &plot
PROC PLOT; PROC PRINT;
PLOT &yvar*&xvar; TITLE "Subset of Data Set &dsn";
RUN; RUN;
)
(程序的第1部分) (程序的第2部分)
SAS实际执行的是如下的语句:
DATA TEMP; SET OLDDATA; IF AGE>=20;
PROC PLOT; PLOT YVAR*XVAR; RUN;
PORC PRINT; TITLE "Subset of Data Set Olddata"; RUN;
2. 宏
宏就是存贮的1个文本,最简单的宏工作起来很象1个宏变量,但是复杂的宏可以做许多宏变量无法完成的事。下面是1个最简单的宏定义。
%MACRO dsn; Newdata %MEND dsn;
宏定义总是由%MACRO语句开始,而且必须包含1个宏名。宏名遵从一般的SAS命名规则,这里的宏名叫DSN。Newdata是宏的内容, %MEND作为宏定义的结束语句;这里的 %MEND语句为了清晰还重复了宏名。要调用1个宏,只需放1个百分号%在宏名的前面,如:
TITLE "Display of Data Set %dsn";
宏处理器执行宏展开,用宏内容去代替%dsn,TITLE语句成为:
TITLE "Display of Data Set Newdata"; 同样,标题必须用双引号括起来。
一个SAS程序可以包含多个宏,一个宏也可以被多次调用,要改变一个象DSN这样简单的宏的内容,只需新重定义宏,对简单的文本代换,用宏变量比定义一个宏效率要高,但是当任务比较复杂的时侯,宏救宏变量要优越得多了。下面的程序创建一个包含整段SAS程序的宏。%MACRO plot; PROC PLOT; PLOT income*age; RUN; %MEND plot;
以后调用宏如下: SAS实际执行的语句是:
DATA temp; DATA TEMP;
SET olddata; SET OLDDATA;
IF age>=20; IF AGE>=20;
%plot PROC PLOT;
PROC PRINT; PLOT INCOME*AGE;
RUN; RUN;
PROC PRINT; RUN;
假设在PROC PLOT中的绘图变量可以改变,你就可以用宏变量引用去替换PLOT语句中的变量名,然后在调用宏之前用%LET语句给宏变量赋值。
%MACRO plot; SAS 实际执行的是如下的语句:
PROC PLOT; DATA TEMP;
PLOT &yvar*&xvar; SET OLDDATA;
RUN; IF AGE>=20;
%MEND plot; PROC PLOT;
DATA temp; PLOT INCOME*AGE;
SET olddata; RUN;
IF age>=20; PROC PLOT;
%LET yvar=income; PLOT INCOME*YRS_EDUC;
%LET xvar=age; RUN;
%plot PROC PRINT;
%LET xvar=yrs_educ; RUN;
%plot
PROC RRINT;
RUN;
宏变量和宏结合起来使用,为灵活设计提供了极大的方便。当程序较长的时侯,就象上面的程序一样,总是使用%LET语句就显得不那么方便了,这时可以使用宏参数, 即定义宏变量作为%MACRO语句的一部分。
%MACRO plot(yvar,xvar); PROC PLOT; PLOT &yvar*&xvar; RUN; %MEND; %MACRO语句的括号中定义的宏变量就是宏参数,当调用宏的时侯,需将值传递给宏参数, 如: %plot(income,age)
宏处理器把第1个值和第1个宏变量匹配,第2个值和第2个宏变量匹配……因此,这些参数叫做位置参数,于是产生下面语句:
PROC PLOT; PLOT INCOME*AGE; RUN;
假设你经常创建数据集TEMP和绘制PLOT图,你可以将数据步做成一个宏,而过程步做成另一个宏,然后在第3个宏里调用它们。注意最好加上注释说明你要做什么。
%MACRO create; %MACRO analyze(yvar,xvar);
Data temp; %*Create the Data Set;
SET olddata; %create
IF age>=20; %*Plot the Variable Selected;
%MEND create; %plot
%MACRO plot; %MEND analyze;
PROC plot;
PLOT &yvar*&xvar;
RUN;
%MEND plot; (右边的程序是左边程序的继续)
在宏ANALYZE中,由“%*”开始的语句是宏注释语句, 它们用于对宏内做的事情作注解。
当宏较长的时侯,是很必要的。当你定义完这些宏后,就可以通过调用宏来运行程序。
%analyze(income,age)
如果你想改变CREATE,重新定义它,则:
%MACRO create;
DATA temp2; SET olddata(obs=100); IF age>=20;
%MEND create;
然后调用ANALYZE,方法是: %analyze(income,age)
当然,你也可以把整个程序段都放入ANALYZE,而不去调用其他宏,但是,定义CREATE和PLOT的优点是你可以在调用ANALYZE之前分别重新定义和调试它们,在规模较大的宏中, 模块结构是很重要的。
当你使用宏ANALYZE时,假设有时你需要执行宏CREATE中的DATA步,而有时又可直接从PROCPLOT开始,重新定义ANALYZE,使用%IF-%THEN语句。
%MACRO analyze(getdata,yvar,xvar);
%IF &getdata=yes %THEN %create; %plot
%MEND analyze;
如下调用ANALYZE,方法为:
%analyze(yes,income,age) %analyze(no,income,yrs_educ) 产生:
DATA TEMP; SET OLDDATA; IF AGE>=20; PROC PLOT; PLOT INCOME*AGE; RUN;
PROC PLOT; PLOT INCOME*YRS_EDUC; RUN;
如果在%THEN从句中有多个语句, 如果我们调用ANALYZE,方法如下:
使用%DO-%END组,如: %analyze(YES,income,age)
%MACRO analyze(getdata,yvar,xvar); 那么产生的语句是:
TITLE; PROC PLOT;
%IF &getdata=yes %THEN %DO; PLOT INCOME*AGE;
%create RUN;
TITLE3 "Data Set Created for This Plot";
%END;
%plot
%MEND analyze;
宏处理器并没有调用宏CREATE去创建一个数据集, 因为宏处理器在把宏展开时并不把小写字母转换为大写。于是%IF语句作的比较是YES=yes,当然为假,事与愿违。为了在调用时不管大小写都正确,应该修改ANALYZE如下:
%MACRO analyze(getdata,yvar,xvar);
%IF %UPCASE(&getdata)=YES %THEN %create;
%plot
%MEND analyze;
其中宏函数%UPCASE把小写转换为大写。
假设你想创建一系列名字用于SAS语句, 如DATA或VAR语句,可以借助宏来实现。
%MACRO names(name,number)
%DO n=1 %TO &number; &name&n %END;
%MEND names;
宏NAMES把参数NAME的值和宏变量的值连接起来形成一系列名字,参数NUMBER的值是宏变量N的终止值。于是在DATA步中调用NAMES,方法如下:
DATA %names(dsn,5); 产生语句: DATA DSN1 DSN2 DSN3 DSN4 DSN5;
3. 数据步中的宏变量
假设你想把一个数据集的观测数放进一个FOOTNOTE语句中,首先定义一个宏,如:
%MACRO create;
DATA temp;
SET olddata END=final; IF age>=20 THEN DO; n+1; OUTPUT; END;
IF final THEN CALL SYMPUT('number',n); RUN;
%MEND create;
以后一旦调用CREATE,则CALL SYMPUT语句将计数变量N的值传递给宏变量NUMBER,NUMBER可以是以前创建的宏变量,也可以是当前的DATA步创建。
这个宏变量以后可以直接引用, 也可以在DATA步用SYMGET函数引用。从而实匣同的数
据步中的值传递,如:
x=SYMGET('number'); 数据步变量X的值就是数据集TEMP中的观测个数。
第5节 SAS实用程序举例
1. 方差分析前如何实现正态性检验和方差齐性检验
用于方差分析的资料要求满足正态性和方差齐性的假设,但用于作方差分析的SAS过程并没有提供这样的检验。下面的程序可用在对数据作方差分析之前, 检查数据是否符合条件。若符合,直接调用方差分析过程;若不符合,看是否能对数据作某种数据变换,使数据近似满足要求;若实在难以办到,只好考虑别的办法,如非参数法。
[SAS程序] [PREANOVA.SAS]
/*Test of normality in groups*/ DATA a2;
PROC SORT DATA=anova; SET a1;
BY a b c d; RUN; v=n-1;
PROC UNIVARIATE DATA=anova NORMAL; _type_=1;
VAR x; BY a b c d; RUN; buf=2*v*LOG(s);
u=1/v;
/* Bartlett's test for H.V. */ RUN;
/* H.V.=homogeneity of variance */ PROC MEANS DATA=a2 NOPRINT;
PROC MEANS DATA=anova NOPRINT; VAR _type_ ss v buf u;
VAR x; BY a b c d; OUTPUT OUT=a3 SUM=k t_ss t_v t1 t2;
OUTPUT OUT=a1 N=n CSS=ss STD=s; RUN;
RUN;
(程序的第1部分) (程序的第2部分)
DATA _NULL_;
FILE PRINT;
SET a3;
chi_sqr=LOG(t_ss/t_v)*t_v-t1;
adjusted=chi_sqr/(1+1/3/(k-1)*(t2-1/t_v));
p=1-PROBCHI(adjusted,k-1);
PUT @12 'Chi Square' @30 'Adjusted Chi Square' @55 'Prob>Chi'
#3 @12 chi_sqr 10.4 @30 adjusted 10.4 @55 p 7.4;
RUN;
(程序的第3部分)
[PREANOVA.SAS修改指导] ①本程序直接调用本章第3节中已建立的SAS数据集ANOVA。②先用SORT过程对该数据集按A,B,C,D进行排序,BY语句中出现的变量,排在前面的变化慢,后面的变化快。例如,假设只有2个变量出现在BY语句中,X在前,Y在后,即'BY X Y;', X有2个值1和2,Y有2个值5和6,则在排完序的数据集中,X=1且Y=5的观测排在最前面,其次是X=1,Y=6,再其次是X=2,Y=5,最后是X=2,Y=6。
③对排完序以后的数据集调用UNIVARIATE过程进行单变量分析,加上选择项NORMAL,表示要进行正态性检验;VAR指定要进行分析的变量;BY语句指明要分组进行分析, 即对A,B,C,D的各种组合都分别进行一次单变量分析。
④[输出结果之一]只是A,B,C,D的24种水平组合产生的24个输出的第1部分,限于篇幅,未全部列出。其中“W:Normal 0.987761 Prob<W 0.9245”是使用了NORMAL选择项的结果,它采用的是W检验,0.9245就是P值,P值越大,一般表明正态性越好。
⑤进行Bartlette方差齐性检验,首先使用了1个MEANS过程。这个过程可以输出许多常用的统计量,如均值、总和、标准差等。NOPRINT表明不产生打印输出, 一般该过程要产生一些输出结果。OUTPUT语句用选择项OUT=规定1个输出数据集名A1,该数据集包含代表观测数的变量N,代表偏差平和的变量SS,代标准差的变量S。在这个OUTPUT语句中,OUT=,N=, CSS=和STD=都是关键词,A1,N,SS和S是自己规定的数据集名轰量名。BY语句指明分组变量。
⑥根据Bartlette法,MEANS过程后紧跟着一个DATA步,里面使用了自然对数函数LOG(),SAS还提供了以10为底的常用对数LOG10()和以2为底的对数函数LOG2()。
⑦DATA步之后再跟1个MEANS过程,输出1个数据集A3,它包含数据集A2中变量_TYPE_,SS,V,BUF和U的总和K,T_SS,T_V,T1和T2。
⑧最后1个DATA步使用空数据集名_NULL_,表明该数据步中不产生数据集, 这样可以节省资源也节省时间,因为在该步中只想产生一些输出。首先用FILE语句指明输出的设备,FILE后面可以跟着1个外部文件名(建立1个外部文件),也可以跟'PRN输出到打印机),或跟1个文件名标记(作用同跟着文件名一样),本例中跟PRINT, 表示要输出到OUTPUT窗口。与FILE语句配合使用的是PUT语句,PUT语句指明要输出的内容,这些内容输出到PUT语句前面最近的1个FILE语句指明的设备。@为列指针,#为行指针,整个PUT语句的意思是说,先把指针移到12列, 输出字符串Chi Square, 接着把指针移到30列,输出字符串Adjusted Chi Square, 再到55列, 输出Prob>Chi,然后把指针移到第3行第12列, 输出变量CHI_SQR的值,输出格式为10.4,再到30列,输出ADJUSTED的值,格式10.4,最后到55列,输出P值,格式7.4。格式W.D的含义是:最多有W位(包括负号和小数点),其中小数部分有D位。
⑨程序还使用了概率函数PROBCHI,本例中,K-1为自由度,ADJUSTED为代表校正的卡方值的变量。由于PROBCHI的返回值为小于卡方值的概率, 但我们需要的却是大于等于该卡方值的概率,故还要经过一次转换。下面给晨分输出结果:
[输出结果之一]
------------------------------ A=1 B=1 C=1 D=1 -------------------------------
Variable=X
Moments
N 4 Sum Wgts 4
Mean 1557.5 Sum 6230
Std Dev 345.1478 Variance 119127
Skewness 0.520517 Kurtosis 0.199882
USS 10060606 CSS 357381
CV 22.16037 Std Mean 172.5739
T:Mean=0 9.025119 Prob>|T| 0.0029
Sgn Rank 5 Prob>|S| 0.1250
Num ^= 0 4
W:Normal 0.987761 Prob<W 0.9245
[输出结果之二]
Chi square Adjusted chi square Prob>chi
40.4205 36.2275 0.0391
[输出结果的解释] 因P=0.0391<0.05,在α=0.05的水平上拒绝方差相等的无效假设,故此需先对此资料作数据变换后方可进行方差分析或采用相应的非参数法分析。
2. 灵活地产生频数表
在实际工作中,频数表的使用是很广泛的。获得一批数据之后,绘制一个频数表, 可以获得对数据的总体的大致印象,对正确选用统计方法进行分析资料是很有好处的。
[SAS程序] [FREQ.SAS]
%LET n=120; m=1;
/*number of observations is 120*/ DO UNTIL (last<=start+m*s);
DATA _NULL_; m=m+1;
ARRAY f(15) _TEMPORARY_; END;
ARRAY sf(15) _TEMPORARY_;
ARRAY percen(15) _TEMPORARY_; n=SYMGET('n')+0;
ARRAY x(&n); DO i=1 TO n;
FILE PRINT; k=m;
INPUT x1-x&n; DO j=1 TO m;
DO i=1 TO 15; IF (start+(j-1)*s)<x(i)<=
f(i)=0; (start+j*s) THEN DO;
sf(i)=0; f(j)=f(j)+1;
END; sf(j)=sf(j)+1;
min_x=MIN(OF x1-x&n); k=j;
max_x=MAX(OF x1-x&n); END;
s=ROUND((max_x-min_x)/10,1); IF j>k THEN sf(j)=sf(j)+1;
start=min_x-1; END;
last=max_x+1; END;
(程序的第1部分) (程序的第2部分)
PUT @2 '______________________________________________________________________'
#2 @4 'Groups' @13 'Frequency' @27 'Cumulative Frequency'
@51 'Cumulative Percent'
#3 @2 '______________________________________________________________________';
DO j=1 TO m;
percen(j)=sf(j)/n*100; buf=start+(j-1)*s;
PUT #(j+3) @2 buf 6. '-' @16 f(j) 3. @37 sf(j) 4. @58 percen(j) 5.1;
END;
PUT #(m+4) @2
'______________________________________________________________________'
#(m+5) @2 'The Largest Value' @19 max_x 6.
@35 'The Smallest Value' @53 min_x 6.;
CARDS;
86 83 77 81 81 80 79 82 82 81 81 87 82 78 80 81 87 81 77 78 77 78 77 77
77 71 95 78 81 79 80 77 76 82 80 82 84 79 90 82 79 82 79 86 76 78 83 75
82 78 73 83 81 81 83 89 81 86 82 82 78 84 84 84 81 81 74 78 78 80 74 78
75 79 85 75 74 71 88 82 76 85 73 78 81 79 77 78 81 87 83 65 64 78 75 82
80 80 77 81 75 83 90 80 85 81 77 78 82 84 85 84 82 85 84 82 85 84 78 78
;
RUN;
(程序的第3部分)
[FREQ.SAS修改指导] ①数据集名用_NULL_,120个观测数据按一维数组的形式读入,SAS系统把它们作为1个观测中的120个变量来处理,这样就可以运用程序设计的许多方法来编程,因为在只有1个观测的情形下,SAS语言吼的计算机语言是非常相似的。
②我们的算法是:准备将数据分成10组左右, 为此先调用MAX和MIN函数计算最大和最小值;然后计算组距S(注:统计书中组距通常用H或i表示);为了形式上的美观,将S取为整数,使用ROUND函数;再对最小值减1,作为第1组的下限值,然后判断实际需要的组数M; 最后计算我们感兴趣的各组的频数、累计频数和累计百分率,以一种比较直观的形式输出。
③ROUND函数的作用是按要求的精度四舍五入,如:
ROUND(10.78,0.1)=10.8; ROUND(10.48,1)=10。
④本程序一般只需在加注释语句的地方(第1个语句处)对数据进行适当的修改,就可用于别的数据。如只有100个样本观测值,只需将120改成100即可。这是因为我们在程序一开始就定义了宏变量N,在后面的程序中多次引用,从这里我们可以看出宏变量的作用。
[输出结果]
表 1.5.3 巧用SAS语言编制的频数表之输出结果
______________________________________________________________________
Groups Frequency Cumulative Frequency Cumulative Percent
______________________________________________________________________
63- 2 2 1.7
66- 0 2 1.7
69- 2 4 3.3
72- 10 14 11.7
75- 29 43 35.8
78- 31 74 61.7
81- 29 103 85.8
84- 12 115 95.8
87- 4 119 99.2
90- 0 119 99.2
93- 1 120 100.0
______________________________________________________________________
The Largest Value 95 The Smallest Value 64
3. R×C表等级资料的Kruskal-Wallis检验、RIDIT分析和CPD分析
Kruskal-Wallis检验(即H检验)是一种经典的非参数统计方法,RIDIT分析为在H检验基础上发展起来,CPD分析为王广仪教授首创,并在他的专著中对3种方法有精彩的论述。下面是一个同时给出Kruskal-Wallis检验、RIDIT 分析和CPD分析结果的程序。
[例1.5.3] 某作者收集了厌色性垂体腺瘤患者手术治疗后视力变化情况资料,见下表。
表1.5.4 厌色性垂体腺瘤患者治疗后视力变化
────────────────────────────────────
例 数
疗 效 ──────────────────────────────
术前视力: 1级 2级 3级 4~5级 合计
────────────────────────────────────
恢复正常 2 3 1 0 6
进 步 51 70 60 45 226
不 变 15 38 27 53 133
恶 化 1 2 2 3 8
────────────────────────────
合计 69 113 90 101 373
────────────────────────────────────
注:术前视力为实验分组指标,疗效为指标分组指标。
[SAS程序] [KRCTEST.SAS]
%let group=4; /* 实验分组(术前视力)的组数为4 */
%let rank=4; /* 指标分组(疗效)的组数为4 */
%let col_row=%eval(&group*&rank); /* 上述2个数需根据具体情况修改 */
(程序的第1部分)
DATA _NULL_; n=n+tcol(j);
row=&rank; END;
col=&group; buf3=n;
ARRAY num(&rank,&group); DO i=1 TO row;
ARRAY cpd(&group); trow(i)=0;
ARRAY ridit(&group); DO j=1 TO col;
ARRAY trow(&rank); trow(i)=trow(i)+num(i,j);
ARRAY tcol(&group); END;
ARRAY y(&rank); buf4=buf4+trow(i)*trow(i)*trow(i);
ARRAY ucpd(&group); buf2=trow(i);
ARRAY uridit(&group); y(i)=buf3-buf1-buf2;
ARRAY t(&group); buf3=y(i);
FILE PRINT; buf1=trow(i);
n=0; buf1=0; END;
buf4=0; buf5=0; DO j=1 TO col;
INPUT num1-num&col_row; cpd(j)=0;
DO j=1 TO col; DO i=1 TO row;
tcol(j)=0; cpd(j)=cpd(j)+num(i,j)*y(i);
DO i=1 TO row; END;
tcol(j)=tcol(j)+num(i,j); ridit(j)=0.5-cpd(j)/(2*tcol(j)*n);
END; t(j)=tcol(j)*(n+1)/2-cpd(j)/2;
(程序的第2部分) (程序的第3部分)
ucpd(j)=cpd(j)/sqrt(tcol(j)*(n-tcol(j))*(n*n*n-buf4)/(3*n*(n-1)));
uridit(j)=(0.5-ridit(j))/sqrt((n*n*n-buf4)/(12*n*n*(n-1)*tcol(j)));
buf5=buf5+t(j)*t(j)/tcol(j); END;
Hc=(12/(n*(n+1))*buf5-3*(n+1))/(1-(buf4-n)/(n*(n*n-1)));
prob_Hc=1-PROBCHI(Hc,col-1);
PUT @27 "Kruskal-Wallis's Test"
#3 @19 'Adjusted Chi-Square(Hc)' @49 'Prob>Hc'
#4 @25 Hc 10.4 @50 prob_Hc 6.4;
PUT #7 @25 'Result of CPD and RIDIT Analysis'
#9 @6 'CPD Values' @20 "CPD's U-Values"
@40 'RIDIT Values' @60 "RIDIT's U-Values";
DO j=1 TO col;
PUT #(j+9) @6 cpd(j) 7. @22 ucpd(j) 8.4
@43 ridit(j) 7.4 @62 uridit(j) 8.4;
END;
buf1=probit(0.05/2/col); buf2=-buf1;
buf3=probit(0.05/2); buf4=-buf3;
PUT #(col+10) @6 "ALPHA= 0.05 CPD's Critical U-Values"
@56 buf1 7.4 @68 buf2 7.4;
PUT #(col+11) @6 " RIDIT's Critical U-Values"
@56 buf3 7.4 @68 buf4 7.4;
buf1=probit(0.01/2/col); buf2=-buf1;
buf3=probit(0.01/2); buf4=-buf3;
PUT #(col+12) @6 "ALPHA= 0.01 CPD's Critical U-values"
@56 buf1 7.4 @68 buf2 7.4;
PUT #(col+13) @6 " RIDIT's Critical U-Values"
@56 buf3 7.4 @68 buf4 7.4;
CARDS;
2 3 1 0
51 70 60 45
15 38 27 53
1 2 2 3
;
RUN; (程序的第4部分)
[KRCTEST.SAS修改指导] ①程序中数据输入方式是:横向为指标分组,纵向为实验分组。
②由于3种方法在本质上是一致的,按其中某种方法计算出结果,通过一些公式的联系,其他方法的结果也就随即得到。本程序以CPD分析的算法进行计算,利用王广仪教授书中的公式算出RIDIT分析和Kruskal-Wallis检验的结果。
③程序中再次定义了宏变量,第1个宏变量的值代表分组数,第2个代表等级数, 对其他资料引用该程序时,只需修改这两个地方即可。
[输出结果] Kruskal-Wallis's Test
Adjusted Chi-Square(Hc) Prob>Hc
21.6427 0.0001
Result of CPD and RIDIT Analysis
CPD Values CPD's U-Values RIDIT Values RIDIT's U-Values 等 级
3894 2.8141 0.4244 2.5405 1
1309 0.7993 0.4845 0.6673 2
1689 1.1077 0.4748 0.9648 3
-6892 -4.3521 0.5915 -3.7165 4
ALPHA= 0.05 CPD's Critical U-Values -2.4977 2.4977
RIDIT's Critical U-Values -1.9600 1.9600
ALPHA= 0.01 CPD's Critical U-Values -3.0233 3.0233
RIDIT's Critical U-Values -2.5758 2.5758
[输出结果的解释和专业结论] 3种方法得到的结果是一致的。Kruskal_Wallis 检验输
出1个卡方值=21.6427和1个P值<0.0001,表明术前视力是影响手术疗效的重要因素, 而RIDIT分析和CPD都对术前视力和疗效的关系作出具体的评价, CPD值从3894到-6892和RIDIT值从0.4244到0.5915分别代表用这两种检验法算得的检验统计量的值, 它们分别反映术前视力从1级到4级的疗效,在α=0.05的水平下,按疗效不同,把术前视力划分为3类, 术前视力为1级的疗效最好(因2.8141>2.4977或2.5405>1.9600), 2级和3级没有明显的区别(因0.7993与 1.1077都介与-2.4977~2.4977之间或0.6673与0.9643都介与-1.9600~1.9600之间),4~5级的疗效最差(因-4.3521<-2.4977或-3.7165<-1.9600)。 |
|