Доброго всем дня, коллеги!
Откопал код функции обратного распределения Стьдента, однако при задании высокого уровня надежности
p она возвращает весьма отдаленное от СТЬЮДЕНТ.ОБР значение. Есть подозрение, что проблема кроется в константах, заданных в подфункции "
NorDev". Если я правильно понял код и Википедию, то
NorDev вроде как должна использовать значения гамма-функции, но так ли это не знаю.
В общем, задача заключается в том, чтобы составить код, возвращающий в точности такие же значения как и СТЬЮДЕНТ.ОБР.
Скрытый текст |
---|
Код |
---|
Function Student(p As Double, f1 As Double) As Double
Dim p4 As Double
Dim x0 As Double
Dim t As Double
Dim sMsg As String
With WorksheetFunction
p4 = p * .Pi()
' точные значения
If (f1 = 1) Then
Student = -Cos(p4) / Sin(p4)
Exit Function
End If
If f1 = 2 Then
Student = Sqr(1# / (2# * p * (1# - p)) - 2#) * Sgn(p - 0.5)
Exit Function
End If
' Модификация Гольдберга и Левина приближения Пизера
If (f1 - 1) * (f1 - 2) <> 0# Then
x0 = NorDev(p)
t = 1# + (1# + x0 * x0) / (4# * f1)
t = t + (3# + 16# * x0 * x0 + 5 * (x0 ^ 4#)) / (96# * f1 * f1)
t = x0 * t
End If
End With
Student = t
Exit Function
End Function
Private Function NorDev(p As Double) As Double
'On Error GoTo ErrH
Dim sMsg As String
'Dim sTit As String
Dim q As Double
Dim X As Double
Dim r As Double
Dim nd As Double
Dim a(4) As Double
Dim b(5) As Double
Dim c(4) As Double
Dim d(3) As Double
a(0) = 2.50662823881
a(1) = -18.6150006252
a(2) = 41.39119773531
a(3) = -25.44106049637
b(0) = 0
b(1) = -8.4735109309
b(2) = 23.08336743743
b(3) = -21.06224101826
b(4) = 3.13082909833
c(0) = -2.78718931138
c(1) = -2.29796479134
c(2) = 4.85014127135
c(3) = 2.32121276858
d(0) = 0
d(1) = 3.54388924762
d(2) = 1.63706781897
q = p - 0.5
If (Abs(q) > 0.42) Then ' хвосты распределения
r = p
If (q > 0) Then r = 1 - p
If (r <= 0) Then ' защита от неправильного ввода
sMsg = "Ошибка данных : P = " & Format(p, "0.000") _
& " и не принадлежит интервалу [0...1]"
Else ' аппроксимация для хвостов распределения
r = Sqr(-Log(r))
X = (((c(3) * r + c(2)) * r + c(1)) * r + c(0))
If (q < 0) Then X = -X
nd = X / ((d(2) * r + d(1)) * r + 1)
End If
Else ' аппроксимация для центральной части распределения
r = q * q
X = q * (((a(3) * r + a(2)) * r + a(1)) * r + a(0))
nd = X / ((((b(4) * r + b(3)) * r + b(2)) * r + b(1)) * r + 1)
End If
NorDev = nd
Exit Function
End Function
|
|
Буду признателен за любую помощь.