macro
IncBeta x p q Bx

#This macro calculates the incomplete beta function
#B(x;p,q) using a series solution from Arfken.

#Variable Definitions:
#x is the upper limit on the beta function integral
#p and q are beta function parameters
#Bx is the value of the incomplete beta function
#lastBx and err are used to track the progress of the summation
#term is the running product term in Arfken's expansion
#n is the loop index

#Example calling statement:
#   mtb> %incbeta 0.1 3 20 k1
#Answer is Bx = 0.379959

#WARNING! The series expansion is asymptotic. You may have to
#change the fractional error that controls the loop to obtain
#the required accuracy.

#References:
#see Arfken, Ed. 3, p292, Exercise 5.2.18
#see AMS55 p. 944, 26.5.1 (Arfken is missing the B(p,q) term)
#B(p,q)=gamma(p)gamma(q)/gamma(p+q) where gamma(p)=(p-1)!
#The F distribution is directly calculable from B(x;p,q). See
#AMS55 p. 945 26.5.28.

#By Paul G. Mathews
#217 Third Street, Fairport Harbor, OH 44077
#pmathews@apk.net
#440-350-0911

#Rev. 1.0 for Minitab 13.3 15 March 2001
#Copyright (c) March 2001 Paul G. Mathews All rights reserved.

mconstant x p q Bx n lastBx term err

note
note This macro is very slow. Please be patient.
note Running macro ...

let Bx=1/p		#the 0th term in Arfken's expansion
let err=1
let term=1
let n=1
while err>0.00001	#set the fractional error here
   let lastBx=Bx
   let term=term*(n-q)/n*x
   let Bx=Bx+term/(p+n)
   let err=abs(1-lastBx/Bx)
   let n=n+1
endwhile
let Bx=gamma(p+q)/gamma(q)/gamma(p)*Bx*x**p
print Bx

endmacro