統計学入門のお供 Julia (1) - 標準得点、偏差値得点 - Qiita
統計学入門のお供 Julia (2) - ピアソンの積率相関係数 - Qiita
統計学入門のお供 Julia (3) - 順位相関係数と自己相関係数 - Qiita
統計学入門のお供 Julia (4) - 最小二乗法と決定係数 - Qiita
統計学入門のお供 Julia (5) - ベイズの定理 - Qiita
統計学入門のお供 Julia (6) - ド・メレの問題 = Oiita
統計学を勉強するため「統計学入門」（東京大学出版会）を読み始めました。練習問題があるのですが、これも勉強中のJuliaを使って計算していきたいと思います。
#１．スタージェスの公式
p22のスタージェスの公式を計算します。
Julia言語のBaseのMathematicsカテゴリのlog10関数を使います。BaseなのでMathematicsカテゴリの関数は特に何もしなくとも利用できます。
Mathematics
https://docs.julialang.org/en/v1/base/math/
k \fallingdotseq 1 + \frac {log_{10} n}{log_{10} 2}
図2.1と表2.1のデータは、n=373なので以下の計算式になります。
k = 1+log10(373)/log10(2)
9.543031820255239
ですからk=10とするのが適当ということになります。
#２．ジニ係数
問2.2(P39)のジニ係数を求めます。ジニ係数の式は以下のようになります。
GI = \sum_i \sum_j \frac {|x_i-x_j|}{2n^2 \bar{x}}
juliaを使って計算します。
標準ライブラリのStatisticsを使います。具体的にはmean関数で平均を求めます。
https://docs.julialang.org/en/v1/stdlib/Statistics/index.html
using Statistics # 標準ライブラリ
a=[0,3,3,5,5,5,5,7,7,10]
b=[0,1,2,3,5,5,7,8,9,10]
c=[3,4,4,5,5,5,5,6,6,7]
function mygini(a)
s = 0.0
n = length(a)
d = 2 * n^2 * mean(a) # mean=平均値
for i in 1:n
for j in 1:n
s += ( abs( (a[i]-a[j]) ) / d )
end
end
s
end
mygini(a)
0.2760000000000002
mygini(b)
0.3760000000000002
mygini(c)
0.12000000000000008
貧富の差はbが最も大きく、cが最も小さくなりました。
##２－１．分布図
データの理解を深めるために分布図を描いてみます。aのデータに関して描きます。
StatsBaseというライブラリを使います。countmap関数で度数を計算していますが、もっと直接的に求める関数があるかもしれません。
http://juliastats.github.io/StatsBase.jl/stable/
import Pkg; Pkg.add("StatsBase")
using StatsBase
cma=StatsBase.countmap(a) # 度数のmap dict
sma=sort(cma) # dictのkeyでソート
sma_x=[key for key in keys(sma)] # x座標=key
sma_y=[sma[key] for key in keys(sma)] # y座標=度数
using Plots
gr()
plot(sma_x, sma_y)
StatsBase.countmapの入力と出力を示します。
a=[0,3,3,5,5,5,5,7,7,10]
cma=StatsBase.countmap(a)
Dict{Int64,Int64} with 5 entries:
0 => 1
7 => 2
10 => 1
3 => 2
5 => 4
sortの出力です。
sma=sort(cma)
OrderedCollections.OrderedDict{Int64,Int64} with 5 entries:
0 => 1
3 => 2
5 => 4
7 => 2
10 => 1
#３．エントロピー
問2.3(p40)のエントロピーを計算します。
Julia言語のBaseのMathematicsカテゴリのlog10関数を使います。BaseなのでMathematicsカテゴリの関数は特に何もしなくとも利用できます。
Mathematics
https://docs.julialang.org/en/v1/base/math/
h1=[32,19,10,24,15]
h2=[28,13,18,29,12]
function myentropy(h)
s = 0.0
n = sum(h)
for i in 1:length(h)
p = h[i]/n
s+=p*log10(p) # log10関数
end
-s
end
myentropy(h1)
0.667724435887455
myentropy(h2)
0.6704368955892825
#４．標準得点、偏差値得点
P40 問2.4を解きます
標準ライブラリのStatisticsの、平均値(mean)や標準偏差（std)を使います。
https://docs.julialang.org/en/v1/stdlib/Statistics/index.html
using Statistics
b=[0,1,2,3,5,5,7,8,9,10]
function myscore(x)
s = 0.0
n = length(x)
m = mean(x) # 平均値
s = std(x) # 標準偏差
z = zeros(n) # サイズnの配列
t = zeros(n) # サイズnの配列
for i in 1:n
z[i] = (x[i]-m)/s # 標準得点
t[i] = 10*z[i]+50 # 偏差値得点
end
println(z)
println(t)
1
end
myscore(b)
[-1.44338, -1.1547, -0.866025, -0.57735, 0.0, 0.0, 0.57735, 0.866025, 1.1547, 1.44338]
[35.5662, 38.453, 41.3397, 44.2265, 50.0, 50.0, 55.7735, 58.6603, 61.547, 64.4338]
ちょっとデータbがあまりにも感覚的に掴みにくいので、ちょっと身近な例も見てみます。６人構成のクラスでの数学、英語、国語の点数を考えます。
b1=[53,55,56,59,62,75] # 数学
b2=[0 ,58,70,75,77,80] # 英吾
b3=[75,78,62,82,84,85] # 国語
myscore(b1)
myscore(b2)
myscore(b3)
[-0.875, -0.625, -0.5, -0.125, 0.25, 1.875] # 数学の標準得点
[41.25, 43.75, 45.0, 48.75, 52.5, 68.75] # 数学の偏差値得点
[-1.97428, -0.0658094, 0.329047, 0.493571, 0.55938, 0.658094] # 英語の標準得点
[30.2572, 49.3419, 53.2905, 54.9357, 55.5938, 56.5809] # 英語の偏差値得点
[-0.311967, 0.0389959, -1.83281, 0.506947, 0.740922, 0.85791] # 国語の標準得点
[46.8803, 50.39, 31.6719, 55.0695, 57.4092, 58.5791] # 国語の偏差値得点
数学と英語、国語で同じ75点でも、偏差値得点で考えれば68.75、54.9357、46.8803の評価になります。数学は良い点数ですが、国語は平均（50点）以下の評価になります。
数学:75->68.75
英語:75->54.9357
国語:75->46.8803
今回は以上です。