ブログ

【Julia入門】JuliaのLinearAlgebraライブラリ① ベクトル・行列の演算

統計や機械学習に関するプログラミングではPythonRが用いられることが多いですが、近年Juliaも注目を集めています。そこで当シリーズではJuliaの基本構文からライブラリの用い方などについて取りまとめます。当記事ではLinearAlgebraライブラリを用いたベクトル・行列の演算について取りまとめを行いました。

Julia入門
https://www.hello-statisticians.com/julia

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

ベクトルの演算

下記のようにLinearAlgebraライブラリを読み込むことで内積を計算するdot関数やクロス積を計算するcross関数を実行することができます。

using LinearAlgebra

println(dot([1, 2, 3], [1, 1, 1]))
println(cross([0, 1, 0], [0, 0, 1]))

・実行結果

6
[1, 0, 0]

LinearAlgebraライブラリではベクトルのノルムを計算するnorm関数やベクトルの正規化を行うnormalize関数などがあり、下記のように実行することができます。

using LinearAlgebra

v = [1, 2, 3]

println(norm(v, 1))
println(norm(v, 2))
println(norm(v, Inf))

println(normalize(v, 1))
println(normalize(v, 2))

・実行結果

6.0
3.7416573867739413
3.0
[0.16666666666666666, 0.3333333333333333, 0.5]
[0.2672612419124244, 0.5345224838248488, 0.8017837257372732]

行列の演算

LinearAlgebraライブラリには行列のトレースを計算するtr関数、行列式を計算するdet関数、逆行列を計算するinv関数が用意されており、それぞれ下記のように実行することができます。

using LinearAlgebra

A = [1 2 3; 4 1 6; 7 8 1]
println(A)
println("===")

println(tr(A))
println(det(A))
println(inv(A))

・実行結果

[1 2 3; 4 1 6; 7 8 1]
===
3
104.0
[-0.4519230769230769 0.2115384615384615 0.08653846153846155; 0.3653846153846154 -0.1923076923076923 0.057692307692307675; 0.24038461538461536 0.057692307692307696 -0.0673076923076923]

参考

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

【Julia入門】Juliaのモジュール② 新しいモジュールの定義

統計や機械学習に関するプログラミングではPythonRが用いられることが多いですが、近年Juliaも注目を集めています。そこで当シリーズではJuliaの基本構文からライブラリの用い方などについて取りまとめます。当記事では新しいモジュールの定義について取りまとめを行いました。

Julia入門
https://www.hello-statisticians.com/julia

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

モジュールの定義

モジュールは下記のように定義することができます。

module Greeting1
hello(name) = println("Hello, $(name).")
end

Greeting1.hello("Julia")

・実行結果

Hello, Julia.

上記ではGreeting1モジュールにhello関数の定義を行いました。Greeting1.helloを実行することでhello関数を実行することができます。

モジュールの相対パス指定

前節ではモジュールの定義のみを行いましたが、Mainモジュールからの相対パスに基づいてモジュールを読み込むことができます。下記ではMainモジュールからの相対パス.Greeting2を用いることでGreeting2の読み込みを行いました。

module Greeting2
export hello
hello(name) = println("Hello, $(name).")
goodbye(name) = println("Goodbye, $(name).")
end

Greeting2.hello("Julia")
Greeting2.goodbye("Julia")

using .Greeting2
hello("Julia")
goodbye("Julia")

・実行結果

Hello, Julia.
Goodbye, Julia.
Hello, Julia.

UndefVarError: goodbye not defined

モジュールを読み込む際にモジュール定義の中でexportで指定した関数はモジュールを読み込むとそのまま使用することができます。上記の例ではexport helloがあるので、helloはそのまま実行できる一方でgoodbyeはそのまま実行できないことが確認できます。

参考

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

Ch.2 「確率分布と期待値」の章末問題の解答例 〜現代数理統計学の基礎(共立出版)〜

当記事は「現代数理統計学の基礎(共立出版)」の読解サポートにあたってChapter.$2$の「確率分布と期待値」の章末問題の解説について行います。
基本的には書籍の購入者向けの解説なので、まだ入手されていない方は購入の上ご確認ください。また、解説はあくまでサイト運営者が独自に作成したものであり、書籍の公式ページではないことにご注意ください。

・解答まとめ
https://www.hello-statisticians.com/answer_textbook_math_stat#basic

章末の演習問題について

問題2.1の解答例

$[1]$
$$
\large
\begin{align}
f(x) = Cx^3, \quad 0 < x < 2
\end{align}
$$

$f(x)$が確率密度関数であるとき、下記が成立する。
$$
\large
\begin{align}
\int_{0}^{2} f(x) dx = C \int_{0}^{2} x^3 dx &= 1 \\
\left[ \frac{C}{4}x^4 \right]_{0}^{2} &= 1 \\
\frac{C}{4} \cdot 2^{4} &= 1 \\
4 C &= 1 \\
C &= \frac{1}{4}
\end{align}
$$

また、累積分布関数$F(x)$は下記のように得られる。
$$
\large
\begin{align}
F(x) &= \int_{0}^{x} \frac{1}{4} u^3 du \\
&= \left[ \frac{1}{16}u^4 \right]_{0}^{x} \\
&= \frac{1}{16}x^4
\end{align}
$$

