SAS-给公众号做一个秩和检验

2019-10-21 17:34:52 浏览数 (1)

日前,小编的公众号的关注人数终于突破了2000了!但还远远没有完成今年的年度计划,于是小编就想看一看,每日增粉数量是否与关注基数有显著的统计学意义~好吧,作为一个不会统计学的菜鸟,写这篇推送的目的是希望各位老师能够指点一下小编,并检阅一下小编给自己布置的作业是否做对了,还请各位老师不吝赐教!

嗯,于是小编从公众号上下载了自2017年11月11日-2018年03月25日的公众号每日增粉相关的数据...接着小编就开始分组了,以500人为区间,分成3个组进行对照研究(group1:<=1000;group2:1000-1500;group3:1501-2000);小编这里想到了临床试验中比较常见的对连续变量进行的描述性统计分析的一个例子,因此,我就套用过来对我公众号每日增粉数量进行分析,并按照临床试验中出三线表的形式,将分析结果进行输出!结果如下:

嗯,看这里的P值是小于0.0001的,得出的结论大概是有显著的统计学意义的,为啥小编这里用秩和检验呢,听说用秩和检验不用考虑数据的正态性,所以就它了。从结果可以看出,均值是逐渐增大的,那么是不是可以说明随着关注量的增大,每日增关量也会增加呢?

接着呢,小编想看一看增粉数量与星期有没有显著性差异,于是小编用了卡方检验,同样得出了P值小于0.0001,是的,有显著性差异!那么还能看出啥呢?小编就不太知道了!

在就是结果中三线表的输出,三线表的输出小编以前虽然写过相关的推送,不过吧,现在水平又提升了一丁点!也正好想找点数据来练手,于是小编就将在本文中在一次的写一写report过程,以及将以前写过的一些基础的专题与推送,简单的连贯起来进行一次实践性质的编程!

首先,我们在写程序的时候,做的一件事是建立文件夹以及逻辑库路径!

嗯,这是我建立的文件夹,数据、文档、日志、宏、输出、程序、源数据都有自己的夹子。如果有兴趣的朋友,可以在公众号对话框回复:文件夹1,既可以下载本文涉及的全部程序以及数据以及macro程序包!解压后程序可以直接运行。

回到正文上来,下载数据 -数据导入:

代码语言:javascript复制
%macro setpath;
%local setup_ runsetup runsetup1 runsetup2 ;
%let setup_= %upcase(%sysget(sas_execfilepath));
%let runsetup=%sysfunc(prxchange(s/(.*)\.*/1/,-1,&setup_));
%let runsetup1=%sysfunc(prxchange(s/(.*)\.*/1/,-1,&runsetup));
%let runsetup2=%sysfunc(prxchange(s/(.*)\.*/1/,-1,&runsetup1));
**run Macro;
%if &SYSVER.=9.4 %then %do;
libname temp "&runsetup2.macros";
%end;
%if &SYSVER.=9.2 %then %do;
libname temp "&runsetup2.macros9.2";
%end;

**run Macro;
Options MAutoSource MStored SASMStore=temp;
%put NOTE:&runsetup2.;
proc datasets library=work kill nolist;
quit;
%mend;

这段代码用到了一个系统宏变量(sysver)来获取你SAS的版本号,小编电脑上装了SAS9.2与SAS9.4,因此将宏分别执行成2个版本宏包(执行后的宏包不能夸版本),所以小编这里利用这个宏变量自动获取当前SAS软件的版本号,然后进行判断进而选择正确的宏包。嗯,下载文件夹后,解压压缩包,不改变压缩包内部文件以及文件夹,在目录下的PGMTableTable.sas程序打开是可以直接运行的。

数据导入后呢,数据集的转换以及建立分组等:

代码语言:javascript复制
proc format ;
value wk 1='星期日' 2='周一' 3='周二' 4='周三' 5='周四' 6='周五' 7='周六';
value gp 1='<=1000' 2='<=1500' 3='<=2000' ;
run;

*建立分组>;
data raw;
set raw.rawdata;
format group gp.;
if '累积关注人数'n <=1000 then group=1;
else if '累积关注人数'n <=1500 then group=2; 
else  group=3; 
weekday=weekday('时间'n);
week=put(weekday,wk.);
run;
ods html close;
ods listing;

接着就到了描述性统计量的计算,这里用mean过程步:

代码语言:javascript复制
proc means data=raw n nmiss mean std median min max;
class group;
var '新关注人数'n;
output out=temp2 std=stdn mean = meann  max=maxn min=minn median=mediann n=nn;
run;
data temp3;
set temp2;
if group=1 then _NAME_="Temp1";
else if group=2 then _NAME_="Temp2";
else if group=3 then _NAME_="Temp3";
else delete;
std=strip(round(stdn,.001));
mean=strip(round(meann,.001));
median=strip(round(mediann,.001));
max=strip(maxn);
min=strip(minn);
n=strip(nn);
run;

计算后的结果:这里为啥小编需要给数据转换成字符型的变量呢,主要原因是这样的,后面小编还要进行数据集的追加,变成字符变量比较好操作!

看到上面的结果:是不是和前面RTF中的排版结构相差很大呢,那是因为没有转置,接下来就来转置一下:

代码语言:javascript复制
proc transpose data=temp3 out=temp3_1 (rename=(_NAME_=Temp0));
var  n  mean std median min max;
ID _NAME_;
IDLABEL GROUP;
run;

结果是这样的:

现在就和前面RTF输出的结构差不多了,好了,然后在进行P值的计算:

代码语言:javascript复制
ods output KruskalWallisTest=temp2 ;
proc npar1way data=raw wilcoxon   ;
class group;
var  '新关注人数'n;
run;
data _null_;
  set temp2;
  if LABEL1='Pr > 卡方' then call symput('frq_pv1',strip(cValue1));
