統計検定1級に合格しました

合格したので感想を。合格した勉強法の1例として参考になれば。
2018年11月下旬に統計検定1級を受験して統計数理と統計応用（理工学を選択）の両方に合格しました。以下、その受験記です。
*1

背景

自分のことを少しだけ書きます。
普段は数学を専攻している大学生です。とは言っても大学で統計学を学んだことはなく、2018年3月ごろから趣味でデータサイエンス周りの勉強を始めました。2018年5月には統計検定2級に合格しています。良ければその受験記も見て下さい。
統計検定2級に合格しました - lazy diary
さて、「2級に合格したなら次は1級でしょ」ということで受験しました。準1級を飛ばしたのは単純に時期が合わなかったからです。

統計検定1級の概要

年1回11月下旬開催の試験です。試験範囲をチェックしましょう。
www.toukei-kentei.jp
試験内容は統計数理と統計応用の2つに分かれています。
この両方に合格すると統計検定1級合格となります。

・統計数理

試験日の午前に行なわれる90分の筆記試験です。60/100点で合格です。
5つの問題から3問を選んで解きます。確率分布に関わる様々な計算がメインです。
確率密度関数の計算、確率母関数を用いた計算、検定論が頻出です。

・統計応用

試験日の午後に行なわれる90分の筆記試験です。60/100点で合格です。
統計応用は「人文科学」「社会科学」「理工学」「医薬生物学」の中から1つの分野を試験申し込み時に選んで解くことになります。
事前に選んだ分野の5つの問題から3問を選んで解きます。

勉強したこと

まず参考書選びです。1級に合格したというブログを調べるとちらほら出てくるのでそこにオススメしてある参考書と、あとは公式サイトにある合格者の声にリストアップされている参考書を書店でチェックしました。
www.toukei-kentei.jp
その中で私が選んだのがこちらの本。

・現代数理統計学の基礎 (共立講座 数学の魅力)

現代数理統計学の基礎 (共立講座 数学の魅力)

この本が無ければ合格は無かったです。それくらい良い本でした。
おそらくこの本は統計検定1級受験を視野に入れて書かれた本です。統計検定1級の範囲となっているが本には掲載されていない事項の一部を著者がWebページに補足資料として載せています。それくらい統計検定に最適化された本です。この本で統計数理のほとんどの試験範囲を抑えることができます。以下、『現代数理統計学の基礎』の目次です。

第1章　確率
第2章　確率分布と期待値
第3章　代表的な確率分布
第4章　多次元確率変数の分布
第5章　標本分布とその近似
第6章　統計的推定
第7章　統計的仮説検定
第8章　統計的区間推定
第9章　線形回帰モデル
第10章　リスク最適性の理論
第11章　計算統計学の方法
第12章　発展的トピック：確率過程
引用：https://www.kyoritsu-pub.co.jp/bookdetail/9784320111660

第1章から第9章までは十分に理解しておくべきでしょう。試験問題の半分くらいは適切に確率密度関数や期待値、分散などの計算ができるかどうかです。この計算がスムーズにできるかどうかが鍵だと思います。
第10章「リスク最適性の理論」は飛ばしても問題ないと思います。
第11章、第12章についても理解しておくべきでしょう。ただし、マルチンゲールといった試験に出ないであろう内容もあるので試験範囲と照らし合わせて必要に応じていくつか飛ばしても良いと思います。

この本を読むのに1ヵ月と少しかかりました。さっそく過去問に取り組みます。2017年の統計数理の問題を解きました。自己採点をしたら合格ライン上の点数だったと思います。それより時間が足りません。明らかに演習不足なのを実感しましたが、同時に解説を読めば理解できるということも分かったので統計数理の勉強はひとまず終了にしました。

次に統計応用の勉強をします。私は「理工学」分野を選択しました。理由は統計数理とさほど範囲が変わらないからです。過去問を見てもらえば分かると思います。「理工学」以外の分野に明るい場合はそれを選べば良いと思います。特にそういうのがないですという場合は、「理工学」分野を選ぶのが無難だと私は思います。
参考書は次の2冊です。

・多変量解析入門――線形から非線形へ

多変量解析入門――線形から非線形へ

１　はじめに
２　線形回帰モデル
３　非線形回帰モデル
４　ロジスティック回帰モデル
５　モデル評価基準
６　判別分析
７　ベイズ判別
８　サポートベクターマシーン
９　主成分分析
10　クラスター分析
引用：https://www.iwanami.co.jp/book/b265441.html

この本はかなり丁寧に記述されています。統計検定云々は関係なく、多変量解析を学びたい人にはオススメしたい1冊でした。
紹介しておいて言うことじゃないかもしれませんが、「試験日まで時間がない」とか「とにかく近道して受かりたい」といった場合はこの本は無くてもおそらく大丈夫ということを強調しておきます。
この本の主な内容は非線形回帰やロジスティック回帰、SVM、部分空間法などです。これまでRやPythonなどで統計モデリングをしてきた方には馴染み深い内容です。これらは準1級試験では出題されることがあるのですが、1級の過去問を見てもこのような内容は1度も出題されていないはずです。
（マーク試験ならまだしも、筆記問題として適切な形で出題できないというのが本音なのでしょうか？？　ただし、試験範囲には入っているので出題されても文句は言えないです）