$[2]$
$$
\large
\begin{align}
f(x) = Ce^{-|x|}, \quad -\infty < x < \infty
\end{align}
$$

$f(x)$が確率密度関数であるとき、下記が成立する。
$$
\large
\begin{align}
\int_{-\infty}^{\infty} f(x) dx = C \int_{-\infty}^{\infty} e^{-|x|} dx &= 1 \\
2C \int_{0}^{\infty} e^{-x} dx &= 1 \\
2C \left[ -e^{-x} \right]_{0}^{\infty} &= 1 \\
2 C &= 1 \\
C &= \frac{1}{2}
\end{align}
$$

また、累積分布関数$F(x)$は下記のように得られる。
・$x \leq 0$のとき
$$
\large
\begin{align}
F(x) &= \int_{-\infty}^{x} \frac{1}{2} e^{u} du \\
&= \left[ \frac{1}{2}e^{u} \right]_{-\infty}^{x} \\
&= \frac{1}{2}e^{x}
\end{align}
$$

・$x > 0$のとき
$$
\large
\begin{align}
F(x) &= \int_{-\infty}^{0} \frac{1}{2} e^{u} du + \int_{0}^{x} \frac{1}{2} e^{-u} du \\
&= \frac{1}{2} + \left[ -\frac{1}{2}e^{-u} \right]_{0}^{x} \\
&= 1 – \frac{1}{2} e^{-x}
\end{align}
$$

問題2.2の解答例

問題2.3の解答例

問題2.4の解答例

問題2.5の解答例

問題2.6の解答例

問題2.7の解答例

問題2.8の解答例

問題2.9の解答例

問題2.10の解答例

問題2.11の解答例

問題2.12の解答例

問題2.13の解答例

問題2.14の解答例

問題2.15の解答例

問題2.16の解答例

問題2.17の解答例

まとめ

【Julia入門】Juliaのモジュール① モジュールの概要と読み込み

統計や機械学習に関するプログラミングではPythonRが用いられることが多いですが、近年Juliaも注目を集めています。そこで当シリーズではJuliaの基本構文からライブラリの用い方などについて取りまとめます。当記事ではモジュールの概要と読み込みについて取りまとめを行いました。

Julia入門
https://www.hello-statisticians.com/julia

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

モジュールの概要

概要

モジュール(module)は何らかの機能を実現するソースコードを取りまとめたものです。当記事では以下、Juliaに標準で用意されたモジュールの活用などについて確認を行います。

Mainモジュール

Juliaを実行する際にデフォルトで使用されるのがMainモジュールです。Mainモジュールでは下記のgcdのように基本的な関数が含まれます。

println(@__MODULE__)

println(gcd(27,72))

・実行結果

Main
9

上記のgcd関数は2つの引数の最大公約数を計算する関数であり、上記では2772の最大公約数である9が出力されたことが確認できます。

using文と既存モジュールの利用

関数を指定した読み込み

using文は下記のように『using モジュール名: 関数名』を実行することでモジュールの読み込みを行います。

using Statistics: mean, std

x = [1, 2, 3, 4, 5]

println(mean(x))
println(std(x))
println(var(x))

・実行結果

3.0
1.5811388300841898

UndefVarError: var not defined

Statisticsモジュールには分散を計算するvar関数が実装されている一方で、上記ではvar関数の読み込みを行わなかったのでUndefVarErrorが出力されました。このように関数名を指定する場合は指定した関数名のみを読み込むことができます。

モジュール全体の読み込みと注意点

モジュール全体を読み込む場合は下記のように『using モジュール名』を実行すれば良いです。

using Statistics

x = [1, 2, 3, 4, 5]

println(mean(x))
println(std(x))
println(var(x))

・実行結果

3.0
1.5811388300841898
2.5

上記のようにusing文を用いてモジュールの読み込みを行う場合、意図しない関数を読みこむ場合が生じやすいです。たとえば上記のStatisticsモジュールのvarなどは他で名称を使用する可能性があり得ます。

このような意図しない関数名の読み込みを行わないようにするにあたっては、前項のusing Statistics: mean, stdように読み込む関数名を指定するか、下記のようにStatisticsモジュールを読み込むなどの方法があります。

using Statistics: Statistics

x = [1, 2, 3, 4, 5]

println(Statistics.mean(x))
println(Statistics.std(x))
println(Statistics.var(x))

・実行結果

3.0
1.5811388300841898
2.5

参考

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

『直感的に理解するTransformer』サポートページ【印刷版】

技術書典$14$で頒布を行った『直感的に理解するTransformer』の印刷版のサポートページです。主に誤植が見つかった場合の正誤表の作成やカラー画像の確認が行えるように作成を行いました。誤植につきましては見つかり次第都度追加いたしますので、お気づきの方は気軽にご指摘ください。

https://hello-stats.booth.pm/items/4781879

追加コンテンツ

Multi Head Attention

直感的に理解するTransformer」の$A.3$の「Multi Head Attention」の処理は正確には論文の内容と一致しないことにご注意ください。具体的には「単に分割」で表したところが「パラメータ行列を用いて線形写像での計算」に対応します。以下に論文の数式などをまとめました。

Transformer論文では下記のようにMulti Head Attentionが表されます。
$$
\large
\begin{align}
\mathrm{MultiHead} (Q,K,V) &= \mathrm{Concat}(\mathrm{head}_{1}, \cdots , \mathrm{head}_{h}) W^{O} \quad (1.1) \\
\mathrm{head}_{i} &= \mathrm{Attention}(QW_{i}^{Q}, KW_{i}^{K}, VW_{i}^{V}) \quad (1.2) \\
\mathrm{Attention}(Q, K, V) &= \mathrm{Softmax} \left( \frac{QK^{\mathrm{T}}}{\sqrt{d_{k}}} \right) \quad (1.3)
\end{align}
$$

Multi Head Attention処理の理解にあたって着目すべきは$(1.2)$式であり、$Q, K, V$にそれぞれ$W_{i}^{Q} \in \mathbb{R}^{d_{model} \times d_{k}}, W_{i}^{K} \in \mathbb{R}^{d_{model} \times d_{k}}, W_{i}^{V} \in \mathbb{R}^{d_{model} \times d_{v}}$をかけることでそれぞれの単語に対応する内部表現を$Q, K$については$d_k$次元、$V$については$d_v$次元にそれぞれ変換した上でDot Product Attention処理を行います。その後に$(1.1)$式で表されるように$h$個のヘッドを連結し、$W_{i}^{O} \in \mathbb{R}^{h d_v \times d_{model}}$を用いて再度パラメータ処理を行います。

Transformer論文ではヘッドの数を表す$h$を$h=8$、各単語の内部表現を$d_{model}=512$で表すのがデフォルトです。また、基本的には$h d_{k} = d_{model}, h d_{v} = d_{model}$であるので$d_{k}=d_{v}=64$がデフォルトになります。詳細の設定は論文のTable.$3$を参照すると良いです。

Transformer論文Table.$3$

上図の「base」がデフォルトのTransformerの構成、「big」がより大きなTransformerの構成にそれぞれ対応します。関連で総パラメータ数の概算についても下記で取り扱ったので合わせてご確認ください。

正誤表

初版第1刷

ページ・行数$\times$
P.57 l.19$P(\mathbf{X})=P(x_0,\cdots,x_7)$の確率が大変小さくなることで$P(\mathbf{X})=P(x_0,\cdots,x_7)$の確率が大変小さくなり
P.75 l.5単語のベクトルを分割してこれまでの処理を「並列で行う」と同義単語 のベクトルを分割してこれまでの処理を「並列で行う」ことに類似
P.83 l.4 教習ラベル教師ラベル

カラー画像の確認

第2章

図$2.1 \,$ 三角関数のグラフ
図$2.3 \,$ 行列の積とMLPの関係性

第3章

図$3.1 \,$ 東京メトロホームページより
図$3.4 \,$ $n=5, d=5$でノードに特徴量を割り当てた場合
図$3.8 \,$ A Comprehensive Survey on Graph Neural Networks. Fig.$2$b

第4章

図$4.1 \,$ Bag of Words
形態素解析の実行例

第5章

図$5.4 \,$ Transformer論文Fig.$1$を改変
図$5.5 \,$ TransformerとMPNN型GNNの各層における処理の概要

Appendix

図A.$7 \,$ 実行結果①
図A.$9 \,$ 実行結果③

Ch.1 「必要な数学知識」の例題・演習問題の解答例 〜言語処理のための機械学習入門〜

当記事は「言語処理のための機械学習入門(コロナ社)」の読解サポートにあたってChapter.$1$の「必要な数学知識」の章末問題の解説について行います。
基本的には書籍の購入者向けの解説なので、まだ入手されていない方は購入の上ご確認ください。また、解説はあくまでサイト運営者が独自に作成したものであり、書籍の公式ページではないことにご注意ください。

例題解答

例題1.1の解答例

$$
\large
\begin{align}
\mathrm{Objective} : \quad & f(\mathbf{x}) = -x_1x_2 \, \rightarrow \, \mathrm{Maximize} \\
\mathrm{Constraint} : \quad & c_1(\mathbf{x}) = x_1-x_2-a = 0
\end{align}
$$

$x_1-x_2-a = 0$より$x_1=x_2+a$を$f(\mathbf{x}) = -x_1x_2$に代入すると下記が得られる。
$$
\large
\begin{align}
f(\mathbf{x}) &= -x_1x_2 \\
&= -x_2(x_2+a) \\
&= – \left( x_2 + \frac{1}{2}a \right)^{2} + \frac{1}{4}a^2
\end{align}
$$

よって$\displaystyle x_1 = \frac{1}{2}a, x_2=-\frac{1}{2}a$のとき$f(\mathbf{x})$は最大値を取る。

例題1.2の解答例

例題1.3の解答例

$$
\large
\begin{align}
f(x) = -x^2
\end{align}
$$

上記のように定義する$f(x)$に関して$0 \leq t \leq 1$である$t$と、$x_1<x_2$である$x_1,x_2$を用いて下記のような式変形が成立する。
$$
\large
\begin{align}
& f(tx_1 + (1-t)x_2) – \left[ tf(x_1) + (1-t)f(x_2) \right] \\
&= -(tx_1 + (1-t)x_2)^2 – \left[ -tx_1^2 – (1-t)x_2^2 \right] \\
&= -t^2x_1^2 – (1-t)^2x_2^2 – 2t(1-t)x_1x_2 + tx_1^2 + (1-t)x_2^2 \\
&= t(1-t)x_1^2 + (1-t)(1-(1-t))x_2^2 – 2t(1-t)x_1x_2 \\
&= t(1-t)x_1^2 + t(1-t)x_2^2 – 2t(1-t)x_1x_2 \\
&= t(1-t)(x_1^2 + x_2^2 – 2x_1x_2) \\
&= t(1-t)(x_1-x_2)^2 \geq 0
\end{align}
$$

上記より、$f(tx_1 + (1-t)x_2) \geq tf(x_1) + (1-t)f(x_2)$が成立するので$f(x)=-x^2$は上に凸の関数である。

例題1.4の解答例

例題1.5の解答例

例題1.6の解答例

例題1.7の解答例

$$
\large
\begin{align}
\mathrm{Objective} : \quad & f(\mathbf{x}) = -(x_1^2+x_2^2) \, \rightarrow \, \mathrm{Maximize} \\
\mathrm{Constraint} : \quad & c_1(\mathbf{x}) = x_1+x_2-1 = 0
\end{align}
$$

上記の制約付き最適化問題のラグランジュ関数を$L(x_1,x_2,\lambda)$とおくと、$L(x_1,x_2,\lambda)$は下記のように表せる。
$$
\large
\begin{align}
L(x_1,x_2,\lambda) &= f(\mathbf{x}) + \lambda c_1(\mathbf{x}) \\
&= -(x_1^2+x_2^2) + \lambda(x_1+x_2-1)
\end{align}
$$

ここでKKT条件より下記が得られる。
$$
\large
\begin{align}
\frac{\partial L}{\partial x_1} &= -2x_1+\lambda = 0 \quad (1) \\
\frac{\partial L}{\partial x_2} &= -2x_2+\lambda = 0 \quad (2) \\
\frac{\partial L}{\partial \lambda} &= x_1+x_2-1 = 0 \quad (3)
\end{align}
$$

$(1),(2)$式より$x_1=x_2$が得られるので、$(3)$式に代入して最適解$\displaystyle (x_1,x_2) = \left( \frac{1}{2},\frac{1}{2} \right)$が得られる。

・考察
目的関数の等高線や制約条件の直線、最適解における勾配ベクトル$\nabla f(\mathbf{x})$と$\lambda \nabla c_1(\mathbf{x})$の釣り合いは下記を実行することでそれぞれ確認することができる。

import numpy as np
import matplotlib.pyplot as plt

x_, y_ = np.arange(-2., 2.01, 0.01), np.arange(-1., 3.01, 0.01)
x, y = np.meshgrid(x_,y_)
z = -(x**2 + y**2)

levs = np.arange(-3.5,-0., 0.5)

plt.scatter(1/2., 1/2., marker="*", s=300, color="blue")
plt.quiver(1/2., 1/2., 1, 1, color="red", angles='xy', scale_units='xy', scale=1)
plt.quiver(1/2., 1/2., -1, -1, color="green", angles='xy', scale_units='xy', scale=1)

plt.plot(x_, -x_+1, "r--")
plt.contour(x,y,z,levels=levs,cmap="viridis")

plt.colorbar()
plt.show()

実行結果

例題1.8の解答例

例題1.9の解答例

例題1.10の解答例

$X^2$の期待値を$E[X^2]$とおくと、$E[X]$は下記のように計算できる。
$$
\large
\begin{align}
E[X^2] &= \sum_{k=1}^{6} k^2 P(X=k) \\
&= \frac{1}{6} \sum_{k=1}^{6} k^2 \\
&= \frac{1}{6} \times \frac{1}{6} \cdot 6 \cdot (6+1) \cdot (12+1) \\
&= \frac{91}{6} \\
& \simeq 15.17
\end{align}
$$

例題1.11の解答例

例題1.12の解答例

標本平均を$\bar{x}$とおくと、標本平均の定義に基づいて下記のように計算を行うことができる。
$$
\large
\begin{align}
\bar{x} &= \frac{1}{8} \sum_{i=1}^{8} x_i \\
&= \frac{1}{8} (3+6+5+5+2+4+5+2) \\
&= 4
\end{align}
$$

例題1.13の解答例

例題1.14の解答例

条件付き確率の定義に基づいて下記が成立する。
$$
\large
\begin{align}
P(x,y) &= P(y|x)P(x) \quad (1) \\
P(x,y) &= P(x|y)P(y) \quad (2)
\end{align}
$$

$(1),(2)$式より下記が成立する。
$$
\large
\begin{align}
P(y|x)P(x) &= P(x|y)P(y) \\
P(y|x) &= \frac{P(x|y)P(y)}{P(x)}
\end{align}
$$

