Julia #2: Bayes Neural Net – Science To Medicine

————————————
ベイズ推論による機械学習入門　須山敦志　著　を読む
————————————
Terminalから目的のフォルダでJulia 0.7を立ち上げて、”]”でパッケージモードへ:

MacBook-Pro-6:~ $ cd Desktop/BayesBook-master
MacBook-Pro-6:~ $ julia

   _       _ _(_)_     |  A fresh approach to technical computing
  (_)     | (_) (_)    |  Documentation: https://docs.julialang.org
   _ _   _| |_  __ _   |  Type "?" for help, "]?" for Pkg help.
  | | | | | | |/ _` |  |
  | | |_| | | | (_| |  |  Version 0.7.0 (2018-08-08 06:46 UTC)
 _/ |\__'_|_|_|\__'_|  |  Official http://julialang.org/ release
|__/                   |  x86_64-apple-darwin14.5.0

(v0.7) pkg>

MacBook-Pro-6:~ $ cd Desktop/BayesBook-master

MacBook-Pro-6:~ $ julia

_ _ _(_)_ | A fresh approach to technical computing

(_) | (_) (_) | Documentation: https://docs.julialang.org

_ _ _| |_ __ _ | Type "?" for help, "]?" for Pkg help.

| | | | | | |/ _` | |

