|
!'標準正規分布
PUBLIC NUMERIC EPS
LET U=0 !'平均
LET M=1 !'分散
LET EPS=1E-13
FOR X=1 TO 400
LET Y=X/100
PRINT Y;":";INTEGRALNORMAL(U,M,Y);INTEGRAL(U,M,0,Y,100)
NEXT X
END
EXTERNAL FUNCTION INTEGRALNORMAL(U,M,X) !'INTEGRAL(EXP(-(X-U)*(X-U)/2/M/M)) 密度関数を項別積分
LET A=X
LET S=X
FOR I=1 TO 1000
LET A=-A*(X-U)*(X-U)/I/2/M/M
LET S=S+A/(2*I+1)
IF ABS(A)<EPS THEN EXIT FOR
NEXT I
LET INTEGRALNORMAL=S/SQR(2*PI)/M
END FUNCTION
EXTERNAL FUNCTION NORMAL(U,M,X) !'正規分布密度関数
LET NORMAL=EXP(-(X-U)*(X-U)/(2*M*M))/SQR(PI*2)/M
END FUNCTION
EXTERNAL FUNCTION INTEGRAL(U,M,A,B,N) !'数値積分
LET N=N*3
LET H=(B-A)/N
FOR K=0 TO N/3-1
LET S=S+3/8*H*NORMAL(U,M,A+H*3*K)+9/8*H*NORMAL(U,M,A+H*(3*K+1))+9/8*H*NORMAL(U,M,A+H*(3*K+2))+3/8*H*NORMAL(U,M,A+H*(3*K+3))
NEXT K
LET INTEGRAL=S
END FUNCTION
--------------------------------------------------------------------------------------------------
ニュートン法を用いてnorminv(z)を求める
PUBLIC NUMERIC EPS
LET U=0 !'平均
LET M=1 !'分散
INPUT PROMPT "norminv(z) z=":Z !' (0~0.5)
LET EPS=1E-13
LET X=0
DO
LET XX=X
LET X=XX-(INTEGRALNORMAL(U,M,XX)-Z)/NORMAL(U,M,XX) !'ニュートン法
LOOP UNTIL ABS(X-XX)<EPS
PRINT X
END
EXTERNAL FUNCTION NORMAL(U,M,X) !'正規分布密度関数
LET S=1
LET A=1
FOR I=1 TO 1000
LET A=-A*(X-U)*(X-U)/I/2/M/M
LET S=S+A
IF ABS(A)<EPS THEN EXIT FOR
NEXT I
LET NORMAL=S/SQR(PI*2)/M
END FUNCTION
EXTERNAL FUNCTION INTEGRALNORMAL(U,M,X) !'INTEGRAL(EXP(-(X-U)*(X-U)/2/M/M)) 密度関数を項別積分
LET A=X
LET S=X
FOR I=1 TO 1000
LET A=-A*(X-U)*(X-U)/I/2/M/M
LET S=S+A/(2*I+1)
IF ABS(A)<EPS THEN EXIT FOR
NEXT I
LET INTEGRALNORMAL=S/SQR(2*PI)/M
END FUNCTION
|
|