MIT OCW, Discrete Stochastic Processes 02日目

Robert Gallager. 6.262 Discrete Stochastic Processes. Spring 2011. Massachusetts Institute of Technology: MIT OpenCourseWare, https://ocw.mit.edu. License: Creative Commons BY-NC-SA.

Lecture videoを要約していく.

• • Lecture 3. Laws of Large Numbers, Convergence
    • Markov, Chebychev, Chernoff bounds
      • Markov inequality
      • Chebyshev inequality
      • Chernoff bound
    • Convergence
      • Weak Law of Large Numbers
        • Definition
    • Central Limit Theorem

Lecture 3. Laws of Large Numbers, Convergence

Markov, Chebychev, Chernoff bounds

Markov inequality

$Y\geq 0, y > 0$なら
$$Pr(Y \geq y) \leq \frac{E[Y]}{y}$$

Chebyshev inequality

$\epsilon > 0$に
$$Pr(|Z-E[Z]|\geq \epsilon) \leq \frac{\sigma_Z^2}{\epsilon^2} = \frac{var[Z]}{\epsilon^2}$$

Chernoff bound

任意の$z>0, r>0$にmoment generationg function $g_Z(r) = E[e^{rZ}]$が定義されているなら,
$$Pr(Z\geq z) \leq g_Z(r)\exp(-rz)$$

proof.

$Y=e^{rZ}$とする. $E[Y]=g_Z(r)$で,$y>0$にMarkov inequalityを使って
$Pr(Y \geq y) \leq g_Z(r)/y$すなわち$Pr(e^{rZ}\geq y) \leq g_Z(r)/y$
$z=(\log y)/r \Leftrightarrow y = e^{rz}$と$z$を導入すると,$Pr(Z\geq z) \leq g_Z(r) \exp(-rz)$がたしかに成立.

Markov, Chebyshev, Chernoffを見比べると,Markovは$y^{-1}$, Chebyshevは$y^{-2}$に従って上限が小さくなるのに対し,Chernoffの上限は指数的に上限が小さくなる.これがChernoff boundの有用性の理由である.

Convergence

Weak Law of Large Numbers

$X_1, ..., X_n$はIIDで,$\bar{X}=E[X_1], \sigma^2 = \sigma_{X_1}^2$とする.
$S_n = X_1 + ...+X_n$とすると$\sigma_{S_n}^2 = n\sigma^2$である.
$S_n/n$の振る舞いを見る.
$var[S_n/n] = n\sigma^2/n^2=\sigma^2/n \rightarrow 0 ( n \rightarrow \infty)$だから
$\lim_{n\rightarrow \infty} E \left[ \left(\frac{S_n}{n}-X\right)^2 \right]=0$. これを$S_n/n$ converges in mean square to $\bar{X}$という. 前者はIIDでないときには成り立たないし,分散が存在しないときにも成り立たないが,このようなときでも実は後者は成立する.

Definition

$Y_1,...,Y_n$が$Y$にconverges in mean square to $Y$
$\Leftrightarrow \lim E[(Y_n-Y)^2]=0$

ここでChebishev’s inequalityを使って
$$Pr(|\frac{S_n}{n} - \bar{X}|\geq \epsilon) \leq \frac{\sigma^2}{n\epsilon^2} \rightarrow 0$$
以上より,weak law of large numbers(WLLN, 大数の弱法則)
$$\forall \epsilon >0. \ \lim_{n\rightarrow \infty} Pr(|\frac{S_n}{n}-\bar{X}|\geq \epsilon)) =0$$
が示せた. (IIDの仮定に注意,varianceは存在しなくても良い)これは$Pr(S_n/n \leq x)$のdistributionが$n$の極限で$x\leq \bar{X}$で0, $x \geq \bar{X}$で$1$をとるステップ関数になるということである

Central Limit Theorem

$$\lim_{n\rightarrow \infty} \left[Pr(\frac{S_n-n\bar{X}}{\sqrt{n}\sigma} \leq y ) \right] = \int^y_{-\infty} \frac{1}{\sqrt{2\pi}}\exp(-\frac{x^2}{2})dx$$

MIT OCW, Discrete Stochastic Processes 01日目

Robert Gallager. 6.262 Discrete Stochastic Processes. Spring 2011. Massachusetts Institute of Technology: MIT OpenCourseWare, https://ocw.mit.edu. License: Creative Commons BY-NC-SA.

Lecture videoを要約していく.

• • Lecture 1
    • モデルを作るときに現れる重要な問題
      • 1. – 完全なモデルは存在しない
      • 2. – その数学的モデルの解が現実で意味を持つか
    • Stochastic Processとは？
      • このコースで学ぶprocessたち
    • Kormogorovの公理
    • Eventの独立性
    • Random Variable
  • Lecture 2. 確率論の復習
    • Expectations
      • example: indicator random variable
    • multiple random variables
      • IID random variables
      • Sample Average
    • The Bernoulli Process
  • Assignment 1. Problem set 1.
    • • Exercise 1.3
      • Exercise 1.9
      • Exercise 1.13
      • Exercise 1.20
      • Exercise 1.26

Lecture 1

well-posed problemを解くのは簡単だが,現実にある現象をモデル化してwell-posed problemに落とし込むのは難しい. このコースでは現実世界での確率と確率の理論を学んだ後discrete processがなんであるかを学び,数あるdiscrete processの内いくつかを学ぶ.
確率論がどこで役に立つか–どこでも役に立つのだが,いくつか例を挙げる
Kormogrovの確率の公理がどのように役に立っているか
確率論の復習

モデルを作るときに現れる重要な問題

1. – 完全なモデルは存在しない

完璧なモデルというのは存在しないが,現実の問題をより詳細に記述するモデル–より複雑なモデルを構築することは出来る. 一方モデルが複雑になるほど理解しづらくなってしまうので,モデルの複雑さと理解のしやすさの間でバランスを取ることが重要になってくる. Whiteheadの警句 “Seek simplicity and distrust it.” は,我々は単純なモデルを正しいと思い込みがちなので,単純なモデルがうまく言っているように見えても,よく検証しなければならないと主張する.

2. – その数学的モデルの解が現実で意味を持つか

確率のモデルの正当性はKormogorvの確率の公理に従っているかで決まる.

Stochastic Processとは？

