FC2ブログ

# 対数正規分布 (逆関数) Logarithmic Normal Distribution (Inverse Function)
# 引数 累積確率 平均 標準偏差 ($X, $Mean, $StandardDeviation)
# 戻り値 対数正規分布 (逆関数) (@InverseFunction)
sub LOGARITHMICNORMALDISTRIBUTIONINVERSE{
my ($X, $Mean, $StandardDeviation) = @_;
my @InverseFunction = ();

# 累積確率 平均 標準偏差の確認
if(($X < 0) || (1 < $X) || ($StandardDeviation <= 0)){
return "Error";
}

if(($X == 0) || ($X == 1)){
return "Inf or 0";
}

# 対数正規分布 (逆関数) Logarithmic Normal Distribution (Inverse Function)
# 下側 Lower
$InverseFunction[0] = &BISECTIONMETHOD($X, $Mean, $StandardDeviation);
# 上側 Upper
$InverseFunction[1] = &BISECTIONMETHOD((1 - $X), $Mean, $StandardDeviation);

return @InverseFunction;
}

# 二分法 Bisection Method
# 引数 累積確率 平均 標準偏差 ($X, $Mean, $StandardDeviation)
# 戻り値 二分法 ($BisectionMethod)
sub BISECTIONMETHOD{
my ($X, $Mean, $StandardDeviation) = @_;
my $BisectionMethod = 0;
my $X1 = 0;
my $X2 = 0;
my $F_m = 0;
my $F_x1 = 0;
my $F_x2 = 0;
my $Middle = 0;
my $PrevMiddle = 0;
my $Limit = 100;
my $Epsilon = 1.0e-20;

# 区間
$X1-- while((&LOGARITHMICNORMALDISTRIBUTION($X1, $Mean, $StandardDeviation) - $X) > 0);
$X2++ while((&LOGARITHMICNORMALDISTRIBUTION($X2, $Mean, $StandardDeviation) - $X) < 0);

# 初期値
$F_x1 = &LOGARITHMICNORMALDISTRIBUTION($X1, $Mean, $StandardDeviation) - $X;
$F_x2 = &LOGARITHMICNORMALDISTRIBUTION($X2, $Mean, $StandardDeviation) - $X;

# 計算
for(my $i = 0; $i < $Limit; $i++){
# 一つ前
$PrevMiddle = $Middle;
# 中間点
$Middle = ($X1 + $X2) / 2;

# f(Middle)
$F_m = &LOGARITHMICNORMALDISTRIBUTION($Middle, $Mean, $StandardDeviation) - $X;

# 置き換え
if(($F_m * $F_x1) > 0){
$X1 = $Middle
}
elsif(($F_m * $F_x2) > 0){
$X2 = $Middle;
}

# 二分法 Bisection Method
$BisectionMethod = $Middle;

# 収束判定
last if(abs($Middle - $PrevMiddle) < $Epsilon);
}

return $BisectionMethod;
}

# 対数正規分布 Logarithmic Normal Distribution
# 引数 変数 平均 標準偏差 ($X, $Mean, $StandardDeviation)
# 戻り値 対数正規分布 ($LogarithmicNormalDistribution)
sub LOGARITHMICNORMALDISTRIBUTION{
my ($X, $Mean, $StandardDeviation) = @_;
my $LogarithmicNormalDistribution = 0;

# 変数の確認
if($X == 0){
return 0;
}

my $LnX = log($X) / log(exp(1));

# 下側累積確率 Lower Probability
$LogarithmicNormalDistribution = 0.5 + (ERRORFUNCTION(($LnX - $Mean) / sqrt(2 * ($StandardDeviation * $StandardDeviation))) / 2);

return $LogarithmicNormalDistribution;
}

# 誤差関数 Error Function
# 引数 値 ($X)
# 戻り値 誤差関数 ($ErrorFunction)
sub ERRORFUNCTION{
my ($X) = @_;
my $ErrorFunction = 0;
my $RombergsMethod = 0;
my @Romberg = ();
my $Sum = 0;
my $a = 0;
my $b = $X;
my $h = 0;
my $N = 10;
my $Power4 = 0;
my $Epsilon = 1.0e-10;
my $Pi = atan2(1, 1) * 4;

if($X == 0){
# 誤差関数 Error Function
$ErrorFunction = 0;

return $ErrorFunction;
}

# ロンバーグ積分 Romberg's Method
# (0, 0)
$Romberg[0] = (($b - $a) / 2) * (&INTEGRAND($a) + &INTEGRAND($b));

for(my $i = 1; $i < $N; $i++){
$Power4 = 1;

for(my $j = $i; $j >= 0; $j--){
if($i == $j){
$Sum = 0;
$h = ($b - $a) / (2 ** $j);

for(my $k = 1; $k <= (2 ** ($i - 1)); $k++){
$Sum += &INTEGRAND($a + (((2 * $k) - 1) * $h));
}

# ($i, 0)
$Romberg[$j] = ($Romberg[$j - 1] / 2) + ($Sum * $h);
}else {
# 累乗
$Power4 = 4 * $Power4;

# リチャードソンの補外 Richardson's Extrapolation
$Romberg[$j] = (($Power4 * $Romberg[$j + 1]) - ($Romberg[$j])) / ($Power4 - 1);

# 収束判定
if(abs($Romberg[$j] - $Romberg[$j + 1]) < $Epsilon){
$Romberg[0] = $Romberg[$j];
$i = $N;
last;
}
}
}
}

# ロンバーグ積分 Romberg's Method
$RombergsMethod = $Romberg[0];

# 誤差関数 Error Function
$ErrorFunction = (2 / sqrt($Pi)) * $RombergsMethod;

return $ErrorFunction;
}

# 被積分関数 Integrand
# 引数 変数 被積分関数 ($X, \@Integrand)
# 戻り値 被積分関数 ($Integrand)
sub INTEGRAND{
my ($X) = @_;
my $Integrand = exp(-($X * $X));

return $Integrand;
}


参考URL
対数正規分布(逆関数) - 高精度計算サイト
オンライン コンパイラ/インタプリタ
テクニカル分析
プロフィール

Author:雨宮
Firefoxを使用しているので気づかなかったけど、IE6でソースコードを上手くコピーできない

5/3
携帯用ならIE6でもソースコードをコピーできる
携帯用

検索フォーム


あわせて読みたいブログパーツ
一寸先は闇 RSS

上記広告は1ヶ月以上更新のないブログに表示されています。新しい記事を書くことで広告を消せます。