Julia に翻訳--183 特殊な指数曲線回帰 b^y = a * x

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

Python による統計処理　
http://aoki2.si.gunma-u.ac.jp/Python/PythonDir.tar.gz

特殊な指数曲線回帰 b**y = a * x
http://aoki2.si.gunma-u.ac.jp/Python/exp3.pdf

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

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

==========#

using Statistics, Printf, Plots

function exp3(x, y; printout=true, graph=true)
    f(x, a, b) = a + b * x
    model = "b^y = a * x"
    n = length(x)
    n >= 2 || error("sample size must be greater than 1")
    all(x .> 0) || error("all x[i] must be positive")
    logx = log.(x)
    intercept, slope = sreg(logx, y)
    b = exp(1 / slope)
    a = exp(intercept / slope)
    #predicted = intercept + slope*logx
    predicted = f.(logx, intercept, slope)
    resid = y - predicted
    println(model)
    println("a = $a, b = $b")
    if printout
        @printf("%12s %12s %12s %12s\n", "x", "y", "pred.", "resid.")
        for i =1:n
            @printf("%12.6f %12.6f %12.6f %12.6f\n", x[i], y[i], predicted[i], resid[i])
        end
    end
    if graph
        plt = scatter(x, y, grid=false, tick_direction=:out, title="\$$model\$",
                      color=:blue, markerstrokecolor=:blue, label="")
        min, max = extrema(x)
        w = 0.05(max - min)
        x2 = range(min - w, max + w, length=1000)
        y2 = f.(log.(x2), intercept, slope)
        plot!(x2, y2, linewidth=0.5, color=:red, label="")
        display(plt)
    end
    Dict(:a => a, :b => b, :intercept => intercept, :slope => slope, :predicted => predicted, :resid => resid, :x => x, :logx => logx, :y => y, :model => model)
end

function sreg(x, y)
    a = cov(x, y) / var(x)
    mean(y) - a * mean(x), a
end

x = collect(1:10);
y = [0.631, 1.262, 1.631, 1.893, 2.096, 2.262, 2.402, 2.524, 2.631, 2.727];
exp3(x, y)
# b^y = a * x
# a = 2.000197950638151, b = 2.999971611554844
#            x            y        pred.       resid.
#     1.000000     0.631000     0.631025    -0.000025
#     2.000000     1.262000     1.261960     0.000040
#     3.000000     1.631000     1.631034    -0.000034
#     4.000000     1.893000     1.892896     0.000104
#     5.000000     2.096000     2.096011    -0.000011
#     6.000000     2.262000     2.261969     0.000031
#     7.000000     2.402000     2.402284    -0.000284
#     8.000000     2.524000     2.523831     0.000169
#     9.000000     2.631000     2.631043    -0.000043
#    10.000000     2.727000     2.726947     0.000053