確率モデルの一種で,sample pointが時間を変数とする関数であるものをstochastic process(確率過程)という. このときあるsample pointは時刻を添字とする確率変数の列の,その時刻における1点の現れと考えることが出来る. この確率変数列の全ての元が離散確率変数であるとき特にdiscrete stochastic processという.

このコースで学ぶprocessたち

• counting process
• Poisson process
• renewal process
• Markov process
• random walkとmartingale

以下はほとんどFundamentals of probabilityでやったが,念の為復習

Kormogorovの公理

$\Omega$はある集合で,$\mathcal{F} \subset 2^{\Omega}$がeventの集合
$\Leftrightarrow$

(1) $\Omega \in \mathcal{F}$
(2) $A_1, A_2, ...\in \mathcal{F} \Rightarrow \cup A_i \in \mathcal{F}$
(3) $A \in \mathcal{F} \Rightarrow A^C \in \mathcal{F}$
つまり$\mathcal{F}$が$\Omega$の$\sigma$-algebraであるということ.

$P$が$\mathcal{F}$に確率を割り当てる($(\Omega, \mathcal{F})$のprobability measureである)
$\Leftrightarrow$

(1) $P(\Omega)=1$
(2) $A \in \mathcal{F} \Rightarrow P(A) \geq 0$
(3) $\{A_i\} \subset \mathcal{F}$が互いに素なら$P(\cup A_i) = \sum P(A_i)$

Eventの独立性

$A_1, A_2 \in \mathcal{F}$がindependent
$\Leftrightarrow P(A_1 \cap A_2) = P(A_1) P(A_2)$

仮に$A_1\in \mathcal{F_1}$が赤いサイコロをふって出る目が1であるというeventで,$A_2 \in \mathcal{F_2}$が白いサイコロをふって出る目が1であるというeventとすると,$A_1 \cap A_2$というのは$(\Omega_1 \times \Omega_2,\mathcal{F_1} \times \mathcal{F}_2, P_1\times P_2)$という新しいprobability spaceがあらわれて,例えば$\Omega_1 \times \Omega_2$は単に$\{1, ... 6\}^2$といデカルト積(とりあえずサイコロの場合はそう).
一方ふるサイコロが両方とも白いとき(区別できないとき)には話は複雑になる.$\Omega_1 \times \Omega_2 = \{(1,1), (1,2), ...,(1,6), (2, 2), (2, 3), ...,(2, 6), (3,3),...,(6,6) \}$となって,$|\Omega_1\times\Omega_2|=21$これはサイコロが区別できるときとは異なる集合である. これは2つのeventを，区別できない状況で組み合わせるときに，組み合わせたeventの確率をどう評価するかという重要な問題の最も簡単な例と言える. この話題はまたあとで扱う.

Random Variable

$X: \Omega \rightarrow \mathbb{R}$が$(\Omega, \mathcal{F}, P)$のrandom variable(r.v., 確率変数)
$\Leftrightarrow \forall r \in \mathbb{R}, \{\omega| X(\omega) \leq r\}\in\mathcal{F} \Leftrightarrow X$は$(\mathcal{F}, \mathcal{B})$-可測$\Leftrightarrow \forall B \in \mathcal{B}. X^{-1}(B)\in \mathcal{F}$

また,$F_X(x) = P(\{\omega| X(\omega) \leq x\})$を$X$のdistributinoという.

Lecture 2. 確率論の復習

Expectations

$X$のexpectation(期待値) $\bar{X} = E[X]$を
$$\begin{aligned} E[X] &= \sum_i a_i p_X(a_i) \ \ \ &\text{for discrete } X \\ E[X] &= \int x f_X (x) dx &\text{for continuous } X \\ E[X] &= \int F^C_X (x) dx &\text{for arbitrary nonegative } X \\ E[X] &= \int^0_{-\infty} F_X(x) dx + \int^\infty_0 F^C_X(x)dx & \text{for arbitrary } X \end{aligned}$$
と定める(上の二つが定義,下の二つが定義から導かれる公式)．
4番目の式の片方の項が$-\infty$,もう一方が$\infty$の場合を考えれば,expectationが定義できないrandom variableが存在することがわかる. 普通$|E[X]|< \infty$のときのみ,$E[X]$が存在するという.
また,$X$のstandard deviation(標準偏差) を$\sigma_X = \sqrt{E[|X-E[X]|^2]}$と定める.
上から3番目の式を,discreteの場合に直感的に正当化する.

figure 1.
fig.1のそれぞれの四角形の面積が$\sum_{i \in \{1,2,3,4\}} a_{i} p_X(a_i)$の各項を表している.$a_0=0$として,$x \in [a_{i-1}, a_{i}]$の$x$軸に垂直な直線と,四角形の重複する部分の長さは$1-F_X(a_i)=F_X^C(A_i)$だから,四角形の面積の和は$\sum_{i \in \{1,2,3,4\}} F^C_X(a_i)$と一致する.

example: indicator random variable

$A \in \mathcal{F}$のindicator random variable $I_A : \Omega \ni \omega \mapsto \begin{cases} 1 \ \ (\omega \in A) \\ 0 \ \ (\omega \notin A) \end{cases}$
について$PI_A (0) = 1 - Pr(A), PI_A(1) = Pr(A)$だから,$E[I_A]=Pr(A)$
$\sigma_{I_A}=\sqrt{Pr(A) (1-Pr(A))^2+(1-Pr(A))(0-Pr(A))^2}=\sqrt{Pr(A)(1-Pr(A))}$

multiple random variables

random variables $X, Y$についてjoint distribution function
$F_{XY} (x, y) = Pr(\{\omega| X(\omega) \leq x\} \cap \{\omega|X(\omega) \leq y\})$が定義できる.これはrandom variableの集合$\{X_i\}$でも動揺に定義できて,$X=(X_1,...,X_n)$がindependentなら
$$F_X(x_1, ..., x_n) = \prod_{m=1}^n F_{X_m}(x_m) \ \ \forall x_1, ...,x_n$$
である.

discrete rv’s $X, Y$について
$$p_{X|Y}(x|y) = \frac{p_{XY}(x,y)}{p_Y(y)} \ \ (\text{if } p_Y(y)\neq0)$$
をconditional probabilityという. $X,Y$がindependentなら

IID random variables

