SAS中文论坛

标题: SAS怎么做RIDIT分析啊？ [打印本页]

作者: shiyiming 时间: 2005-7-4 17:13
标题: SAS怎么做RIDIT分析啊？
这个可以做吗？怎么实现啊？
谢谢！

作者: shiyiming 时间: 2005-7-5 10:01
标题: Re: SAS怎么做RIDIT分析啊？
[b:88781]one example[/b:88781]

[quote:88781]%let group=3;  /* 三种药物 */
%let rank=4; /* 每种药物的疗效有四个等级 */
%let col_row=%eval(&group*&rank);
DATA _NULL_;
row=&rank;
col=&group;
ARRAY num(&rank,&group);
ARRAY cpd(&group);
ARRAY ridit(&group);
ARRAY trow(&rank);
ARRAY tcol(&group);
ARRAY y(&rank);
ARRAY ucpd(&group);
ARRAY uridit(&group);
ARRAY t(&group);
FILE PRINT;
n=0;    buf1=0;
buf4=0;  buf5=0;
INPUT num1-num&col_row;
DO j=1 TO col;
  tcol(j)=0;
  DO i=1 TO row;
tcol(j)=tcol(j)+num(i,j);
  END;
n=n+tcol(j);
END;
buf3=n;
DO i=1 TO row;
trow(i)=0;
DO j=1 TO col;
  trow(i)=trow(i)+num(i,j);
END;
buf4=buf4+trow(i)*trow(i)*trow(i);
buf2=trow(i);
y(i)=buf3-buf1-buf2;
buf3=y(i);
buf1=trow(i);
END;
DO j=1 TO col;
cpd(j)=0;
DO i=1 TO row;
  cpd(j)=cpd(j)+num(i,j)*y(i);
END;
ridit(j)=0.5-cpd(j)/(2*tcol(j)*n);
t(j)=tcol(j)*(n+1)/2-cpd(j)/2;
  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 @5 "秩和检验结果"
   #3 @5 '秩和检验统计量值' @30 '秩和检验P值'
   #4 @5 Hc 10.4 @30 prob_Hc 6.4;
PUT #7 @5 'Ridit分析结果'
   #9 @5 'Ridit值' @15 "检验Ridit所对应的U统计量的值";
DO j=1 TO col;
  PUT #(j+9) @5 ridit(j) 7.4  @15 uridit(j) 8.4;
END;
buf1=probit(0.05/2/col);  buf2=-buf1;
buf3=probit(0.05/2);    buf4=-buf3;
PUT #(col+12) @5 "U 统计量值的95％置信区间"
            @30 buf3 7.4 @40 buf4 7.4;
buf1=probit(0.01/2/col);  buf2=-buf1;
buf3=probit(0.01/2);    buf4=-buf3;
PUT #(col+14) @5 "U 统计量值的99％置信区间"
            @30 buf3 7.4 @40 buf4 7.4;
CARDS;
15  4  1
49  9 15
31 50 45
5 22 24
;
RUN;[/quote:88781]

作者: shiyiming 时间: 2005-7-5 18:42
标题: Re: SAS怎么做RIDIT分析啊？
不过要是只要2种药物的话不知道程序哪里需要改一下啊？
呵呵
我刚接触这个
不懂啊
再次谢谢！

作者: shiyiming 时间: 2005-7-6 07:20
标题: %let group=3; /* 三种药物 */
%let group=3; /* 三种药物 */

欢迎光临 SAS中文论坛 (http://mysas.net/forum/)