Edited at

Juliaで学ぶ統計学入門(第2章:1次元のデータ)

More than 5 years have passed since last update.


はじめに

Facebookグループ『東大の統計学教科書を読み進める会』で、

統計学入門(通称:赤本)の読み進めを行っています。

読み進める中で、「ここの処理をJuliaではどうやるの?」が気になったのでコードを書き連ねました。


統計学入門(通称:赤本)とは?

字義通り。「東京大学教養学部統計学教室 編」の統計学入門書です。

参考:統計学入門


Juliaとは?

ポストR。ベターMATLAB。処理が超早い。

参考:Juliaの公式サイト


調査対象の章

第2章「1次元のデータ」


データ

データをJuliaへ読み込ませる。


DataRead.jl

# [変数定義]

##整数
x = 1
##小数
x = 1.0
##文字列
x = "abc"
##配列
x = [1,2,3]
##連想配列
x = {:id:1, name:"nezuq", address:["japan","tokyo","nakameguro"]}
##行列
x = [1 2 3 4 5; 6 7 8 9 10; 11 12 13 14 15]
##データフレーム
using DataFrames
df = DataFrame(縦軸=["a", "a", "c", "b", "b"], 横軸=[10, 20, 30, 40, 50])

# [ファイルの読み込み]
##ファイルのオープン
txt = open("abc.txt")
##ファイルの読み込み
###テキストファイルの場合
x = readlines(txt)
###CSVファイルの場合
x = readcsv(txt)
##ファイルのクローズ
close(txt)

# [Rデータセットの呼び出し]
##Rデータセット用パッケージのインストール
using RDatasets
##irisデータを呼び出す
x = dataset("datasets", "iris")

# [SQLiteの読み込み]
##SQLite管理用パッケージのインストール
using SQLite
##SQLiteファイルへの接続
con = SQLite.connect("abc.db")
##クエリの発行
x = query("select * from m_user", con)
##ファイルのクローズ
close(con)



基本統計量

基本統計量を算出する。


Summary.jl

data = reduce((sum,a) -> vcat(sum,a), [repeat([d[1]],inner=[d[2]]) for d = [(1,43), (2,179), (3,238), (4,38), (5,2)]])

#[度数]
##Baseパッケージ
x = hist(data, int(ceil(1 + log2(length(data))))) # => (0.0:1.0:5.0,[43,179,238,38,2])
##StatsBaseパッケージ
Pkg.add("StatsBase")
using StatsBase
x = counts(data, minimum(data):maximum(data)) # => [43, 179, 238, 38, 2]

#[相対度数]
using StatsBase
x = proportions(data, minimum(data):maximum(data)) # => [0.086, 0.358, 0.476, 0.076, 0.004]

#[平均値]
##算術平均
x = mean(data) # => 2.554
##幾何平均
using StatsBase
x = geomean([1.218, 1.305, 1.536, 1.500, 1.129]) # => 1.3282714663720916
##調和平均
using StatsBase
x = harmmean([25, 15]) # => 18.75

#[中央値]
x = median(data) # => 3.0

#[最頻値]
x = mode(data) # => 3

#[ミッド・レンジ]
x = (maximum(data) - minimum(data)) / 2 # => 2.0

#[分位点]
x = quantile(data, 1/4) # => 2.0

#[分散]
x = var(data) # => 0.5962765531062114

#[標準偏差]
x = std(data) # => 0.7721894541537144



グラフ描画

データからグラフを描画する。


Plot.jl

data = reduce((sum,a) -> vcat(sum,a), [repeat([d[1]],inner=[d[2]]) for d = [(1,43),(2,179),(3,238),(4,38),(5,2)]])

x = hist(data, int(ceil(1 + log2(length(data)))))

#[棒グラフ]
Pkg.add("Gadfly")
using Gadfly
plot(x = collect(x[1])[2:end], y = x[2], Geom.bar)

#[線グラフ]
using Gadfly
plot(x = collect(x[1])[2:end], y = x[2], Geom.smooth)



練習問題の解答


prac02.jl

#[代表値の間隔]

#参照:統計局ホームページ/家計調査報告(貯蓄・負債編)-平成25年(2013年)平均結果速報-(二人以上の世帯)
# http://www.stat.go.jp/data/sav/sokuhou/nen/
k = (1023 - 50) / (1739 - 1023) # => 1.3589385474860336

#[ジニ係数]
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]
##平均差
diffmean(nary) = reduce((sum1,x1) -> sum1 + reduce((sum2,x2) -> sum2 + abs(x1 - x2), 0, nary), 0, nary) / (length(nary) ^ 2)
diffmean(A) # => 2.76
diffmean(B) # => 3.76
diffmean(C) # => 1.2
##ジニ係数
geni(nary) = reduce((sum1,x1) -> sum1 + reduce((sum2,x2) -> sum2 + abs(x1 - x2), 0, nary), 0, nary) / 2(length(nary) ^ 2)mean(nary)
geni(A) # => 0.276
geni(B) # => 0.376
geni(C) # => 0.12

#[エントロピー]
cnts_now = [32, 19, 10, 24, 15]
cnts_old = [28, 13, 18, 29, 12]
Pkg.add("StatsBase")
using StatsBase
entropy(cnts_now / sum(cnts_now)) # => 1.537492332302312(logの底がeで計算)
entropy(cnts_old / sum(cnts_old)) # => 1.5437380015770872(logの底がeで計算)

#[諸得点]
##標準得点
using StatsBase
b_now = zscore([19], mean(cnts_now), std(cnts_now)) # => -0.118262
b_old = zscore([13], mean(cnts_old), std(cnts_old)) # => -0.864923
##偏差値得点
10 * b_now + 50 # => 48.8174
10 * b_old + 50 # => 41.3508



感想

Juliaで統計処理を行う為のコマンド(方法)を集めたJLJPWikiが欲しい。