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) || 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

.pdf(t, v) ⇒ Object 



5
6
7
 
# 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 



66
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
 
# File 'lib/distribution/t/ruby.rb', line 66

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

  return ptsub(q, n) if (n <= 5)
  return ptsub(q, n) if q <= 5.0e-3 && n <= 13

  f1 = 4.0 * (f = n.to_f)
  f5 = (f4 = (f3 = (f2 = f * f) * f) * f) * f
  f2 *= 96.0
  f3 *= 384.0
  f4 *= 92_160.0
  f5 *= 368_640.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 + 17_955.0 * u) / f5
  u + w
end

.ptsub(q, n) ⇒ Object

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



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

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 = 10_000.0
  w = 0.0
  loop do
    w += s
    return w if (s <= eps)
    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

.quantile(y, n) ⇒ Object Also known as: p_value

Returns the P-value of tdist().



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

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