・日本統計学会公式認定 統計検定1級対応 統計学

日本統計学会公式認定 統計検定1級対応 統計学

第I部　統計数理
第1章　確率と確率変数
第2章　種々の確率分布
第3章　統計的推定
第4章　仮説検定
第5章　主なデータ解析手法

第II部統計応用
第6章　統計応用共通手法
第7章　人文科学分野キーワード
第8章　社会科学分野キーワード
第9章　理工学分野キーワード
第10章　医学生物学分野キーワード
引用：http://www.tokyo-tosho.co.jp/books/ISBN978-4-489-02150-3.html

次は公式の参考書です。上2冊で学習できていない分散分析や直交表などを補いつつ、これまで学んだ事を一気に復習します。この本は某レビューで酷評されていますが、私は使い方の問題だと思っていて、この本でイチから学ぼうとするとおそらく失敗します。この本は復習に使うと良いです。良い意味でも悪い意味でもコンパクトに試験範囲がまとまっています。

さて、これで統計応用の範囲も全てチェックしました。
私は比較的時間にゆとりがあったのでこの間にも別の本をいくつかやっていますが、これだけで十分だと思います。後はひたすら演習です。時間をきっちり計って過去問を解きまくりましょう。
試験日まで時間がある人はこのあたりで気づくと思うのですが、年1回しか試験がないので過去問が少ないです。おいおい全然演習が足りないよと思うと思います。そんな人にはこれがオススメです。
Exam past papers
王立統計学会（RSS）が実施している試験の過去問です。これも公式サイトに載っている合格者の声で知りました。このうち「Graduate Diploma」の過去問を解いてみるのが良いと思います。難易度は1級よりやや簡単ですが、問題に慣れるのには十分でしょう。実際に似たような問題が1級試験でも出題されています。

受験当日

受験票に写真を貼ることと、電卓、時計を忘れないようにしましょう。
会場に行きます。受験者は男性ばかりでしたね。最初にアンケートとWeb合格発表を希望するかどうかもあるのでしっかりマークしておきます。
さて試験開始です。統計数理はどこかで見たことがある問題を見つけ1問サクッと解答。残り2問は最後がダメでしたがそれ以外は一通り解けたと思います。そしてお昼休憩を挟んで、午後から統計応用の試験です。
計算量が多そうな正規分布周りの計算を含む問題はスルーして、とにかく分かる問題を落とさずしっかり得点するという気持ちで解いていきました。

12月半ばにWeb合格発表が行われ無事合格しました。成績上位者に贈られる優秀賞は残念ながら逃しました。

最後に

統計検定1級を通して、いわゆる古典的な統計は一通り学べたと思っています。今は〇〇netを筆頭にDeep Learningが流行っていますが、検定だったり確率分布の操作だったりはそれらを理解するための基本だと思うので、その理解度をチェックする上で統計検定を受験してみてはいかがでしょうか。

*1:統計検定の公式キャラクターらしいです http://www.toukei-kentei.jp/info/toukeiakari/

非線形SVMによる分類

参考書はこちら

多変量解析入門――線形から非線形へ

ハードマージンSVM、ソフトマージンSVMと続いて今回でSVMは最後です。
今までは全て線形SVMでしたが、今回は非線形のSVMとなります。非線形SVMでは特徴空間次元の内積計算が多数必要になり、それを避けるためにカーネル法という手法が用いられます。カーネル法にも種類は様々あるのですが今回は最も一般的なガウスカーネルを用います。データはソフトマージンSVMのときと同じデータです。

from matplotlib import pyplot as plt
import numpy as np
import pandas as pd
import cvxopt
from cvxopt import matrix

