Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Ruby script to calculate an adjusted coefficient of determination.
#! /usr/local/bin/ruby
#*******************************************************************
# Ruby script to calculate an adjusted coefficient of determination.
#*******************************************************************
#
class AdjustedCoefficientOfDetermination
# 説明変数と目的変数
X = [
[17.5, 17.0, 18.5, 16.0, 19.0, 19.5, 16.0, 18.0, 19.0, 19.5],
[30 , 25 , 20 , 30 , 45 , 35 , 25 , 35 , 35 , 40 ]
]
Y = [45, 38, 41, 34, 59, 47, 35, 43, 54, 52]
def exec
puts "説明変数 X1 = [#{X[0].join(', ')}]"
puts "説明変数 X2 = [#{X[1].join(', ')}]"
puts "目的変数 Y = [#{Y.join(', ')}]"
puts "---"
# 重回帰式算出
reg = X.reg_multi(Y)
puts " a = %20.16f" % reg[0]
puts " b = %20.16f" % reg[1]
puts " c = %20.16f" % reg[2]
# 推定値配列
y_e = calc_estimations(X, Y, reg)
# 標本値 Y (目的変数)の平均
y_b = Y.inject(0) { |s, a| s += a } / Y.size.to_f
# 自由度調整済み決定係数
p, n = X.size, Y.size
r_2 = calc_s_e(Y, y_e) / (n - p - 1)
r_2 /= calc_s_r(Y, y_b) / (n - 1)
r_2 = 1 - r_2
puts "---"
puts "自由度調整済み決定係数 R2f = %20.16f" % r_2
rescue => e
$stderr.puts "[#{e.class}] #{e.message}"
e.backtrace.each{ |tr| $stderr.puts "\t#{tr}" }
exit 1
end
private
# 推定値
#
# @param xs: 説明変数配列(2個)
# @param y: 目的変数配列
# @param reg: 重回帰式の値 a, b, c
# @return y_e: 推定値配列
def calc_estimations(xs, y, reg)
y_e = Array.new
begin
a, b, c = reg
y.zip(xs[0], xs[1]).each do |y, x1, x2|
y_e << a + b * x1 + c * x2
end
return y_e
rescue => e
raise
end
end
# 標本値の変動
#
# @param y_s: 標本値(目的変数)配列
# @param y_b: 標本値(目的変数)の平均
# @return s_r: 標本値の変動
def calc_s_r(y_s, y_b)
s_r = 0.0
begin
y_s.each do |a|
v = a - y_b
s_r += v * v
end
return s_r
rescue => e
raise
end
end
# 残差の変動
#
# @param y_s: 標本値(目的変数)配列
# @param y_e: 推定値配列
# @return s_e: 残差の変動
def calc_s_e(y_s, y_e)
s_e = 0.0
begin
y_s.zip(y_e).each do |a, b|
v = a - b
s_e += v * v
end
return s_e
rescue => e
raise
end
end
end
class Array
def reg_multi(y)
# 元の数、変量内のサンプル数
e_size, s_size = self.size, y.size
# 以下の場合は例外スロー
# - 引数の配列が Array クラスでない
# - 自身配列が空
# - 配列サイズが異なれば例外
raise "Argument is not a Array class!" unless y.class == Array
raise "Self array is nil!" if e_size == 0
raise "Argument array size is invalid!" unless self[0].size == s_size
1.upto(e_size - 1) do |i|
raise "Argument array size is invalid!" unless self[0].size == self[i].size
end
x1, x2 = self
sum_x1 = x1.sum
sum_x2 = x2.sum
sum_y = y.sum
sum_x1x1 = x1.map { |a| a * a }.sum
sum_x1x2 = x1.zip(x2).map { |a, b| a * b }.sum
sum_x1y = x1.zip(y).map { |a, b| a * b }.sum
sum_x2x1 = sum_x1x2
sum_x2x2 = x2.map { |a| a * a }.sum
sum_x2y = x2.zip(y).map { |a, b| a * b }.sum
mtx = [
[s_size, sum_x1 , sum_x2 , sum_y ],
[sum_x1, sum_x1x1, sum_x1x2, sum_x1y],
[sum_x2, sum_x2x1, sum_x2x2, sum_x2y]
]
# 連立方程式を解く (ガウスの消去法)
return gauss_e(mtx)
end
private
# ガウスの消去法
def gauss_e(ary)
# 行数
n = ary.size
# 前進消去
0.upto(n - 2) do |k|
(k + 1).upto(n - 1) do |i|
if ary[k][k] == 0
puts "解けない!"
exit 1
end
d = ary[i][k] / ary[k][k].to_f
(k + 1).upto(n) do |j|
ary[i][j] -= ary[k][j] * d
end
end
end
# 後退代入
(n - 1).downto(0) do |i|
if ary[i][i] == 0
puts "解けない!"
exit 1
end
d = ary[i][n]
(i + 1).upto(n - 1) do |j|
d -= ary[i][j] * ary[j][n]
end
ary[i][n] = d / ary[i][i].to_f
end
return ary.map { |a| a[-1] }
end
end
AdjustedCoefficientOfDetermination.new.exec if __FILE__ == $0
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment