SAS中文论坛

 找回密码
 立即注册

扫一扫,访问微社区

查看: 857|回复: 3
打印 上一主题 下一主题

关于一个循环的SAS程序

[复制链接]

49

主题

76

帖子

1462

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
1462
楼主
 楼主| 发表于 2004-11-18 09:24:59 | 只看该作者

关于一个循环的SAS程序

我这里有一个SAS程序,有个问题想请教各位高手。

data  eq ;
do n=1 to    5          ;

   do aa=0  to  n       ;
     an=n-aa;
     do ab=0  to an     ;
      cn=an-ab ;
       do ac=0 to cn     ;
         ad=n-aa-ab-ac;
am1=aa+ab;an1=aa+ac;
total_p=0;
s_bc=ab+ac;
binomial=probbnml(0.5,s_bc,min(ab,ac));
doub_bin=2*binomial;
if binomial=0.5 then doub_bin=1;

  pi=(am1+an1)/(2*n) ;

  apa=am1/n;
  apb=an1/n;
  diff_bc=abs(ab-ac);
   p_diff=(am1-an1)/n;
q=1-pi;



  do m1=0 to n ;      m2=n-m1;
   do a=0 to m1 ;      b=m1-a;
    do c=0 to m2 ;      d=m2-c;
                         n1=a+c;
                           n2=n-n1;
        bc=b+c;
        b_c=b-c;
        b_c_abs=abs(b_c);
        mn1=m1+n1;
        mn2=m2+n2;
        rate_a=m1/n;
        rate_b=n1/n;
        rate_a_b=abs(rate_a-rate_b);
  if rate_a_b>=p_diff then do ;

   if abs(b_c)>=diff_bc  then do;
   if m1>0 then
     pm1=probbnml(pi,n,m1)-probbnml(pi,n,(m1-1));
        if m1=0 then pm1=probbnml(pi,n,m1) ;

   if a>0 then
     pa =probbnml(pi,m1,a)-probbnml(pi,m1,(a -1));
        if a=0 and m1>0 then pa=probbnml(pi,m1,a) ;
   if m1=0 then pa=1;

   if c>0 then
     pc  =probbnml(pi,m2,c)-probbnml(pi,m2,(c-1));
        if c=0 and m2>0 then pc=probbnml(pi,m2,c)  ;
        if m2=0 then pc=1;

     prob=pm1*pa*pc   ;

     total_p=total_p+prob;


      if (n*(m1+n1)*(m2+n2)) ne 0 then
      x2_p=2*(m1*n2-m2*n1)**2/(n*(m1+n1)*(m2+n2));
                        P_x2_p=1-probchi(x2_p,1);
      if (b+c) ne 0 then
      x2_m=(b-c)**2/(b+c);
                       P_x2_m=1-probchi(x2_m,1);

       two_tail=2*probbnml(0.5,bc,b);
            if two_tail>1 then two_tail=1;
      if (b+c) ne 0 then
      x2_mc=(abs(b-c)-1)**2/(b+c);
                        P_x2_mc=1-probchi(x2_mc,1);

if ((a+d)*(b+c)+4*b*c) ne 0 then
      x2_tj=((b-c)**2*(n-1))/((a+d)*(b+c)+4*b*c);
                        P_x2_tj=1-probchi(x2_tj,1);


         r1=(p_x2_m-two_tail)/two_tail;
         r2=(p_x2_mc-two_tail)/two_tail;
         r3=(p_x2_tj-two_tail)/two_tail;output;

      end;
      end;
     end;
    end;
  end;end;end;end;end;
keep a b c d pm1 pa pc prob bc b_c b_c_abs x2_p x2_m x2_mc
      p_x2_p p_x2_m two_tail p_x2_mc x2_tj p_x2_tj r1-r3  total_p  n aa ab ac ad doub_bin two_tail;
run;



proc print;
var n aa ab ac ad total_p doub_bin p_x2_p p_x2_m p_x2_mc p_x2_tj two_tail;
run;
问题:为什么循环以后在结果中aa, ab, ac , ad出现几个相同的结果?如aa, ab, ac, ad分别为0, 0, 0 , 1的情况出现了四次,请问出错在哪里?
回复 支持 反对

使用道具 举报

49

主题

76

帖子

1462

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
1462
沙发
 楼主| 发表于 2004-11-18 14:14:14 | 只看该作者

呵呵简化一下吧

天啦
这么大的一堆程序
谁也看不进去
你在模拟什么贝努利实验吗?
如果你能分成小块来问.
可能我们能帮上忙.
<!-- s:oops: --><img src="{SMILIES_PATH}/icon_redface.gif" alt=":oops:" title="Embarassed" /><!-- s:oops: -->
回复 支持 反对

使用道具 举报

49

主题

76

帖子

1462

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
1462
板凳
 楼主| 发表于 2004-11-22 20:08:06 | 只看该作者

解释一下

那,前20行是程序的第一部分,是一个循环,也是一个贝努里实验,求二项分布的概率。
后面是程序的第二部分,也是一个循环,prob=pm1*pa*pb, total_p是prob的累积和。
我的目的是让循环中不同情况下的total_p与doub_bin显示在同一个结果窗口中,即在同一行上,
回复 支持 反对

使用道具 举报

49

主题

76

帖子

1462

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
1462
地板
 楼主| 发表于 2004-12-9 14:35:44 | 只看该作者

可能缺少OUTPUT语句吧,多用几个看看!

可能缺少OUTPUT语句吧,多用几个看看!
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|小黑屋|手机版|Archiver|SAS中文论坛  

GMT+8, 2025-6-8 21:13 , Processed in 0.146151 second(s), 20 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.

快速回复 返回顶部 返回列表