class NonLinearSoftMarginSVM():
  
  def __init__(self, X, y, lam, sigma):
    self.X = X
    self.y = y
    self.lam = lam
    self.sigma = sigma
    
  #ガウスカーネル
  def kernel(self, x, y):
    return np.exp(np.dot(x - y, x - y) / (-2 * self.sigma))
    
  def fit(self):
    num_of_data = self.X.shape[0]
    label = self.y * 2.0 - 1.0 #ラベルを-1と1に変換する
    
    #2次計画問題を解く
    alpha = np.zeros((num_of_data, 1))
    q = matrix(np.ones((num_of_data ,1)))
    P = np.zeros((num_of_data, num_of_data))
    for i in range(num_of_data):
      for j in range(num_of_data):
        P[i,j] = label.iloc[i] * label.iloc[j] * self.kernel(self.X.iloc[i], self.X.iloc[j])
    P = matrix(P)
    A = matrix(np.array(label), (1, num_of_data))
    b = matrix(0.0, (1,1))
    #条件緩和によって変化した部分
    G_1 = (-1) * np.identity(num_of_data)
    G_2 = np.identity(num_of_data)
    G = matrix(np.vstack((G_1, G_2)))
    h_1 = np.zeros((num_of_data, 1))
    h_2 = np.ones((num_of_data, 1)) * self.lam
    h = matrix(np.vstack((h_1, h_2)))
    sol = cvxopt.solvers.qp(P, -q, G, h, A, b)
    alpha = sol["x"]
    #2次計画問題終了
    
    #パラメータを求める
    w = 0
    for i in range(num_of_data):
      w = w + alpha[i] * label.iloc[i] * self.X.iloc[i]
    
    #切片bを求めるためのサポートベクター集合S
    epsilon = 0.00001 #これより小さいものを0とみなす
    Sup_index = []
    for i in range(num_of_data):
      if alpha[i] < epsilon or alpha[i] > self.lam: #条件緩和によって変化した部分
        continue
      Sup_index.append(i)
    b = 0
    temp = 0
    print(Sup_index)
    for i in range(len(Sup_index)): #分離超平面より上と下の個数を数える
      if label.iloc[Sup_index[i]] == 1:
        temp = temp + 1
      for j in range(num_of_data):
        b = b + alpha[j] * label.iloc[j] * self.kernel(self.X.iloc[j], self.X.iloc[Sup_index[i]])
    b = ((temp - (len(Sup_index) - temp)) - b) / len(Sup_index)
    
    return [w, b, Sup_index, alpha]

前回のソフトマージンSVMと比べて変化しているのは、内積計算が全てカーネルによる計算に置き換えられているという点だけです。続いて、ここにデータを入れてみます。

lam = 1.0
sigma = 1.0
X = data[["x", "y"]]
y = data["c"]
label = data["c"] * 2 - 1
svm = NonLinearSoftMarginSVM(X , y, lam, sigma)
[w, b, Sup_index, alpha] = svm.fit()

ここまでは実はさほど計算時間はかかりません。決定境界をプロットしてみましょう。

def f(x, alpha, Sup_index, b, X, label):
  result = b
  for i in range(len(Sup_index)):
    j = Sup_index[i]
    result = result + alpha[j] * label[j]  * svm.kernel(X.iloc[j], x)
  return result

#色の設定
col = []
for i in range(len(data["c"])):
  if data["c"].iloc[i] == 1:
    col.append("red")
  else:
    col.append("blue")

#メッシュグリッドを作成して決定境界をプロット
seq = np.linspace(0, 10, 50)
X1, X2 = np.meshgrid(seq, seq)
Z = []
x1 = X1.reshape(-1, 1)
x2 = X2.reshape(-1, 1)
for i in range(X1.size):
  Z.append(f(np.array([x1[i][0], x2[i][0]]), alpha, Sup_index, b, X, label))

plt.title("Non Lnear Soft Margin SVM")
plt.xlabel("x")
plt.ylabel("y")
Z = np.array(Z).reshape(X1.shape)
plt.contour(X1, X2, Z, colors = "black", levels = [0.0])
plt.scatter(data["x"], data["y"], c = col)

このように非線形の決定境界となります。SVM自体の話はここまでで、あとはどうやってモデルを良いものにしていくかが難しいところです。例えば今回はソフトマージンなので正則化係数、さらにはガウスカーネルなので分散を事前に設定する必要があります。他には多項式カーネルであれば2つのハイパーパラメータを設定しなければなりません。クロスバリデーションなどを適用しながらこれらを適切に設定することが求められます。

参考

非線形の決定境界を書く際に参考になりました。
python で等高線を描くなら meshgrid して contour | コード７区

ソフトマージンSVMによる分類

参考書はこちら

多変量解析入門――線形から非線形へ

前回はハードマージンSVMを実装しました。今回はそれの拡張のようなもので、ソフトマージンSVMを行ないます。ハードマージンSVMではデータは完全に線形分離できるものにしてくださいと書きました。
ハードマージンSVMによる分類 - lazy diary
しかし現実問題そんな都合の良いデータばかりではないでしょう。なんらかの理由でラベルの異なるグループに入ってしまうデータもあるはずです。ソフトマージンSVMでは、そのようなデータが入ることを許すよう条件を緩和します。しかし、どの程度まで許すかということが問題になります。そこで条件を緩和しすぎず、なおかつマージンを大きくすることを目指したのがソフトマージンSVMです。
実際に解析するにあたってその程度を決める指標として正則化係数を私たちは与える必要があります。これについてはクロスバリデーションで決める方法がおそらく最もシンプルだと思われます。

コード自体はハードマージンSVMとほとんど変わりません。二次計画問題で、条件が増える点と、サポートベクトル集合にも条件がさらに加わることだけ注意します。

from matplotlib import pyplot as plt
import numpy as np
import pandas as pd
import cvxopt
from cvxopt import matrix

