用SAS实现单变量统计分析

愿太阳照亮我

用SAS实现单变量统计分析
　　　　对一组单指标试验数据进行分析常采取2种方法,即图示法和计算几种简单的描述统计量
　　。图示法包括 绘制资料的干叶图(A STEM-AND-LEAF PLOT)、箱图(A BOX PLOT)和正态概率图(A NORMAL PROBABILITY PLOT)、简单的描述统计量包括矩(MOMENTS)、分位数(QUANTILES)、极端值(EXTREMES)以及资料的频数表,详见例子和输出结 果。
　　　　[例1.6.1]　68例健康成年人血清NSE酶(ng/ml)含量的测定结果参见[UNIV1.PRG], 试用SAS中UNIVARIATE过程进行单变量统计分析。
　　　　[SAS程序]　　　　　　　　　　　　　[UNIV1.PRG]
　　DATA c1;　　　　　　　　　　　　　　　 6.5 7.6 3.5 5.7 5.4 2.3 7.4
　　INPUT x @@;　　　　　　　　　　　　　 2.7 4.2 6.4 6.9 7.2 6.7 6.5
　　CARDS;　　　　　　　　　　　　　　　　4.0 7.3 1.1 4.9 2.5 2.9 1.9
　　6.5 5.1 6.1 3.9 3.5 7.7 2.1　　　　　　3.6 1.4 2.5 4.4 2.5
　　1.9 9.6 7.9 7.6 7.8 4.6 6.1　　　　　　;
　　6.4 2.8 7.6 2.5 4.6 8.1 4.8　　　　　　PROC UNIVARIATE plot normal
　　6.9 5.1 2.0 6.4 6.0 4.5 8.0　　　　　　　　　　　　　　def=4 freq;
　　8.0 8.6 6.4 4.9 6.4 6.8 4.7　　　　　　 VAR x;
　　3.4 7.7 1.2 2.8 0.5 2.6 3.2　　　　　　RUN;
　　(程序的第1部分)　　　　　　　　　　　　(程序的第2部分)
　 　　　[UNIV1.PRG修改指导]　过程步中各选择项的意义是 : PLOT绘制干叶图、箱图和正态概率图;NORMAL对资料进行正态性检验; DEF=4从5种计算分位数的定义式中挑选了第4种(DEF=4是该选择项的缺省值);FREQ输出资料的频数表,此表不实用, 因为分组过细。本篇第4章第5节中已给出了1个很实用的编制频数表的SAS程序, 后面将用1个简便的方法实现这一目的。
　　　　[输出结果及其解释]　　　　UNIVARIATE PROCEDURE
　　　　　　　　　　　　　　　　　Moments(矩或称动差)
　　样本含量(N)→　　　N　　　　　　　　68　Sum Wgts　　　　 68　　 ←观测值权众之和
　　算术均值→　　　　 Mean　　　 5.019118　Sum　　　　　 341.3　 ←样本观测值之总和
　　样本标准差→　　　 Std Dev　　2.238753　Variance　 5.012017　　　　　 ←样本方差
　　偏度系数→　　　　 Skewness　 -0.12124　Kurtosis　　-1.0954　　　 ←样本峰度系数
　　未校正的平和→　 USS　　　　 2048.83　CSS　　　　335.8051　　　 ←校正的平和
　　样本变异系数→　　 CV　　　　 44.60452　Std Mean　 0.271489 ←样本均值的标准误差
　　t统计量的值→　　　T:Mean=0　 18.48739　Prob>|T|　　 0.0001　←与t检验对应的概率
　　符号秩统计量的值→ Sgn Rank　　　 1173　Prob>|S|　　 0.0001 ←与符号秩检验对应的
　　非零数据个数　　　 Num ^= 0　　　　 68　　　　　　　　　　　　　　　　　　　概率
　　W统计量的值→　　　W:Normal　 0.948119　Prob　 　　　[说明] USS─即各观测值平之和、CSS─即各观测值的离均差平之和、t统计量的值t=18.48739, 其检验假设为H0:总体均值μ= 0 、中心化符号秩统计量的值R=1173,其检验假设为H0:μ=0、W统计量的值W=0.948119,其检验假设为H0:资料服从正态。
　　　　本例t=18.48739, P<0.0001, 说明资料并非抽自均值为零的总体; 又因W=0.948119, P=0.0147,说明该资料不服从正态。如果N>2000,SAS系统将采用D检验法进行正态性检验。
　　　　　　　　　　　　　　　　Quantiles(Def=4)(分位数)
　　第100百分位数(较大值)→　　100% Max　　9.6　　　 99%　　9.6　　　 ←第99百分位数
　　第75百分位数(第3四分位数)→ 75% Q3　　 6.875　　 95%　　8.055　　 ←第95百分位数
　　第50百分位数(中位数)→　　　50% Med　　5　　　　 90%　　7.81　　　←第90百分位数
　　第25百分位数(第1四分位数)→ 25% Q1　　 2.825　　 10%　　1.99　　　←第10百分位数
　　第0百分位数(最小值)→　　　　0% Min　　0.5　　　　5%　　1.29　　　←第5百分位数
　　　　　　　　　　　　　　　　　　　　　　　　　　　1%　　0.5　　　 ←第1百分位数
　　极差(即Max-Min)→　　　　　　 Range　　9.1
　　四分位数间距→　　　　　　　　Q3-Q1　　4.05
　　众数(出现次数最多的值)→　　　 Mode　　6.4
　　　　　　　　　　　　　　　　　　Extremes(极端值)
　　　　　　　　　　　　　Lowest　　Obs　　 Highest　　Obs
　　这是资料中　　　　　　0.5(　　　40)　　　　8(　　　28)　　　　　　 这是资料中
　　最小的5 个　　　　　　1.1(　　　59)　　　　8(　　　29)　　　　　　 较大的5 个
　　数值及其在　　　　　　1.2(　　　38)　　　8.1(　　　20)　　　　　　 数值及其在
　　资料中的顺　　　　　　1.4(　　　65)　　　8.6(　　　30)　　　　　　 资料中的顺
　　序号　　　　　　　　　1.9(　　　63)　　　9.6(　　　 9)　　　　　　 序号
　　　　　　　　　　　　　Stem Leaf(干叶图)　　　　　　　　　#　　 Boxplot(箱图)
　　有5个数≥8　　　　　　　 8 00166　　　　　　　　　　　　 5　　　　|
　　有25个数≥6且<8　　　　　6 0114444455578992346667789　　25　　 +-----+
　　有14个数≥4且<6　　　　　4 02456678991147　　　　　　　 14　　 *--+--*
　　有18个数≥2且<4　　　　　2 013555567889245569　　　　　 18　　 +-----+
　　有6个数≥0且<2　　　　　 0 512499　　　　　　　　　　　　6　　　　|
　　水平标尺→　　　　　　　　 ----+----+----+----+----+
　　　　[说明]　干叶图(Stem　Leaf　Plot)与箱图(Box　 plot)直观地反映了68个原始数据主
　　要集中在[2,8)之间;干叶图中最下面一行“0　512499　　6”的含意是:原始数据≥0且小于
　　2的数据个数有6个,它们的个位数分别为5,1,2,4,9,9;同理,可解释其他行的含意。
　　　　　　　　　　　　　　　　Normal Probability Plot(正态概率图)
　　　　　　　　　　9+　　　　　　　　　　　　　　　　　　　 +**+*+*++ *
　　　　　　　　　　 |　　　　　　　　　　　　　 **********+**
　　　　　　　　　　5+　　　　　　　　　　 ******+++
　　　　　　　　　　 |　　　　　　****+*****
　　　　　　　　　　1+ * ++*+*+**+
　　　　　　　　　　　+----+----+----+----+----+----+----+----+----+----+
　　　　　　　　　　　　　-2　　　　-1　　　　 0　　　　+1　　　　+2
　　　　　　　　　　　　　　　　 Frequency Table(频数表,略)
　　　　正态概率图显示各散点并非呈直线趋势而呈阶梯形着,说明资料偏离正态较远。
　　　　第2节　用SAS语言编程求小样本率的置信区间
　　　　[例1.6.2]　从某地随机抽取15人,经检查,其中有7人具有某种阳性体征, 试求该阳性体征在该地人群中的发生率的95％置信区间。
　　　　[分析与解答]　已知:n=15,X=7,α=0.05,即α/2=0.025,所需的SAS程序如下:
　　[SAS程序]　　　　　　　　[ERXFBCI.PRG]
　　DATA a;
　　x=7;　 n=15;　　　　　　　　　　 /* The second line */
　　p=x/n;
　　DO pp=0.001 TO p BY 0.0001;
　　　　p1=1-PROBBNML(pp,n,x-1);
　　　 IF ABS(p1-0.025)<0.00005 THEN GOTO ok1;　 /* line 6 */
　　END;
　　ok1: PUT x n p pp p1;
　　DO qq=p TO 0.999 BY 0.0001;
　　　　p2=PROBBNML(qq,n,x);
　　　 IF ABS(p2-0.025)<0.00005 THEN GOTO ok2;　 /* line 11 */
　　END;
　　ok2: PUT x n p qq p2;
　　PROC PRINT;
　　VAR x n p p1 p2 pp qq;
　　RUN;
　　　　[ERXFBCI.PRG修改指导]　SAS程序中字母大写的词表明是SAS语句固定的写法,小写者表明是用户定义的变量名等,下同。从运行的角度看,大小写没有任何差别。
　 　　　第2行中的数字7和15需根据用户的具体问题作相应修改;第6行与第11行是求95％置信区间的关键语句,若要求99％置信区间, 只需将这2句中的0.025改成0.005即可。要想使估计的结果更较精确, 可把第6行与第11行中的0.00005改成更小的正数。本程序是根据二项原理编制的。一般统计书上给出的“百分率的95%可信区间”表就是这样算得的。
　　　　[输出结果]　X　　 n　　　 p　　　　p1　　　　p2　　　　 pp　　　 qq
　　　　　　　　　　7　　15　 0.46667　0.024959　 0.025040　　0.2126　 0.73407
　　　　[输出结果的解释]　P=0.46667为观测所得的阳性率,[0.213, 0.734]为总体阳性率π的
　　95％置信区间,P1,P2分别为左侧和右侧的近似概率,即当P=0.734,n=15时, X≤7的概率为P1=
　　0.024959≈0.025;当P=0.213,n=15时,X≥7的概率为P2=0.025040≈0.025。欲使P1,P2的值更接近于0.025,需进一步提高误差的精度,即把SAS程序中的0.00005改成更小的值,比如0.000000005。
　　　　[专业结论]　该阳性体征在该地人群中的发生率的95％置信区间为[0.213,0.734]。
　　　　第3节　用SAS语言编程实现相对比误差传递公式的计算
　　　　在第1篇第3章第9节中介绍了２变量之商(相对比)的误差传递公式, 并给出了应用实例,当时,用计算器给出了计算结果。下面仍沿用该例, 并用SAS语言编程实现相对比误差传递公式的计算。
　　　　[例1.6.3]　求均数宏准差用[MNSTD.PRG];求存活率及其标准误差用[SDWCCD.PRG]。
　　　　[SAS程序]　　　　　　　　　　　　[MNSTD.PRG]
　　DATA chl1;　　　　　　　　　　　　　 PROC MEANS NOPRINT;
　　DO a=0 TO 6; DO n=1 TO 3;　　　　　　CLASS a;
　　　　INPUT x @@;　OUTPUT;　　　　　　　VAR x;
　　END; END;　CARDS;　　　　　　　　　　OUTPUT OUT=aaa MEAN=m STD=std;
　　CARDS;　　　　　　　　　　　　　　　 RUN;
　　611 608 622　　　　　　　　　　　　　DATA chl2;
　　358 329 341　　　　　　　　　　　　　 SET aaa;
　　253 239 284　　　　　　　　　　　　　 IF a='.' THEN DELETE;
　　111 124 137　　　　　　　　　　　　　PROC PRINT;
　　68　41　55　　　　　　　　　　　　　 ID a;
　　35　24　11　　　　　　　　　　　　　RUN;
　　　0　 0　 0
　　;
　　(程序的第1部分)　　　　　　　　　　　(程序的第2部分)
　　　　[MNSTD.PRG修改指导]　只需修改A(实验组数)和N(各组重复实验次数)的取值以及实验数据。A=0代表的是对照组。如果各组重复实验次数不等, 请参看本书第2篇第2章SAS程序─
　　[YBTTEST1.PRG],修改本程序的数据步。
　　　　[SAS程序]　　　　　　　　　　　　　[SDWCCD.PRG]
　　DATA abc;　　　　　　　　　　　　　　　342.667　　14.5717
　　m0=613.667;　　　　　　　　　　　　　 258.667　　23.0290
　　std0=7.3711;　　　　　　　　　　　　　124.000　　13.0000
　　INPUT M STD;　　　　　　　　　　　　　 54.667　　13.5031
　　P=m/m0;　　　　　　　　　　　　　　　　23.333　　12.0139
　　sp=SQRT((std/m0)**2+　　　　　　　　　　0.000　　 0.0000
　　　　(m*std0/m0/m0)**2);　　　　　　　 ;
　　p=p*100;　　　　　　　　　　　　　　 PROC PRINT;
　　sp=sp*100;　　　　　　　　　　　　　　VAR p sp;
　　CARDS;　　　　　　　　　　　　　　　　RUN;
　　(程序的第1部分)　　　　　　　　　　　(程序的第2部分)
　　　　[SDWCCD.PRG修改指导]　将[MNSTD.PRG]的输出结果(各组的均值宏准差)作为输入数据,M0、STD0分别为对照组的均值与标准差。
　　[MNSTD.PRG的输出结果]　　　　　　　　 [SDWCCD.PRG的输出结果]
　　　A　 _FREQ_　　 M　　　　　STD　　　　OBS　　　　P　　　　　 SP
　　组别　频数　　 均值　　　标准差　　　 组别　 存活率(％)　P的标准误差(％)
　　　0　　　3　　613.667　　　7.371　　　　
　　　1　　　3　　342.667　　 14.572　　　　 1　　　55.8392　　 2.46744
　　　2　　　3　　258.667　　 23.029　　　　 2　　　42.1510　　 3.78669
　　　3　　　3　　124.000　　 13.000　　　　 3　　　20.2064　　 2.13227
　　　4　　　3　　 54.667　　 13.503　　　　 4　　　 8.9083　　 2.20300
　　　5　　　3　　 23.333　　 12.014　　　　 5　　　 3.8022　　 1.95826
　　　6　　　3　　　0.000　　　0.000　　　　 6　　　 0.0000　　 0.00000
　　过程步与右边相同。