run;
%put NOTE: p值:&frq_pv1.; 
data temp1;
length temp0 Temp4 $400.;
temp0="新关注人数^R/RTF'bsuperfs22'[3]" ;Temp4="&frq_pv1.";
run;

data task1;
set temp1 temp3_1(in=a);
order=1;
if a then temp0=upcase(temp0);
i=_n_;
if a then temp0=strip("^R/RTF'u32?u32?u32?u32?'")||strip(temp0);
run;

运行完后的结果是这样的: 下面首列的一串字符是啥意思呢,这是RTF标记语言,实现一些排版缩进等操作!

这里一部分算是完成了,接下来在来看第二部分。卡方检验结果的输出!其实和上面的也是很类似,都整到数据集中,在进行转置啊等几步数据结构的操作,这样就可以实现想要输出结构的排版,这里就不细说,直接贴代码了!

代码语言:javascript复制
ods listing;
Ods Output CrossTabFreqs=temp1  Chisq=Chisq    ;
proc freq data=raw ; 
table  group *week/ Chisq ;
weight '新关注人数'n;
run;


data _null_;
set chisq;
if STATISTIC='卡方' then call symput('frq_pv1',strip(vvalue(PROB)));
run;
%put NOTE: p值:&frq_pv1.; 
data temp2;
set temp1;
value=strip(FREQUENCY)||strip('(')||strip(round(ROWPERCENT,.1))||strip('%)');
if missing(week) then delete;
if missing(group) then delete;
if group=1 then _NAME_="Temp1";
else if group=2 then _NAME_="Temp2";
else if group=3 then _NAME_="Temp3";
else  delete;;
rename group=_LABEL_;
run;

proc sort data=temp2  out=temp2 ;by week  _NAME_  ;;quit;
proc transpose data=temp2 out=temp3(drop=_NAME_) ;
by week;
var value;
run;

data temp1;
length  Temp4 $400.;
week="星期^R/RTF'bsuperfs22'[4]" ;Temp4="&frq_pv1.";
run;

data Task2;
length  week temp1-temp4 $400.;
format  week temp1-temp4 $400.;
informat week temp1-temp4 $400.;
set temp1 temp3(in=a);
order=2;
i=_n_;
rename week=temp0;
if a then week=strip("^R/RTF'u32?u32?u32?u32?'")||strip(week);
run;
/*数据集合并*/
data final;
length   temp0-temp4 $500.;
set Task1-Task2;
page=1;
run;

proc datasets library=work  nolist memtype=data;
save final Task1-Task2;
quit;

最后在来看看要输出的final数据集:

做到这一步,整个数据的操作基本上是完成了,接下来就是数据集的输出,输出到RTF中。这里小编就用了直接写好的宏进行输出。

首先,我先导入输出RTF中的标题、备注等信息,前面可以看到小编输出的rtf里面是带有备注等信息的,其实就是通过这里控制。顺便看一看titles的Excel模板是啥样的。为啥要做Excel中呢,放在外部,其实也便于修改与管理,以及用Macro来实现自动添加这些信息等等!第二个sheet有一些常用的RTF标记的例子,方便忘记的时候随时查找!

代码语言:javascript复制
/*import title*/
%macro excel2sas(path,sheet,outds);
proc import out=&outds.
datafile= "&path." 
dbms=excel replace;
range="&sheet.$"; 
getnames=yes;
mixed=no;
scantext=yes;
usedate=yes;
scantime=yes;
run;
%mend;
%excel2sas(path=&Macpath.titles.xlsx,sheet=sheet1,outds=titles);

紧接着,就是一些Macro的调用了,这里就不过多的说了,有兴趣的朋友可以私下找我讨论。

代码语言:javascript复制
/*import template*/
%rtf_ods_temp(font=%str(Courier New),size=8pt);
%rtf_ods_title(pgmname=Table.sas,tablename=微信公众号:xiaocgn,inds=titles,tableid=Table_1,ftyn=N)
*pgmname :title显示程序名称  tablename :title显示表格名称  inds:title的数据  tableid:对饮表格编号 ftyn:title 是否在body中 ;
/*Set output file name and path*/
%rtf_ods_st(path=&Table.,report=输出的结果,style=style_tb);

ods rtf exclude none;

然后就到了proc report过程步了...

代码语言:javascript复制
ods rtf exclude none;
proc report data=final nowd headskip headline split='|' missing nocenter
;
  column page order temp0 -temp4   ;
  define page / order order=internal noprint;
  define order /  order order=internal noprint;
  define temp0 / display " " style(header)=[just=left] style(column)=[just=left cellwidth=38% asis=on] flow;
  define temp1 / display "Group1|(<=1000)"   style(header)=[just=center] style(column)=[cellwidth=15% just=center asis=on] flow;
  define temp2 / display "Group2|(1000-1500)"   style(header)=[just=center] style(column)=[cellwidth=15% just=center asis=on] flow;
  define temp3 / display "Group3|(1501-2000)"   style(header)=[just=center] style(column)=[cellwidth=15% just=center asis=on] flow;
  define temp4 / display "P值"   style(header)=[just=center] style(column)=[cellwidth=15% just=center asis=on] flow;
  break after page/page;

 compute before order;
    line " ";
  endcomp;;

  %rtf_ods_on;
run;
%rtf_ods_end;

到这里,基本上整个程序是写完了!如对程序有兴趣,可以在可以在公众号对话框回复:文件夹1。

小编是初学者,各位老师如果觉得文中涉及统计存在问题,欢迎指正!以免其他初学者被小编带偏了!

今天就这么多了,后续内容,敬请期待~

0 人点赞