class SoftMarginSVM():
  
  def __init__(self, X, y, lam):
    self.X = X
    self.y = y
    self.lam = lam
    
  def fit(self):
    num_of_data = self.X.shape[0]
    label = self.y * 2.0 - 1.0 #ラベルを-1と1に変換する
    
    #2次計画問題を解く
    alpha = np.zeros((num_of_data, 1))
    q = matrix(np.ones((num_of_data ,1)))
    P = np.zeros((num_of_data, num_of_data))
    for i in range(num_of_data):
      for j in range(num_of_data):
        P[i,j] = label.iloc[i] * label.iloc[j] * np.dot(self.X.iloc[i], self.X.iloc[j])
    P = matrix(P)
    A = matrix(np.array(label), (1, num_of_data))
    b = matrix(0.0, (1,1))
    #条件緩和によって変化した部分
    G_1 = (-1) * np.identity(num_of_data)
    G_2 = np.identity(num_of_data)
    G = matrix(np.vstack((G_1, G_2)))
    h_1 = np.zeros((num_of_data, 1))
    h_2 = np.ones((num_of_data, 1)) * self.lam
    h = matrix(np.vstack((h_1, h_2)))
    sol = cvxopt.solvers.qp(P, -q, G, h, A, b)
    alpha = sol["x"]
    #2次計画問題終了
    
    #パラメータを求める
    w = 0
    for i in range(num_of_data):
      w = w + alpha[i] * label.iloc[i] * self.X.iloc[i]
    
    #切片bを求めるためのサポートベクター集合S
    epsilon = 0.00001 #これより小さいものを0とみなす
    Sup_index = []
    for i in range(num_of_data):
      if alpha[i] < epsilon or alpha[i] > self.lam: #条件緩和によって変化した部分
        continue
      Sup_index.append(i)
    b = 0
    temp = 0
    print(Sup_index)
    for i in range(len(Sup_index)): #分離超平面より上と下の個数を数える
      if label.iloc[Sup_index[i]] == 1:
        temp = temp + 1
      b = b + np.dot(w, self.X.iloc[Sup_index[i]])
    b = ((temp - (len(Sup_index) - temp)) - b) / len(Sup_index)
    
    return [w, b, Sup_index]

データを作ります。ハードマージンSVMと違い多少は入り乱れたデータでも平気です。

#create data
x_1 = 4.2
y_1 = 7
x_2 = 7
y_2 = 4.2
x1 = x_1 + np.random.randn(50)
y1 = y_1 + np.random.randn(50)
x2 = x_2 + np.random.randn(50)
y2 = y_2 + np.random.randn(50)
#上を元にデータをpandas.dataframeで作成
x = pd.Series(np.hstack([x1, x2]))
y = pd.Series(np.hstack([y1, y2]))
c = pd.Series(np.repeat([1, 0], 50))
data = pd.DataFrame({"x":x, "y":y, "c":c})

これをSVMにかけてみましょう。

lam = 10.0
svm = SoftMarginSVM(data[["x", "y"]], data["c"], lam)
[w, b, Sup_index] = svm.fit()

col = []
for i in range(len(data["c"])):
  if data["c"].iloc[i] == 1:
    col.append("red")
  else:
    col.append("blue")
seq = np.arange(np.min(data["x"]) - 1, np.max(data["x"]), 0.1)
plt.scatter(data["x"], data["y"], c = col)
plt.plot(seq, (-w[0] * seq - b) / (w[1]), c = "black")
plt.title("Soft Margin SVM")
plt.xlabel("x")
plt.ylabel("y")

このように、ハードマージンSVMの厳しすぎる条件を緩和することでソフトマージンSVMでは分類を達成します。

ハードマージンSVMによる分類

参考書はこちら

多変量解析入門――線形から非線形へ

今回はハードマージンSVMというものを実装します。
ソフトでもない上に非線形でもないという最も基本的なSVMです。ハードマージンSVMさえ理解してしまえば他のSVMはその拡張みたいなものです。

必要なものをimportしておきます。

from matplotlib import pyplot as plt
import numpy as np
import pandas as pd
import cvxopt
from cvxopt import matrix

このcvxoptというのは二次計画問題のソルバーです。詳しくはこの記事の最後にあるブログを参考にしてください。google colabにもデフォルトではインストールされていないようなので、

pip install cvxopt

で導入しておきましょう。

データは相変わらずいい加減に作ったものです。とは言いつつも今回はハードマージンSVMですので、線形分離できないデータは困ります。かならず超平面で分離できるデータをチョイスしてください。2次元なら直線で完全に分離できるものです。

#create data
x_1 = 3
y_1 = 7
x_2 = 7
y_2 = 2
x1 = x_1 + np.random.randn(50)
y1 = y_1 + np.random.randn(50)
x2 = x_2 + np.random.randn(50)
y2 = y_2 + np.random.randn(50)
#上を元にデータをpandas.dataframeで作成
x = pd.Series(np.hstack([x1, x2]))
y = pd.Series(np.hstack([y1, y2]))
c = pd.Series(np.repeat([1, 0], 50))
data = pd.DataFrame({"x":x, "y":y, "c":c})

ここからが本編です。流れとしては
1. データから二次計画問題を解く
2. パラメータwを求める
3. サポートベクトル集合を求める
4. サポートベクトル集合から切片bを求める
となっています。