例題1.15の解答例

$P(x_1,x_2,x_3)$は条件付き確率の定義式を元に下記のように表すことができる。
$$
\large
\begin{align}
P(x_1,x_2,x_3) &= P(x_1|x_2,x_3)P(x_2,x_3) \\
&= P(x_1|x_2,x_3)P(x_2|x_3)P(x_3) \quad (1)
\end{align}
$$

また、下記のようにも表すことができる。
$$
\large
\begin{align}
P(x_1,x_2,x_3) &= P(x_2|x_1,x_3)P(x_1,x_3) \\
&= P(x_2|x_1,x_3)P(x_1|x_3)P(x_3) \quad (2)
\end{align}
$$

$(1), (2)$式より下記が成立する。
$$
\large
\begin{align}
P(x_2|x_1,x_3)P(x_1|x_3) \cancel{P(x_3)} &= P(x_1|x_2,x_3)P(x_2|x_3) \cancel{P(x_3)} \\
P(x_2|x_1,x_3)P(x_1|x_3) &= P(x_1|x_2,x_3)P(x_2|x_3) \\
P(x_2|x_1,x_3) &= \frac{P(x_2|x_3)P(x_1|x_2,x_3)}{P(x_1|x_3)}
\end{align}
$$

例題1.16の解答例

$E[X]$は下記のように計算できる。
$$
\large
\begin{align}
E[X] &= 1 \cdot p + 0 \cdot (1-p) \\
&= p
\end{align}
$$

例題1.17の解答例

$X \sim \mathrm{Bin}(5,0.8)$である確率変数$X$に対し、下記のように$P(X=3)$を計算すれば良い。
$$
\large
\begin{align}
P(X=3) &= {}_{5} C_{3} \cdot 0.8^{3} \cdot 0.2^{2} \\
&= 10 \cdot 0.8^{3} \cdot 0.2^{2} \\
&= 0.2048
\end{align}
$$

例題1.18の解答例

観測値を$\mathbf{x}=(x_1,\cdots,x_6)$、$P$氏の確率を$P(\mathbf{x}|P)$、$N$氏の確率を$P(\mathbf{x}|N)$とおく。このとき$P(\mathbf{x}|P)$と$P(\mathbf{x}|N)$はそれぞれ下記のように計算できる。
$$
\large
\begin{align}
P(\mathbf{x}|P) &= \frac{6!}{2!1!1!2!} 0.4^{2} 0.1^{1} 0.3^{1} 0.2^{2} = 0.03456 \\
P(\mathbf{x}|N) &= \frac{6!}{2!1!1!2!} 0.1^{2} 0.4^{1} 0.2^{1} 0.3^{2} = 0.01296
\end{align}
$$

例題1.19の解答例

パラメータ$\mu$に関する尤度を$L(\mu)$とおくと、$L(\mu)$は下記のように得られる。
$$
\large
\begin{align}
L(\mu) &= P(D|\mu) \\
&= \prod_{i=1}^{N} \frac{\mu^{x^{(i)} e^{-\mu}}}{x^{(i)}!}
\end{align}
$$

上記に対し、対数尤度$\log{L(\mu)}$は下記のように得られる。
$$
\large
\begin{align}
\log{L(\mu)} &= \log{\left[ \prod_{i=1}^{N} \frac{\mu^{x^{(i)}} e^{-\mu}}{x^{(i)}!} \right]} \\
&= \sum_{i=1}^{N} \log{\left[ \frac{\mu^{x^{(i)}} e^{-\mu}}{x^{(i)}!} \right]} \\
&= \sum_{i=1}^{N} (x^{(i)}\log{\mu} – \mu -\log{x^{(i)}!})
\end{align}
$$

上記の$\mu$に関する微分は標本平均$\bar{x}$を用いて下記のように表せる。
$$
\large
\begin{align}
\frac{\partial \log{L(\mu)}}{\partial \mu} &= \sum_{i=1}^{N} \left[ \frac{x^{(i)}}{\mu} – 1 \right] \\
&= \frac{N \bar{x}}{\mu} – N \\
&= \frac{N(\bar{x}-\mu)}{\mu}
\end{align}
$$

上記より、対数尤度$\log{L(\mu)}$と尤度$L(\mu)$を最大にする$\mu$は$\mu=\bar{x}$であることが得られる。

例題1.20の解答例

$N=n_{\mathrm{good}}+n_{\mathrm{bad}}+n_{\mathrm{exciting}}+n_{\mathrm{boring}}=10$であるので、それぞれ下記が成立する。
$$
\large
\begin{align}
p_{\mathrm{good}} &= \frac{n_{\mathrm{good}}}{N} \\
&= \frac{5}{10} = 0.5 \\
p_{\mathrm{bad}} &= \frac{n_{\mathrm{bad}}}{N} \\
&= \frac{1}{10} = 0.1 \\
p_{\mathrm{exciting}} &= \frac{n_{\mathrm{exciting}}}{N} \\
&= \frac{4}{10} = 0.4 \\
p_{\mathrm{boring}} &= \frac{n_{\mathrm{boring}}}{N} \\
&= \frac{0}{10} = 0
\end{align}
$$