| | |_| | | | (_| | | Version 0.7.0 (2018-08-08 06:46 UTC)

_/ |\__'_|_|_|\__'_| | Official http://julialang.org/ release

|__/ | x86_64-apple-darwin14.5.0

(v0.7) pkg>

generate 命令でバッケージのテンプレートを作成し、add命令で使用する追加モジュールDistributionsをtomlファイルへ登録する。

(v0.7) pkg> generate BayesNeuralNet
Generating project BayesNeuralNet:
    BayesNeuralNet/Project.toml
    BayesNeuralNet/src/BayesNeuralNet.jl

(v0.7) pkg> generate BayesNeuralNet

Generating project BayesNeuralNet:

BayesNeuralNet/Project.toml

BayesNeuralNet/src/BayesNeuralNet.jl

一度Juliaを終了し、作成したBayesNeuralNetディレクトリに入ります。(これ省略すると、一つ上位のフォルダにパッケージpoml等が更新されてしまうので要注意！
Juliaを起動して]を押してパッケージモードにして、パッケージの中で利用するパッケージ、ここではDistributionsを以下のようにしてManifest.tomlに登録。

(v0.7) pkg> activate .
(BayesNeuralNet) pkg> add Distributions
  Updating git-repo `https://github.com/JuliaRegistries/General.git`
 Resolving package versions...
  Updating `Project.toml`
  [31c24e10] + Distributions v0.21.0
  Updating `Manifest.toml`
..................

(v0.7) pkg> activate .

(BayesNeuralNet) pkg> add Distributions

Updating git-repo `https://github.com/JuliaRegistries/General.git`

Resolving package versions...

Updating `Project.toml`

[31c24e10] + Distributions v0.21.0

Updating `Manifest.toml`

..................

BayesNeuralNetのフォルダにはProject.tomlが出来上がっているとともに、Manifest.tomlの内容が更新されている。
delキーを押してパッケージモードを終了した後、

julia> using BayesNeuralNet
[ Info: Precompiling BayesNeuralNet [db08b2e8-c615-11e9-225a-adac7eced69b]

1 2	julia> using BayesNeuralNet [ Info: Precompiling BayesNeuralNet [db08b2e8-c615-11e9-225a-adac7eced69b]

Jupyter Notebookを立ち上げる：

julia> using Pkg

julia> Pkg.add("IJulia")
  Updating registry at `~/.julia/registries/General`
  Updating git-repo `https://github.com/JuliaRegistries/General.git`
 Resolving package versions...
 Installed Colors ─ v0.9.6
  Updating `~/.julia/environments/v0.7/Project.toml`
 [no changes]
  Updating `~/.julia/environments/v0.7/Manifest.toml`
  [5ae59095] ↑ Colors v0.9.5 ⇒ v0.9.6

julia> using IJulia

julia> IJulia.notebook()
[ Info: running setenv(`/Users/********/.julia/conda/3/bin/jupyter notebook`,["XP
............

julia> using Pkg

julia> Pkg.add("IJulia")

Updating registry at `~/.julia/registries/General`

Updating git-repo `https://github.com/JuliaRegistries/General.git`

Resolving package versions...

Installed Colors ─ v0.9.6

Updating `~/.julia/environments/v0.7/Project.toml`

[no changes]

Updating `~/.julia/environments/v0.7/Manifest.toml`

[5ae59095] ↑ Colors v0.9.5 ⇒ v0.9.6

julia> using IJulia

julia> IJulia.notebook()

[ Info: running setenv(`/Users/********/.julia/conda/3/bin/jupyter notebook`,["XP

............

立ち上がったJupyterから, Julia 0.7を選択してdemoプロうグラムを入力：

####################################
## Demo script for Bayesian neural network.

using PyPlot, PyCall

####################################

## Demo script for Bayesian neural network.

using PyPlot, PyCall

作成した自作パッケージBayesNeuralNetへのパスを間違えないように、

push!(LOAD_PATH, "/Users/*****/Desktop/BayesBook-master/src/BayesNeuralNet")
import BayesNeuralNet

1 2	push!(LOAD_PATH, "/Users/*****/Desktop/BayesBook-master/src/BayesNeuralNet") import BayesNeuralNet

入力して、エラーがでなければBayesNeuralNetパッケージ導入成功！

"""
Sample neural nets from prior.
"""
function sample_test()
    # model parameters
    D = 1 # output
    K = 3 # hidden
    M = 2 # input
    sigma2_w = 10.0
    sigma2_y = 0.1
    
    xmin = -5
    xmax = 5
    N_lin = 1000
    X_lin = ones(M, N_lin)
    X_lin[1,:] = linspace(xmin, xmax, N_lin)
    X_lin[2,:] = 1 # bias

    # visualize
    num_samples = 5
    figure("Function samples")
    clf()
    for i in 1 : num_samples
        _, Y_true, _, _ = BayesNeuralNet.sample_data_from_prior(X_lin, sigma2_w, sigma2_y, D, K)
        plot(X_lin[1,:], Y_true)
        xlim([xmin, xmax])
    end
    ratey = (ylim()[2] - ylim()[1]) * 0.1
    ratex = (xlim()[2] - xlim()[1]) * 0.1
    text(xlim()[1] + ratex, ylim()[2] - ratey, @sprintf("K=%d", K), fontsize=18)
    show()
end

"""

Sample neural nets from prior.

"""

function sample_test()

# model parameters

D = 1 # output

K = 3 # hidden

M = 2 # input

sigma2_w = 10.0

sigma2_y = 0.1

xmin = -5

xmax = 5

N_lin = 1000

X_lin = ones(M, N_lin)

X_lin[1,:] = linspace(xmin, xmax, N_lin)

X_lin[2,:] = 1 # bias

# visualize

num_samples = 5

figure("Function samples")

clf()

for i in 1 : num_samples

_, Y_true, _, _ = BayesNeuralNet.sample_data_from_prior(X_lin, sigma2_w, sigma2_y, D, K)

plot(X_lin[1,:], Y_true)

xlim([xmin, xmax])

end

ratey = (ylim()[2] - ylim()[1]) * 0.1

ratex = (xlim()[2] - xlim()[1]) * 0.1

text(xlim()[1] + ratex, ylim()[2] - ratey, @sprintf("K=%d", K), fontsize=18)

show()

end

Out: sample_test

"""
Run a test script of variational inference for Bayesian neural net.
"""
function test()
    #################
    # prepara data
    
    # data size
    D = 1 # output
    M = 2 # input
        
    # function setting
    xmin = -2
    xmax = 4
    N_lin = 1000
    X_lin = ones(M, N_lin)
    X_lin[1,:] = linspace(xmin, xmax, N_lin)
    X_lin[2,:] = 1 # bias

    # training data
    N = 50 # data size
    X = 2*rand(M, N) - 0.0 # input
    X[2,:] = 1.0 # bias
    Y = 0.5*sin.(2*pi * X[1,:]/3) + 0.05 * randn(N)    
    
    # model parameters
    K = 5
    sigma2_w = 10.0
    sigma2_y = 0.01
    
    ################
    # inference
    alpha = 1.0e-5
    max_iter = 100000
    mu1, rho1, mu2, rho2 = BayesNeuralNet.VI(Y, X, sigma2_w, sigma2_y, K, alpha, max_iter)
    Y_mean = [mu2'* tanh.(mu1'X_lin[:,n]) for n in 1 : N_lin]

    ################
    # visualize        
    figure("result")
    clf()
    Y_list = []
    num_samples = 100
    for i in 1 : num_samples
        Y_est, _ = BayesNeuralNet.sample_data_from_posterior(X_lin, mu1, rho1, mu2, rho2, sigma2_y, D)
        push!(Y_list, Y_est)
        plot(X_lin[1,:], Y_est, "-c", alpha=0.25)
    end
    plot(X[1,:], Y, "ok")
    plot(X_lin[1,:], Y_mean, "b-")
    xlim([xmin, xmax])
    xlabel("x")
    ylabel("y")
    show()
end

"""

Run a test script of variational inference for Bayesian neural net.

"""

function test()

#################

# prepara data

# data size

D = 1 # output

M = 2 # input

# function setting

xmin = -2

xmax = 4

N_lin = 1000

X_lin = ones(M, N_lin)

X_lin[1,:] = linspace(xmin, xmax, N_lin)

X_lin[2,:] = 1 # bias

# training data

N = 50 # data size

X = 2*rand(M, N) - 0.0 # input

X[2,:] = 1.0 # bias

Y = 0.5*sin.(2*pi * X[1,:]/3) + 0.05 * randn(N)

# model parameters

K = 5

sigma2_w = 10.0

sigma2_y = 0.01

################

# inference

alpha = 1.0e-5

max_iter = 100000

mu1, rho1, mu2, rho2 = BayesNeuralNet.VI(Y, X, sigma2_w, sigma2_y, K, alpha, max_iter)

Y_mean = [mu2'* tanh.(mu1'X_lin[:,n]) for n in 1 : N_lin]

################

# visualize

figure("result")

clf()

Y_list = []

num_samples = 100

for i in 1 : num_samples

Y_est, _ = BayesNeuralNet.sample_data_from_posterior(X_lin, mu1, rho1, mu2, rho2, sigma2_y, D)

push!(Y_list, Y_est)

plot(X_lin[1,:], Y_est, "-c", alpha=0.25)

end

plot(X[1,:], Y, "ok")

plot(X_lin[1,:], Y_mean, "b-")

xlim([xmin, xmax])

xlabel("x")

ylabel("y")

show()

end

Out: test

#sample_test()
test()

1 2	#sample_test() test()

いっぱいwarningが出力されるので不安になるが、using LinearAlgebraということなので、カーネル中断して、

using LinearAlgebra

1	using LinearAlgebra

で再度、

#sample_test()
test()

1 2	#sample_test() test()

Science To Medicine

Just My Daily Study Note by ts.anesth.kpum