class HardMarginSVM():
  
  def __init__(self, X, y):
    self.X = X
    self.y = y
    
  def fit(self):
    num_of_data = self.X.shape[0]
    label = self.y * 2.0 - 1.0 #ラベルを-1と1に変換する
    
    #2次計画問題を解く
    alpha = np.zeros((num_of_data, 1))
    q = matrix(np.ones((num_of_data ,1)))
    P = np.zeros((num_of_data, num_of_data))
    for i in range(num_of_data):
      for j in range(num_of_data):
        P[i,j] = label.iloc[i] * label.iloc[j] * np.dot(self.X.iloc[i], self.X.iloc[j])
    P = matrix(P)
    A = matrix(np.array(label), (1, num_of_data))
    b = matrix(0.0, (1,1))
    G = matrix((-1) * np.identity(num_of_data))
    h = matrix(np.zeros((num_of_data, 1)))
    sol = cvxopt.solvers.qp(P, -q, G, h, A, b)
    alpha = sol["x"]
    
    #パラメータを求める
    w = 0
    for i in range(num_of_data):
      w = w + alpha[i] * label.iloc[i] * self.X.iloc[i]
    
    #切片bを求めるためのサポートベクター集合S
    epsilon = 0.00001 #これより小さいものを0とみなす
    Sup_index = []
    for i in range(num_of_data):
      if alpha[i] < epsilon:
        continue
      Sup_index.append(i)
    b = 0
    temp = 0
    for i in range(len(Sup_index)): #分離超平面より上と下の個数を数える
      if label.iloc[Sup_index[i]] == 1:
        temp = temp + 1
      b = b + np.dot(w, self.X.iloc[Sup_index[i]])
    b = ((temp - (len(Sup_index) - temp)) - b) / len(Sup_index)
    
    return [w, b, Sup_index]

これを上のデータに当てはめてみましょう。

svm = HardMarginSVM(data[["x", "y"]], data["c"])
[w, b, Sup_index] = svm.fit()

col = []
for i in range(len(data["c"])):
  if data["c"].iloc[i] == 1:
    col.append("red")
  else:
    col.append("blue")
seq = np.arange(np.min(data["x"]) - 1, np.max(data["x"]), 0.1)
plt.scatter(data["x"], data["y"], c = col)
plt.plot(seq, (-w[0] * seq - b) / (w[1]), c = "black")
plt.scatter(data["x"].iloc[Sup_index], data["y"].iloc[Sup_index], c = "yellow", marker = "D")
plt.title("Hard Margin SVM")
plt.xlabel("x")
plt.ylabel("y")

黄色になっている部分はサポートベクトルです。SVMの特徴として実質的にサポートベクトルにしか依存しないので、次のような外れ値が入っても影響がありません。

data.loc[100] = [0, 20, 1]
svm = HardMarginSVM(data[["x", "y"]], data["c"])
[w, b, Sup_index] = svm.fit()

col = []
for i in range(len(data["c"])):
  if data["c"].iloc[i] == 1:
    col.append("red")
  else:
    col.append("blue")
seq = np.arange(np.min(data["x"]) - 1, np.max(data["x"]), 0.1)
plt.scatter(data["x"], data["y"], c = col)
plt.plot(seq, (-w[0] * seq - b) / (w[1]), c = "black")
plt.scatter(data["x"].iloc[Sup_index], data["y"].iloc[Sup_index], c = "yellow", marker = "D")
plt.title("Hard Margin SVM")
plt.xlabel("x")
plt.ylabel("y")

参考

今回、二次計画問題を解くにあたって参考にしました。SMOアルゴリズムについてはまだ理解できていないため、今回はcvxoptというソルバーに頼りました。
二次計画法(Quadratic Programming)の概要とPythonソルバーcvxoptの使い方 - MyEnigma

XORのようなデータを非線形ロジスティック回帰で分類する

参考書は相変わらずこちら

多変量解析入門――線形から非線形へ

今回は、基底関数展開を利用した非線形のロジスティック回帰を行ないます。基底関数を求めて、データを変換、最後にロジスティック回帰を行なうだけです。今回は基底関数としてガウスカーネルを使いました。
非線形の恩恵を確かめるために、通常の線形ロジスティック回帰では分類できないXORのようなデータを今回は使います。

import numpy as np
import pandas as pd
from matplotlib import pyplot as plt
from sklearn.cluster import KMeans
from sklearn import linear_model

from scipy.stats import multivariate_normal

cov = 0.5 * np.identity(2)
data1 = pd.DataFrame(multivariate_normal.rvs(mean = (-2.0, -2.0), cov = cov, size = 25))
data2 = pd.DataFrame(multivariate_normal.rvs(mean = (2.0, -2.0), cov = cov, size = 25))
data3 = pd.DataFrame(multivariate_normal.rvs(mean = (2.0, 2.0), cov = cov, size = 25))
data4 = pd.DataFrame(multivariate_normal.rvs(mean = (-2.0, 2.0), cov = cov, size = 25))

data1["c"] = 0
data2["c"] = 1
data3["c"] = 0
data4["c"] = 1

data = pd.concat([data1, data2, data3, data4], join = "inner")
data = data.reset_index(drop = True)
data = data.rename(columns = {0:"x", 1:"y", "c":"c"})