例題1.21の解答例

“Yes”の発言回数を確率変数$X$で表すと、$X \sim \mathrm{Bin}(5,p)$が成立する。このとき尤度$L(p)=P(X_1=2,X_2=3,X_3=3)$は下記のように表せる。
$$
\large
\begin{align}
L(p) &= P(X_1=2, X_2=3, X_3=3) \\
&= P(X_1=2)P(X_2=3)P(X_3=3) \\
&= {}_{5} C_{2} ({}_{5} C_{3})^{2} p^{8} (1-p)^{7}
\end{align}
$$

上記に対し、対数尤度$\log{L(p)}$は下記のように得られる。
$$
\large
\begin{align}
\log{L(p)} = 8 \log{p} + 7 \log{(1-p)}
\end{align}
$$

上記の$p$に関する微分は下記のように得られる。
$$
\large
\begin{align}
\frac{\partial \log{L(p)}}{\partial p} &= \frac{8}{p} – \frac{7}{1-p} \\
&= \frac{8(1-p)-7p}{p(1-p)} \\
&= \frac{8-15p}{p(1-p)}
\end{align}
$$

上記より最尤推定解は$\displaystyle p=\frac{8}{15}$である。

例題1.22の解答例

例題1.23の解答例

$P_{\mathrm{Bin}}(X=0), P_{\mathrm{Bin}}(X=1), P_{\mathrm{Bin}}(X=2), P_{\mathrm{Bin}}(X=3)$はそれぞれ下記のように表せる。
$$
\large
\begin{align}
P_{\mathrm{Bin}}(X=0) &= {}_{3} C_{0} \times 0.8^0 \times 0.2^3 = 0.2^3 \\
P_{\mathrm{Bin}}(X=1) &= {}_{3} C_{1} \times 0.8^1 \times 0.2^2 = 3 \times 0.8 \times 0.2^2 \\
P_{\mathrm{Bin}}(X=2) &= {}_{3} C_{2} \times 0.8^2 \times 0.2^1 = 3 \times 0.8^2 \times 0.2 \\
P_{\mathrm{Bin}}(X=3) &= {}_{3} C_{3} \times 0.8^3 \times 0.2^0 = 0.8^3
\end{align}
$$

上記を元にエントロピーの定義より$H(P_{\mathrm{Bin}})$は下記のように計算できる。
$$
\large
\begin{align}
H(P_{\mathrm{Bin}}) &= -\sum_{k=0}^{3} P_{\mathrm{Bin}}(X=k) \log{P_{\mathrm{Bin}}(X=k)} \\
&= -0.2^3 \log{0.2^3} – (3 \times 0.8 \times 0.2^2) \log{(3 \times 0.8 \times 0.2^2)} – (3 \times 0.8^2 \times 0.2) \log{(3 \times 0.8^2 \times 0.2)} – 0.8^3 \log{0.8^3} \\
&=
\end{align}
$$

例題1.24の解答例

例題1.25の解答例

演習問題解答

問題1.1の解答例

問題1.2の解答例

問題1.3の解答例

問題1.4の解答例

問題1.5の解答例

問題1.6の解答例

問題1.7の解答例

問題1.8の解答例

期待値を$E[X]$とおくと$E[X]$は下記のように計算できる。
$$
\large
\begin{align}
E[X] &= \frac{1}{6}(1+2+3+4+5+6) \\
&= \frac{21}{6} \\
&= \frac{7}{2}
\end{align}
$$

問題1.9の解答例

$P(x_1,x_2,x_3,x_4)$は条件付き確率の定義に基づいて下記のように表せる。
$$
\large
\begin{align}
P(x_1,x_2,x_3,x_4) = P(x_1,x_2|x_3,x_4)P(x_3,x_4) \quad (1)
\end{align}
$$

同様に条件付き確率の定義を用いることで下記のように表すこともできる。
$$
\large
\begin{align}
P(x_1,x_2,x_3,x_4) &= P(x_1|x_2,x_3,x_4)P(x_2,x_3,x_4) \\
&= P(x_1|x_2,x_3,x_4)P(x_2|x_3,x_4)P(x_3,x_4) \quad (2)
\end{align}
$$

$(1), (2)$式より下記が成立する。
$$
\large
\begin{align}
P(x_1,x_2|x_3,x_4) \cancel{P(x_3,x_4)} &= P(x_1|x_2,x_3,x_4) P(x_2|x_3,x_4) \cancel{P(x_3,x_4)} \\
P(x_1,x_2|x_3,x_4) &= P(x_1|x_2,x_3,x_4) P(x_2|x_3,x_4)
\end{align}
$$

上記より与式が示される。

問題1.10の解答例

$X \sim \mathrm{Bin}(5,0.8)$である確率変数$X$に対し、下記のように$P(X=3)+P(X=4)+P(X=5)$を計算すれば良い。
$$
\begin{align}
P(X=3) + P(X=4) + P(X=5) &= {}_{5} C_{3} \cdot 0.8^{3} \cdot 0.2^{2} + {}_{5} C_{4} \cdot 0.8^{4} \cdot 0.2^{1} + {}_{5} C_{5} \cdot 0.8^{5} \cdot 0.2^{0} \\
&= 10 \cdot 0.8^{3} \cdot 0.2^{2} + 5 \cdot 0.8^{4} \cdot 0.2^{1} + 0.8^{5} \cdot 0.2^{0} \\
&= 0.942 \cdots
\end{align}
$$

