Module: Distribution::T::Ruby_

Defined in:
lib/distribution/t/ruby.rb

Class Method Summary collapse

Class Method Details

.cdf(t, n) ⇒ Object

Returns the integral of t-distribution with n degrees of freedom over (-Infty, x].


10
11
12
# File 'lib/distribution/t/ruby.rb', line 10

def cdf(t, n)
  p_t(n, t)
end

.p_t(df, t) ⇒ Object

t-distribution ([1]) (-infty, x]


16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
# File 'lib/distribution/t/ruby.rb', line 16

def p_t(df, t)
  if df.to_i!=df
    x=(t+Math.sqrt(t**2+df)) / (2*Math.sqrt(t**2+df))
    return Math.regularized_beta(x,df/2.0,df/2.0)
  end
  df=df.to_i
  c2 = df.to_f / (df + t * t);
  s = Math.sqrt(1.0 - c2)
  s = -s if t < 0.0
  p = 0.0;
  i = df % 2 + 2
  while i <= df
    p += s
    s *= (i - 1) * c2 / i
    i += 2
  end
  
  if df.is_a? Float or df & 1 != 0
    0.5+(p*Math.sqrt(c2)+Math.atan(t/Math.sqrt(df))) / Math::PI
  else
    (1.0 + p) / 2.0
  end
end

.p_value(y, n) ⇒ Object

Returns the P-value of tdist().


99
100
101
102
103
104
105
# File 'lib/distribution/t/ruby.rb', line 99

def p_value(y,n)
  if y > 0.5
    pt(2.0 - y*2.0, n)
  else
    - pt(y*2.0, n)
  end
end

.pdf(t, v) ⇒ Object


5
6
7
8
# File 'lib/distribution/t/ruby.rb', line 5

def pdf(t,v)
  ((Math.gamma((v+1) / 2.0)) / (Math.sqrt(v*Math::PI)*Math.gamma(v/2.0))) * ((1+(t**2 / v.to_f))**(-(v+1) / 2.0))
  
end

.pt(q, n) ⇒ Object


67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
# File 'lib/distribution/t/ruby.rb', line 67

def pt(q, n)
  q = q.to_f
  if(q < 1.0e-5 || q > 1.0 || n < 1)
  $stderr.printf("Error : Illegal parameter in pt()!\n")
  return 0.0
  end
  
  if(n <= 5) then return ptsub(q, n) end
  if(q <= 5.0e-3 && n <= 13) then return ptsub(q, n) end
  
  f1 = 4.0 * (f = n.to_f)
  f5 = (f4 = (f3 = (f2 = f * f) * f) * f) * f
  f2 *= 96.0
  f3 *= 384.0
  f4 *= 92160.0
  f5 *= 368640.0
  u = Normal.p_value(1.0 - q / 2.0)
  
  w0 = (u2 = u * u) * u
  w1 = w0 * u2
  w2 = w1 * u2
  w3 = w2 * u2
  w4 = w3 * u2
  w = (w0 + u) / f1
  w += (5.0 * w1 + 16.0 * w0 + 3.0 * u) / f2
  w += (3.0 * w2 + 19.0 * w1 + 17.0 * w0 - 15.0 * u) / f3
  w += (79.0 * w3 + 776.0 * w2 + 1482.0 * w1 - 1920.0 * w0 - 9450.0 * u) / f4
  w += (27.0 * w4 + 339.0 * w3 + 930.0 * w2 - 1782.0 * w1 - 765.0 * w0 + 17955.0 * u) / f5
  u + w
end

.ptsub(q, n) ⇒ Object

inverse of t-distribution ([2]) (-infty, -q/2] + [q/2, infty)


43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
# File 'lib/distribution/t/ruby.rb', line 43

def ptsub(q, n)
q = q.to_f
if(n == 1 && 0.001 < q && q < 0.01)
eps = 1.0e-4
elsif (n == 2 && q < 0.0001)
eps = 1.0e-4
elsif (n == 1 && q < 0.001)
eps = 1.0e-2
else
eps = 1.0e-5
end
s = 10000.0
w = 0.0
loop do
w += s
if(s <= eps) then return w end
if((qe = 2.0 - p_t(n, w)*2.0 - q) == 0.0) then return w end
if(qe < 0.0)
  w -= s
  s /= 10.0 #/
end
end
end