x = data["x"]
y = data["y"]
clusters = data["c"]

plt.scatter(x, y, c = clusters)

色付けをまともにしてないので見づらいですが、たしかにXORです。ここから前回と同じようにk-meansからの基底関数の構成およびデータの変換を行ないます。

#k-meansなどでクラスタリングする
dist = []
for i in range(1,30):
  km = KMeans(n_clusters = i)
  km.fit(data[["x", "y"]])
  dist.append(km.inertia_)
  
plt.plot(range(1,30), dist, marker = ".")

エルボー則に従って4くらいかなと考えられます。

from sklearn import linear_model
#クラスタ数を4とする

#ガウス型基底関数
def gauss_func(x, mu, sigma):
  return np.exp((np.sum(np.power(x - mu, 2))) / (-2 * sigma))

m = 4
km = KMeans(n_clusters = m)
km.fit(data[["x", "y"]])
cluster_center = km.cluster_centers_
dist = np.power(km.transform(data[["x", "y"]]), 2)
cluster_num = km.predict(data[["x", "y"]])
least_dist = np.zeros(shape = (m, 1))
for i in range(len(data)):
  for j in range(m):
    if cluster_num[i] == j:
      least_dist[j] += dist[i][j]
for i in range(m):
  least_dist[i] = least_dist[i] / np.sum(cluster_num == i)

transformed_X = pd.DataFrame(np.zeros(shape = (len(data), m)))
for i in range(len(data)):
  for j in range(m):
    transformed_X.iloc[i, j] = gauss_func(data[["x", "y"]].iloc[i,:], cluster_center[j], least_dist[j])

さて、カーネルパラメータの推定とカーネルを使ってデータの変換まで終わったのでロジスティック回帰を実行します。

#ロジスティック回帰の実行
log_reg = linear_model.LogisticRegression(penalty = "l2", solver = "saga", C = 1000.0, max_iter = 1000).fit(transformed_X, clusters)
print(log_reg.score(transformed_X, clusters)) #->1.0

とりあえず分類できました。決定境界をプロットしてみます。

seq = np.linspace(data["x"].min(), data["x"].max(), 50)
X1, X2 = np.meshgrid(seq, seq)
Z = []
x1 = X1.reshape(-1, 1)
x2 = X2.reshape(-1, 1)
testdata = np.hstack((x1, x2))
transformed_testdata = pd.DataFrame(np.zeros(shape = (len(x1), m)))
for i in range(len(x1)):
  for j in range(m):
    transformed_testdata.iloc[i, j] = gauss_func(testdata[i,:], cluster_center[j], least_dist[j])
Z = log_reg.predict_proba(transformed_testdata)

plt.title("Non Lnear Soft Margin SVM")
plt.xlabel("x")
plt.ylabel("y")
Z = np.array(Z[:,0]).reshape(X1.shape)
plt.contour(X1, X2, Z, colors = "black", levels = [0.5])
plt.scatter(data["x"], data["y"], c = clusters)

ガウスカーネルっぽい丸みを帯びた決定境界です。

僕の図書館（2019/02/16更新）

ちゃんと読んだもの一覧

Python

・Python学習講座~PythonエンジニアによるPython3学習支援サイト~

これで十分。あとはリファレンスを読む。
PythonエンジニアによるPython3学習サイト

・Aidemy

有名なやつです。Python使いならぜひ。
（2018/05/07追記：いつからか有料になっていました。当初は全講座が無料でした。また、当初よりコンテンツがかなり増えているようです。）
aidemy.net

・【ゼロから始めるデータ分析】 ビジネスケースで学ぶPythonデータサイエンス入門

Udemyの講座です。分析の一通りの手順を教えてくれます。
実際にSignateに出てくる問題を扱っているので、実践形式で学べるところが良し。
https://www.udemy.com/optworks_1/

R

・全人類がわかる統計学

Rの基礎をこのサイトで覚えました。あとは必要に応じて調べるだけで今のところ十分のような気がしています。
to-kei.net

機械学習

・Coursera：Machine Learning

この世界なら知らない人はいないであろうAndrew Ng先生の授業です。
内容は機械学習の基礎で、大半は日本語字幕があります。WEEK1-11まであり、各WEEKの最後にはクイズやプログラミングの宿題が課されます。具体的な授業内容は別の記事に書いたので、良ければ参考に。
Coursera：Machine Learningを受講しました - lazy diary
www.coursera.org

・ゼロから作るDeep Learning――Pythonで学ぶディープラーニングの理論と実装

NNをPythonでスクラッチして理論を学んで最終的にCNNまでやってみようという本です。分かりやすいと人気の高い本ですが、第5章の計算グラフを使ったback propagationの説明は好みが分かれるような気がします。しかしながら第5章以外、特に第6章の学習に関するテクニック（最適化の方法）は知らないこともあったので勉強になりました。

ゼロから作るDeep Learning ―Pythonで学ぶディープラーニングの理論と実装

伝統的な統計学

古典的な統計学。いわゆる〇〇検定とかの本です。

・統計学入門 (基礎統計学Ⅰ)

