24小时热门版块排行榜    

Znn3bq.jpeg
查看: 1236  |  回复: 4
本帖产生 1 个 程序强帖 ,点击这里进行查看
当前只显示满足指定条件的回帖,点击这里查看本话题的所有回帖

ian_zhangty

木虫 (著名写手)

【答案】应助回帖

★ ★ ★ ★ ★ ★
jjdg(金币+1): 感谢分享你的解决方案 2011-08-06 15:28:45
dubo(金币+5, 程序强帖+1): good 2011-08-21 14:01:09
已经解决,方案如下,谢谢

FUNCTION betai(a,b,x)
REAL betai,a,b,x
!USES betacf,gammln
REAL bt,betacf,gammln
if(x<0..or.x>1.) pause 'bad argument x in betai'
if(x==0..or.x==1.) then
  bt=0.
else
  bt=exp(gammln(a+b)-gammln(a)-gammln(b)+a*log(x)&
     +b*log(1.-x))
endif
if(x<(a+1.)/(a+b+2.)) then
  betai=bt*betacf(a,b,x)/a
  return
else
  betai=1.-bt*betacf(b,a,1.-x)/b
  return
endif
END FUNCTION betai

FUNCTION betacf(a,b,x)
INTEGER maxit
REAL betacf,a,b,x,EPS,fpmin
PARAMETER (maxit=100,EPS=3.e-7,fpmin=1.e-30)
INTEGER m,m2
REAL aa,c,d,del,h,qab,qam,qap
qab=a+b
qap=a+1.
qam=a-1.
c=1.
d=1.-qab*x/qap
if(abs(d) d=1./d
h=d
do m=1,maxit
  m2=2*m
  aa=m*(b-m)*x/((qam+m2)*(a+m2))
  d=1.+aa*d
  if(abs(d)   c=1.+aa/c
  if(abs(c)   d=1./d
  h=h*d*c
  aa=-(a+m)*(qab+m)*x/((a+m2)*(qap+m2))
  d=1.+aa*d
  if(abs(d)   c=1.+aa/c
  if(abs(c)   d=1./d
  del=d*c
  h=h*del
  if(abs(del-1.)     betacf=h
    return
  end if
end do
pause 'a or b too big, or maxit too small in betacf'
END FUNCTION betacf

FUNCTION gammln(xx)
REAL gammln,xx
INTEGER j
DOUBLE PRECISION ser,stp,tmp,x,y,cof(6)
SAVE cof,stp
DATA cof,stp/76.18009172947146d0,-86.50532032941677d0,&
     24.01409824083091d0,-1.231739572450155d0,&
         .1208650973866179d-2,-.5395239384953d-5,&
         2.5066282746310005d0/
x=xx
y=x
tmp=x+5.5d0
tmp=(x+0.5d0)*log(tmp)-tmp
ser=1.000000000190015d0
do j=1,6
  y=y+1.d0
  ser=ser+cof(j)/y
end do
gammln=tmp+log(stp*ser/x)
END FUNCTION gammln



program pvalue
        REAL df, r, t, prob
        df=17-2   ! 17 is sample size

        r=0.188446645   ! pearson coefficient
        t=(abs(r)*sqrt(df))/sqrt(1-r**2)
        prob=betai(0.5*df,0.5,df/(df+t**2))
        write (*,*) prob





end program pvalue
5楼2011-08-06 10:41:54
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

智能机器人

Robot (super robot)

我们都爱小木虫

相关版块跳转 我要订阅楼主 瘦竹云淡风轻 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[基金申请] 国自然面上和省基金B类撒花 +9 花田半亩~白 2026-04-21 9/450 2026-04-22 16:50 by zhouyanli11
[教师之家] 又一批高校组建人工智能学院 师资行吗 不是骗人吗 +5 yexuqing 2026-04-19 5/250 2026-04-22 10:01 by easeheart
[论文投稿] 急需审稿人!!! +3 陆小果画大饼 2026-04-21 3/150 2026-04-21 23:54 by jzy_123456
[考博] 申博/考博 +4 啃面包的小书虫 2026-04-17 8/400 2026-04-21 16:26 by 啃面包的小书虫
[考研] 295分求调剂 +6 ?要上岸? 2026-04-17 6/300 2026-04-21 08:18 by Equinoxhua
[考研] 一志愿中科大材料与化工,353分还有调剂学校吗 +11 否极泰来2026 2026-04-15 13/650 2026-04-20 22:31 by Equinoxhua
[考研] 085600材料与化工调剂 5+3 孜孜不倦2002 2026-04-19 6/300 2026-04-20 21:25 by babero
[论文投稿] 有没有接收比较快的sci期刊呀,最好在一个月之内的,研三孩子求毕业 20+4 之护着 2026-04-16 7/350 2026-04-20 15:45 by 豆豆7758
[考研] 337求调剂 +3 jyz04 2026-04-18 3/150 2026-04-20 12:24 by 研可安
[考博] 申博 +3 Xyyx. 2026-04-18 3/150 2026-04-20 10:44 by YuY66
[考研] 求计算机方向调剂 +3 Toffee2 2026-04-16 6/300 2026-04-19 22:37 by ll叶
[考研] 294求调剂 +8 淡然654321 2026-04-17 9/450 2026-04-19 19:51 by Equinoxhua
[考研] 304求调剂 +8 castLight 2026-04-16 8/400 2026-04-19 17:14 by 中豫男
[考研] 300求调剂 +12 橙a777 2026-04-15 12/600 2026-04-18 23:51 by 路病情
[考研] 接受任何调剂 +6 也就是栗子 2026-04-17 7/350 2026-04-18 17:20 by 涵竹刘
[考研] 收到复试调剂但是去不了 +8 小蜗牛* 2026-04-16 8/400 2026-04-18 11:15 by zixin2025
[考研] 260求调剂 +4 Zyt1314520.. 2026-04-17 5/250 2026-04-18 08:28 by babysonlkd
[考研] 急需调剂 +9 绝不放弃22 2026-04-15 10/500 2026-04-18 08:09 by chixmc
[有机交流] 二苯甲酮酸类衍生物 50+3 小白爱主人 2026-04-17 6/300 2026-04-17 18:47 by kf2781974
[考研] 322求调剂 +6 tekuzu 2026-04-17 6/300 2026-04-17 13:48 by Espannnnnol
信息提示
请填处理意见