Julia に翻訳--125 ウィリアムズ法，平均値の多重比較

#==========
Julia の修行をするときに，いろいろなプログラムを書き換えるのは有効な方法だ。
以下のプログラムを Julia に翻訳してみる。

ウィリアムズの方法による平均値の多重比較
http://aoki2.si.gunma-u.ac.jp/R/Williams.html

ファイル名: williams.jl　　関数名: williams

翻訳するときに書いたメモ

p 値は求めない

==========#

using Statistics, Rmath

function williams(data, group; method = "up") # or method = "down"
    func = method == "down" ? minimum : maximum
    index, ni = table(group)
    sumi = [sum(data[group .== g]) for g in index]
    meani = [mean(data[group .== g]) for g in index]
    vi = [var(data[group .== g]) for g in index]
    a = length(ni)
    phie = sum(ni) - a
    ve = sum((ni .- 1) .* vi) / phie
    t = [(func(cumsum(reverse(sumi[2:p])) ./
          cumsum(reverse(ni[2:p]))) .- meani[1]) /
          sqrt(ve * (1 / ni[1] + 1 / ni[p])) for p = a:-1:2]
    println("phi.e = $phie")
    println("t = $t")
    Dict(:phie => phie, :t => method == "down" ? -t : t)
end

function table(x) # indices が少ないとき
    indices = sort(unique(x))
    counts = zeros(Int, length(indices))
    for i in indexin(x, indices)
        counts[i] += 1
    end
    return indices, counts
end

data = [
    415, 380, 391, 413, 372, 359, 401,
    387, 378, 359, 391, 362, 351, 348,
    357, 379, 401, 412, 392, 356, 366,
    361, 351, 378, 332, 318, 344, 315,
    299, 308, 323, 351, 311, 285, 297
];
group = [repeat([i], 7) for i = 1:5]
group = vcat(group...)

williams(data, group, method = "down")
# phi.e = 30
# t = [-7.076037814702816, -4.217674245029326, -1.4164779467493078, -1.9690949035528496]