統計の赤本の名で有名な統計学の入門書です。線形代数だと行列の基本的な扱い、解析だと重積分あたりが前提知識かと思います。日常でも聞く平均などの統計量の説明から始まり、そこに確率変数を導入して分布の説明、次に標本の話と標本からパラメータを推定したり仮説検定を行なったりする話、最後に少しですが回帰分析の話が載っています。
私はこの本を統計検定2級を受験する際の教科書として利用しました。統計検定2級の受験記も良ければ見て下さい。
統計検定2級に合格しました - lazy diary

統計学入門 (基礎統計学?)

・現代数理統計学の基礎

分布収束といった収束の概念を使って細かい説明がなされています（※測度論の話は数か所で出てくる程度です）。この本の良い点は誤植や演習問題の解答、補足等を著者がWEBサイトで公開している点だと思います。数学書でこれほどありがたいことはないです。
私はこの本を統計検定1級を受験する際の教科書として利用しました。統計検定1級の受験記も良ければ見て下さい。
統計検定1級に合格しました - lazy diary

現代数理統計学の基礎 (共立講座 数学の魅力)

統計モデリング

モデリングの話。ベイズベイズしていない本。

・一般化線形モデル入門

一般化線形モデルを数理統計の視点で解説しています。この本を読めば一般化線形モデルの数理の側面は一通り理解できると思います。

一般化線形モデル入門 原著第2版

・多変量解析入門――線形から非線形へ

線形のモデルを理解したら次は非線形と思ったので読みました。この本は線形モデルの知識を前提とはしていませんが、知っていると捗ると思います。基底関数展開や判別分析、SVM、主成分分析を丁寧に数式を展開しながら説明しています。最後に教師なし学習が少しだけ説明されています。丁寧な数式展開、誤植が見当たらないなど数学書としてはこの上ない素晴らしい本だと思います。以下、本で紹介されているものをいくつか私がPythonで実装したものです。
XORのようなデータを非線形ロジスティック回帰で分類する - lazy diary
ハードマージンSVMによる分類 - lazy diary
ソフトマージンSVMによる分類 - lazy diary
非線形SVMによる分類 - lazy diary

多変量解析入門――線形から非線形へ

・経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

沖本本の名で有名な時系列解析の本です。1通りまとめてみたものの、GARCHモデルらへんはまだあんまり理解できた自信がないです。

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

・続・わかりやすいパターン認識―教師なし学習入門―

前半はベイズ推論の入口です（ベイズの定理、事前分布とは？）。そのあと教師あり学習との対比をしたうえで教師なし学習で頻出のEMアルゴリズム、HMM（Viterbiアルゴリズム、Baum-Welchアルゴリズム）があります。次に、k-means法、凸クラスタリング、最後にノンパラベイズ（ディリクレ過程混合モデル）の説明があります。丁寧に記述されていてとても読みやすいのですが、ノンパラベイズは使いこなせる自信がないです。

続・わかりやすいパターン認識―教師なし学習入門―

ベイズモデリング

ベイズモデリングについての本です。

・データ解析のための統計モデリング入門――一般化線形モデル・階層ベイズモデル・MCMC (確率と情報の科学)

みどり本の名で有名な統計モデリングの本です。モデルの説明もありますがそれよりモデルを作っていくときの考え方を紹介した本だと私は思います。下記のPyStanに移植したというありがたい記事がに参考になりました。
データ解析のための統計モデリング入門（緑本）の読書メモ（PythonとStanで） - Qiita

データ解析のための統計モデリング入門――一般化線形モデル・階層ベイズモデル・MCMC (確率と情報の科学)

・パターン認識と機械学習

PRMLの名で有名な本です。ベイジアンなら読まないとモグリだと思うので読みました。
やっぱりこの本は難しいですね。難解な説明でPRMLを読むことが目的になっているような気がしないでもないです。これを理解して使いこなせるようになるまではまだまだ時間がかかるような気がします。

パターン認識と機械学習 上

パターン認識と機械学習 下 (ベイズ理論による統計的予測)

・ベイズ統計の理論と方法

渡辺ベイズ本とかベイズの赤本とか言われるやつです。ベイズ統計学の基礎からWAICの導出がメインだと思われます。ベイズ推論における理論のエッセンスが詰まっています。「質問と回答」を読むだけでも十分に価値があると思っています。ベイズモデリングやるならこの本は必携だと思いました。

ベイズ統計の理論と方法

・StanとRでベイズ統計モデリング (Wonderful R)

アヒル本と呼ばれるやつです。Stanの使い方と緑本や一般化線形モデル入門で扱われているモデリングの良い復習になりました。この本はRで書かれていますが、私は一部を除いて全てPyhonで実装しました。そのためPythonによる可視化の練習にもなりました。以下の記事が実装で困ったときの参考になりました。
StanとRでベイズ統計モデリング(アヒル本)をPythonにしてみる - 目次 - Qiita

StanとRでベイズ統計モデリング (Wonderful R)

・ベイズ統計で実践モデリング: 認知モデルのトレーニング（現在、取り組み中）