問題1.11の解答例

問題1.12の解答例

問題1.13の解答例

問題1.14の解答例

問題1.15の解答例

問題1.16の解答例

問題1.17の解答例

まとめ

【Julia入門】Juliaの多次元配列⑤ サブ配列(subarray)

統計や機械学習に関するプログラミングではPythonRが用いられることが多いですが、近年Juliaも注目を集めています。そこで当シリーズではJuliaの基本構文からライブラリの用い方などについて取りまとめます。当記事ではJuliaのサブ配列(subarray)について取りまとめを行いました。

Julia入門
https://www.hello-statisticians.com/julia

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

サブ配列の概要

概要

サブ配列(subarray)は配列の一部を表すオブジェクトです。view関数を用いることで配列からサブ配列を作成することができます。

サブ配列の使用例

サブ配列はview関数を用いて下記のように作成することができます。

using Random
Random.seed!(0)

X = rand(3,3)

println(X)
println(view(X, 1:2, 1))
println(view(X, 1:2, 1:2))

・実行結果

[0.4552384158732863 0.9405848223512736 0.7468008914093891; 0.5476424498276177 0.02964765308691042 0.9766699015845924; 0.7733535276924052 0.74694291453392 0.08694684883050086]
[0.4552384158732863, 0.5476424498276177]
[0.4552384158732863 0.9405848223512736; 0.5476424498276177 0.02964765308691042]

参考

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

【Julia入門】Juliaの多次元配列④ map、reduce、filterを用いた操作

統計や機械学習に関するプログラミングではPythonRが用いられることが多いですが、近年Juliaも注目を集めています。そこで当シリーズではJuliaの基本構文からライブラリの用い方などについて取りまとめます。当記事ではJuliamapreducefilterを用いた多次元配列の操作について取りまとめを行いました。

Julia入門
https://www.hello-statisticians.com/julia

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

map・reduce・filterの概要

配列の各要素に関数を使用した値の変換や集約、フィルタリングにあたってはそれぞれmapreducefilterが有用です。

使用例

map

mapは配列の各要素に関数を適用して値を変換する際などに用いられる関数です。下記のように実行することができます。

using Random
Random.seed!(1)

X = randn(2,2)

println(X)
println(map(x -> x+2.0, X))
println(map(x-> x*2.0-3.0, X))

・実行結果

[-0.07058313895389791 -0.806852326006714; 0.5314767537831963 2.456991333983293]
[1.929416861046102 1.193147673993286; 2.5314767537831964 4.456991333983293]
[-3.141166277907796 -4.613704652013428; -1.9370464924336075 1.913982667966586]

reduce

reduceは配列の要約などに用いられる関数です。下記のように実行することができます。

using Random
Random.seed!(1)

X = randn(5,5)

println(X)
println(reduce(+,X))
println(reduce(+,X)/25.0)

・実行結果

[0.06193274031408013 -1.5765649225859841 0.21787878613277867 -1.5242735305753605 -0.5841980481085709; 0.2784058141640002 0.1759399913010747 -0.6559505611509957 -0.7829327668434358 -0.29942199646186524; -0.5958244153640522 0.8653808054093252 0.26751404691066 -0.5370139742021761 -2.3830630673003137; 0.04665938957338174 -2.790281005549307 0.0073860677115008865 -2.0074061601841917 -1.0026135242092884; 1.0857940215432762 -1.8920155582259128 1.0607244834908975 0.16218731766404215 0.972024394360624]
-11.429731672185813
-0.4571892668874325

filter

filterはフィルタリングなどに用いられる関数です。下記のように実行することができます。

using Random
Random.seed!(0)

X = randn(5,5)

println(X)
println(filter(x -> x > 1, X))

・実行結果

[-0.23190906957695134 1.1722448627021573 -0.26988492037793593 0.13457571047366457 -1.0745769949749375; 0.940390210248594 -1.696813019281842 0.5780036022342506 0.16554255091922548 -1.3819748518823003; 0.5967616711335215 -2.1161459685604527 1.1603442420756227 -1.0720974154467635 0.9518527642027133; 1.9978240930577937 0.5583656057835097 0.28788799558526296 -1.2601642953877346 -0.2922912657969977; -0.05156184220322609 -0.8647295755028067 -0.44119282303422414 0.196408003082906 0.23759838081412232]
[1.9978240930577937, 1.1722448627021573, 1.1603442420756227]

参考

Julia 1.8 Documentation
https://docs.julialang.org/en/v1/

Python・matplotlib-vennを用いたベン図(Venn diagrams)の描画

集合と要素や集合に基づいて必要条件・十分条件・必要十分条件などを取り扱う際に用いると良いのがベン図(Venn diagrams)です。当記事ではPythonを用いてベン図を描画するにあたって、matplotlib-vennの用法などの取りまとめを行いました。