$X_1, ..., X_n$がindependent and identically distributed(IID)
$\Leftrightarrow \forall x_1,...,x_n F_X(x_1,...,x_n)=\prod_{k=1}^n F_{X_1}(x_k)$かつ$\forall x.F_{X_i}(x)=F_{X_1}(x)$
$\Omega=\mathbb{R}$で,その上のr.v. $X$があるとき,i.i.d.として並べて$\mathbb{R}^n$上のr.v. $(X_1,...,X_n)$を考えることが出来る(extended modelという).
($\Omega$がなんであっても$\Omega^n$を考えればいいような気がするが・・・・)

Sample Average

$S_n / n =(X_1 +...+X_n)/n$はある実数に収束して,その極限が,extended modelが現実世界における試行の繰り返しであるとするなら,その結果の平均であるというのが,大数の法則の主張である.
しかし,いかに述べる問題から,正しいモデルが作れないこともある.
1. 現実での試行の列というのは,それぞれが十分似ていなかったり，独立でなかったりして,i.i.d.でモデル化できないかもしれない.
2. もとのモデルが間違っているかもしれない.例えばコイントスが表が出る確率0.5としたが,実際には0.45かもしれない

実験によって得られる結果というのはsample pointであって,確率ではない. extended modelが現実と合致しているときにのみ大数の法則やその関連を使って分布を考えることが出来る.
また，$S_n-n\bar{X}$は平均$0$,分散$n\sigma^2$のr.v.で,$(S_n-n\bar{X})/\sqrt{n}\sigma$は平均0, 分散1. これが$n\rightarrow \infty$で$N(0,1)$に分布収束するというのがcentral limit theorem(CLT, 中心極限定理)の主張である.(characteristic functionの各点収束を示して証明した)

The Bernoulli Process

$p_Y(1)=p>0, p_Y(0)=1-p=q >0$がi.i.d.に並んだ$Y_1, ..., Y_n$があるとき,
$S_n = Y_1 + ... +Y_n$の分布を調べる.この節では$p < 1/2$と約束する.
$(Y_1,...,Y_n) = (\underline{1,...,1}_\text{k個},\underline{0,...,0}_{\text{n-k個}})$となる確率は$p^kq^{n-k}$.$k=0$で最大となり,$k$の増加とともに減少する.
また,$(Y_1,...,Y_n)$のうちちょうど$k$個が$1$,ほかが$0$である場合の数は$\left(\begin{array}{} n \\k \end{array} \right)$個で,それぞれの確率は$p^kq^{n-k}$だから,確率の和は$\left(\begin{array}{} n \\k \end{array} \right)p^kq^{n-k}$.よって
$$p_{S_n}(k) = \left(\begin{array}{} n \\k \end{array} \right)p^kq^{n-k}$$
$k$による増減は,
$$\begin{aligned} \frac{p_{S_n}(k+1)}{p_{S_n}(k)} &= \frac{n!}{(k+1)!(n-k-1)!}\frac{k!(n-k)!}{n!}\frac{p^{k+1}q^{n-k-1}}{p^kq^{n-k}} \\ &= \frac{n-k}{k+1}\frac{p}{q} \end{aligned}$$
だから,$k$の増加とともに狭義単調減少する.
さらに,
$$\frac{p_{S_n}(k+1)}{p_{S_n}(k)} = \begin{cases} < 1 \ \ &\text{for } k \geq pn \\ \sim 1 &\text{for } k < pn < k+1 \\ >1 &\text{for } k+1 < pn \end{cases}$$

(以下CLTの成立の証明が続く)

Assignment 1. Problem set 1.

Exercise 1.3

$A_1, A_2, ... \in \mathcal{F}$はdisjointで,$Pr(A_n)=2^{-n-1},\Omega=\cup A_n$とする.
(a) この過程が確率の公理に反することを示せ
(b) $Pr(\cup_{i\geq 1} A_i)=\sum_{i\geq 1}Pr(A_i)$を$Pr(\cup_{i=1}^n A_i)= \sum_{i=1}^n Pr(A_i)$で置き換えると,上の過程がその公理系を満たすことを示せ
この結果から,countable additivityがfinite additivityよりも強い概念であることがわかる.

答案.

(a)
$$\sum_{i \geq 1} Pr(A_i) = \sum_{n\geq 1} 2^{-n-1} = 1/2 \neq Pr(\cup_{i\geq 1} A_i)=Pr(\Omega)=1$$
よって示せた.
(b)
$\sum_{i=1}^n Pr(A_i)= (1/2)(1-2^{-n})<1$. これと$Pr(\Omega)=1$は矛盾しない.

Exercise 1.9

$X$はr.v.で,distributionは$F_X$とする.以下のr.v.たちのdistributionを与えよ
(a) $X$にIIDなr.v.たちの最大値のr.v.
(b) $X$にIIDなr.v.だちの最小値のr.v.
(c) (a)と(b)の差のr.v. ただし$X$はPDF $f_X$を持つとする.

答案.

(a)
$\max(X_1,...,X_n) \leq x \Leftrightarrow X_1 \leq x \land ...\land X_n \leq x$から
$Pr(\max(X_1,...,X_n)\leq x ) = Pr(X_1 \leq x, ....,X_n\leq x) = Pr(X_1\leq x)...Pr(X_n\leq x)$(独立性)
よって$F_{\max{(X_1,...,X_n)}}(x) = \prod_{i=1}^n F_{X_i}(x) = F_X(x)^n$
(b)
(a)とほとんど同様に,
$F_(\min{(X_1,...,X_n)})=(1-F_X(x))^n$
(c) (模範解答)
求めるr.v. を$R$とする.
$MAX= \max(X_1,...,X_n), MIN=\min(X_1,...,X_n)$とする.
$MAX\leq x \land MIN > y \Leftrightarrow \forall i. y<X_i \leq x$から
$Pr(MAX \leq x, MIN>y) = [F_X(x)-F_Y(y)]^n$.
$$Pr(R\leq r) = \int \frac{\partial Pr(MAX \leq x, MIN >y)}{\partial x} |_{y=x-r}dx$$
$X$がPDF $f_X$を持つことから,これは$\int n f_X(x)[F_X(x)-F_X(x-r)]^{n-1}dx$である.

Exercise 1.13

$X_1,...$はPDF $f_X(x)$をもつr.v.のIIDな無限列とする. $n\geq 2$について,$X_n$をrecord-to-date,すなわち$\forall i< n. X_n > X_i$と定める. IIDの対称性を使って以下の問いに答えよ.
(a) $X_2$がrecord-to-dateである確率を求めよ.
(b) $X_n$がrecord-to-dateである確率を$n$の関数として求めよ
(c) 任意の$m$について,最初の$m$回の試行におけるrecord-to-dateの個数の期待値を求めよ. $m\rightarrow \infty$によって期待値が$\infty$に発散すると示せ.