コワい本と呼ばれるやつです。

ベイズ統計で実践モデリング: 認知モデルのトレーニング

その他

・入門 機械学習による異常検知―Rによる実践ガイド

異常検知に興味が湧いたので読みました。式はあまり行間が飛ぶことなく丁寧に書かれており、著者の親切さ（のようなもの？）をひしひしと感じました。異常検知の流れを掴むことができました。

入門 機械学習による異常検知―Rによる実践ガイド

・赤池情報量規準AIC―モデリング・予測・知識発見―

本書はその主役である赤池博士の業績および、博士によって切り拓かれた情報量に基づく情報処理の最先端の成果を紹介するものとして企画されたみたいです。AICを考えるようになったきっかけから始まりAICとその周辺の事柄がまとめられています。また、モデリングに対する姿勢も描かれています。この本はあくまで読み物として手に取るのが良いと思われます。

赤池情報量規準AIC―モデリング・予測・知識発見

統計検定2級に合格しました

合格したので感想を。
*1

統計検定とは

www.toukei-kentei.jp
そもそも統計検定って何？って方もいると思うのでザックリと説明を。

「統計検定」とは、統計に関する知識や活用力を評価する全国統一試験です。
データに基づいて客観的に判断し、科学的に問題を解決する能力は、仕事や研究をするための21世紀型スキルとして国際社会で広く認められています。
日本統計学会は、中高生・大学生・職業人を対象に、各レベルに応じて体系的に国際通用性のある統計活用能力評価システムを研究開発し、統計検定として実施します。

うむ、よく分からんということで具体的な検定内容は次のようになっています。

試験の種別	試験内容
統計検定1級	実社会の様々な分野でのデータ解析を遂行する統計専門力
統計検定準1級	統計学の活用力 ─ データサイエンスの基礎
統計検定2級	大学基礎統計学の知識と問題解決力
統計検定3級	データの分析において重要な概念を身に付け、身近な問題に活かす力
統計検定4級	データや表・グラフ、確率に関する基本的な知識と具体的な文脈の中での活用力
統計調査士	統計に関する基本的知識と利活用
専門統計調査士	調査全般に関わる高度な専門的知識と利活用手法

2/3級については、CBT方式というもので受験が可能です。この方式で申し込みをすれば指定された全国の会場でほぼ毎日受検することが出来ます。
一方で、1級と準1級については受験日が決まっています。1級は11月、準1級は6月に試験が開催されます。全て合格点は60/100点だったと思います。
統計調査士、専門統計調査士についてはよく知らないので略。

この統計検定始まったのは2011年とまだ最近です。

勉強したこと

統計検定2級の範囲をまずチェックしましょう。
私は数学畑の人間とは言っても統計学については素人です。確率分布って何？という状態です。そこで合格した人のブログをいくつか調べて、この本で学ぶことにしました。

統計学入門 (基礎統計学?)

検定理論あたりで少し苦戦しましたがそれ以外はとても読みやすい本でした。この本をノートにまとめながら読むこと約2週間（※練習問題は飛ばしました）。
これだけでどれくらい点数が取れるのだろうかと過去問を解くことにしました。

日本統計学会公式認定 統計検定 2級 公式問題集[2014〜2016年]

ギリギリですが合格点を取ることが出来ました。私はこのとき初めてを2級の範囲を見たのですが（最初にチェックしておけよって話ですね）、『統計学入門』にいくつか載っていない項目があったのでそこをググりながら補うことにしました。これに1週間ほどかけました。その後、過去問の再度解き直しに1週間をかけ、安定して合格できるくらいまで知識を確認しました。

受験当日

事前に日程をこちらが決めてCBT方式で受験しました。会場に行ってパソコンの前に着席してから説明を聞いて、すぐに試験開始です。メモ用紙が貰えるでそこでガリガリ計算していきました。偶然かもしれませんが、過去問よりいくぶんか簡単でした。その場で合否が分かるので合格発表をドキドキして待つこともありません。

最後に

ここまで言及していなかったのですが、忘れてはいけないのが電卓です。この試験をサクサク乗り切るためには電卓をうまく使うことが結構重要なのではと過去問を解きながら思っていました。私は電卓といわれるとスマホの電卓しか長らく使っていませんでしたが、スマホの電卓はNGなので今回この試験のために買いました。

カシオ 電卓 時間・税計算 手帳タイプ 8桁 SL-300A-N

あまり使う機会ないですしお金が無い学生にはこれで十分でした。

参考書を選んで学習、最後に過去問で演習という流れで合格できました。参考書は範囲をおおよそカバーしているのものの中から気に入ったのを選べば良いと思います。*2『統計学入門』もその１つだと思います。

次は統計検定1級に挑戦してみようと思います。頑張ります。
【追記】
2018年度の統計検定1級試験を受験して合格しました。1級の受験記はこちらです。
mt19.hatenablog.com

*1:統計検定の公式キャラクターらしいです http://www.toukei-kentei.jp/info/toukeiakari/

*2:公式の参考書もありますがAmazonではなかなかの酷評です。私は公式のものを見ていないので評価できません。