・プログラミングまとめ
https://www.hello-statisticians.com/program

matplotlib-vennライブラリ

matplotlib-vennの概要

matplotlib-vennmatplotlibを用いて事象が2もしくは3のベン図の描画を行うにあたっての関数を実装したライブラリです。

matplotlib-vennのインストール

下記を実行することでmatplotlib-vennのインストールを行うことが可能です。

> pip install matplotlib-venn

matplotlib-vennの主な関数

matplotlib-vennの主な関数はvenn2venn2_circlesvenn3venn3_circles4つでそれぞれ下記のような描画を行うことができます。

関数 概要
venn22つの事象に基づくベン図を作成する
venn2_circles2つの事象に基づく色・ラベルのないベン図を作成する
venn33つの事象に基づくベン図を作成する
venn3_circles3つの事象に基づく色・ラベルのないベン図を作成する

ベン図の描画

基本的な実行例

PyPIのドキュメントでは下記のような実行例が紹介されています。

from matplotlib import pyplot as plt
from matplotlib_venn import venn2, venn2_circles, venn3, venn3_circles

figure, axes = plt.subplots(2, 2)
venn2(subsets={'10': 1, '01': 1, '11': 1}, set_labels = ('A', 'B'), ax=axes[0][0])
venn2_circles((1, 2, 3), ax=axes[0][1])
venn3(subsets=(1, 1, 1, 1, 1, 1, 1), set_labels = ('A', 'B', 'C'), ax=axes[1][0])
venn3_circles({'001': 10, '100': 20, '010': 21, '110': 13, '011': 14}, ax=axes[1][1])

plt.show()

・実行結果

上記より、venn2venn2_circlesvenn3venn3_circlesの主な4つの関数の出力結果が確認できます。

以下では、venn2関数について確認します。

from matplotlib import pyplot as plt
from matplotlib_venn import venn2

venn2(subsets = (3, 2, 1))

plt.show()

・実行結果

上図では引数のsubsetsでそれぞれの領域の大きさを指定しましたが、下記のようにset関数を元に集合を引数に与えることも可能です。

from matplotlib import pyplot as plt
from matplotlib_venn import venn2

venn2([set(['A', 'B', 'C', 'D']), set(['D', 'E', 'F'])])

plt.show()

・実行結果

和集合の描画

和集合は下記のように描画を行うことができます。

from matplotlib import pyplot as plt
from matplotlib_venn import venn2

v = venn2(subsets = (1, 1, 1))

v.get_patch_by_id('10').set_color('limegreen')
v.get_patch_by_id('11').set_color('limegreen')
v.get_patch_by_id('01').set_color('limegreen')

v.get_patch_by_id('01').set_alpha(1.0)
v.get_patch_by_id('11').set_alpha(1.0)
v.get_patch_by_id('10').set_alpha(1.0)

v.get_patch_by_id('11').set_edgecolor('white')

v.get_label_by_id('10').set_text('')
v.get_label_by_id('11').set_text('')
v.get_label_by_id('01').set_text('')

plt.title("Union of two sets")
plt.show()

・実行結果

積集合の描画

積集合は下記のように描画を行うことができます。

from matplotlib import pyplot as plt
from matplotlib_venn import venn2

v = venn2(subsets = (1, 1, 1))

v.get_patch_by_id('10').set_color('limegreen')
v.get_patch_by_id('11').set_color('limegreen')
v.get_patch_by_id('01').set_color('limegreen')

v.get_patch_by_id('10').set_alpha(0.1)
v.get_patch_by_id('11').set_alpha(1.0)
v.get_patch_by_id('01').set_alpha(0.1)

v.get_label_by_id('10').set_text('')
v.get_label_by_id('11').set_text('')
v.get_label_by_id('01').set_text('')

plt.show()

・実行結果

参考

https://pypi.org/project/matplotlib-venn/
https://github.com/konstantint/matplotlib-venn

統計検定準1級 問題解説 ~2016年6月実施 選択問題及び部分記述問題 問6~

過去問

過去問題は統計検定公式が問題と解答例を公開しています。こちらを参照してください。

解説

[1] 解答

$\boxed{\mathsf{6}}$ : $①$

データの個数を $n$ ,説明変数の個数を $p$ とすると,自由度調整済み決定係数 $\bar{R^2}$ は,
$\bar{R^2} = 1-\dfrac{n-1}{n-p-1}(1-R^2)$ と表されるので, $\bar{R^2} = 1-\dfrac{21-1}{21-1-1}(1-0.0372) = -0.01347$
となる.

[2] 解答

$\boxed{\mathsf{7}}$ : $⑤$

ダービン・ワトソン統計量 $DW$ と自己相関係数 $\hat{\rho}$ には,$DW \approx 2(1-\hat{\rho})$ の関係がある.いま,$DW = 2.98$ より,$\rho \approx -0.49$ .
このことから,強い負の相関があることがわかる.(強い)負の相関があるとき残差は上下に(激しく)動く.(なぜなら,$1$ つ前が大きいときにはその次のときは小さくなり,逆に1つ前が小さいときは次のときに大きくなる)

したがって,⑤である.