答案.

(a) $Pr(X_1<X_2)= Pr(X_1>X_1)=1/2$
(b) $X_1 > \max(\cup\{X_i\}\backslash\{X_1\}), X_2 > \max(\cup\{X_i\}\backslash\{X_2\},..., X_n > \max(\cup \{X_i\}\backslash\{X_n\})$
のどれか一つが確立1で成立し,対称性から$Pr(X_n>X_1,...,X_{n-1})=1/n$.
(c) $\underline{1}_{\text{個数}}\cdot (\underline{1/2}_{(2)} + \underline{1/3}_{(3)} + \cdots \underline{1/m}_{m})$ ($(i)$は$X_i$がrecord-to-dateである確率)
$\sum_{i \geq 2} 1/i$が発散することは有名である.
(模範解答では$X_1$は必ずrecord-to-dateになっていた.全件否定からたしかに$X_1$はrecod-to-date)

Exercise 1.20

(a) $X, Y, Z$は確率1/2で1,確率1/2で0をとるbinary rv’sとする. $\{X,Y,Z\}$がdependentだが,それぞれ2つを選ぶとindependentになる例を挙げ,PMF $p_{XYZ}(x,y,z)$を求めよ (hint: もっとも単純な例では,ただ4つのjoint probabilityが正となる)
(b) pairwise independenceは
$$E[\prod_{i=1}^n X_i] = \prod_{i=1}^n E[X_i]$$
の十分条件か.

模範解答.

(a) $Z= (X+Y)_{mod2}$
(b) この例で, $Pr(XYZ=0)=1$.よって十分条件でない.

Exercise 1.26

r.v. $X$はcontinuousで,distributionは$F_X(x)$とする. $Y=F_X(X)$を新たなr.v.として考える.すなわち,$\omega \in \Omega$に,$X(\omega)=x$ならば$Y(\omega)=F_X(x)$ということである. $Y$は$[0, 1]$上でuniformly distributedであることを示せ.

答案.

$F_Y(y) = y\ \ (0\leq y \leq 1)$を言えば良い.
$F_Y(y) = Pr(Y\leq y) = Pr(\{\omega|Y(\omega)\leq y \})=Pr(\{\omega| F_X(X(\omega)) \leq y \})$
$F_X$の単調性から,ある$r$があって,$x \leq r \Rightarrow F_X(x)\leq y$.$r$の最大値(連続性から存在)を$F^{-1}(y)$とすると$F_Y(y) = Pr(X\leq F^{-1}_X(y)) = F_X(F^{-1}_X(y))=y$.

MIT OCW, Machine Learning 08日目 カーネル

Rohit Singh, Tommi Jaakkola, and Ali Mohammad. 6.867 Machine Learning. Fall 2006. Massachusetts Institute of Technology: MIT OpenCourseWare, https://ocw.mit.edu. License: Creative Commons BY-NC-SA.

• • Lecture 6.
    • Active Learning (cont)
    • Non-linear Predictions, Kernels
  • Lecture 7.
    • Linear Regression and Kernels
    • Kernels

Lecture 6.

Active Learning (cont)

$y = \theta^{*}\mathbf{x} + \theta^*_0 + \epsilon, \ \ \ \epsilon \sim N(0,\sigma^2)$というlinear modelについて,最尤法で推測されるパラメータ$\hat{\theta}, \hat{\theta_0}$のMSEは
$$E\left[\left\|\left[\begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right] - \left[\begin{array}{} {\theta}^{*} \\ {\theta_0}^* \end{array} \right] \right\|^2 | \mathbf{X} \right] = \sigma^{*2} Tr[(\mathbf{X^TX})^{-1}]$$
となることから,$\mathbf{X}$をうまく設計することで少ないexampleからよりよい推測を行うことをactive learningといった. $\mathbf{X}$の設計で最も単純な方法は,$\mathbf{x_1}, ..., \mathbf{x_k}$があるときに,$Tr[\mathbf{X^TX}]$が最少になるように$\mathbf{x}_{k+1}$を選ぶという操作を繰り返すというのがある. すでに$\mathbf{X}$があって,$\mathbf{A}=(\mathbf{X^TX})^{-1}$とする. $[\mathbf{x}^T, 1]$を$\mathbf{X}$の行に新たに加えることを考える.
$$\left[\begin{array}{} \mathbf{X} \\ \mathbf{x}^T 1 \end{array} \right]^T\left[\begin{array}{} \mathbf{X} \\ \mathbf{x}^T 1 \end{array} \right] = (\mathbf{X^TX}) + \left[\begin{array}{} \mathbf{x} \\ 1 \end{array} \right]\left[\begin{array}{} \mathbf{x} \\ 1 \end{array} \right]^T = \mathbf{A}^{-1} + \mathbf{vv^T} \ \ \ (\mathbf{v}=[\mathbf{x}^T, 1]^T)$$
$$Tr[(\mathbf{A}^{-1} + \mathbf{vv^T})^{-1}]$$
を最小化する$\mathbf{v}$を考える.
$$(\mathbf{A}^{-1} + \mathbf{vv^T})^{-1} = \mathbf{A} - \frac{1}{1 + \mathbf{v^TAv}} \mathbf{Avv^TA}$$
であって,$Tr(A+B) = Tr(A)+Tr(B), Tr(AB)=Tr(BA)$を考えれば
$$Tr[(\mathbf{A}^{-1}+\mathbf{vv^T})^{-1}] = Tr[A] - \frac{\mathbf{v^T AAv}}{1 + \mathbf{v^TAv}}$$
が成立する. ($\mathbf{v^TAAv}$は実数で,traceはその実数そのもの)
任意の$\mathbf{v}$に$\frac{\mathbf{v^TAAv}}{1+\mathbf{v^TAv}}>$だから,どのような$\mathbf{x}$を加えたとしてもMSEは減少するが,減少量が最大であるような$\mathbf{x}$を求めたい.
$$\frac{\mathbf{v^TAAv}}{1+\mathbf{v^TAv}}$$
の大きさは$A$の最大の固有値が上限である. 言い換えると,新しいexampleによってパラメータ空間からせいぜい1つだけ自由度を減じることが出来る. $\mathbf{x}$に制限がなければ,$\mathbf{A}$の最大の固有値に対応する固有ベクトルに平行な長さ無限のベクトルを$\mathbf{v}$とするのだが,$\|\mathbf{v}\|\leq c$という制限が有る場合には,最大固有値に対応する固有ベクトルと平行でながさ$c$のベクトルを$\mathbf{v}$とする. ほかにも$\mathbf{x}$に制限が有る場合には,$\mathbf{v}$もその制限を考慮することになる.

これまでMSEを推定量の良さの基準としてきたが,今度はvarianceを考える.
$$\begin{aligned} var[y|\mathbf{x, X}] &= E\left[(\hat{\theta}^T\mathbf{x}+\hat{\theta}_0 - \theta^{*T}\mathbf{x} -\theta^{*}_0)^2 | \mathbf{x, X} \right]\\ &=E\left[\left[\begin{array}{} \mathbf{x} \\ 1 \end{array} \right]^T \left(\left[ \begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right] - \left[ \begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right]\right) \left(\left[ \begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right] - \left[ \begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right]\right)^T\left[\begin{array}{} \mathbf{x} \\ 1 \end{array} \right]| \mathbf{x, X}\right]\\ &= \left[\begin{array}{} \mathbf{x} \\ 1 \end{array} \right]^T \sigma^{*2}(\mathbf{X^TX})^{-1}\left[\begin{array}{} \mathbf{x} \\ 1 \end{array} \right] \\ &= \sigma^{*2}\cdot \mathbf{v^TAv}\end{aligned}$$

よって,$\mathbf{v^TAv}$が最大になるような$\mathbf{v}$がよいが,これはMSEを小さくするような$\mathbf{v}$と同じである.
(MSEを小さくしつつvarianceを小さくすることが対立することを言いたいのかと思ったら,varianceを大きくしたいらしい・・・)

Non-linear Predictions, Kernels

$\mathbf{x}$の非線形な写像に対する像$\phi(\mathbf{x})$に対してこれまで議論してきた方法が使える.例えば$y = \theta x + \theta_0 + \epsilon, \epsilon \sim N(0, \sigma^2)$というlinear modelが有るとき,$x$を$x^2$を含む高次元のベクトルに写像してquadratic(二次) modelが得られ,$x^3$を含む高次元のベクトルに写像するとthird order modelが得られる.
$\phi(x) = [1, \sqrt{2}x, x^2]^T, \phi(x)=[1, \sqrt{3}x, \sqrt{3}x^2, x^3]^T$のような感じである.$\sqrt{2}$や$\sqrt{3}$の意味は後で見る.
新しいpolynomial regression modelは
$$y = \theta^T \phi(x) + \theta_0 + \epsilon, \ \ \epsilon \sim N(0, \sigma^2)$$
となる. 高次元空間に写像してから線形回帰するわけだが,このときregularizationを行わないとoverfittingが起きることが多い.(figure 2)
!

$\mathbf{x}$が多次元の場合も,
$$\mathbf{x}=[x_1, x_2]^T\mapsto^{\phi} [1, x_1, x_2, \sqrt{2}x_1x_2, x_1^2,x_2^2]^T = \phi(\mathbf{x})$$
というふうにしてより高次元な空間に写像できる.
高次元な空間への変換は計算コストが膨大になることが有るが,$\phi$を直接計算せずとも,例えば
$$\begin{aligned} \phi(x) &= [1, \sqrt{3}x, \sqrt{3}x^2, x^3]^T \\ \phi(x') &= [1, \sqrt{3}x', \sqrt{3}x^{'2}, x^{'3}]^T \\ \phi(x)^T\phi(x') &= 1 + 3xx' + 3(xx')^2 + (xx')^3 = (1+xx')^3 \end{aligned}$$
のように,$\phi(x)^T\phi(x')=k(x,x')$と,$\phi$を暗黙に表現する計算が簡単な$K$が存在することが有る(存在するように$\phi$を定めたのである). $\phi$ではなく計算が簡単な$K$を使うように問題を書き換えることを考える.

Lecture 7.

Linear Regression and Kernels

$\theta_0$を外したモデル$y = \theta^T \phi(\mathbf{x}) + \epsilon$はの推測は
$$J(\theta) = \sum_{t=1}^n (y_t-\theta^T \phi(\mathbf{x}_t))^2 + \lambda\|\theta\|^2$$
の最適化問題である. 前節で述べたとおり,$\phi$ではなく$k$でこの最適化問題を表現する.
regularizationによって$\theta$は$0$に圧縮され,$\theta$のtraining feature vectorと関係ない次元は$0$になる. よってこの問題の解は$\{\phi(\mathbf{x}_t)\}$の張る空間の元である.
proof.

局地の条件を考えると
$$\frac{dJ}{d\theta} = -2 \sum_{t=1}^n \underline{(y_t-\theta^T\phi(\mathbf{x}_t))}_{\alpha_t}\phi(\mathbf{x}_t) + 2\lambda \theta=0$$
$$\theta = \frac{1}{\lambda} \sum_{t=1}^n \alpha_t \phi(\mathbf{x}_t)$$
は$\frac{dJ}{d\theta}=0$を満たして,最適解である.

$$\alpha_t = y_t - \theta^T \phi(\mathbf{x}_t)=y_t - \frac{1}{\lambda}\sum_{t'=1}^n \alpha_{t'} \phi(\mathbf{x_{t'}})^T \phi(\mathbf{x}_t)$$
が成立するから,$\alpha_t$は$y_t$と$\phi(\mathbf{x}), \phi(\mathbf{x'})$だけで決まる.
Gram行列
$$\mathbf{K} = \left[\begin{array}{} \phi(\mathbf{x}_1)^T\phi(\mathbf{x}_1) & \phi(\mathbf{x}_1)^T\phi(\mathbf{x}_2) & \cdots & \phi(\mathbf{x}_1)^T\phi(\mathbf{x}_n) \\ \cdots & \cdots & \cdots & \cdots \\ \phi(\mathbf{x}_n)^T\phi(\mathbf{x}_1) & \cdots & \cdots &\phi(\mathbf{x}_n)^T\phi(\mathbf{x}_n) \end{array} \right]$$
によってベクトルで書くと
$$\begin{aligned} \mathbf{a} &= [\alpha_1, ..., \alpha_n]^T \\ \mathbf{y} &= [y_1, ..., y_n]^T \\ \mathbf{a} &= \mathbf{y} - \frac{1}{\lambda} \mathbf{Ka} \end{aligned}$$
そして解は
$$\hat{\mathbf{a}} = \lambda(\lambda \mathbf{I} + K)^{-1} \mathbf{y}$$
$\hat{\alpha}_t$が得られたら,
$$y = \hat{\theta}^T \phi(\mathbf{x}) = \sum_{t=1}^n (\hat{\alpha_t}/\lambda)\phi(\mathbf{x}_{t'})^T\phi(\mathbf{x})=\sum_{t=1}^n\hat{\alpha}_tK(\mathbf{x_{t'}, x})$$
によって,新しいexample $\mathbf{x}$に対してresponseの推測$y$が計算できる.ここで$K(\mathbf{x_{t'}, x})$はkernel functionという.

Kernels

以上で, regularized linear regressionをkernel formに変形できた. kernel function $K$を変えることで,例えば任意の次数のpolynomial expansionが実現できるし,polynomial expansion以外の$\mathbf{x}$を高次元に写した像を使ったlinear regressionも実現できる.
実現される高次元への写像の種類によって$K$を分類することが有る.例えば
- Polynomial kernel

$$K(\mathbf{x', x})=(1 + \mathbf{x^Tx'})^p, \ \ p = 1,2,...$$
- Radial basis kernel
$$K(\mathbf{x', x}) = \exp \left(-\frac{\beta}{2}\|\mathbf{x}-\mathbf{x}'\|^2 \right), \ \ \beta>0$$

polynomial kernelは,$\mathbf{x}=[x_1,...,x_n]^T$を$(x_1+\cdots +x_n)^p$を二項展開したときの各項へと写す写像$\phi$を考えたときのkernel functionで, radial basis kernelは無限次元空間への写像のkernel functionである. radial basis functionは$\mathbf{x}$と$\mathbf{x'}$の近さを表していると考えることが出来る.

MIT OCW, Fundamentals of Probability 24日目 Markov Chain III

David Gamarnik, and John Tsitsiklis. 6.436J Fundamentals of Probability. Fall 2008. Massachusetts Institute of Technology: MIT OpenCourseWare, https://ocw.mit.edu. License: Creative Commons BY-NC-SA.

• • Lecture 25. Markov Chain III. Periodicity, Mixing, Absorption
    • 25.1 Periodicity
      • Lemma 25-1
      • Lemma 25-2 (証明略)
      • Theorem 25-3 (証明略)
      • Theorem 25-4 (証明略)
    • 25.2 Absorption Probabilities and Expected Time to Absorption
      • Example: Gambbler’s Ruin

Lecture 25. Markov Chain III. Periodicity, Mixing, Absorption

25.1 Periodicity

$x \in \mathcal{X}$がrecurrentであるとき,$x$がそれ自身からaccessibleである時刻,すなわち$I_x = \{n\geq 1: p_{xx}^{(n)} > 0\}$を考える. $I_x$は和に対して閉じている(i.e. $m,n\in I_x \Rightarrow m+n\in I_x$). $d_x$を$I_x$の元の最大公約数として,$x$のperiodという. periodの諸性質を論じる.

Lemma 25-1

$x, y$が同じrecurrentにあるとき,$d_x = d_y$である.

proof.

$p^{(m)}_{xy}, p^{(n)}_{yx} > 0$である$m, n$を選ぶ(同じrecurrentだから存在する). $p_{yy}^{(m+n)} \geq p_{xy}^{(m)}p_{yx}^{(n)} > 0$だから$d_y$は$m+n$を割り切る. また$l$を$p_{xx}^{(l)} > 0$なる$l \in \mathbb{N}$とすると,$p_{yy}^{m+n+l} \geq p_{yx}^{(n)}p_{xx}^{(l)}p_{xy}^{(m)}>0$だから,$d_y$は$m+n+l$を割り切り,故に$l$を割り切る. したがって$d_y$は$d_x$を割り切る. 同じ論法で$d_x$が$d_y$を割り切ることも言えて,以上より$d_x=d_y$

$d > 1$であるようなrecurrent classをperiodicといい,$d=1$であるときにはaperidicという. periodicityは$p_{xy}^{(n)}$が$\pi_y$に収束することを妨げている. $y$がperiodicなrecurrent classの元とすると,$p_{yy}^{(n)}=0$が,$n$が$d$の倍数でない限り成立するが,$\pi_y>0$である. 一方$d=1$(aperiodic)であれば,十分大きな全ての$n$に,markov chainが$y$に戻ってくる確率が正になる.

Lemma 25-2 (証明略)

$d_y=1$であれば,$N \geq 1$があって,$n \geq N\Rightarrow p_{yy}^{(n)}>0$である.

Markov chainがただ一つのrecurrent classをもち(irreducible),かつaperiodicであるとき,steady stateの振る舞いはstationary distributionによって与えられる.この事実をmixingという.

Theorem 25-3 (証明略)

irreduibleかつaperiodicなMarkov chainがあるとき,任意のstateの組$x, y$について$\lim_{n \rightarrow \infty} p_{xy}^{(n)}=\pi_y$

periodicな場合には,$p_{xy}^{(n)}$の部分列の収束に関する定理が有るがここでは扱わない.
$\pi_x p_{xy} = \pi_y p_{yx}$が任意の$x, y \in \mathcal{X}$に成り立つとき,そのMarkov chainはreversibleであるという. Theorem 25-3の仮定にreversible性を加ええ場合の重要な定理が知られている.

Theorem 25-4 (証明略)

irreducible, aperiodic, reversibleなMarkov chainについて,任意の$x, y \in \mathcal{X}$に$|p_{xy}^{(n)} - \pi_y| \leq C|\lambda_2|^n$が成り立つような定数$C$が存在する.ただし$\lambda_2$は$P$の二番目に絶対値が大きいeigenvalueとする.

$|\lambda_2|<1$だから,これは$p_{xy}^{(n)}$の$pi_y$への収束の速さを与える.

25.2 Absorption Probabilities and Expected Time to Absorption

Markov chainの長期的な振る舞いを見てきたが,今度は短期的な振る舞いを議論する. 特にtransientなstateから始まったchainの振る舞いを考える. 簡単のため,recurrent state $i$はabsorbingであるとする.すなわち$p_{ii}=1$である. これから考察するMarkov chainはtransient classのほかは全てabsorbingとする.
absorbing state $i$がただ一つであるときには$\pi_i = 1$であって,必ず$i$に到達する. 一方absorbing stateが複数存在するときには,どのabsorbing stateに至るかは確率的に決まる.
$$a_{ki} = P(X_n \text{eventually equals i} | X_0 = k)$$
をabsorbing probabilityという. $j$がabsorbingなら$a_{jj}=0, a_{ji}=0$である.
$k$がtransientなら
$$\begin{aligned} a_{ki} &=P(\exists n: X_n=i|X_0=k) \\ &=\sum_{j=1}^N P(\exists n: X_n=i|X_1=j)p_{kj} \\ &= \sum_{j=1}^N a_{ji}p_{kj} \end{aligned}$$
だから,この線形連立方程式を解くことでabsorption probabilityを得られる.

Example: Gambbler’s Ruin

あるギャンブラーが一回の勝負ごとに$p$の確率で1ドルを得て,$1-p$の確率で1ドルを失うとする. それぞれの勝負は独立であるとする. ギャンブラーは$m$ドルを稼ぐか金を全て失うまで勝負を続ける. 彼が全財産を失う確率を求めよ

$i$はギャンブラーの持つ金額として,Markov chain $\mathcal{X} = \{0, 1, ..., m\}$を考える. $i=0$なるとき,彼は全財産を失うったということであり,$i=m$となるとき,彼は目的を達成したということである. $0, m$はabsorbing stateであると言える.
transition probabilityは$p_{i, i+1} = p, p_{i, i-1}=1-p$が全ての$i=1, ..., m-1$で成立する.また$p_{00}=p_{mm}=1$である. $i=0$のabsorbing probabilityは
$$\begin{aligned} a_{00} &= 1 \\ a_{m0} &= 0 \\ a_{mm}&=1\\ a_{i0} &=(1-p)a_{i-1, 0} + p(a_{i+1}, 0) \ \ \ \ \ \text{for } i=1, ..., m-1\end{aligned}$$
によって計算できる.$b_i = a_{i0} - a_{i+1, 0}, \rho = (1-p)/p$とすると,上の方程式から
$$\begin{aligned} (1-p)(a_{i-1, 0} -a_{i, 0}) &= p(a_{i0} -a_{i+1, 0}) \\ b_i &= \rho b_{i-1} \end{aligned}$$
であって,故に$b_i = \rho^i b_0$である. $b_0+b_1+\cdots +b_{m-1} = a_{00} - a_{m0} = 1$であって, $(1+\rho + ...+rho^{m-1})b_0 = 1$であって,
$$b_i = \begin{cases} \frac{\rho^i (1-\rho)}{1-\rho^m} \ \ \ &if \rho \neq 1 \\ \frac{1}{m} & otherwise\end{cases}$$
さらに$a_{i,0}$は$\rho \neq 1$ならば$i=1,...,m-1$について
$$a_{i0} = a_{00} -b_{i-1} - ... -b_0 = \frac{\rho^i - \rho^m}{1-\rho^m}$$
$\rho = 1$ならば
$$a_{i0} = \frac{m-i}{m}$$
したがって$i$がいかなる値でも$m$が大きくなると全財産を失う確率が$1$に近づく.

Discrete Stochastic Processesに続く

MIT OCW, Machine Learning 07日目 リッジ回帰

Rohit Singh, Tommi Jaakkola, and Ali Mohammad. 6.867 Machine Learning. Fall 2006. Massachusetts Institute of Technology: MIT OpenCourseWare, https://ocw.mit.edu. License: Creative Commons BY-NC-SA.

Penalized Log-Likelihood and Ridge Regression

training dataが,その各exampleの次元$d$に対して十分に大きくないときには,パラメータをregularizeすることが多い. prior distributionを$P(\theta, \theta_0)$にassign することで,どのようにregularizeすればよいかを見る. prior distributionは, パラメータの推測値の絶対値を小さくするために導入する.
prior distributionを平均0のnormal distributionとする.つまり
$$P(\theta, \theta_0, \sigma^{'2})=N(\theta_0;0, \sigma^{'2})$$
をlikelihood $L$に追加すると
$$\begin{aligned}l'(\theta, \theta_0, \sigma^2) &= \sum_{t=1}^n \log \left[\frac{1}{\sqrt{2\pi \sigma^2}} \exp \left(-\frac{1}{2\sigma^2} (y_t-\theta^T \mathbf{x}_t - \theta_0)^2 \right) \right] + \log P(\theta, \theta_0; \sigma^{`2}) \\ &= \text{const.} - \frac{n}{2} \log \sigma^2 - \frac{1}{2\sigma^2} \sum_{t=1}^n (y_t - \theta^T \mathbf{x}_t - \theta_0)^2 \\ &-\frac{1}{2\sigma^{'2}}(\theta_0^2 + \sum_{j=1}^d \theta_j^2) - \frac{d+1}{2} \log \theta^{'2} \end{aligned}$$
また,$\sigma^{'2}=\sigma^2/ \lambda$とすることも多い. $\sigma^2$が小さいときにはoverfittingのおそれがあるので,よりpenallityを大きくしてパラメータを0に近づけるのである. training dataが小さいときには$\sigma^2$が小さくなりなちなので,この節のはじめに説明したregularizationをする動機と合目的である.
$\sigma^{'2}=\sigma^2/\lambda$を$l'$に代入すると
$$\begin{aligned} l'(\theta, \theta_0, \sigma^2) &= \text{const.} - \frac{n}{2}\log \sigma^2 - \frac{1}{2\sigma^2}\sum_{t=1}^n (y_t - \theta^T \mathbf{x}_t - \theta_0)^2 \\ &-\frac{\lambda}{2\sigma^2}(\theta^2_0 + \sum_{j=1}^d \theta_j^2) - \frac{d+1}{2} \log (\sigma^2/ \lambda) \\ &= \text{const.} - \frac{n+d+1}{2} \log\sigma^2 + \frac{d+1}{2} \log \lambda \\ &-\frac{1}{2\sigma^2} \left[ \sum_{t=1}^n (y_t-\theta^T \mathbf{x}_t - \theta_0)^2 + \lambda (\theta^2_0 + \sum_{j=1}^d \theta^2_j)\right]\end{aligned}$$
このregularization problemの解を求めることをRidge regressionという.
その解$\hat{\theta}, \hat{\theta_0}$は,
$$\left[\begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right] = (\lambda \mathbf{I} + \mathbf{X^TX})^{-1} \mathbf{Xy}$$
で与えられる.
$$\begin{aligned}E \left[\left[\begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right]|\mathbf{X} \right] &=(\lambda \mathbf{I}+\mathbf{X^TX})^{-1}\mathbf{X^TX} \left[\begin{array}{}{\theta^*} \\ {\theta_0^*} \end{array} \right] \\&=(\lambda \mathbf{I} +\mathbf{X^TX})^{-1} (\mathbf{X^TX + \lambda I -\lambda I})\left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] \\ &= \left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] - \lambda(\lambda \mathbf{I} + \mathbf{X^TX})^{-1} \left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] \\ &= (\mathbf{I} - \lambda(\lambda \mathbf{I} + \mathbf{X^TX})^{-1}) \left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] \end{aligned}$$
だから,$\hat{\theta}, \hat{\theta_0}$はbiasedな推測である. また$(\mathbf{I} - \lambda(\lambda \mathbf{I} + \mathbf{X^TX})^{-1})$は固有値が1未満の正定値行列で,$\lambda$が大きくなるとともに$\theta, \theta_0$は$0$へと近づいていく. 以前やったのと同じ方法で MSEを計算すると,
$$\begin{aligned} E\left[ \left\| \left[\begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right] - \left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] \right\| | \mathbf{X} \right] &= \sigma^{*2} \cdot Tr[(\lambda \mathbf{I} + \mathbf{X^TX})^{-1} - \lambda(\lambda \mathbf{I} + \mathbf{X^TX})^{-2}] \\& +\lambda^2 \left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right]^T (\lambda \mathbf{I} + \mathbf{X^TX})^{-2}\left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] \end{aligned}$$
であって,これはregularizationを考えない場合のMSE $\sigma^{*2} \cdot Tr[(\mathbf{X^TX})^{-1}]$よりも小さく出来る.

Active Learning

training data $\{\mathbf{x}_1, ..., \mathbf{x}_n\}$を能動的に選んでestimation errorを小さくすることを,active learning問題という. 例えば画像の分類で,すでにたくさんのtraining dataのもととなるlabelなしの画像データが有るが,そこからできるだけ少なくデータを選んでラベル付けし(ときにラベル付は画像そのものの収集よりコストがかかる),training dataとする状況を考える. 推測の正確性を犠牲にせずに,できるだけ選ぶ画像データを少なくする方法を考えるのである.
この問題を考察するため,regularizationの無い場合のestimation errorを再掲する.
$$\begin{aligned} E\left[ \left\| \left[\begin{array}{} \hat{\theta} \\ \hat{\theta_0} \end{array} \right] - \left[\begin{array}{} {\theta}^* \\ {\theta_0}^* \end{array} \right] \right\| | \mathbf{X} \right] =\sigma^{*2} Tr[(\mathbf{X^TX})^{-1}] \end{aligned}$$
$\sigma^{*2}$はtraining dataの選び方によらないので,$Tr[(\mathbf{X^TX})^{-1}]$が小さくなるようにすれば良い. ただし,この方法はexampleと推定値の写像の線形性を仮定しているから,そうでない場合には使えない.

MIT OCW, Machine Learning 06日目 宿題2

Rohit Singh, Tommi Jaakkola, and Ali Mohammad. 6.867 Machine Learning. Fall 2006. Massachusetts Institute of Technology: MIT OpenCourseWare, https://ocw.mit.edu. License: Creative Commons BY-NC-SA.

Assignment

Problem set 1 Section B

問

1.

パーセプトロンの実装.pythonを使う.

import numpy as np
import matplotlib.pyplot as plt

def sign(r):
    if r >= 1:
        return 1
    else:
        return -1

class PerceptronClassifier:
    def get_params(self):
        gamma = np.argmin(np.abs([np.dot(self.theta, x) for x in self.train_X]))
        gamma_geom = gamma / np.linalg.norm(self.theta)
        print("theta is: {0}, k till the convergence is {1}".format(self.theta, self.k))
        print("The angle with [1, 0] and theta is: {0} (rad)".format(np.arccos(self.theta[0]/np.linalg.norm(self.theta))))
        print("The geometric margin is {0}".format(gamma_geom))

    def perceptron_train(self, X, y):
        self.train_X = X
        self.train_y = y
        self.theta = np.zeros(len(X[0]))
        self.k = 0
        while True:
            cnt = 0
            for i in range(len(self.train_y)):
                if y[i] != sign(np.dot(self.theta, self.train_X[i])):
                    self.k += 1
                    self.theta = self.theta + y[i]*X[i]
                    cnt += 1
            if cnt == 0:
                break

    def perceptron_test(self, test_X, test_y):
        self.test_X = test_X
        self.test_y = test_y
        errors = 0
        for i in range(len(test_y)):
            if test_y[i] != sign(np.dot(self.theta, test_X[i])):
                errors += 1
        print("The error ratio is: {0}".format(errors/len(test_y)))

    def draw_graph(self, train=True):
        if train:
            X = self.train_X
            y = self.train_y
        else:
            X = self.test_X
            y = self.test_y

        plus = np.array([x for x in X if np.dot(self.theta, x)>=0 ])
        minus = np.array([x for x in X if np.dot(self.theta, x) < 0])

        plt.scatter(plus[:, 0], plus[:, 1], color='red', s=2)
        plt.scatter(minus[:, 0], minus[:, 1], color='blue', s=2)
        plt.show()

模範解答とはことなった結果を示すが,収束までの更新回数$k$や$\gamma_{geom}$は更新を行う前に定義する$\theta$の初期値にわりと鋭敏に反応するので,深く考えなくてもいいかもしれない(MATLABとPythonの精度も関係しているかも？).

X_a = np.loadtxt('p1_a_X.dat' )
y_a = np.loadtxt('p1_a_y.dat')
X_b = np.loadtxt('p1_b_X.dat')
y_b = np.loadtxt('p1_b_y.dat')

per_cla = PerceptronClassifier()
per_cla.perceptron_train(X_a, y_a)
per_cla.perceptron_test(X_a, y_a)
per_cla.draw_graph('Dataset A')
per_cla.get_params()

per_cla = PerceptronClassifier()
per_cla.perceptron_train(X_b, y_b)
per_cla.perceptron_test(X_b, y_b)
per_cla.draw_graph('Dataset B')
per_cla.get_params()

2, 3.

SVMの実装. quadratic programを解く関数を使っていいらしいがpythonだとpipにも入ってないから飛ばす