Class: Statsample::MLE::Normal

Inherits:
BaseMLE
  • Object
show all
Defined in:
lib/statsample/mle/normal.rb

Overview

Normal Distribution MLE estimation. Usage:

mle=Statsample::MLE::Normal.new
mle.newton_raphson(x,y)
beta=mle.parameters
likehood=mle.likehood(x,y,beta)
iterations=mle.iterations

Constant Summary

Constants inherited from BaseMLE

BaseMLE::ITERATIONS, BaseMLE::MIN_DIFF, BaseMLE::MIN_DIFF_PARAMETERS

Instance Attribute Summary

Attributes inherited from BaseMLE

#iterations, #output, #parameters, #stop_criteria, #var_cov_matrix, #verbose

Instance Method Summary collapse

Methods inherited from BaseMLE

#initialize, #likehood, #newton_raphson, #set_default_parameters

Constructor Details

This class inherits a constructor from Statsample::MLE::BaseMLE

Instance Method Details

#first_derivative(x, y, p) ⇒ Object

First derivative for Normal Model. p should be [k+1,1], because the last parameter is sigma^2



24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
 
# File 'lib/statsample/mle/normal.rb', line 24

def first_derivative(x,y,p)
    raise "x.rows!=y.rows" if x.row_size!=y.row_size
    raise "x.columns+1!=p.rows" if x.column_size+1!=p.row_size
    
    n = x.row_size
    k = x.column_size
    b = Array.new(k)
    k.times{|i| b[i]=[p[i,0]]}
    beta = Matrix.rows(b)
    sigma2 = p[k,0]
    sigma4=sigma2*sigma2
    e = y-(x*(beta))
    xte = x.transpose*(e)
    ete = e.transpose*(e)
    #rows of the Jacobian
    rows = Array.new(k+1)
    k.times{|i| rows[i] = [xte[i,0] / sigma2]}
    rows[k] = [ete[0,0] / (2*sigma4) - n / (2*sigma2)]
    fd = Matrix.rows(rows, true)
end

#log_likehood(x, y, b) ⇒ Object

Total MLE for given X, Y and B matrices



14
15
16
17
18
19
20
21
 
# File 'lib/statsample/mle/normal.rb', line 14

def log_likehood(x,y,b)
    n=x.row_size.to_f
    sigma2=b[b.row_size-1,0]
    betas=Matrix.columns([b.column(0). to_a[0...b.row_size-1]])
    e=y-(x*betas)
    last=(1 / (2*sigma2))*e.t*e
    (-(n / 2.0) * Math::log(2*Math::PI))-((n / 2.0)*Math::log(sigma2)) - last[0,0]
end

#second_derivative(x, y, p) ⇒ Object

p should be [k+1,1], because the last parameter is sigma^2



47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
 
# File 'lib/statsample/mle/normal.rb', line 47

def second_derivative(x,y,p)
    raise "x.rows!=y.rows" if x.row_size!=y.row_size
    raise "x.columns+1!=p.rows" if x.column_size+1!=p.row_size

    n = x.row_size
    k = x.column_size
    b = Array.new(k)
    k.times{|i| b[i]=[p[i,0]]}
    beta = Matrix.rows(b)
    sigma2 = p[k,0]
    sigma4=sigma2*sigma2
    sigma6 = sigma2*sigma2*sigma2
    e = y-(x*(beta))
    xtx = x.transpose*(x)
    xte = x.transpose*(e)
    ete = e.transpose*(e)
    #rows of the Hessian
    rows = Array.new(k+1)
    k.times do |i|
        row = Array.new(k+1)
        k.times do |j|
            row[j] = -xtx[i,j] / sigma2
        end
        row[k] = -xte[i,0] / sigma4
        rows[i] = row
    end
    last_row = Array.new(k+1)
    k.times do |i|
        last_row[i] = -xte[i,0] / sigma4
    end
    last_row[k] = 2*sigma4 - ete[0,0] / sigma6
    rows[k] = last_row
    sd = Matrix.rows(rows, true)
end